This page contains the Classical Package documentation.
Classical AMG (Ruge-Stuben AMG)
Create a multilevel solver using Classical AMG (Ruge-Stuben AMG)
Parameters : | A : csr_matrix
strength : [‘symmetric’, ‘classical’, ‘evolution’, None]
CF : {string}
presmoother : {string or dict}
postsmoother : {string or dict}
max_levels: {integer} : default 10
max_coarse: {integer} : default 500
keep: {bool} : default False
|
---|---|
Returns : | ml : multilevel_solver
|
See also
aggregation.smoothed_aggregation_solver, multilevel_solver, 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.
References
[R10] | 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)
Compatible Relaxation
Use Compatible Relaxation to compute a C/F splitting
Parameters : | S : csr_matrix
method : {‘habituated’,’concurrent’}
maxiter : int
|
---|---|
Returns : | splitting : array
|
References
[R11] | Livne, O.E., “Coarsening by compatible relaxation.” Numer. Linear Algebra Appl. 11, No. 2-3, 205-227 (2004). |
Examples
>>> from pyamg.gallery import poisson
>>> from pyamg.classical.cr import CR
>>> A = poisson((20,20),format='csr')
>>> splitting = CR(A)
Binormalize matrix A. Attempt to create unit l_1 norm rows.
Parameters : | A : csr_matrix
tol : float
x : array
maxiter : int
|
---|---|
Returns : | C : csr_matrix
|
Notes
References
[R12] | 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 |
Examples
>>> from pyamg.gallery import poisson
>>> from pyamg.classical import binormalize
>>> A = poisson((10,),format='csr')
>>> C = binormalize(A)
Classical AMG Interpolation methods
Create prolongator using direct interpolation
Parameters : | A : {csr_matrix}
C : {csr_matrix}
splitting : array
|
---|---|
Returns : | P : {csr_matrix}
|
Examples
>>> 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. ]]
Functions to compute C/F splittings for use in Classical AMG
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.
C/F splitting is represented by an array with ones for all the C-nodes and zeros for the F-nodes.
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
[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. |
Compute a C/F splitting using Ruge-Stuben coarsening
Parameters : | S : csr_matrix
|
---|---|
Returns : | splitting : ndarray
|
See also
amg_core.rs_cf_splitting
References
[R13] | 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. |
Examples
>>> from pyamg.gallery import poisson
>>> from pyamg.classical import RS
>>> S = poisson((7,), format='csr') # 1D mesh with 7 vertices
>>> splitting = RS(S)
C/F splitting using the Parallel Modified Independent Set method
Parameters : | S : csr_matrix
|
---|---|
Returns : | splitting : ndarray
|
See also
References
[R14] | 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. |
Examples
>>> from pyamg.gallery import poisson
>>> from pyamg.classical import PMIS
>>> S = poisson((7,), format='csr') # 1D mesh with 7 vertices
>>> splitting = PMIS(S)
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.
Parameters : | S : csr_matrix
method : string
|
---|---|
Returns : | splitting : array
|
See also
References
[R15] | David M. Alber and Luke N. Olson “Parallel coarse-grid selection” Numerical Linear Algebra with Applications 2007; 14:611-643. |
Examples
>>> from pyamg.gallery import poisson
>>> from pyamg.classical import PMISc
>>> S = poisson((7,), format='csr') # 1D mesh with 7 vertices
>>> splitting = PMISc(S)
Compute a maximal independent set of a graph in parallel
Parameters : | G : csr_matrix
weights : ndarray
maxiter : int
|
---|---|
Returns : | mis : array
|
See also
fn
Examples
>>> 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)