classical Package¶
classical Package¶
Classical AMG
- pyamg.classical.CR(S, method='habituated', maxiter=20)¶
Use Compatible Relaxation to compute a C/F splitting
- S : csr_matrix
- sparse matrix (n x n) usually matrix A of Ax=b
- method : {‘habituated’,’concurrent’}
- Method used during relaxation:
- concurrent: GS relaxation on F-points, leaving e_c = 0
- habituated: full relaxation, setting e_c = 0
- maxiter : int
- maximum number of outer iterations (lambda)
- splitting : array
- C/F list of 1’s (coarse pt) and 0’s (fine pt) (n x 1)
[1] Livne, O.E., “Coarsening by compatible relaxation.” Numer. Linear Algebra Appl. 11, No. 2-3, 205-227 (2004). >>> from pyamg.gallery import poisson >>> from pyamg.classical.cr import CR >>> A = poisson((20,20),format='csr') >>> splitting = CR(A)
- pyamg.classical.MIS(G, weights, maxiter=None)¶
Compute a maximal independent set of a graph in parallel
- G : csr_matrix
- Matrix graph, G[i,j] != 0 indicates an edge
- weights : ndarray
- Array of weights for each vertex in the graph G
- maxiter : int
- Maximum number of iterations (default: None)
- mis : array
- Array of length of G of zeros/ones indicating the independent set
>>> from pyamg.gallery import poisson >>> from pyamg.classical import MIS >>> import numpy >>> G = poisson((7,), format='csr') # 1D mesh with 7 vertices >>> w = numpy.ones((G.shape[0],1)).ravel() >>> mis = MIS(G,w)
fn = amg_core.maximal_independent_set_parallel
- pyamg.classical.PMIS(S)¶
C/F splitting using the Parallel Modified Independent Set method
- S : csr_matrix
- Strength of connection matrix indicating the strength between nodes i and j (S_ij)
- splitting : ndarray
- Array of length of S of ones (coarse) and zeros (fine)
>>> from pyamg.gallery import poisson >>> from pyamg.classical import PMIS >>> S = poisson((7,), format='csr') # 1D mesh with 7 vertices >>> splitting = PMIS(S)
MIS
[1] Hans De Sterck, Ulrike M Yang, and Jeffrey J Heys “Reducing complexity in parallel algebraic multigrid preconditioners” SIAM Journal on Matrix Analysis and Applications 2006; 27:1019-1039.
- pyamg.classical.PMISc(S, method='JP')¶
C/F splitting using Parallel Modified Independent Set (in color)
PMIS-c, or PMIS in color, improves PMIS by perturbing the initial random weights with weights determined by a vertex coloring.
- S : csr_matrix
- Strength of connection matrix indicating the strength between nodes i and j (S_ij)
- method : string
- Algorithm used to compute the initial vertex coloring:
- ‘MIS’ - Maximal Independent Set
- ‘JP’ - Jones-Plassmann (parallel)
- ‘LDF’ - Largest-Degree-First (parallel)
- splitting : array
- Array of length of S of ones (coarse) and zeros (fine)
>>> from pyamg.gallery import poisson >>> from pyamg.classical import PMISc >>> S = poisson((7,), format='csr') # 1D mesh with 7 vertices >>> splitting = PMISc(S)
MIS
[1] David M. Alber and Luke N. Olson “Parallel coarse-grid selection” Numerical Linear Algebra with Applications 2007; 14:611-643.
- pyamg.classical.RS(S)¶
Compute a C/F splitting using Ruge-Stuben coarsening
- S : csr_matrix
- Strength of connection matrix indicating the strength between nodes i and j (S_ij)
- splitting : ndarray
- Array of length of S of ones (coarse) and zeros (fine)
>>> from pyamg.gallery import poisson >>> from pyamg.classical import RS >>> S = poisson((7,), format='csr') # 1D mesh with 7 vertices >>> splitting = RS(S)
amg_core.rs_cf_splitting
[1] Ruge JW, Stuben K. “Algebraic multigrid (AMG)” In Multigrid Methods, McCormick SF (ed.), Frontiers in Applied Mathematics, vol. 3. SIAM: Philadelphia, PA, 1987; 73-130.
- pyamg.classical.binormalize(A, tol=1e-05, maxiter=10)¶
Binormalize matrix A. Attempt to create unit l_1 norm rows.
- A : csr_matrix
- sparse matrix (n x n)
- tol : float
- tolerance
- x : array
- guess at the diagonal
- maxiter : int
- maximum number of iterations to try
- C : csr_matrix
- diagonally scaled A, C=DAD
- Goal: Scale A so that l_1 norm of the rows are equal to 1:
- B = DAD
- want row sum of B = 1
- easily done with tol=0 if B=DA, but this is not symmetric
- algorithm is O(N log (1.0/tol))
>>> from pyamg.gallery import poisson >>> from pyamg.classical import binormalize >>> A = poisson((10,),format='csr') >>> C = binormalize(A)
[1] Livne, Golub, “Scaling by Binormalization” Tech Report SCCM-03-12, SCCM, Stanford, 2003 http://citeseerx.ist.psu.edu/viewdoc/summary?doi=10.1.1.3.1679
- pyamg.classical.direct_interpolation(A, C, splitting)¶
Create prolongator using direct interpolation
- A : {csr_matrix}
- NxN matrix in CSR format
- C : {csr_matrix}
- Strength-of-Connection matrix Must have zero diagonal
- splitting : array
- C/F splitting stored in an array of length N
- P : {csr_matrix}
- Prolongator using direct interpolation
>>> from pyamg.gallery import poisson >>> from pyamg.classical import direct_interpolation >>> import numpy >>> A = poisson((5,),format='csr') >>> splitting = numpy.array([1,0,1,0,1], dtype='intc') >>> P = direct_interpolation(A, A, splitting) >>> print P.todense() [[ 1. 0. 0. ] [ 0.5 0.5 0. ] [ 0. 1. 0. ] [ 0. 0.5 0.5] [ 0. 0. 1. ]]
- pyamg.classical.ruge_stuben_solver(A, strength=('classical', {'theta': 0.25}), CF='RS', presmoother=('gauss_seidel', {'sweep': 'symmetric'}), postsmoother=('gauss_seidel', {'sweep': 'symmetric'}), max_levels=10, max_coarse=500, keep=False, **kwargs)¶
Create a multilevel solver using Classical AMG (Ruge-Stuben AMG)
- A : csr_matrix
- Square matrix in CSR format
- strength : [‘symmetric’, ‘classical’, ‘evolution’, 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’
- Method used for coarse grid selection (C/F splitting) Supported methods are RS, PMIS, PMISc, CLJP, and CLJPc
- presmoother : {string 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.
- postsmoother : {string or dict}
- Postsmoothing method with the same usage as presmoother
- max_levels: {integer} : default 10
- Maximum number of levels to be used in the multilevel solver.
- max_coarse: {integer} : default 500
- Maximum number of variables permitted on the coarse grid.
- keep: {bool} : default False
- Flag to indicate keeping extra operators in the hierarchy for diagnostics. For example, if True, then strength of connection (C) and tentative prolongation (T) are kept.
- ml : multilevel_solver
- Multigrid hierarchy of matrices and prolongation operators
>>> from pyamg.gallery import poisson >>> from pyamg import ruge_stuben_solver >>> A = poisson((10,),format='csr') >>> ml = ruge_stuben_solver(A,max_coarse=3)
“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.
[1] Trottenberg, U., Oosterlee, C. W., and Schuller, A., “Multigrid” San Diego: Academic Press, 2001. Appendix A aggregation.smoothed_aggregation_solver, multilevel_solver, aggregation.rootnode_solver
classical Module¶
Classical AMG (Ruge-Stuben AMG)
- pyamg.classical.classical.ruge_stuben_solver(A, strength=('classical', {'theta': 0.25}), CF='RS', presmoother=('gauss_seidel', {'sweep': 'symmetric'}), postsmoother=('gauss_seidel', {'sweep': 'symmetric'}), max_levels=10, max_coarse=500, keep=False, **kwargs)[source]¶
Create a multilevel solver using Classical AMG (Ruge-Stuben AMG)
- A : csr_matrix
- Square matrix in CSR format
- strength : [‘symmetric’, ‘classical’, ‘evolution’, 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’
- Method used for coarse grid selection (C/F splitting) Supported methods are RS, PMIS, PMISc, CLJP, and CLJPc
- presmoother : {string 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.
- postsmoother : {string or dict}
- Postsmoothing method with the same usage as presmoother
- max_levels: {integer} : default 10
- Maximum number of levels to be used in the multilevel solver.
- max_coarse: {integer} : default 500
- Maximum number of variables permitted on the coarse grid.
- keep: {bool} : default False
- Flag to indicate keeping extra operators in the hierarchy for diagnostics. For example, if True, then strength of connection (C) and tentative prolongation (T) are kept.
- ml : multilevel_solver
- Multigrid hierarchy of matrices and prolongation operators
>>> from pyamg.gallery import poisson >>> from pyamg import ruge_stuben_solver >>> A = poisson((10,),format='csr') >>> ml = ruge_stuben_solver(A,max_coarse=3)
“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.
[1] Trottenberg, U., Oosterlee, C. W., and Schuller, A., “Multigrid” San Diego: Academic Press, 2001. Appendix A aggregation.smoothed_aggregation_solver, multilevel_solver, aggregation.rootnode_solver
cr Module¶
Compatible Relaxation
- pyamg.classical.cr.CR(S, method='habituated', maxiter=20)[source]¶
Use Compatible Relaxation to compute a C/F splitting
- S : csr_matrix
- sparse matrix (n x n) usually matrix A of Ax=b
- method : {‘habituated’,’concurrent’}
- Method used during relaxation:
- concurrent: GS relaxation on F-points, leaving e_c = 0
- habituated: full relaxation, setting e_c = 0
- maxiter : int
- maximum number of outer iterations (lambda)
- splitting : array
- C/F list of 1’s (coarse pt) and 0’s (fine pt) (n x 1)
[1] Livne, O.E., “Coarsening by compatible relaxation.” Numer. Linear Algebra Appl. 11, No. 2-3, 205-227 (2004). >>> from pyamg.gallery import poisson >>> from pyamg.classical.cr import CR >>> A = poisson((20,20),format='csr') >>> splitting = CR(A)
- pyamg.classical.cr.binormalize(A, tol=1e-05, maxiter=10)[source]¶
Binormalize matrix A. Attempt to create unit l_1 norm rows.
- A : csr_matrix
- sparse matrix (n x n)
- tol : float
- tolerance
- x : array
- guess at the diagonal
- maxiter : int
- maximum number of iterations to try
- C : csr_matrix
- diagonally scaled A, C=DAD
- Goal: Scale A so that l_1 norm of the rows are equal to 1:
- B = DAD
- want row sum of B = 1
- easily done with tol=0 if B=DA, but this is not symmetric
- algorithm is O(N log (1.0/tol))
>>> from pyamg.gallery import poisson >>> from pyamg.classical import binormalize >>> A = poisson((10,),format='csr') >>> C = binormalize(A)
[1] Livne, Golub, “Scaling by Binormalization” Tech Report SCCM-03-12, SCCM, Stanford, 2003 http://citeseerx.ist.psu.edu/viewdoc/summary?doi=10.1.1.3.1679
interpolate Module¶
Classical AMG Interpolation methods
- pyamg.classical.interpolate.direct_interpolation(A, C, splitting)[source]¶
Create prolongator using direct interpolation
- A : {csr_matrix}
- NxN matrix in CSR format
- C : {csr_matrix}
- Strength-of-Connection matrix Must have zero diagonal
- splitting : array
- C/F splitting stored in an array of length N
- P : {csr_matrix}
- Prolongator using direct interpolation
>>> from pyamg.gallery import poisson >>> from pyamg.classical import direct_interpolation >>> import numpy >>> A = poisson((5,),format='csr') >>> splitting = numpy.array([1,0,1,0,1], dtype='intc') >>> P = direct_interpolation(A, A, splitting) >>> print P.todense() [[ 1. 0. 0. ] [ 0.5 0.5 0. ] [ 0. 1. 0. ] [ 0. 0.5 0.5] [ 0. 0. 1. ]]
split Module¶
Functions to compute C/F splittings for use in Classical AMG
Overview¶
A C/F splitting is a partitioning of the nodes in the graph of as connection matrix (denoted S for strength) into sets of C (coarse) and F (fine) nodes. The C-nodes are promoted to the coarser grid while the F-nodes are retained on the finer grid. Ideally, the C-nodes, which represent the coarse-level unknowns, should be far fewer in number than the F-nodes. Furthermore, algebraically smooth error must be well-approximated by the coarse level degrees of freedom.
Representation¶
C/F splitting is represented by an array with ones for all the C-nodes and zeros for the F-nodes.
C/F Splitting Methods¶
- RS : Original Ruge-Stuben method
- Produces good C/F splittings but is inherently serial.
- May produce AMG hierarchies with relatively high operator complexities.
- See References [1,4]
- PMIS: Parallel Modified Independent Set
- Very fast construction with low operator complexity.
- Convergence can deteriorate with increasing problem size on structured meshes.
- Uses method similar to Luby’s Maximal Independent Set algorithm.
- See References [1,3]
- PMISc: Parallel Modified Independent Set in Color
- Fast construction with low operator complexity.
- Better scalability than PMIS on structured meshes.
- Augments random weights with a (graph) vertex coloring
- See References [1]
- CLJP: Clearly-Luby-Jones-Plassmann
- Parallel method with cost and complexity comparable to Ruge-Stuben.
- Convergence can deteriorate with increasing problem size on structured meshes.
- See References [1,2]
- CLJP-c: Clearly-Luby-Jones-Plassmann in Color
- Parallel method with cost and complexity comparable to Ruge-Stuben.
- Better scalability than CLJP on structured meshes.
- See References [1]
Summary¶
In general, methods that use a graph coloring perform better on structured meshes [1]. Unstructured meshes do not appear to benefit substantially from coloring.
method parallel in color cost RS no no moderate PMIS yes no very low PMISc yes yes low CLJP yes no moderate CLJPc yes yes moderate
References¶
[1] | David M. Alber and Luke N. Olson “Parallel coarse-grid selection” Numerical Linear Algebra with Applications 2007; 14:611-643. |
[2] | Cleary AJ, Falgout RD, Henson VE, Jones JE. “Coarse-grid selection for parallel algebraic multigrid” Proceedings of the 5th International Symposium on Solving Irregularly Structured Problems in Parallel. Springer: Berlin, 1998; 104-115. |
[3] | Hans De Sterck, Ulrike M Yang, and Jeffrey J Heys “Reducing complexity in parallel algebraic multigrid preconditioners” SIAM Journal on Matrix Analysis and Applications 2006; 27:1019-1039. |
[4] | Ruge JW, Stuben K. “Algebraic multigrid (AMG)” In Multigrid Methods, McCormick SF (ed.), Frontiers in Applied Mathematics, vol. 3. SIAM: Philadelphia, PA, 1987; 73-130. |
- pyamg.classical.split.RS(S)[source]¶
Compute a C/F splitting using Ruge-Stuben coarsening
- S : csr_matrix
- Strength of connection matrix indicating the strength between nodes i and j (S_ij)
- splitting : ndarray
- Array of length of S of ones (coarse) and zeros (fine)
>>> from pyamg.gallery import poisson >>> from pyamg.classical import RS >>> S = poisson((7,), format='csr') # 1D mesh with 7 vertices >>> splitting = RS(S)
amg_core.rs_cf_splitting
[1] Ruge JW, Stuben K. “Algebraic multigrid (AMG)” In Multigrid Methods, McCormick SF (ed.), Frontiers in Applied Mathematics, vol. 3. SIAM: Philadelphia, PA, 1987; 73-130.
- pyamg.classical.split.PMIS(S)[source]¶
C/F splitting using the Parallel Modified Independent Set method
- S : csr_matrix
- Strength of connection matrix indicating the strength between nodes i and j (S_ij)
- splitting : ndarray
- Array of length of S of ones (coarse) and zeros (fine)
>>> from pyamg.gallery import poisson >>> from pyamg.classical import PMIS >>> S = poisson((7,), format='csr') # 1D mesh with 7 vertices >>> splitting = PMIS(S)
MIS
[1] Hans De Sterck, Ulrike M Yang, and Jeffrey J Heys “Reducing complexity in parallel algebraic multigrid preconditioners” SIAM Journal on Matrix Analysis and Applications 2006; 27:1019-1039.
- pyamg.classical.split.PMISc(S, method='JP')[source]¶
C/F splitting using Parallel Modified Independent Set (in color)
PMIS-c, or PMIS in color, improves PMIS by perturbing the initial random weights with weights determined by a vertex coloring.
- S : csr_matrix
- Strength of connection matrix indicating the strength between nodes i and j (S_ij)
- method : string
- Algorithm used to compute the initial vertex coloring:
- ‘MIS’ - Maximal Independent Set
- ‘JP’ - Jones-Plassmann (parallel)
- ‘LDF’ - Largest-Degree-First (parallel)
- splitting : array
- Array of length of S of ones (coarse) and zeros (fine)
>>> from pyamg.gallery import poisson >>> from pyamg.classical import PMISc >>> S = poisson((7,), format='csr') # 1D mesh with 7 vertices >>> splitting = PMISc(S)
MIS
[1] David M. Alber and Luke N. Olson “Parallel coarse-grid selection” Numerical Linear Algebra with Applications 2007; 14:611-643.
- pyamg.classical.split.MIS(G, weights, maxiter=None)[source]¶
Compute a maximal independent set of a graph in parallel
- G : csr_matrix
- Matrix graph, G[i,j] != 0 indicates an edge
- weights : ndarray
- Array of weights for each vertex in the graph G
- maxiter : int
- Maximum number of iterations (default: None)
- mis : array
- Array of length of G of zeros/ones indicating the independent set
>>> from pyamg.gallery import poisson >>> from pyamg.classical import MIS >>> import numpy >>> G = poisson((7,), format='csr') # 1D mesh with 7 vertices >>> w = numpy.ones((G.shape[0],1)).ravel() >>> mis = MIS(G,w)
fn = amg_core.maximal_independent_set_parallel