 |
Block-Structured AMR Software Framework
|
|
Block-Structured AMR Software Framework
|
Algebraic multigrid solver for SpMatrix / AlgVector systems. More...
#include <AMReX_AlgMG.H>
Public Types | |
| using | TS = std::conditional_t<std::is_same_v<T,double>, float, T> |
| using | BottomSolver = AlgMGBottomSolver |
| using | KrylovSolver = AlgMGKrylovSolver |
| using | InterpType = AlgMGInterpType |
| using | Smoother = AlgMGSmoother |
Public Member Functions | |
| AlgMG ()=default | |
| AlgMG (SpMatrix< T > &a_mat) | |
| void | define (SpMatrix< T > &a_mat) |
| Define the solver with a square matrix. Must be called before solve. | |
| bool | updateMatrix (SpMatrix< T > &a_mat) |
| Replace the matrix, keeping the hierarchy. | |
| void | solve (AlgVector< T > &a_sol, AlgVector< T > const &a_rhs) |
| Solve A x = b. a_sol holds the initial guess on entry. | |
| void | setVerbose (int v) |
| void | setPrintIndentation (std::string s) |
| Prefix of the lines the solver prints, e.g. to nest them in MLMG output. | |
| void | setMaxIter (int n) |
| void | setFixedIter (int n) |
| If positive, run exactly this many V-cycles regardless of tolerance. | |
| void | setRelTol (T t) |
| void | setAbsTol (T t) |
| Converged when the residual 2-norm is below max(reltol*bnorm, abstol). | |
| void | setSmoother (Smoother sm) |
| Takes effect at the next solve; the hierarchy is kept. | |
| void | setRelaxWeight (T w) |
| Relaxation weight of the smoother; negative restores the default. Next solve. | |
| void | setChebyshevDegree (int d) |
| Chebyshev polynomial degree. Next solve. | |
| void | setChebyshevRatio (T r) |
| Chebyshev eigenvalue ratio lambda_max / lambda_min. Next solve. | |
| void | setPreSmooth (int nu) |
| Number of pre-smoothing sweeps; negative restores the smoother's default. | |
| void | setPostSmooth (int nu) |
| Number of post-smoothing sweeps; negative restores the smoother's default. | |
| void | setBottomSmooth (int nu) |
| void | setStrongThreshold (T theta) |
| void | setMaxLevels (int n) |
| void | setMaxCoarseSize (Long n) |
| Stop coarsening when a level has at most this many rows. Next setup. | |
| void | setAggressiveNumLevels (int n) |
| Number of levels with aggressive coarsening, starting at the finest. Next setup. | |
| void | setAggressiveDirectInterp (bool b) |
| void | setSingular (bool b) |
| void | setInterpType (InterpType it) |
| Takes effect at the next setup. | |
| void | setPMaxElmts (int n) |
| void | setTruncFactor (T f) |
| Drop interpolation weights below this fraction of the row maximum. Next setup. | |
| void | setKrylovSolver (KrylovSolver a) |
| Krylov solver with the V-cycle as preconditioner (see KrylovSolver). Default None. | |
| void | precond (AlgVector< T > &x, AlgVector< T > const &b) |
| void | setBottomSolver (BottomSolver bs) |
| void | setBottomTol (T t) |
| void | setBottomMaxIter (int n) |
| void | setBottomVerbose (int v) |
| void | setThrowException (bool b) |
| int | getNumIters () const |
| int | getNumSetups () const |
| Number of setups run so far (for checks that the setup is reused). | |
| T | getResidualNorm () const |
| int | numLevels () const |
| Number of levels; 0 before setup. | |
| std::size_t | getSetupPeakBytes () const |
| Peak bytes in use in The_Arena during setup (0 if not available). | |
| Long | numGlobalRows (int lev) const |
| Global number of rows on a level. Valid after setup. | |
| SpMatrix< T > const & | getMatrix (int lev) const |
| SpMatrix< T > const & | getInterp (int lev) const |
| void | setup () |
| bool | coarsen_level (int lev, SpMatrix< T > &A, SpMatrix< T > &P, AlgPartition &cpart, bool &aggressive, Long &ncoarse, bool &truncated, Long &nnz_untruncated) |
| SpMatrix< TS > | create_soc (SpMatrix< T > const &A, T scale) const |
| void | pmis (SpMatrix< TS > &S, SpMatrix< TS > &ST, Gpu::DeviceVector< int > &state, bool isolated_as_coarse=false) const |
| void | project_out (AlgVector< T > &x) |
| SpMatrix< T > | interpolation (SpMatrix< T > &A, SpMatrix< TS > &S, SpMatrix< TS > const &ST, Gpu::DeviceVector< int > const &state, Gpu::DeviceVector< Long > const &cidx, AlgPartition const &cpart, T scale) const |
| Interpolation of the selected type from the C-points in state. | |
Static Public Member Functions | |
| static T | strength_scale (SpMatrix< T > const &A) |
| Largest |a_ii|, the unit of the strength matrix entries. | |
| static AlgPartition | coarse_numbering (Gpu::DeviceVector< int > const &state, Gpu::DeviceVector< Long > &cidx, int which=PointType::coarse) |
| Local index among the rows of type which and their partition. | |
| static SpMatrix< T > | build_interp (SpMatrix< T > &A, SpMatrix< TS > const &S, Gpu::DeviceVector< int > const &state, Gpu::DeviceVector< Long > const &cidx, AlgPartition const &cpart) |
| Classical direct interpolation. | |
| static SpMatrix< T > | truncate_interp (SpMatrix< T > &P, AlgPartition const &cpart, int max_elmts, T trunc_factor) |
| static void | normalize_rows (SpMatrix< T > &P) |
| Scales the rows of P to sum to one; aborts on an empty row. | |
| static SpMatrix< TS > | second_pass_strength (SpMatrix< TS > &S, Gpu::DeviceVector< int > const &state, Gpu::DeviceVector< Long > const &cidx, AlgPartition const &cpart) |
| static SpMatrix< T > | select_coarse_rows (SpMatrix< T > const &P, Gpu::DeviceVector< int > const &state, Gpu::DeviceVector< Long > const &cidx, AlgPartition const &cpart) |
| Rows of P at the C-points, in C-point numbering. | |
| static SpMatrix< T > | build_interp_mm (SpMatrix< T > &A, SpMatrix< TS > &S, SpMatrix< TS > const &ST, Gpu::DeviceVector< int > const &state, Gpu::DeviceVector< Long > const &cidx, AlgPartition const &cpart, T scale, bool plus_i, int const *row_state=nullptr, int max_elmts=0, T trunc_factor=T(0), Long *nnz_untruncated=nullptr) |
Static Public Attributes | |
| static constexpr int | max_p_elmts = 32 |
| static constexpr int | default_p_max_elmts = 4 |
Algebraic multigrid solver for SpMatrix / AlgVector systems.
The setup selects coarse points with PMIS [2] on the classical strength of connection [1], builds the interpolation P (extended+i [3] in the matrix-matrix form of [4] by default, or direct [1]) with at most p_max_elmts entries per row [3], and forms A_c = P^T A P until the coarsest level has at most max_coarse_size rows. The first aggressive_num_levels levels use aggressive coarsening [5]. The solve runs V-cycles with Chebyshev [7], l1-Jacobi [6], weighted Jacobi or, on CPUs, l1 hybrid Gauss-Seidel [6] smoothing, on their own or as the preconditioner of a Krylov solver (setKrylovSolver).
Singular matrices with the constant null vector (setSingular) are handled by removing the mean of the right-hand side and of the coarsest-level residuals.
The matrix must be square with rows sorted by column index (see SpMatrix::sortCSR), and it must stay alive while the solver is in use. With R = P^T, the coarsening assumes a symmetric or nearly symmetric matrix. Parameters are set per object with the setters below; the Linear Solvers chapter of the documentation describes their effect.
References (full citations in the documentation): [1] Ruge & Stüben, Algebraic multigrid, SIAM 1987. doi:10.1137/1.9781611971057.ch4 [2] De Sterck, Yang & Heys, SIAM J. Matrix Anal. Appl. 27 (2006). doi:10.1137/040615729 [3] De Sterck, Falgout, Nolting & Yang, Numer. Linear Algebra Appl. 15 (2008). doi:10.1002/nla.559 [4] Li, Sjögreen & Yang, SIAM J. Sci. Comput. 43 (2021). doi:10.1137/20M134931X [5] Yang, Numer. Linear Algebra Appl. 17 (2010). doi:10.1002/nla.689 [6] Baker, Falgout, Kolev & Yang, SIAM J. Sci. Comput. 33 (2011). doi:10.1137/100798806 [7] Adams, Brezina, Hu & Tuminaro, J. Comput. Phys. 188 (2003). doi:10.1016/S0021-9991(03)00194-3
| using amrex::AlgMG< T >::BottomSolver = AlgMGBottomSolver |
| using amrex::AlgMG< T >::InterpType = AlgMGInterpType |
| using amrex::AlgMG< T >::KrylovSolver = AlgMGKrylovSolver |
| using amrex::AlgMG< T >::Smoother = AlgMGSmoother |
| using amrex::AlgMG< T >::TS = std::conditional_t<std::is_same_v<T,double>, float, T> |
Value type of the strength matrices: float when T is double. Their entries are the matrix entries divided by strength_scale.
|
default |
|
explicit |
|
static |
Classical direct interpolation.
|
static |
Extended (plus_i = false) or extended+i interpolation in matrix-matrix form. ST is the transpose of S, scale its unit. With row_state, only its C-points get a row of P; the other rows are empty. In CPU builds, P is also truncated as in truncate_interp if max_elmts or trunc_factor is positive, and nnz_untruncated (if not null) gets the number of weights before.
|
static |
Local index among the rows of type which and their partition.
| bool amrex::AlgMG< T >::coarsen_level | ( | int | lev, |
| SpMatrix< T > & | A, | ||
| SpMatrix< T > & | P, | ||
| AlgPartition & | cpart, | ||
| bool & | aggressive, | ||
| Long & | ncoarse, | ||
| bool & | truncated, | ||
| Long & | nnz_untruncated ) |
Selects the C-points of level lev and builds its interpolation P (rows partitioned by cpart). Returns false when coarsening stops. The strength matrices are freed on return, before the triple product. If truncated, P is already truncated and nnz_untruncated holds its number of weights before.
|
nodiscard |
Strong part of A: S_ij = a_ij / scale if i strongly depends on j, else absent. scale is strength_scale(A).
| void amrex::AlgMG< T >::define | ( | SpMatrix< T > & | a_mat | ) |
Define the solver with a square matrix. Must be called before solve.
|
inlinenodiscard |
|
inlinenodiscard |
|
inlinenodiscard |
|
inlinenodiscard |
Number of setups run so far (for checks that the setup is reused).
|
inlinenodiscard |
|
inlinenodiscard |
Peak bytes in use in The_Arena during setup (0 if not available).
| SpMatrix< T > amrex::AlgMG< T >::interpolation | ( | SpMatrix< T > & | A, |
| SpMatrix< TS > & | S, | ||
| SpMatrix< TS > const & | ST, | ||
| Gpu::DeviceVector< int > const & | state, | ||
| Gpu::DeviceVector< Long > const & | cidx, | ||
| AlgPartition const & | cpart, | ||
| T | scale ) const |
Interpolation of the selected type from the C-points in state.
|
static |
Scales the rows of P to sum to one; aborts on an empty row.
|
inlinenodiscard |
Global number of rows on a level. Valid after setup.
|
inlinenodiscard |
Number of levels; 0 before setup.
| void amrex::AlgMG< T >::pmis | ( | SpMatrix< TS > & | S, |
| SpMatrix< TS > & | ST, | ||
| Gpu::DeviceVector< int > & | state, | ||
| bool | isolated_as_coarse = false ) const |
PMIS C/F selection. state gets PointType values. Points without any strong connection become F, or C if isolated_as_coarse.
| void amrex::AlgMG< T >::precond | ( | AlgVector< T > & | x, |
| AlgVector< T > const & | b ) |
Apply one V-cycle to b with zero initial guess: x = M^{-1} b. This is the preconditioner used by the Krylov solver; usable stand-alone after the first solve or setup.
|
inline |
Singular systems: remove the mean of x (the constant null vector is the same on every level since interpolation preserves constants).
|
static |
Strength graph between C-points connected through a strong path of length two, in C-point numbering (rows partitioned by cpart).
|
static |
Rows of P at the C-points, in C-point numbering.
|
inline |
Converged when the residual 2-norm is below max(reltol*bnorm, abstol).
|
inline |
Use classical direct interpolation for the first stage of the two-stage aggressive interpolation (cheaper setup, more cycles).
|
inline |
Number of levels with aggressive coarsening, starting at the finest. Next setup.
|
inline |
|
inline |
Number of smoother sweeps on the coarsest level with the Jacobi bottom solver; negative restores the default.
|
inline |
Coarsest-level solver (see BottomSolver). Default direct, which uses smoother sweeps instead beyond 1024 rows.
|
inline |
|
inline |
|
inline |
Chebyshev polynomial degree. Next solve.
|
inline |
Chebyshev eigenvalue ratio lambda_max / lambda_min. Next solve.
|
inline |
If positive, run exactly this many V-cycles regardless of tolerance.
|
inline |
Takes effect at the next setup.
|
inline |
Krylov solver with the V-cycle as preconditioner (see KrylovSolver). Default None.
|
inline |
Stop coarsening when a level has at most this many rows. Next setup.
|
inline |
|
inline |
Maximum number of levels. Two means a two-level method. Takes effect at the next setup, i.e., the next solve after define.
|
inline |
Maximum number of entries per row of P (0: no limit; at most max_p_elmts; negative restores the default). Next setup.
|
inline |
Number of post-smoothing sweeps; negative restores the smoother's default.
|
inline |
Number of pre-smoothing sweeps; negative restores the smoother's default.
|
inline |
Prefix of the lines the solver prints, e.g. to nest them in MLMG output.
|
inline |
Relaxation weight of the smoother; negative restores the default. Next solve.
|
inline |
|
inline |
Singular matrix with the constant vector in its null space (e.g. Poisson with periodic or Neumann boundaries). The mean is removed from the right-hand side and the coarsest residuals, and the solution is returned with zero mean. Next setup.
|
inline |
Takes effect at the next solve; the hierarchy is kept.
|
inline |
Strength threshold, typically 0.1 to 0.5 (0.25): connections at least this fraction of the strongest one in the row are strong. A row whose sum exceeds 0.9 times its diagonal has no strong connection. Takes effect at the next setup, i.e., the next solve after define.
|
inline |
|
inline |
Drop interpolation weights below this fraction of the row maximum. Next setup.
| void amrex::AlgMG< T >::setup | ( | ) |
|
inline |
| void amrex::AlgMG< T >::solve | ( | AlgVector< T > & | a_sol, |
| AlgVector< T > const & | a_rhs ) |
Solve A x = b. a_sol holds the initial guess on entry.
|
static |
Largest |a_ii|, the unit of the strength matrix entries.
|
static |
Keep the max_elmts largest weights of each row of P (0: all) and drop weights below trunc_factor times the row maximum. Kept weights are rescaled to preserve the row's positive and negative sums.
| bool amrex::AlgMG< T >::updateMatrix | ( | SpMatrix< T > & | a_mat | ) |
Replace the matrix, keeping the hierarchy.
The interpolation and coarse matrices of the last setup are kept, and only the smoothers (and a direct bottom solve on the finest level) are rebuilt, so the multigrid cycle then only approximates the new matrix. Cheaper than define() when the matrix changes little between solves. The matrix must have the same row partition. If no hierarchy has been set up since the last define() (no solve yet), this is define().
a_mat
|
staticconstexpr |
|
staticconstexpr |
1.15.0