Sparse LU decomposition (Gaussian elimination) is used by default to
solve linear systems of equations in FEniCS programs. This is a very
robust and simple method. It is the recommended method for systems
with up to a few thousand unknowns and may hence be the method of
choice for many 2D and smaller 3D problems. However, sparse LU
decomposition becomes slow and one quickly runs out of memory for
larger problems. For large problems, we instead need to use *iterative
methods* which are faster and require much less memory. We will now
look at how to take advantage of state-of-the-art iterative solution
methods in FEniCS.

Preconditioned Krylov solvers is a type of popular iterative methods that are easily accessible in FEniCS programs. The Poisson equation results in a symmetric, positive definite system matrix, for which the optimal Krylov solver is the Conjugate Gradient (CG) method. For non-symmetric problems, a Krylov solver for non-symmetric systems, such as GMRES, is a better choice. Incomplete LU factorization (ILU) is a popular and robust all-round preconditioner, so let us try the GMRES-ILU pair:

```
solve(a == L, u, bc,
solver_parameters={'linear_solver': 'gmres',
'preconditioner': 'ilu'})
# Alternative syntax
solve(a == L, u, bc,
solver_parameters=dict(linear_solver='gmres',
preconditioner='ilu'))
```

the section List of linear solver methods and preconditioners lists the most popular choices of Krylov solvers and preconditioners available in FEniCS.

The actual GMRES and ILU implementations that are brought into action
depend on the choice of linear algebra package. FEniCS interfaces
several linear algebra packages, called *linear algebra backends* in
FEniCS terminology. PETSc is the default choice if FEniCS is compiled
with PETSc. If PETSc is not available, then FEniCS falls back to using
the Eigen backend. The linear algebra backend in FEniCS can be set
using the following command:

```
parameters.linear_algebra_backend = backendname
```

where `backendname`

is a string. To see which linear algebra backends
are available, you can call the FEniCS function
`list_linear_algebra_backends`

. Similarly, one may check which
linear algebra backend is currently being used by the following
command:

```
print(parameters.linear_algebra_backend)
```

We will normally want to control the tolerance in the stopping
criterion and the maximum number of iterations when running an
iterative method. Such parameters can be controlled at both a *global*
and a *local* level. We will start by looking at how to set global
parameters. For more advanced programs, one may want to use a number
of different linear solvers and set different tolerances and other
parameters. Then it becomes important to control the parameters at a
*local* level. We will return to this issue in the section Linear variational problem and solver objects.

Changing a parameter in the global FEniCS parameter database affects
all linear solvers (created *after* the parameter has been set).
The global FEniCS parameter database is simply called `parameters`

and
it behaves as a nested dictionary. Write

```
info(parameters, verbose=True)
```

to list all parameters and their default values in the database.
The nesting of parameter sets is indicated through indentation in the
output from `info`

.
According to this output, the relevant parameter set is
named `'krylov_solver'`

, and the parameters are set like this:

```
prm = parameters.krylov_solver # short form
prm.absolute_tolerance = 1E-10
prm.relative_tolerance = 1E-6
prm.maximum_iterations = 1000
```

Stopping criteria for Krylov solvers usually involve some norm of
the residual, which must be smaller than the absolute tolerance
parameter *or* smaller than the relative tolerance parameter times
the initial residual.

We remark that default values for the global parameter database can be defined in an XML file. To generate such a file from the current set of parameters in a program, run

```
File('parameters.xml') << parameters
```

If a `dolfin_parameters.xml`

file is found in the directory where a
FEniCS program is run, this file is read and used to initialize the
`parameters`

object. Otherwise, the file
`.config/fenics/dolfin_parameters.xml`

in the user's home directory is
read, if it exists. Another alternative is to load the XML file (with any
name) manually in the program:

```
File('parameters.xml') >> parameters
```

The XML file can also be in gzip'ed form with the extension `.xml.gz`

.

We may extend the previous solver function from
`ft12_poisson_solver.py`
in the section A more general solver function
such that it also offers the GMRES+ILU
preconditioned Krylov solver:

This new `solver`

function, found in the file
`ft10_poisson_extended.py`,
replaces the one in
`ft12_poisson_solver.py`.
It has all the functionality of the previous `solver`

function, but
can also solve the linear system with iterative methods.

Regarding verification of the new `solver`

function in terms of unit
tests, it turns out that unit testing for a problem where the
approximation error vanishes gets more complicated when we use
iterative methods. The problem is to keep the error due to iterative
solution smaller than the tolerance used in the verification
tests. First of all, this means that the tolerances used in the Krylov
solvers must be smaller than the tolerance used in the `assert`

test,
but this is no guarantee to keep the linear solver error this small.
For linear elements and small meshes, a tolerance of \( 10^{-11} \) works
well in the case of Krylov solvers too (using a tolerance \( 10^{-12} \)
in those solvers). The interested reader is referred to the
`demo_solvers`

function in
`ft10_poisson_extended.py`
for details:
this function tests the numerical solution for direct and iterative
linear solvers, for different meshes, and different degrees of the
polynomials in the finite element basis functions.

Which linear solvers and preconditioners that are available in FEniCS depends on how FEniCS has been configured and which linear algebra backend is currently active. The following table shows an example of which linear solvers that can be available through FEniCS when the PETSc backend is active:

Name | Method |

`'bicgstab'` | Biconjugate gradient stabilized method |

`'cg'` | Conjugate gradient method |

`'gmres'` | Generalized minimal residual method |

`'minres'` | Minimal residual method |

`'petsc'` | PETSc built in LU solver |

`'richardson'` | Richardson method |

`'superlu_dist'` | Parallel SuperLU |

`'tfqmr'` | Transpose-free quasi-minimal residual method |

`'umfpack'` | UMFPACK |

The set of available preconditioners also depends on configuration and linear algebra backend. The following table shows an example of which preconditioners may be available:

Name | Method |

`'icc'` | Incomplete Cholesky factorization |

`'ilu'` | Incomplete LU factorization |

`'petsc_amg'` | PETSc algebraic multigrid |

`'sor'` | Successive over-relaxation |

An up-to-date list of the available solvers and preconditioners for your FEniCS installation can be produced by

```
list_linear_solver_methods()
list_krylov_solver_preconditioners()
```