Merge pull request #777 from KulaginVladimir/UpdDocs

Documentation: stepsize and Newton Solver

docs/source/userguide/index.rst

@@ -17,5 +17,6 @@ User's Guide
    traps
    settings
    stepsize
+   newton_solver
    export_post_processing
    troubleshooting

docs/source/userguide/newton_solver.rst

@@ -0,0 +1,108 @@
+=============
+Newton solver
+=============
+
+For advanced simulations, the parameters of the Newton solver can be adapted depending on a specific problem. 
+
+-----------------
+Built-in options
+-----------------
+A limited set of the solver parameters can be accessed via the built-in attributes of classes. 
+
+The parameters of the Newton solver for :class:`festim.HTransportProblem` can be chosen in :class:`festim.Settings` (see :ref:`settings_ug`). Absolute and relative tolerances of the Newton solver
+are defined with ``absolute_tolerance`` and ``relative_tolerance`` attributes, respectively. The maximum number of the solver iterations can be set using 
+the ``maximum_iterations`` parameter. Additionally, there is an option to choose linear solver and preconditioning methods that may be more suitable for particular problems.
+
+The linear solver method can be set with the ``linear_solver`` attribute. The list of available linear solvers can be viewed with: ``print(fenics.list_linear_solver_methods())``.
+
+.. dropdown:: Linear solver methods
+
+    * "bicgstab" - Biconjugate gradient stabilized method
+    * "cg" - Conjugate gradient method
+    * "gmres" - Generalized minimal residual method
+    * "minres" - Minimal residual method
+    * "mumps" - MUMPS (MUltifrontal Massively Parallel Sparse direct Solver)
+    * "petsc" - PETSc built in LU solver
+    * "richardson" - Richardson method 
+    * "superlu" - SuperLU
+    * "superlu_dist" - Parallel SuperLU
+    * "tfqmr" - Transpose-free quasi-minimal residual method
+    * "umfpack" - UMFPACK (Unsymmetric MultiFrontal sparse LU factorization)
+
+The preconditioner can be set with the ``preconditioner`` attribute. The list of available preconditioners can be viewed with: ``print(fenics.list_krylov_solver_preconditioners())``.
+
+.. dropdown:: Preconditioners
+
+    * "amg" - Algebraic multigrid
+    * "hypre_amg" - Hypre algebraic multigrid (BoomerAMG)
+    * "hypre_euclid" - Hypre parallel incomplete LU factorization
+    * "hypre_parasails" - Hypre parallel sparse approximate inverse
+    * "icc" - Incomplete Cholesky factorization
+    * "ilu" - Incomplete LU factorization
+    * "jacobi" - Jacobi iteration 
+    * "petsc_amg" - PETSc algebraic multigrid
+    * "sor" - Successive over-relaxation
+
+Similarly, the Newton solver parameters of :class:`festim.HeatTransferProblem`, :class:`festim.ExtrinsicTrap`, or :class:`festim.NeutronInducedTrap` 
+can be defined if needed. Here is an example for the heat transfer problem:
+
+.. code-block:: python
+
+    model.T = HeatTransferProblem(
+        transient=True,
+        initial_condition=300,
+        absolute_tolerance=1.0,
+        relative_tolerance=1e-10,
+        maximum_iterations=50,
+        linear_solver="gmres",
+        preconditioner="icc")
+
+
+--------------
+Custom solver
+--------------
+
+For a finer control, the built-in Newton solver can be overwritten with a custom solver based on the ``fenics.NewtonSolver`` class.
+
+.. warning::
+
+    Defining a custom Newton solver will override the solver parameters given with the built-in settings.
+
+A user-defined Newton solver can be provided after :class:`festim.Simulation.initialise()`. Here is a simple example for the H transport problem:
+
+.. code-block:: python
+
+    custom_solver = fenics.NewtonSolver()
+    custom_solver.parameters["error_on_nonconvergence"] = False
+    custom_solver.parameters["absolute_tolerance"] = 1e10
+    custom_solver.parameters["relative_tolerance"] = 1e-10
+    custom_solver.parameters["maximum_iterations"] = 100
+    custom_solver.parameters["linear_solver"] = "gmres"
+    custom_solver.parameters["preconditioner"] = "ilu"
+
+    model.initialise()
+
+    model.h_transport_problem.newton_solver = custom_solver
+
+    model.run()
+
+.. warning::
+
+    For a stationary heat transfer problem, a custom Newton solver has to be provided before the simulation initialisation! 
+
+To extend the functionality, the `NewtonSolver <https://bitbucket.org/fenics-project/dolfin/src/master/dolfin/nls/NewtonSolver.cpp>`_ class 
+can be overwritten: 
+
+.. code-block:: python
+
+    class CustomSolver(f.NewtonSolver):
+        def __init__(self):
+            super().__init__()
+
+        def converged(self, r, problem, iteration):
+            if iteration == 0:
+                self.r0 = r.norm("l2")
+            print(f"Iteration {iteration}, relative residual {r.norm('l2')/self.r0}")
+            return super().converged(r, problem, iteration)
+
+In this example, the relative residual will be printed after each Newton solver iteration.

docs/source/userguide/settings.rst

@@ -1,3 +1,5 @@
+.. _settings_ug:
+
 ========
 Settings
 ========
@@ -16,19 +18,21 @@ The settings of a FESTIM simulation are defined with a :class:`festim.Settings`
     )
 
 
-Here you define:
+Here you define with:
 
-* wether the simulation is transient or steady-state
-* the final time of the simulation
-* wether to run the simulation with conservation of chemical potential at interfaces (only useful for multi-materials)
-* wether to turn the Soret effect on
-* the absolute and relative tolerance of the Newton solver
-* the maximum iterations of the Newton solver
+* ``transient``: wether the simulation is transient or steady-state
+* ``final_time``: the final time of the simulation
+* ``chemical_pot``: wether to run the simulation with conservation of chemical potential at interfaces (only useful for multi-materials)
+* ``soret``: wether to turn the Soret effect on or not
+* ``absolute_tolerance``: the absolute tolerance of the Newton solver
+* ``relative_tolerance``: the relative tolerance of the Newton solver
+* ``maximum_iterations``: the maximum iterations of the Newton solver
 
 More advanced settings are also available:
 
-* the type of finite elements for traps (DG elements can be useful to account for discontinuities)
-* Wether to update the jacobian at each iteration or not
-* the linear solver
+* ``traps_element_type``: the type of finite elements for traps (DG elements can be useful to account for discontinuities)
+* ``update_jacobian``: wether to update the jacobian at each iteration or not
+* ``linear_solver``: linear solver method for the Newton solver
+* ``preconditioner``: preconditioning method for the Newton solver
 
 See :ref:`settings_api` for more details.

docs/source/userguide/stepsize.rst

@@ -1,3 +1,4 @@
+========
 Stepsize
 ========
 
@@ -11,7 +12,7 @@ Here is an example creating a stepsize of 1.2 seconds:
     import festim as F
     my_stepsize = F.Stepsize(initial_value=1.2)
 
-To use of the adaptive time stepping implemented in FESTIM, the arguments ``stepsize_change_ratio`` needs to be set to a value above 1.
+To use the adaptive time stepping implemented in FESTIM, the arguments ``stepsize_change_ratio`` needs to be set to a value above 1.
 
 .. code-block:: python
 
@@ -26,15 +27,35 @@ To cap the stepsize after some time, the parameters ``t_stop`` and ``stepsize_st
 
 .. code-block:: python
 
-    my_stepsize = F.Stepsize(initial_value=1.2, stepsize_change_ratio=1.5, dt_min=1e-6, t_stop=10, stepsize_stop_max=1.5)
+    my_stepsize = F.Stepsize(
+        initial_value=1.2, 
+        stepsize_change_ratio=1.5, 
+        dt_min=1e-6, 
+        t_stop=10, 
+        stepsize_stop_max=1.5)
 
 .. warning::
 
-    Please note that the parameters ``t_stop`` and ``stepsize_stop_max`` will be deprecated in a future release. To set the maximal value of the stepsize, consider using the ``max_stepsize`` parameter:
-
-    .. code-block:: python
+    Please note that parameters ``t_stop`` and ``stepsize_stop_max`` will be deprecated in a future release. 
 
-        my_stepsize = F.Stepsize(initial_value=1.2, stepsize_change_ratio=1.5, dt_min=1e-6, max_stepsize=lambda t: 1 if t < 1 else 2)
+Another option for controlling the stepsize is to use the ``max_stepsize`` parameter. This parameter defines the maximal value of the stepsize during simulations, 
+and it can be set as a constant or a callable function of time:
+
+.. code-block:: python
+
+    def max_stepsize(t):
+        if t <= 5:
+            return 1.5
+        elif t > 5 and t < 10:
+            return 2.5
+        else:
+            return None
+
+    my_stepsize = F.Stepsize(
+        initial_value=1.2, 
+        stepsize_change_ratio=1.5, 
+        dt_min=1e-6, 
+        max_stepsize=max_stepsize)
 
 The ``milestones`` argument can be used to make sure the simulation passes through specific times.
 This will modify the stepsize as needed.
@@ -45,7 +66,6 @@ This will modify the stepsize as needed.
         initial_value=1.2,
         stepsize_change_ratio=1.5,
         dt_min=1e-6,
-        t_stop=10,
-        stepsize_stop_max=1.5,
+        max_stepsize=5,
         milestones=[1, 5, 6, 10]
         )