pyamg.classical#
Classical AMG.
- pyamg.classical.air_solver(A, strength=('classical', {'norm': 'min', 'theta': 0.3}), CF=('RS', {'second_pass': True}), interpolation='one_point', restrict=('air', {'degree': 2, 'theta': 0.05}), presmoother=None, postsmoother=('fc_jacobi', {'c_iterations': 1, 'f_iterations': 2, 'iterations': 1, 'omega': 1.0, 'withrho': False}), filter_operator=None, max_levels=20, max_coarse=20, keep=False, **kwargs)[source]#
Create a multilevel solver using approximate ideal restriction (AIR) AMG.
- Parameters:
- Acsr_array
Square (non)symmetric matrix in CSR format
- strength[‘symmetric’, ‘classical’, ‘evolution’, ‘distance’,
‘algebraic_distance’,’affinity’, ‘energy_based’, None]
Method used to determine the strength of connection between unknowns of the linear system. Method-specific parameters may be passed in using a tuple, e.g. strength=(‘symmetric’,{‘theta’ : 0.25 }). If strength=None, all nonzero entries of the matrix are considered strong.
- CF{string}default ‘RS’ with second pass
Method used for coarse grid selection (C/F splitting) Supported methods are RS, PMIS, PMISc, CLJP, CLJPc, and CR.
- interpolation{string}default ‘one_point’
Options include ‘direct’, ‘classical’, ‘inject’ and ‘one-point’.
- restrict{string}default distance-2 AIR, with theta = 0.05.
Option is ‘air’ for local approximate ideal restriction (lAIR), with inner options specifying degree, strength tolerance, etc.
- presmoother{string or dict}default None
Method used for presmoothing at each level. Method-specific parameters may be passed in using a tuple.
- postsmoother{string or dict}
Postsmoothing method with the same usage as presmoother. postsmoother=(‘fc_jacobi’, … ) with 2 F-sweeps, 1 C-sweep is default.
- filter_operator(bool, tol)default None
Remove small entries in operators on each level if True. Entries are considered “small” if |a_ij| < tol |a_ii|.
- max_levels: {integer}default 20
Maximum number of levels to be used in the multilevel solver.
- max_coarse: {integer}default 20
Maximum number of variables permitted on the coarse grid.
- keep: {bool}default False
Flag to indicate keeping strength of connection matrix (C) in hierarchy.
- kwargsdict
Extra keywords passed to the Multilevel class
- Returns:
- mlmultilevel_solver
Multigrid hierarchy of matrices and prolongation operators
See also
aggregation.smoothed_aggregation_solver
,multilevel_solver
aggregation.rootnode_solver
,ruge_stuben_solver
Notes
“coarse_solver” is an optional argument and is the solver used at the coarsest grid. The default is a pseudo-inverse. Most simply, coarse_solver can be one of [‘splu’, ‘lu’, ‘cholesky, ‘pinv’, ‘gauss_seidel’, … ]. Additionally, coarse_solver may be a tuple (fn, args), where fn is a string such as [‘splu’, ‘lu’, …] or a callable function, and args is a dictionary of arguments to be passed to fn. See [2001TrOoSc] for additional details.
References
[1] Manteuffel, T. A., Münzenmaier, S., Ruge, J., & Southworth, B. S. (2019). Nonsymmetric reduction-based algebraic multigrid. SIAM Journal on Scientific Computing, 41(5), S242-S268.
[2] Manteuffel, T. A., Ruge, J., & Southworth, B. S. (2018). Nonsymmetric algebraic multigrid based on local approximate ideal restriction (lAIR). SIAM Journal on Scientific Computing, 40(6), A4105-A4130.
Examples
>>> from pyamg.gallery import poisson >>> from pyamg import air_solver >>> A = poisson((10,),format='csr') >>> ml = air_solver(A,max_coarse=3)
- pyamg.classical.ruge_stuben_solver(A, strength=('classical', {'theta': 0.25}), CF=('RS', {'second_pass': False}), interpolation='classical', presmoother=('gauss_seidel', {'sweep': 'symmetric'}), postsmoother=('gauss_seidel', {'sweep': 'symmetric'}), max_levels=30, max_coarse=10, keep=False, **kwargs)[source]#
Create a multilevel solver using Classical AMG (Ruge-Stuben AMG).
- Parameters:
- Acsr_array
Square matrix in CSR format
- strengthstr
Valid strings are [‘symmetric’, ‘classical’, ‘evolution’, ‘distance’, ‘algebraic_distance’,’affinity’, ‘energy_based’, None]. Method used to determine the strength of connection between unknowns of the linear system. Method-specific parameters may be passed in using a tuple, e.g. strength=(‘symmetric’,{‘theta’: 0.25 }). If strength=None, all nonzero entries of the matrix are considered strong.
- CFstr or tuple, default ‘RS’
Method used for coarse grid selection (C/F splitting) Supported methods are RS, PMIS, PMISc, CLJP, CLJPc, and CR.
- interpolationstr, default ‘classical’
Method for interpolation. Options include ‘direct’, ‘classical’.
- presmootherstr or dict
Method used for presmoothing at each level. Method-specific parameters may be passed in using a tuple, e.g. presmoother=(‘gauss_seidel’,{‘sweep’:’symmetric}), the default.
- postsmootherstr or dict
Postsmoothing method with the same usage as presmoother
- max_levelsint, default 30
Maximum number of levels to be used in the multilevel solver.
- max_coarseint, default 20
Maximum number of variables permitted on the coarse grid.
- keepbool, default False
Flag to indicate keeping strength of connection (C) in the hierarchy for diagnostics.
- kwargsdict
Extra keywords passed to Multilevel class
- Returns:
- mlMultilevelSolver
Multigrid hierarchy of matrices and prolongation operators
See also
aggregation.smoothed_aggregation_solver
,MultilevelSolver
aggregation.rootnode_solver
Notes
“coarse_solver” is an optional argument and is the solver used at the coarsest grid. The default is a pseudo-inverse. Most simply, coarse_solver can be one of [‘splu’, ‘lu’, ‘cholesky, ‘pinv’, ‘gauss_seidel’, … ]. Additionally, coarse_solver may be a tuple (fn, args), where fn is a string such as [‘splu’, ‘lu’, …] or a callable function, and args is a dictionary of arguments to be passed to fn. See [2001TrOoSc] for additional details.
References
[2001TrOoSc]Trottenberg, U., Oosterlee, C. W., and Schuller, A., “Multigrid” San Diego: Academic Press, 2001. Appendix A
Examples
>>> from pyamg.gallery import poisson >>> from pyamg import ruge_stuben_solver >>> A = poisson((10,),format='csr') >>> ml = ruge_stuben_solver(A,max_coarse=3)