Iterative solvers

AMGCL provides several iterative solvers, but it should be easy to use AMGCL preconditioners with a user-provided solver as well. Each solver in AMGCL is a class template. Its single template parameter specifies the backend to use. This allows to preallocate necessary resources at class construction. Obviously, the solver backend has to coincide with the preconditioner backend.

Each of the solvers in AMGCL provides two overloads for the operator():

boost::tuple<size_t, scalar_type> operator()(const Matrix &A, const Precond &P, const Vec1 &rhs, Vec2 &&x) const

Computes the solution for the given system matrix A and the right-hand side rhs. Returns the number of iterations made and the achieved residual as a boost::tuple. The solution vector x provides initial approximation on input and holds the computed solution on output.

The system matrix may differ from the matrix used during initialization. This may be used for the solution of non-stationary problems with slowly changing coefficients. There is a strong chance that a preconditioner built for a time step will act as a reasonably good preconditioner for several subsequent time steps [DeSh12].

boost::tuple<size_t, scalar_type> operator()(const Precond &P, const Vec1 &rhs, Vec2 &&x) const

Computes the solution for the given right-hand side rhs. The system matrix is the same that was used for the setup of the preconditioner P. Returns the number of iterations made and the achieved residual as a boost::tuple. The solution vector x provides initial approximation on input and holds the computed solution on output.

Conjugate Gradient

#include <amgcl/solver/cg.hpp>

template <class Backend, class InnerProduct = detail::default_inner_product>
class amgcl::solver::cg

Conjugate Gradients method.

An effective method for symmetric positive definite systems [Barr94].

Public Functions

cg(size_t n, const params &prm = params (), const backend_params &backend_prm = backend_params(), const InnerProduct &inner_product = InnerProduct())

Preallocates necessary data structures for the system of size n.

struct params

Solver parameters.

Public Members

template<>
size_t maxiter

Maximum number of iterations.

template<>
scalar_type tol

Target relative residual error.

template<>
scalar_type abstol

Target absolute residual error.

BiConjugate Gradient Stabilized (BiCGSTAB)

#include <amgcl/solver/bicgstab.hpp>

template <class Backend, class InnerProduct = detail::default_inner_product>
class amgcl::solver::bicgstab

BiConjugate Gradient Stabilized (BiCGSTAB) method.

The BiConjugate Gradient Stabilized method (Bi-CGSTAB) was developed to solve nonsymmetric linear systems while avoiding the often irregular convergence patterns of the Conjugate Gradient [Barr94].

Public Functions

bicgstab(size_t n, const params &prm = params (), const backend_params &backend_prm = backend_params(), const InnerProduct &inner_product = InnerProduct())

Preallocates necessary data structures for the system of size n.

struct params

Solver parameters.

Public Members

template<>
preconditioner::side::type pside

Preconditioning kind (left/right).

template<>
size_t maxiter

Maximum number of iterations.

template<>
scalar_type tol

Target relative residual error.

template<>
scalar_type abstol

Target absolute residual error.

BiCGSTAB(L)

#include <amgcl/solver/bicgstabl.hpp>

template <class Backend, class InnerProduct = detail::default_inner_product>
class amgcl::solver::bicgstabl

BiCGStab(L) method.

Generalization of BiCGStab method [SlDi93].

Public Functions

bicgstabl(size_t n, const params &prm = params (), const backend_params &backend_prm = backend_params(), const InnerProduct &inner_product = InnerProduct())

Preallocates necessary data structures for the system of size n.

struct params

Solver parameters.

GMRES

#include <amgcl/solver/gmres.hpp>

template <class Backend, class InnerProduct = detail::default_inner_product>
class amgcl::solver::gmres

Generalized Minimal Residual (GMRES) method.

The Generalized Minimal Residual method is an extension of MINRES (which is only applicable to symmetric systems) to unsymmetric systems [Barr94].

Public Functions

gmres(size_t n, const params &prm = params (), const backend_params &backend_prm = backend_params(), const InnerProduct &inner_product = InnerProduct())

Preallocates necessary data structures for the system of size n.

struct params

Solver parameters.

Public Members

template<>
unsigned M

Number of iterations before restart.

template<>
preconditioner::side::type pside

Preconditioning kind (left/right).

template<>
unsigned maxiter

Maximum number of iterations.

template<>
scalar_type tol

Target relative residual error.

template<>
scalar_type abstol

Target absolute residual error.

LGMRES

#include <amgcl/solver/lgmres.hpp>

template <class Backend, class InnerProduct = detail::default_inner_product>
class amgcl::solver::lgmres

“Loose” GMRES.

The LGMRES algorithm [BaJM05] is designed to avoid some problems in the convergence in restarted GMRES, and often converges in fewer iterations.

Public Functions

lgmres(size_t n, const params &prm = params (), const backend_params &bprm = backend_params(), const InnerProduct &inner_product = InnerProduct())

Preallocates necessary data structures for the system of size n.

struct params

Solver parameters.

Public Members

template<>
unsigned M

Number of inner GMRES iterations per each outer iteration.

template<>
unsigned K

Number of vectors to carry between inner GMRES iterations.

According to [BaJM05], good values are in the range of 1...3. However, note that if you want to use the additional vectors to accelerate solving multiple similar problems, larger values may be beneficial.

template<>
bool always_reset

Reset augmented vectors between solves.

If the solver is used to repeatedly solve similar problems, then keeping the augmented vectors between solves may speed up subsequent solves. This flag, when set, resets the augmented vectors at the beginning of each solve.

template<>
bool store_Av

Whether LGMRES should store also A*v in addition to vectors v.

template<>
preconditioner::side::type pside

Preconditioning kind (left/right).

template<>
size_t maxiter

Maximum number of iterations.

template<>
scalar_type tol

Target relative residual error.

template<>
scalar_type abstol

Target absolute residual error.

FGMRES

#include <amgcl/solver/fgmres.hpp>

template <class Backend, class InnerProduct = detail::default_inner_product>
class amgcl::solver::fgmres

Flexible GMRES method.

Flexible version of the GMRES method [Saad03].

Public Functions

fgmres(size_t n, const params &prm = params (), const backend_params &bprm = backend_params(), const InnerProduct &inner_product = InnerProduct())

Preallocates necessary data structures for the system of size n.

struct params

Solver parameters.

Public Members

template<>
unsigned M

Number of inner GMRES iterations per each outer iteration.

template<>
unsigned maxiter

Maximum number of iterations.

template<>
scalar_type tol

Target relative residual error.

template<>
scalar_type abstol

Target absolute residual error.