17. Linear solver

The linear solver is used for solving the pressure Poisson equations and when quasi-implicit discretization is used for diffusion, the solver is also used for solving the intermediate solution of the momentum equations.

The linear solver can be selected by setting numerical/linear_solver/@solver_choice to one of the identifiers listed below. ComFLOW supports the following solver-preconditioner combinations (in order of recommendation):

name description OpenMP parallel MPI parallel
"bicgstab_ilu" BiCGSTAB + ILU(0) preconditioner with special drop tolerance yes yes
"sor" SOR with improved search algorithm yes (outdated) no
"sor_legacy" SOR with legacy search algorithm yes (outdated) no
"cg_milu" BiCGSTAB with MILU(0) preconditioner (Cartesian grids only) no no
"1" CG with MILU(0) preconditioner (Cartesian grids only) no no

Only the first of the listed solvers (@solver_choice="bicgstab_ilu") is compatible with all simulations settings of ComFLOW. The preconditioner in this solver takes care of the special matrix structure near the GABC boundaries and is able to cope with large density variations in two-phase flow simulations. Typically the best convergence rate is obtained by selecting this solver. Note that parallelization is available for a limited number of solver, in particular where it concerns MPI parallelization.

For most one-phase simulations the SOR solver (@solver_choice="sor") provides a satisfactory alternative. The SOR solver cannot be used for two-phase simulations or simulations including a GABC boundary condition. In case of an absorbing boundary condition (GABC) the BiCGSTAB solver should be used (@solver_choice="bicgstab_ilu"). There are two options for the SOR solver. When @solver_choice="sor_legacy", the legacy SOR solver will be used, where the automatically adjusted optimal relaxation parameter follows from a robust, but slower, search procedure. When @solver_choice="sor" (recommended choice when using SOR), an alternative search procedure is used to determine the optimal relaxation parameter, which in general leads to a better convergence rate.

The other solvers that are available in ComFLOW do not support local grid refinement and have not been parallelized.

Section 17.4 outlines the recommended settings for the various solvers.

17.1. Stopping criterion

The settings for the stopping criterion can be specified in the section numerical/linear_solver/stopping_criterion/. The number of solver iterations is determined by two mechanisms.

Firstly, the number of iterations is bounded by @iter_max. This value should be chosen large enough to allow the solver to reach the desired precision.

Secondly, the solver iterations are controlled by a tolerance on the residual. All solvers that are available in ComFLOW apply an inner iteration loop and an outer iteration loop. In the inner loop the convergence criterion is based on a 2-norm of the residual, which can be computed more efficiently than a maximum norm. In the outer iteration loop the convergence criterion is based on a more expensive maximum norm of the residual. More precisely, the stopping criterion in the outer iteration loop is based on a cell-wise maximum norm of the residual of the Poisson equation, row-wise normalized by the solution magnitude, roughly something like: \(|\text{diag}[1./x](Ax-b)\|_\infty\), where matrix \(A\) and vector \(b\) describe the pressure Poisson problem, \(x\) is the pressure solution obtained so far, \(\text{diag}[*]\) is a matrix with * on the diagonal and the division is applied point-wise. For different simulations a different number can be chosen. In general, smaller values result in more iterations, but improve the accuracy of the solution. The convergence behavior (error w.r.t. iteration count) depends on the chosen solver. A value of "1.0E-6" is generally a good value, but a useful value in case of wave simulations is "1E-8".

The tolerance level for the outer iteration loop determines the accuracy of the final answer and can be specified in the attributes @outer_p and @outer_uvw for the pressure Poisson and quasi-implicit diffusion scheme, respectively. The tolerance level for the inner iteration loop is specified by attributes @inner_p and @inner_uvw and should be set to a slightly larger value than the outer tolerance level. Every time when the inner iteration loop converges, the solution is tested for the outer (and final) stopping criterion. If necessary the inner iteration loop is automatically continued using a smaller inner tolerance level.

17.2. The ILU preconditioner

The BiCGSTAB solver (@solver_choice="bicgstab_ilu") employs an ILU preconditioner which calculates an incomplete LU decomposition based on a drop tolerance approach. In the section numerical/linear_solver/ilu_preconditioner/ two drop tolerances can be specified. In most of the computational cells a general drop tolerance is applied, which can be set with attribute @droptol. In boundary cells containing a (G)ABC boundary conditions a separate drop tolerance is applied, which can be set with attribute @droptolbc. Here it is recommended to use a smaller dropping criterion because otherwise too much information is lost and the convergence rate is negatively affected.

A (default) value of @droptol="1.0e-2" is usually suitable for the general drop tolerance combines with a stricter drop tolerance of @droptolbc="1.0e-6" for the (G)ABC-related equations. The memory requirements for storage of the ILU preconditioner increase when using lowering the drop tolerance. On the other hand, the number of iterations decreases with a better preconditioner, obtained with a smaller value of the drop-tolerance.

17.3. The SOR search algorithm

The initial relaxation factor for the SOR solver can be set in numerical/linear_solver/relaxation/@omega Preferably use a value of "1.0" for safety. Sometimes little speed up may be achieved for values of "1.3". During the solver iterations the relaxation factor is automatically adjusted.

17.5. Diagnostics

For diagnostics output the attribute numerical/linear_solver/@verbosity can be set to a value larger than 0. ComFLOW distinguishes the following verbosity levels:

  1. no output
  2. basic statistics
  3. detailed statistics per solver iteration

Footnotes

[1](1, 2) Note: the convergence of the SOR solver can slow down in case of complicated geometries.