Relaxation Documentation

This page contains the Relaxation Package documentation.

The chebyshev Module

Compute coefficients for polynomial smoothers

pyamg.relaxation.chebyshev.chebyshev_polynomial_coefficients(a, b, degree)[source]

Chebyshev polynomial coefficients for the interval [a,b]

Parameters :

a,b : float

The left and right endpoints of the interval.

degree : int

Degree of desired Chebyshev polynomial

Returns :

Coefficients of the Chebyshev polynomial C(t) with minimum :

magnitude on the interval [a,b] such that C(0) = 1.0. :

The coefficients are returned in descending order. :

Notes

a,b typically represent the interval of the spectrum for some matrix that you wish to damp with a Chebyshev smoother.

Examples

>>> from pyamg.relaxation.chebyshev import chebyshev_polynomial_coefficients
>>> print chebyshev_polynomial_coefficients(1.0,2.0, 3)
[-0.32323232  1.45454545 -2.12121212  1.        ]

The relaxation Module

Relaxation methods for linear systems

pyamg.relaxation.relaxation.sor(A, x, b, omega, iterations=1, sweep='forward')[source]

Perform SOR iteration on the linear system Ax=b

Parameters :

A : {csr_matrix, bsr_matrix}

Sparse NxN matrix

x : ndarray

Approximate solution (length N)

b : ndarray

Right-hand side (length N)

omega : scalar

Damping parameter

iterations : int

Number of iterations to perform

sweep : {‘forward’,’backward’,’symmetric’}

Direction of sweep

Returns :

Nothing, x will be modified in place. :

Notes

When omega=1.0, SOR is equivalent to Gauss-Seidel.

Examples

>>> ## Use SOR as stand-along solver
>>> from pyamg.relaxation import sor
>>> from pyamg.gallery import poisson
>>> from pyamg.util.linalg import norm
>>> import numpy
>>> A = poisson((10,10), format='csr')
>>> x0 = numpy.zeros((A.shape[0],1))
>>> b = numpy.ones((A.shape[0],1))
>>> sor(A, x0, b, 1.33, iterations=10)
>>> print norm(b-A*x0)
3.03888724811
>>> #
>>> ## Use SOR as the multigrid smoother 
>>> from pyamg import smoothed_aggregation_solver
>>> sa = smoothed_aggregation_solver(A, B=numpy.ones((A.shape[0],1)),
...         coarse_solver='pinv2', max_coarse=50,
...         presmoother=('sor', {'sweep':'symmetric', 'omega' : 1.33}), 
...         postsmoother=('sor', {'sweep':'symmetric', 'omega' : 1.33}))
>>> x0 = numpy.zeros((A.shape[0],1))
>>> residuals=[]
>>> x = sa.solve(b, x0=x0, tol=1e-8, residuals=residuals)
pyamg.relaxation.relaxation.gauss_seidel(A, x, b, iterations=1, sweep='forward')[source]

Perform Gauss-Seidel iteration on the linear system Ax=b

Parameters :

A : {csr_matrix, bsr_matrix}

Sparse NxN matrix

x : ndarray

Approximate solution (length N)

b : ndarray

Right-hand side (length N)

iterations : int

Number of iterations to perform

sweep : {‘forward’,’backward’,’symmetric’}

Direction of sweep

Returns :

Nothing, x will be modified in place. :

Examples

>>> ## Use Gauss-Seidel as a Stand-Alone Solver
>>> from pyamg.relaxation import *
>>> from pyamg.gallery import poisson
>>> from pyamg.util.linalg import norm
>>> import numpy
>>> A = poisson((10,10), format='csr')
>>> x0 = numpy.zeros((A.shape[0],1))
>>> b = numpy.ones((A.shape[0],1))
>>> gauss_seidel(A, x0, b, iterations=10)
>>> print norm(b-A*x0)
4.00733716236
>>> #
>>> ## Use Gauss-Seidel as the Multigrid Smoother
>>> from pyamg import smoothed_aggregation_solver
>>> sa = smoothed_aggregation_solver(A, B=numpy.ones((A.shape[0],1)),
...         coarse_solver='pinv2', max_coarse=50,
...         presmoother=('gauss_seidel', {'sweep':'symmetric'}), 
...         postsmoother=('gauss_seidel', {'sweep':'symmetric'}))
>>> x0=numpy.zeros((A.shape[0],1))
>>> residuals=[]
>>> x = sa.solve(b, x0=x0, tol=1e-8, residuals=residuals)
pyamg.relaxation.relaxation.jacobi(A, x, b, iterations=1, omega=1.0)[source]

Perform Jacobi iteration on the linear system Ax=b

Parameters :

A : csr_matrix

Sparse NxN matrix

x : ndarray

Approximate solution (length N)

b : ndarray

Right-hand side (length N)

iterations : int

Number of iterations to perform

omega : scalar

Damping parameter

Returns :

Nothing, x will be modified in place. :

Examples

>>> ## Use Jacobi as a Stand-Alone Solver
>>> from pyamg.relaxation import *
>>> from pyamg.gallery import poisson
>>> from pyamg.util.linalg import norm
>>> import numpy
>>> A = poisson((10,10), format='csr')
>>> x0 = numpy.zeros((A.shape[0],1))
>>> b = numpy.ones((A.shape[0],1))
>>> jacobi(A, x0, b, iterations=10, omega=1.0)
>>> print norm(b-A*x0)
5.83475132751
>>> #
>>> ## Use Jacobi as the Multigrid Smoother
>>> from pyamg import smoothed_aggregation_solver
>>> sa = smoothed_aggregation_solver(A, B=numpy.ones((A.shape[0],1)),
...         coarse_solver='pinv2', max_coarse=50,
...         presmoother=('jacobi', {'omega': 4.0/3.0, 'iterations' : 2}), 
...         postsmoother=('jacobi', {'omega': 4.0/3.0, 'iterations' : 2}))
>>> x0=numpy.zeros((A.shape[0],1))
>>> residuals=[]
>>> x = sa.solve(b, x0=x0, tol=1e-8, residuals=residuals)
pyamg.relaxation.relaxation.polynomial(A, x, b, coefficients, iterations=1)[source]

Apply a polynomial smoother to the system Ax=b

Parameters :

A : sparse matrix

Sparse NxN matrix

x : ndarray

Approximate solution (length N)

b : ndarray

Right-hand side (length N)

coefficients : {array_like}

Coefficients of the polynomial. See Notes section for details.

iterations : int

Number of iterations to perform

Returns :

Nothing, x will be modified in place. :

Notes

The smoother has the form x[:] = x + p(A) (b - A*x) where p(A) is a polynomial in A whose scalar coefficients are specified (in descending order) by argument ‘coefficients’.

  • Richardson iteration p(A) = c_0:

    polynomial_smoother(A, x, b, [c_0])

  • Linear smoother p(A) = c_1*A + c_0:

    polynomial_smoother(A, x, b, [c_1, c_0])

  • Quadratic smoother p(A) = c_2*A^2 + c_1*A + c_0:

    polynomial_smoother(A, x, b, [c_2, c_1, c_0])

Here, Horner’s Rule is applied to avoid computing A^k directly.

For efficience, the method detects the case x = 0 one matrix-vector product is avoided (since (b - A*x) is b).

Examples

>>> ## The polynomial smoother is not currently used directly 
>>> ## in PyAMG.  It is only used by the chebyshev smoothing option,
>>> ## which automatically calculates the correct coefficients.
>>> from pyamg.gallery import poisson
>>> from pyamg.util.linalg import norm
>>> import numpy
>>> from pyamg.aggregation import smoothed_aggregation_solver
>>> A = poisson((10,10), format='csr')
>>> b = numpy.ones((A.shape[0],1))
>>> sa = smoothed_aggregation_solver(A, B=numpy.ones((A.shape[0],1)),
...         coarse_solver='pinv2', max_coarse=50,
...         presmoother=('chebyshev', {'degree':3, 'iterations':1}), 
...         postsmoother=('chebyshev', {'degree':3, 'iterations':1}))
>>> x0=numpy.zeros((A.shape[0],1))
>>> residuals=[]
>>> x = sa.solve(b, x0=x0, tol=1e-8, residuals=residuals)
pyamg.relaxation.relaxation.schwarz(A, x, b, iterations=1, subdomain=None, subdomain_ptr=None, inv_subblock=None, inv_subblock_ptr=None, sweep='forward')[source]

Perform Overlapping multiplicative Schwarz on the linear system Ax=b

Parameters :

A : {csr_matrix, bsr_matrix}

Sparse NxN matrix

x : ndarray

Approximate solution (length N)

b : ndarray

Right-hand side (length N)

iterations : int

Number of iterations to perform

subdomain : {int array}

Linear array containing each subdomain’s elements

subdomain_ptr : {int array}

Pointer in subdomain, such that subdomain[subdomain_ptr[i]:subdomain_ptr[i+1]]] contains the _sorted_ indices in subdomain i

inv_subblock : {int_array}

Linear array containing each subdomain’s inverted diagonal block of A

inv_subblock_ptr : {int array}

Pointer in inv_subblock, such that inv_subblock[inv_subblock_ptr[i]:inv_subblock_ptr[i+1]]] contains the inverted diagonal block of A for the i-th subdomain in _row_ major order

sweep : {‘forward’,’backward’,’symmetric’}

Direction of sweep

Returns :

Nothing, x will be modified in place. :

Notes

If subdomains is None, then a point-wise iteration takes place, with the overlapping region defined by each degree-of-freedom’s neighbors in the matrix graph.

If subdomains is not None, but subblocks is, then the subblocks are formed internally.

Currently only supports CSR matrices

Examples

>>> ## Use Overlapping Schwarz as a Stand-Alone Solver
>>> from pyamg.relaxation import *
>>> from pyamg.gallery import poisson
>>> from pyamg.util.linalg import norm
>>> import numpy
>>> A = poisson((10,10), format='csr')
>>> x0 = numpy.zeros((A.shape[0],1))
>>> b = numpy.ones((A.shape[0],1))
>>> schwarz(A, x0, b, iterations=10)
>>> print norm(b-A*x0)
0.126326160522
>>> #
>>> ## Schwarz as the Multigrid Smoother
>>> from pyamg import smoothed_aggregation_solver
>>> sa = smoothed_aggregation_solver(A, B=numpy.ones((A.shape[0],1)),
...         coarse_solver='pinv2', max_coarse=50,
...         presmoother='schwarz', 
...         postsmoother='schwarz')
>>> x0=numpy.zeros((A.shape[0],1))
>>> residuals=[]
>>> x = sa.solve(b, x0=x0, tol=1e-8, residuals=residuals)
pyamg.relaxation.relaxation.jacobi_ne(A, x, b, iterations=1, omega=1.0)[source]

Perform Jacobi iterations on the linear system A A.H x = A.H b (Also known as Cimmino relaxation)

Parameters :

A : csr_matrix

Sparse NxN matrix

x : ndarray

Approximate solution (length N)

b : ndarray

Right-hand side (length N)

iterations : int

Number of iterations to perform

omega : scalar

Damping parameter

Returns :

Nothing, x will be modified in place. :

References

[R28]Brandt, Ta’asan. “Multigrid Method For Nearly Singular And Slightly Indefinite Problems.” 1985. NASA Technical Report Numbers: ICASE-85-57; NAS 1.26:178026; NASA-CR-178026;
[R29]Kaczmarz. Angenaeherte Aufloesung von Systemen Linearer Gleichungen. Bull. Acad. Polon. Sci. Lett. A 35, 355-57. 1937
[R30]Cimmino. La ricerca scientifica ser. II 1. Pubbliz. dell’Inst. pre le Appl. del Calculo 34, 326-333, 1938.

Examples

>>> ## Use NE Jacobi as a Stand-Alone Solver
>>> from pyamg.relaxation import *
>>> from pyamg.gallery import poisson
>>> from pyamg.util.linalg import norm
>>> import numpy
>>> A = poisson((50,50), format='csr')
>>> x0 = numpy.zeros((A.shape[0],1))
>>> b = numpy.ones((A.shape[0],1))
>>> jacobi_ne(A, x0, b, iterations=10, omega=2.0/3.0)
>>> print norm(b-A*x0)
49.3886046066
>>> #
>>> ## Use NE Jacobi as the Multigrid Smoother
>>> from pyamg import smoothed_aggregation_solver
>>> sa = smoothed_aggregation_solver(A, B=numpy.ones((A.shape[0],1)),
...         coarse_solver='pinv2', max_coarse=50,
...         presmoother=('jacobi_ne', {'iterations' : 2, 'omega' : 4.0/3.0}), 
...         postsmoother=('jacobi_ne', {'iterations' : 2, 'omega' : 4.0/3.0}))
>>> x0=numpy.zeros((A.shape[0],1))
>>> residuals=[]
>>> x = sa.solve(b, x0=x0, tol=1e-8, residuals=residuals)
pyamg.relaxation.relaxation.gauss_seidel_ne(A, x, b, iterations=1, sweep='forward', omega=1.0, Dinv=None)[source]

Perform Gauss-Seidel iterations on the linear system A A.H x = b (Also known as Kaczmarz relaxation)

Parameters :

A : csr_matrix

Sparse NxN matrix

x : { ndarray }

Approximate solution (length N)

b : { ndarray }

Right-hand side (length N)

iterations : { int }

Number of iterations to perform

sweep : {‘forward’,’backward’,’symmetric’}

Direction of sweep

omega : { float}

Relaxation parameter typically in (0, 2) if omega != 1.0, then algorithm becomes SOR on A A.H

Dinv : { ndarray}

Inverse of diag(A A.H), (length N)

Returns :

Nothing, x will be modified in place. :

References

[R31]Brandt, Ta’asan. “Multigrid Method For Nearly Singular And Slightly Indefinite Problems.” 1985. NASA Technical Report Numbers: ICASE-85-57; NAS 1.26:178026; NASA-CR-178026;
[R32]Kaczmarz. Angenaeherte Aufloesung von Systemen Linearer Gleichungen. Bull. Acad. Polon. Sci. Lett. A 35, 355-57. 1937

Examples

>>> ## Use NE Gauss-Seidel as a Stand-Alone Solver
>>> from pyamg.relaxation import *
>>> from pyamg.gallery import poisson
>>> from pyamg.util.linalg import norm
>>> import numpy
>>> A = poisson((10,10), format='csr')
>>> x0 = numpy.zeros((A.shape[0],1))
>>> b = numpy.ones((A.shape[0],1))
>>> gauss_seidel_ne(A, x0, b, iterations=10, sweep='symmetric')
>>> print norm(b-A*x0)
8.47576806771
>>> #
>>> ## Use NE Gauss-Seidel as the Multigrid Smoother
>>> from pyamg import smoothed_aggregation_solver
>>> sa = smoothed_aggregation_solver(A, B=numpy.ones((A.shape[0],1)),
...         coarse_solver='pinv2', max_coarse=50,
...         presmoother=('gauss_seidel_ne', {'sweep' : 'symmetric'}), 
...         postsmoother=('gauss_seidel_ne', {'sweep' : 'symmetric'}))
>>> x0=numpy.zeros((A.shape[0],1))
>>> residuals=[]
>>> x = sa.solve(b, x0=x0, tol=1e-8, residuals=residuals)
pyamg.relaxation.relaxation.gauss_seidel_nr(A, x, b, iterations=1, sweep='forward', omega=1.0, Dinv=None)[source]

Perform Gauss-Seidel iterations on the linear system A.H A x = A.H b

Parameters :

A : csr_matrix

Sparse NxN matrix

x : { ndarray }

Approximate solution (length N)

b : { ndarray }

Right-hand side (length N)

iterations : { int }

Number of iterations to perform

sweep : {‘forward’,’backward’,’symmetric’}

Direction of sweep

omega : { float}

Relaxation parameter typically in (0, 2) if omega != 1.0, then algorithm becomes SOR on A.H A

Dinv : { ndarray}

Inverse of diag(A.H A), (length N)

Returns :

Nothing, x will be modified in place. :

References

[R33]Yousef Saad, “Iterative Methods for Sparse Linear Systems, Second Edition”, SIAM, pp. 247-9, 2003 http://www-users.cs.umn.edu/~saad/books.html

Examples

>>> ## Use NR Gauss-Seidel as a Stand-Alone Solver
>>> from pyamg.relaxation import *
>>> from pyamg.gallery import poisson
>>> from pyamg.util.linalg import norm
>>> import numpy
>>> A = poisson((10,10), format='csr')
>>> x0 = numpy.zeros((A.shape[0],1))
>>> b = numpy.ones((A.shape[0],1))
>>> gauss_seidel_nr(A, x0, b, iterations=10, sweep='symmetric')
>>> print norm(b-A*x0)
8.45044864352
>>> #
>>> ## Use NR Gauss-Seidel as the Multigrid Smoother
>>> from pyamg import smoothed_aggregation_solver
>>> sa = smoothed_aggregation_solver(A, B=numpy.ones((A.shape[0],1)),
...      coarse_solver='pinv2', max_coarse=50,
...      presmoother=('gauss_seidel_nr', {'sweep' : 'symmetric'}), 
...      postsmoother=('gauss_seidel_nr', {'sweep' : 'symmetric'}))
>>> x0=numpy.zeros((A.shape[0],1))
>>> residuals=[]
>>> x = sa.solve(b, x0=x0, tol=1e-8, residuals=residuals)
pyamg.relaxation.relaxation.gauss_seidel_indexed(A, x, b, indices, iterations=1, sweep='forward')[source]

Perform indexed Gauss-Seidel iteration on the linear system Ax=b

In indexed Gauss-Seidel, the sequence in which unknowns are relaxed is specified explicitly. In contrast, the standard Gauss-Seidel method always performs complete sweeps of all variables in increasing or decreasing order. The indexed method may be used to implement specialized smoothers, like F-smoothing in Classical AMG.

Parameters :

A : csr_matrix

Sparse NxN matrix

x : ndarray

Approximate solution (length N)

b : ndarray

Right-hand side (length N)

indices : ndarray

Row indices to relax.

iterations : int

Number of iterations to perform

sweep : {‘forward’,’backward’,’symmetric’}

Direction of sweep

Returns :

Nothing, x will be modified in place. :

Examples

>>> from pyamg.gallery import poisson
>>> import numpy
>>> A = poisson((4,), format='csr')
>>> x = numpy.array([0.0, 0.0, 0.0, 0.0])
>>> b = numpy.array([0.0, 1.0, 2.0, 3.0])
>>> gauss_seidel_indexed(A, x, b, [0,1,2,3])                #relax all four rows, in order
>>> gauss_seidel_indexed(A, x, b, [0,1])                    #relax first two rows
>>> gauss_seidel_indexed(A, x, b, [2,0])                    #relax row 2, then row 0
>>> gauss_seidel_indexed(A, x, b, [2,3], sweep='backward')  #relax row 3, then row 2
>>> gauss_seidel_indexed(A, x, b, [2,0,2])                  #relax row 2, then 0, then 2 again
pyamg.relaxation.relaxation.block_jacobi(A, x, b, Dinv=None, blocksize=1, iterations=1, omega=1.0)[source]

Perform block Jacobi iteration on the linear system Ax=b

Parameters :

A : csr_matrix or bsr_matrix

Sparse NxN matrix

x : ndarray

Approximate solution (length N)

b : ndarray

Right-hand side (length N)

Dinv : array

Array holding block diagonal inverses of A size (N/blocksize, blocksize, blocksize)

blocksize : int

Desired dimension of blocks

iterations : int

Number of iterations to perform

omega : scalar

Damping parameter

Returns :

Nothing, x will be modified in place. :

Examples

>>> ## Use block Jacobi as a Stand-Alone Solver
>>> from pyamg.relaxation import *
>>> from pyamg.gallery import poisson
>>> from pyamg.util.linalg import norm
>>> import numpy
>>> A = poisson((10,10), format='csr')
>>> x0 = numpy.zeros((A.shape[0],1))
>>> b = numpy.ones((A.shape[0],1))
>>> block_jacobi(A, x0, b, blocksize=4, iterations=10, omega=1.0)
>>> print norm(b-A*x0)
4.66474230129
>>> #
>>> ## Use block Jacobi as the Multigrid Smoother
>>> from pyamg import smoothed_aggregation_solver
>>> sa = smoothed_aggregation_solver(A, B=numpy.ones((A.shape[0],1)),
...        coarse_solver='pinv2', max_coarse=50,
...        presmoother=('block_jacobi', {'omega': 4.0/3.0, 'iterations' : 2, 'blocksize' : 4}), 
...        postsmoother=('block_jacobi', {'omega': 4.0/3.0, 'iterations' : 2, 'blocksize' : 4}))
>>> x0=numpy.zeros((A.shape[0],1))
>>> residuals=[]
>>> x = sa.solve(b, x0=x0, tol=1e-8, residuals=residuals)
pyamg.relaxation.relaxation.block_gauss_seidel(A, x, b, iterations=1, sweep='forward', blocksize=1, Dinv=None)[source]

Perform block Gauss-Seidel iteration on the linear system Ax=b

Parameters :

A : {csr_matrix, bsr_matrix}

Sparse NxN matrix

x : ndarray

Approximate solution (length N)

b : ndarray

Right-hand side (length N)

iterations : int

Number of iterations to perform

sweep : {‘forward’,’backward’,’symmetric’}

Direction of sweep

Dinv : array

Array holding block diagonal inverses of A size (N/blocksize, blocksize, blocksize)

blocksize : int

Desired dimension of blocks

Returns :

Nothing, x will be modified in place. :

Examples

>>> ## Use Gauss-Seidel as a Stand-Alone Solver
>>> from pyamg.relaxation import *
>>> from pyamg.gallery import poisson
>>> from pyamg.util.linalg import norm
>>> import numpy
>>> A = poisson((10,10), format='csr')
>>> x0 = numpy.zeros((A.shape[0],1))
>>> b = numpy.ones((A.shape[0],1))
>>> block_gauss_seidel(A, x0, b, iterations=10, blocksize=4, sweep='symmetric')
>>> print norm(b-A*x0)
0.958333817624
>>> #
>>> ## Use Gauss-Seidel as the Multigrid Smoother
>>> from pyamg import smoothed_aggregation_solver
>>> sa = smoothed_aggregation_solver(A, B=numpy.ones((A.shape[0],1)),
...        coarse_solver='pinv2', max_coarse=50,
...        presmoother=('block_gauss_seidel', {'sweep':'symmetric', 'blocksize' : 4}), 
...        postsmoother=('block_gauss_seidel', {'sweep':'symmetric', 'blocksize' : 4}))
>>> x0=numpy.zeros((A.shape[0],1))
>>> residuals=[]
>>> x = sa.solve(b, x0=x0, tol=1e-8, residuals=residuals)

The setup Module

pyamg.relaxation.setup.configuration(parent_package='', top_path=None)[source]

The smoothing Module

Method to create pre- and post-smoothers on the levels of a multilevel_solver

pyamg.relaxation.smoothing.change_smoothers(ml, presmoother, postsmoother)[source]

Initialize pre- and post- smoothers throughout a multilevel_solver, with the option of having different smoothers at different levels

For each level of the multilevel_solver ‘ml’ (except the coarsest level), initialize the .presmoother() and .postsmoother() methods used in the multigrid cycle.

Parameters :

ml : {pyamg multilevel hierarchy}

Data structure that stores the multigrid hierarchy.

presmoother : {None, string, tuple, list}

presmoother can be (1) the name of a supported smoother, e.g. “gauss_seidel”, (2) a tuple of the form (‘method’,’opts’) where ‘method’ is the name of a supported smoother and ‘opts’ a dict of keyword arguments to the smoother, or (3) a list of instances of options 1 or 2. See the Examples section for illustrations of the format.

If presmoother is a list, presmoother[i] determines the smoothing strategy for level i. Else, presmoother defines the same strategy for all levels.

If len(presmoother) < len(ml.levels), then presmoother[-1] is used for all remaining levels

If len(presmoother) > len(ml.levels), then the remaining smoothing strategies are ignored

postsmoother : {string, tuple, list}

Defines postsmoother in identical fashion to presmoother

Returns :

ml changed in place :

ml.level[i].presmoother <=== presmoother[i] :

ml.level[i].postsmoother <=== postsmoother[i] :

Notes

  • Parameter ‘omega’ of the Jacobi, Richardson, and jacobi_ne methods is scaled by the spectral radius of the matrix on each level. Therefore ‘omega’ should be in the interval (0,2).

  • Parameter ‘withrho’ (default: True) controls whether the omega is rescaled by the spectral radius in jacobi, block_jacobi, and jacobi_ne

  • By initializing the smoothers after the hierarchy has been setup, allows for “algebraically” directed relaxation, such as strength_based_schwarz, which uses only the strong connections of a degree-of-freedom to define overlapping regions

  • Available smoother methods:

    gauss_seidel
    block_gauss_seidel
    jacobi
    block_jacobi
    richardson
    sor
    chebyshev
    gauss_seidel_nr
    gauss_seidel_ne
    jacobi_ne
    cg
    gmres
    cgne
    cgnr
    schwarz
    strength_based_schwarz
    None
    

Examples

>>> from pyamg.gallery import poisson
>>> from pyamg.aggregation import smoothed_aggregation_solver
>>> from pyamg.relaxation.smoothing import change_smoothers
>>> from pyamg.util.linalg import norm
>>> from scipy import rand, array, mean
>>> A = poisson((10,10), format='csr')
>>> b = rand(A.shape[0],)
>>> ml = smoothed_aggregation_solver(A, max_coarse=10)
>>> #
>>> ## Set all levels to use gauss_seidel's defaults
>>> smoothers = 'gauss_seidel'
>>> change_smoothers(ml, presmoother=smoothers, postsmoother=smoothers)
>>> residuals=[]
>>> x = ml.solve(b, tol=1e-8, residuals=residuals)
>>> #
>>> ## Set all levels to use three iterations of gauss_seidel's defaults
>>> smoothers = ('gauss_seidel', {'iterations' : 3})
>>> change_smoothers(ml, presmoother=smoothers, postsmoother=None)
>>> residuals=[]
>>> x = ml.solve(b, tol=1e-8, residuals=residuals)
>>> #
>>> ## Set level 0 to use gauss_seidel's defaults, and all  
>>> ## subsequent levels to use 5 iterations of cgnr
>>> smoothers = ['gauss_seidel', ('cgnr', {'maxiter' : 5})]
>>> change_smoothers(ml, presmoother=smoothers, postsmoother=smoothers)
>>> residuals=[]
>>> x = ml.solve(b, tol=1e-8, residuals=residuals)

Table Of Contents

Previous topic

Krylov Documentation

Next topic

Testing Documentation

This Page