Coarsening Strategies

A coarsening strategy defines various options for creating coarse systems in the AMG hierarchy. A coarsening strategy takes the system matrix \(A\) at the current level, and returns prolongation operator \(P\) and the corresponding restriction operator \(R\).

Ruge-Stuben

template<class Backend>
class amgcl::coarsening::ruge_stuben

Include <amgcl/coarsening/ruge_stuben>

The classic Ruge-Stuben coarsening with direct interpolation [Stue99].

class params
float eps_strong = 0.25

Parameter \(\varepsilon_{str}\) defining strong couplings.

Variable \(i\) is defined to be strongly negatively coupled to another variable, \(j\), if

\[-a_{ij} \geq \varepsilon_{str}\max\limits_{a_{ik}<0}|a_{ik}|\quad \text{with fixed} \quad 0 < \varepsilon_{str} < 1.\]

In practice, a value of \(\varepsilon_{str}=0.25\) is usually taken.

bool do_trunc = true

Prolongation operator truncation. Interpolation operators, and, hence coarse operators may increase substantially towards coarser levels. Without truncation, this may become too costly. Truncation ignores all interpolatory connections which are smaller (in absolute value) than the largest one by a factor of \(\varepsilon_{tr}\). The remaining weights are rescaled so that the total sum remains unchanged. In practice, a value of \(\varepsilon_{tr}=0.2\) is usually taken.

float eps_trunc = 0.2

Truncation parameter \(\varepsilon_{tr}\).

Aggregation-based coarsening

The aggregation-base class of coarsening methods split the nodes at the fine grid into disjoint sets of nodes, the so-called aggregates that act as nodes on the coarse grid. The prolongation operators are then built by first constructing a tentative prolongator using the knowledge of zero energy modes of the principal part of the differential operator with natural boundary conditions (e.g., near null-space vectors, or rigid body modes for elasticity). In case of smoothed aggregation the prolongation operator is then smoothed by a carefully selected iteration.

All of the aggregation based methods take the aggregation and the nullspace parameters:

class amgcl::coarsening::pointwise_aggregates

Pointwise aggregation. When the system matrix has a block structure, it is converted to a poinwise matrix (single value per block), and the aggregates are created for this reduced matrix instead.

class params

The aggregation parameters.

float eps_strong = 0.08

Parameter \(\varepsilon_{strong}\) defining strong couplings. Connectivity is defined in a symmetric way, that is, two variables \(i\) and \(j\) are considered to be connected to each other if \(\frac{a_{ij}^2}{a_{ii}a_{jj}} > \varepsilon_{strong}\) with fixed \(0 < \varepsilon_{strong} < 1\).

int block_size = 1

The block size in case the system matrix has a block structure.

class amgcl::coarsening::nullspace_params

The nullspace parameters.

int cols = 0

The number of near nullspace vectors.

std::vector<double> B

The near nullspace vectors. The vectors are represented as columns of a 2D matrix stored in row-major order.

Aggregation

template<class Backend>
class amgcl::coarsening::aggregation

Include <amgcl/coarsening/aggregation.hpp>

The non-smoothed aggregation coarsening [Stue99].

class params

The aggregation coarsening parameters

amgcl::coarsening::pointwise_aggregates::params aggr;

The aggregation parameters

amgcl::coarsening::nullspace_params nullspace

The near nullspace parameters

float over_interp = 1.5

Over-interpolation factor \(\alpha\) [Stue99]. In case of aggregation coarsening, coarse-grid correction of smooth error, and by this the overall convergence, can often be substantially improved by using “over-interpolation”, that is, by multiplying the actual correction (corresponding to piecewise constant interpolation) by some factor \(\alpha > 1\). Equivalently, this means that the coarse-level Galerkin operator is re-scaled by \(1/\alpha\):

\[I_h^HA_hI_H^h \to \frac{1}{\alpha}I_h^HA_hI_H^h.\]

Smoothed Aggregation

template<class Backend>
class amgcl::coarsening::smoothed_aggregation

Include <amgcl/coarsening/smoothed_aggregation.hpp>

The smoothed aggregation coarsening [VaMB96].

class params

The smoothed aggregation coarsening parameters

amgcl::coarsening::pointwise_aggregates::params aggr;

The aggregation parameters

amgcl::coarsening::nullspace_params nullspace

The near nullspace parameters

float relax = 1.0

The relaxation factor \(r\). Used as a scaling for the damping factor \(\omega\). When estimate_spectral_radius is set, then

\[\omega = r * (4/3) / \rho.\]

where \(\rho\) is the spectral radius of the system matrix. Otherwise

\[\omega = r * (2/3).\]

The tentative prolongation \(\tilde P\) from the non-smoothed aggregation is improved by smoothing to get the final prolongation matrix \(P\). Simple Jacobi smoother is used here, giving the prolongation matrix

\[P = \left( I - \omega D^{-1} A^F \right) \tilde P.\]

Here \(A^F = (a_{ij}^F)\) is the filtered matrix given by

\[\begin{split}\begin{aligned} a_{ij}^F &= \begin{cases} a_{ij} \quad \text{if} \; j \in N_i\\ 0 \quad \text{otherwise} \end{cases}, \quad \text{if}\; i \neq j, \\ \quad a_{ii}^F &= a_{ii} - \sum\limits_{j=1,j\neq i}^n \left(a_{ij} - a_{ij}^F \right), \end{aligned}\end{split}\]

where \(N_i\) is the set of variables strongly coupled to variable \(i\), and \(D\) denotes the diagonal of \(A^F\).

bool estimate_spectral_radius = false

Estimate the matrix spectral radius. This usually improves the convergence rate and results in faster solves, but costs some time during setup.

int power_iters = 0

The number of power iterations to apply for the spectral radius estimation. Use Gershgorin disk theorem when power_iters = 0.

Smoothed Aggregation with Energy Minimization

template<class Backend>
class amgcl::coarsening::smoothed_aggr_emin

Include <amgcl/coarsening/smoothed_aggr_emin.hpp>

The smoothed aggregation with energy minimization coarsening [SaTu08].

class params

The smoothed aggregation with energy minimization coarsening parameters

amgcl::coarsening::pointwise_aggregates::params aggr;

The aggregation parameters

amgcl::coarsening::nullspace_params nullspace

The near nullspace parameters

Block to Scalar convertor

template<template<class> class Coarsening>
class amgcl::coarsening::as_scalar

Include <amgcl/coarsening/as_scalar.hpp>

Wrapper for the specified coarsening. Converts the input matrix from block to scalar format, applies the base coarsening, converts the resulting transfer operators back to block format. See the Using near null-space vectors tutorial.

template <class Backend>
class type

The resulting coarsening class.