Preconditioners¶

Aside from the AMG, AMGCL implements preconditioners for some common problem types. For example, there is a Schur complement pressure correction preconditioner for Navie-Stokes type problems, or CPR preconditioner for reservoir simulations. Also, it is possible to use single level relaxation method as a preconditioner.

General preconditioners¶

These preconditioners do not take the origin or structure of matrix into account, and may be useful both on their own, as well as building blocks for the composite preconditioners.

AMG¶

template<class Backend, template<class> class Coarsening, template<class> class Relax> class amgcl::amg¶

Include <amgcl/amg.hpp>

AMG is one the most effective methods for the solution of large sparse unstructured systems of equations, arising, for example, from discretization of PDEs on unstructured grids [TrOS01]. The method may be used as a black-box solver for various computational problems, since it does not require any information about the underlying geometry.

The three template parameters allow the user to select the exact components of the method:

The Backend to transfer the constructed hierarchy to,
The Coarsening strategy for the hierarchy construction, and
The Relaxation scheme (smoother to use during the solution phase).

typedef typename Backend::params backend_params¶: The backend parameters

typedef typename Backend::value_type value_type¶: The value type of the system matrix

typedef typename amgcl::math::scalar_of<value_type>::type scalar_type¶: The scalar type corresponding to the value type. For example, when the value type is std::complex<double>, then the scalar type is double.

typedef Coarsening<Backend> coarsening_type¶: The coarsening class instantiated for the specified backend

typedef Relax<Backend> relax_type¶: The relaxation class instantiated for the specified backend

class params¶

The AMG parameters. The coarsening and the relaxation parameters are encapsulated as part of the AMG parameters.

typedef typename coarsening_type::params coarsening_params¶: The type of the coarsening parameters

typedef typename relax_type::params relax_params¶: The type of the relaxation parameters

coarsening_params coarsening¶: The coarsening parameters

relax_params relax¶: The relaxation parameters

unsigned coarse_enough = Backend::direct_solver::coarse_enough()¶: Specifies when a hierarchy level is coarse enough to be solved directly. If the number of variables at the level is lower than this threshold, then the hierarchy construction is stopped and the linear system is solved directly at this level. The default value comes from the direct solver class defined as part of the backend.

bool direct_coarse = true¶: Use direct solver at the coarsest level. When set, the coarsest level is solved with a direct solver. Otherwise a smoother is used as a solver. This may be useful when the system is singular, but is still possible to solve with an iterative solver.

unsigned max_levels = std::numeric_limits<unsigned>::max()¶: The maximum number of levels. If this number is reached even when the size of the last level is greater that coarse_enough, then the hierarchy construction is stopped. The coarsest level will not be solved directly, but will use a smoother.

unsigned npre = 1¶: The number of pre-relaxations.

unsigned npost = 1¶: The number of post-relaxations.

unsigned ncycle = 1¶: The shape of AMG cycle (1 for V-cycle, 2 for W-cycle, etc).

unsigned pre_cycles = 1¶: The number of cycles to make as part of preconditioning.

bool allow_rebuild = false¶: Keep transfer operator matrices (\(P\) and \(R\)) in internal format to allow for quick rebuild of the hierarchy. See amgcl::amg::rebuild().

template<class Matrix> void rebuild(const Matrix &A, const backend_params &bprm = backend_params())¶: Rebuild the hierarchy using the new system matrix. Requires allow_rebuild to be set in parameters. The transfer operators created during the initial setup are reused in order to create the coarse system matrices: \(A_c^{\text{new}} = R A^{new} P\).

Single-level relaxation¶

template<class Backend, template<class> class Relax> class amgcl::relaxation::as_preconditioner¶

Include <amgcl/relaxation/as_preconditioner.hpp>

Allows to use a relaxation method as a standalone preconditioner.

Relax<backend> smoother;: The relaxation class instantiated for the specified backend

typedef typename smoother::params params¶: The relaxation params are inherited as the parameters for the preconditioner

Dummy¶

template<class Backend> class amgcl::preconditioner::dummy¶

Include <amgcl/preconditioner/dummy.hpp>

The dummy preconditioner, equivalent to an identity matrix. May be used to test the convergence of unpreconditioned iterative solvers.

class params¶: There are no parameters

Composite preconditioners¶

The preconditioners in this section take the into account the block structure of the system and properties of the individual blocks. Most often the preconditioners are used for the solution of saddle point or Stokes-like systems, where the system matrix may be represented in the following form:

(1)¶\[\begin{split} \begin{pmatrix} A & B_1^T \\ B_2 & C \end{pmatrix} \begin{pmatrix} u \\ p \end{pmatrix} = \begin{pmatrix} b_u \\ b_p \end{pmatrix}\end{split}\]

CPR¶

template<class PPrecond, class SPrecond> class amgcl::preconditioner::cpr¶

Include <amgcl/preconditioner/cpr.hpp>

The Constrained Pressure Residual (CPR) preconditioner [Stue07]. The CPR preconditioners are based on the idea that coupled system solutions are mainly determined by the solution of their elliptic components (i.e., pressure). Thus, the procedure consists of extracting and accurately solving pressure subsystems. Residuals associated with this solution are corrected with an additional preconditioning step that recovers part of the global information contained in the original system.

The template parameters PPrecond and SPrecond for the CPR preconditioner specify which preconditioner to use with the pressure subblock (the \(C\) matrix in (1)), and with the complete system.

The system matrix should be ordered by grid nodes, so that the pressure and suturation/concentration unknowns belonging to the same grid node are compactly located together in the vector of unknowns. The pressure should be the first unknown in the block of unknowns associated with a grid node.

class params¶

The CPR preconditioner parameters

typedef typename SPrecond::value_type value_type¶: The value type of the system matrix

typedef typename PPrecond::params pprecond_params¶: The type of the pressure preconditioner parameters

typedef typename SPrecond::params sprecond_params¶: The type of the global preconditioner parameters

pprecond_params pprecond¶: The pressure preconditioner parameters

sprecond_params sprecond¶: The global preconditioner parameters

int block_size = 2¶: The number of unknowns associated with each grid node. The default value is 2 when the system matrix has scalar value type. Otherwise, the block size of the system matrix value type is used.

size_t active_rows = 0¶: When non-zero, only unknowns below this number are considered to be pressure. May be used when a system matrix contains unstructured tail block (for example, the unknowns associated with wells).

CPR (DRS)¶

template<class PPrecond, class SPrecond> class amgcl::preconditioner::cpr_drs¶

Include <amgcl/preconditioner/cpr.hpp>

The Constrained Pressure Residual (CPR) preconditioner with weighted dynamic row sum (WDRS) [Grie14], [BrCC15].

The template parameters PPrecond and SPrecond for the CPR WDRS preconditioner specify which preconditioner to use with the pressure subblock (the \(C\) matrix in (1)), and with the complete system.

The system matrix should be ordered by grid nodes, so that the pressure and suturation/concentration unknowns belonging to the same grid node are compactly located together in the vector of unknowns. The pressure should be the first unknown in the block of unknowns associated with a grid node.

class params¶

The CPR preconditioner parameters

typedef typename SPrecond::value_type value_type¶: The value type of the system matrix

typedef typename PPrecond::params pprecond_params¶: The type of the pressure preconditioner parameters

typedef typename SPrecond::params sprecond_params¶: The type of the global preconditioner parameters

pprecond_params pprecond¶: The pressure preconditioner parameters

sprecond_params sprecond¶: The global preconditioner parameters

int block_size = 2¶: The number of unknowns associated with each grid node. The default value is 2 when the system matrix has scalar value type. Otherwise, the block size of the system matrix value type is used.

size_t active_rows = 0¶: When non-zero, only unknowns below this number are considered to be pressure. May be used when a system matrix contains unstructured tail block (for example, the unknowns associated with wells).

double eps_dd = 0.2¶: Controls the severity of the violation of diagonal dominance. See [Grie14] for more details.

double eps_ps = 0.02¶: Controls the pressure/saturation coupling. See [Grie14] for more details.

std::vector<double> weights¶: The weights for the weighted DRS method. See [BrCC15] for more details.

Schur Pressure Correction¶

template<class USolver, class PSolver> class amgcl::preconditioner::schur_pressure_correction¶

Include <amgcl/preconditioner/schur_pressure_correction.hpp>

The system (1) may be rewritten as

\[\begin{split}\begin{pmatrix} A & B_1^T \\ 0 & S \end{pmatrix} \begin{pmatrix} u \\ p \end{pmatrix} = \begin{pmatrix} b_u \\ b_p - B_2 A^{-1} b_u \end{pmatrix}\end{split}\]

where \(S = C - B_2 A^{-1} B_1^T\) is the Schur complement. The Schur complement pressure correction preconditioner uses this representation and an approximation to the Schur complement matrix in order to decouple the pressure and the velocity parts of the system [ElHS08].

The two template parameters for the method, USolver and PSolver, specify the preconditioned solvers for the \(A\) and \(S\) blocks.

class params¶

The parameters for the Schur pressure correction preconditioner

typedef typename USolver::params usolver_params¶: The type of the USolver parameters

typedef typename PSolver::params psolver_params¶: The type of the PSolver parameters

usolver_params usolver¶: The USolver parameters

psolver_params psolver¶: The PSolver parameters

std::vector<char> pmask¶: The indicator vector, containing 1 for pressure unknowns, and 0 otherwise.

int type = 1¶

The variant of the block preconditioner to use.

When type = 1:

\[\begin{split}\begin{aligned} S p &= b_p - B_2 A^{-1} b_u \\ A u &= b_u - B_1^T p \end{aligned}\end{split}\]
When type = 2:

\[\begin{split}\begin{aligned} S p &= b_p \\ A u &= b_u - B_1^T p \end{aligned}\end{split}\]

bool approx_schur = false¶: When set, approximate \(A^{-1}\) as \(\mathrm{diag}(A)^{-1}\) during computation of the matrix-less Schur complement when solving the \(Sp=b_p\) system. Otherwise, the full solve using USolver is used.

int adjust_p = 1¶

Adjust the matrix used to construct the preconditioner for the Schur complement system.

When adjust_p = 0, use the unmodified \(C\) matrix;
When adjust_p = 1, use \(C - \mathrm{diag}(B_2 \mathrm{diag}(A)^{-1} B_1^T)\);
When adjust_p = 2, use \(C - B_2 \mathrm{diag}(A)^{-1} B_1^T\).

bool simplec_dia = true¶: When set, use \(\frac{1}{\sum_j \|A_{i,j}\|}\) instead of \(\mathrm{diag}(A)^{-1}\) as the approximation for \(A^{-1}\) (similar to the SIMPLEC algorithm).

int verbose = 0¶

When verbose >= 1, show the number of iterations and the relative residual achieved after each nested solve.
When verbose >= 2, save the \(A\) and \(C\) submatrices as Kuu.mtx and Kpp.mtx.