pyamg.krylov

Krylov solvers.

This module contains several Krylov subspace methods, in addition to two simple iterations, to solve linear systems iteratively. These methods often use multigrid as a preconditioner to accelerate convergence to the solution. See [1] and [2].

Functions

• gmres

• fgmres

• cgne

• cgnr

• cg

• bicgstab

• steepest descent, (simple iteration)

• minimal residual (MR), (simple iteration)

References

[1]

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

[2]

Richard Barrett et al. “Templates for the Solution of Linear Systems: Building Blocks for Iterative Methods, 2nd Edition”, SIAM http://www.netlib.org/linalg/html_templates/Templates.html http://www.netlib.org/templates/

pyamg.krylov.bicgstab(A, b, x0=None, tol=1e-05, criteria='rr', maxiter=None, M=None, callback=None, residuals=None)

Biconjugate Gradient Algorithm with Stabilization.

Solves the linear system Ax = b. Left preconditioning is supported.

Parameters:

Aarray, matrix, sparse matrix, LinearOperator

Linear system of size (n,n) to solve.

barray, matrix

Right hand side of size (n,) or (n,1).

x0array, matrix

Initial guess, default is a vector of zeros.

tolfloat

Tolerance for stopping criteria.

criteriastr

Stopping criteria, let r=r_k, x=x_k

‘rr’: ||r|| < tol ||b|| ‘rr+’: ||r|| < tol (||b|| + ||A||_F ||x||)

if ||b||=0, then set ||b||=1 for these tests.

maxiterint

Maximum number of iterations allowed.

Marray, matrix, sparse matrix, LinearOperator

Inverted preconditioner of size (n,n), i.e. solve M A x = M b.

callbackfunction

User-supplied function is called after each iteration as callback(xk), where xk is the current solution vector.

residualslist

Residual history in the 2-norm, including the initial residual.

Returns:

array

Updated guess after k iterations to the solution of Ax = b.

int

Halting status

0	successful exit
>0	convergence to tolerance not achieved, return iteration count instead.
<0	numerical breakdown, or illegal input

Notes

The LinearOperator class is in scipy.sparse.linalg. Use this class if you prefer to define A or M as a mat-vec routine as opposed to explicitly constructing the matrix.

References

[1]

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

Examples

>>> from pyamg.krylov import bicgstab
>>> from pyamg.util.linalg import norm
>>> import numpy as np
>>> from pyamg.gallery import poisson
>>> A = poisson((10,10))
>>> b = np.ones((A.shape[0],))
>>> (x,flag) = bicgstab(A,b, maxiter=2, tol=1e-8)
>>> print(f'{norm(b - A@x):.6}')
4.68163

pyamg.krylov.cg(A, b, x0=None, tol=1e-05, criteria='rr', maxiter=None, M=None, callback=None, residuals=None)

Conjugate Gradient algorithm.

Solves the linear system Ax = b. Left preconditioning is supported.

Parameters:

Aarray, matrix, sparse matrix, LinearOperator

Linear system of size (n,n) to solve.

barray, matrix

Right hand side of size (n,) or (n,1).

x0array, matrix

Initial guess, default is a vector of zeros.

tolfloat

Tolerance for stopping criteria.

criteriastr

Stopping criteria, let r=r_k, x=x_k

‘rr’: ||r|| < tol ||b|| ‘rr+’: ||r|| < tol (||b|| + ||A||_F ||x||) ‘MrMr’: ||M r|| < tol ||M b|| ‘rMr’: <r, Mr>^1/2 < tol

if ||b||=0, then set ||b||=1 for these tests.

maxiterint

Maximum number of iterations allowed.

Marray, matrix, sparse matrix, LinearOperator

Inverse preconditioner of size (n,n), i.e. solve M A x = M b.

callbackfunction

User-supplied function is called after each iteration as callback(xk), where xk is the current solution vector.

residualslist

Residual history in the 2-norm, including the initial residual.

Returns:

array

Updated guess after k iterations to the solution of Ax = b.

int

Halting status

0	successful exit
>0	convergence to tolerance not achieved, return iteration count instead.
<0	numerical breakdown, or illegal input

Notes

The LinearOperator class is in scipy.sparse.linalg. Use this class if you prefer to define A or M as a mat-vec routine as opposed to explicitly constructing the matrix.

References

[1]

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

Examples

>>> from pyamg.krylov import cg
>>> from pyamg.util.linalg import norm
>>> import numpy as np
>>> from pyamg.gallery import poisson
>>> A = poisson((10,10))
>>> b = np.ones((A.shape[0],))
>>> (x,flag) = cg(A,b, maxiter=2, tol=1e-8)
>>> print(f'{norm(b - A@x):.6}')
10.9371

pyamg.krylov.cgne(A, b, x0=None, tol=1e-05, criteria='rr', maxiter=None, M=None, callback=None, residuals=None)

Conjugate Gradient, Normal Error algorithm.

Applies CG to the normal equations, A A.H x = b. Left preconditioning is supported. Note that unless A is well-conditioned, the use of CGNE is inadvisable.

Parameters:

Aarray, matrix, sparse matrix, LinearOperator

Linear system of size (n,n) to solve.

barray, matrix

Right hand side of size (n,) or (n,1).

x0array, matrix

Initial guess, default is a vector of zeros.

tolfloat

Tolerance for stopping criteria.

criteriastr

Stopping criteria, let r=r_k, x=x_k

‘rr’: ||r|| < tol ||b|| ‘rr+’: ||r|| < tol (||b|| + ||A||_F ||x||) ‘MrMr’: ||M r|| < tol ||M b|| ‘rMr’: <r, Mr>^1/2 < tol

if ||b||=0, then set ||b||=1 for these tests.

maxiterint

Maximum number of iterations allowed.

Marray, matrix, sparse matrix, LinearOperator

Inverted preconditioner of size (n,n), i.e. solve M A A.H x = M b.

callbackfunction

User-supplied function is called after each iteration as callback(xk), where xk is the current solution vector.

residualslist

Residual history in the 2-norm, including the initial residual.

Returns:

array

Updated guess after k iterations to the solution of Ax = b.

int

Halting status

0	successful exit
>0	convergence to tolerance not achieved, return iteration count instead.
<0	numerical breakdown, or illegal input

Notes

The LinearOperator class is in scipy.sparse.linalg. Use this class if you prefer to define A or M as a mat-vec routine as opposed to explicitly constructing the matrix.

References

[1]

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

Examples

>>> from pyamg.krylov import cgne
>>> from pyamg.util.linalg import norm
>>> import numpy as np
>>> from pyamg.gallery import poisson
>>> A = poisson((10,10))
>>> b = np.ones((A.shape[0],))
>>> (x,flag) = cgne(A,b, maxiter=2, tol=1e-8)
>>> print(f'{norm(b - A@x):.6}')
46.1547

pyamg.krylov.cgnr(A, b, x0=None, tol=1e-05, criteria='rr', maxiter=None, M=None, callback=None, residuals=None)

Conjugate Gradient, Normal Residual algorithm.

Applies CG to the normal equations, A.H A x = A.H b. Left preconditioning is supported. Note that unless A is well-conditioned, the use of CGNR is inadvisable

Parameters:

Aarray, matrix, sparse matrix, LinearOperator

Linear system of size (n,n) to solve.

barray, matrix

Right hand side of size (n,) or (n,1).

x0array, matrix

Initial guess, default is a vector of zeros.

tolfloat

Tolerance for stopping criteria.

criteriastr

Stopping criteria, let r=r_k, x=x_k

‘rr’: ||r|| < tol ||b|| ‘rr+’: ||r|| < tol (||b|| + ||A||_F ||x||) ‘MrMr’: ||M r|| < tol ||M b|| ‘rMr’: <r, Mr>^1/2 < tol

if ||b||=0, then set ||b||=1 for these tests.

maxiterint

Maximum number of iterations allowed.

Marray, matrix, sparse matrix, LinearOperator

Inverted preconditioner of size (n,n), i.e. solve M A.H A x = M A.H b.

callbackfunction

User-supplied function is called after each iteration as callback(xk), where xk is the current solution vector.

residualslist

Residual history in the 2-norm, including the initial residual.

Returns:

array

Updated guess after k iterations to the solution of Ax = b.

int

Halting status

0	successful exit
>0	convergence to tolerance not achieved, return iteration count instead.
<0	numerical breakdown, or illegal input

Notes

The LinearOperator class is in scipy.sparse.linalg. Use this class if you prefer to define A or M as a mat-vec routine as opposed to explicitly constructing the matrix.

References

[1]

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

Examples

>>> from pyamg.krylov import cgnr
>>> from pyamg.util.linalg import norm
>>> import numpy as np
>>> from pyamg.gallery import poisson
>>> A = poisson((10,10))
>>> b = np.ones((A.shape[0],))
>>> (x,flag) = cgnr(A,b, maxiter=2, tol=1e-8)
>>> print(f'{norm(b - A@x):.6}')
9.39102

pyamg.krylov.cr(A, b, x0=None, tol=1e-05, criteria='rr', maxiter=None, M=None, callback=None, residuals=None)

Conjugate Residual algorithm.

Solves the linear system Ax = b. Left preconditioning is supported. The matrix A must be Hermitian symmetric (but not necessarily definite).

Parameters:

Aarray, matrix, sparse matrix, LinearOperator

Linear system of size to solve.

barray, matrix

Right hand side of size (n,) or (n,1).

x0array, matrix

Initial guess, default is a vector of zeros.

tolfloat

Tolerance for stopping criteria.

criteriastr

Stopping criteria, let r=r_k, x=x_k

‘rr’: ||r|| < tol ||b|| ‘rr+’: ||r|| < tol (||b|| + ||A||_F ||x||) ‘MrMr’: ||M r|| < tol ||M b||

if ||b||=0, then set ||b||=1 for these tests.

maxiterint

Maximum number of iterations allowed.

Marray, matrix, sparse matrix, LinearOperator

Inverted preconditioner of size (n,n), i.e. solve M A x = M b.

callbackfunction

User-supplied function is called after each iteration as callback(xk), where xk is the current solution vector.

residualslist

Residual history in the 2-norm, including the initial residual.

Returns:

array

Updated guess after k iterations to the solution of Ax = b.

int

Halting status

0	successful exit
>0	convergence to tolerance not achieved, return iteration count instead.
<0	numerical breakdown, or illegal input

Notes

The LinearOperator class is in scipy.sparse.linalg. Use this class if you prefer to define A or M as a mat-vec routine as opposed to explicitly constructing the matrix.

References

[1]

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

Examples

>>> from pyamg.krylov import cr
>>> from pyamg.util.linalg import norm
>>> import numpy as np
>>> from pyamg.gallery import poisson
>>> A = poisson((10,10))
>>> b = np.ones((A.shape[0],))
>>> (x,flag) = cr(A,b, maxiter=2, tol=1e-8)
>>> print(f'{norm(b - A@x):.6}')
6.54282

pyamg.krylov.fgmres(A, b, x0=None, tol=1e-05, restart=None, maxiter=None, M=None, callback=None, residuals=None, restrt=None)

Flexible Generalized Minimum Residual Method (fGMRES).

fGMRES iteratively refines the initial solution guess to the system Ax = b. fGMRES is flexible in the sense that the right preconditioner (M) can vary from iteration to iteration.

Parameters:

Aarray, matrix, sparse matrix, LinearOperator

Linear system of size (n,n) to solve.

barray, matrix

Right hand side of size (n,) or (n,1).

x0array, matrix

Initial guess, default is a vector of zeros.

tolfloat

Tolerance for stopping criteria, let r=r_k

||r|| < tol ||b||

if ||b||=0, then set ||b||=1 for these tests.

restartNone, int

• If int, restart is max number of inner iterations and maxiter is the max number of outer iterations.

• If None, do not restart GMRES, and max number of inner iterations is maxiter.

maxiterNone, int

• If restart is None, maxiter is the max number of inner iterations and GMRES does not restart.

• If restart is int, maxiter is the max number of outer iterations, and restart is the max number of inner iterations.

• Defaults to min(n,40) if restart=None.

Marray, matrix, sparse matrix, LinearOperator

Inverted preconditioner of size (n,n), i.e. solve M A x = M b. M need not be stationary for fgmres.

callbackfunction

User-supplied function is called after each iteration as callback(xk), where xk is the current solution vector.

residualslist

Residual history in the 2-norm, including the initial residual.

restrtNone, int

Deprecated. See restart.

Returns:

array

Updated guess after k iterations to the solution of Ax = b.

int

Halting status

0	successful exit
>0	convergence to tolerance not achieved, return iteration count instead.
<0	numerical breakdown, or illegal input

Notes

The LinearOperator class is in scipy.sparse.linalg. Use this class if you prefer to define A or M as a mat-vec routine as opposed to explicitly constructing the matrix.

fGMRES allows for non-stationary preconditioners, as opposed to GMRES.

For robustness, Householder reflections are used to orthonormalize the Krylov Space. Givens Rotations are used to provide the residual norm each iteration. Flexibility implies that the right preconditioner, M, can vary from iteration to iteration.

References

[1]

Yousef Saad, “Iterative Methods for Sparse Linear Systems, Second Edition”, SIAM, pp. 151-172, pp. 272-275, 2003 http://www-users.cs.umn.edu/~saad/books.html

Examples

>>> from pyamg.krylov import fgmres
>>> from pyamg.util.linalg import norm
>>> import numpy as np
>>> from pyamg.gallery import poisson
>>> A = poisson((10,10))
>>> b = np.ones((A.shape[0],))
>>> (x,flag) = fgmres(A,b, maxiter=2, tol=1e-8)
>>> print(f'{norm(b - A@x):.6}')
6.54282

pyamg.krylov.gmres(A, b, x0=None, tol=1e-05, restart=None, maxiter=None, M=None, callback=None, residuals=None, orthog='householder', restrt=None, **kwargs)

Generalized Minimum Residual Method (GMRES).

GMRES iteratively refines the initial solution guess to the system Ax = b. Left preconditioned. Residuals are preconditioned residuals.

Parameters:

Aarray, matrix, sparse matrix, LinearOperator

Linear system of size (n,n) to solve.

barray, matrix

Right hand side of size (n,) or (n,1).

x0array, matrix

Initial guess, default is a vector of zeros.

tolfloat

Tolerance for stopping criteria, let r=r_k

||M r|| < tol ||M b||

if ||b||=0, then set ||M b||=1 for these tests.

restartNone, int

• If int, restart is max number of inner iterations and maxiter is the max number of outer iterations.

• If None, do not restart GMRES, and max number of inner iterations is maxiter.

maxiterNone, int

• If restart is None, maxiter is the max number of inner iterations and GMRES does not restart.

• If restart is int, maxiter is the max number of outer iterations, and restart is the max number of inner iterations.

• Defaults to min(n,40) if restart=None.

Marray, matrix, sparse matrix, LinearOperator

Inverted preconditioner of size (n,n), i.e. solve M A x = M b.

callbackfunction

User-supplied function is called after each iteration as callback(xk), where xk is the current solution vector.

residualslist

Preconditioned residual history in the 2-norm, including the initial residual.

orthogstr

• ‘householder’ calls _gmres_householder which uses Householder reflections to find the orthogonal basis for the Krylov space.

• ‘mgs’ calls _gmres_mgs which uses modified Gram-Schmidt to find the orthogonal basis for the Krylov space.

restrtNone, int

Deprecated. See restart.

**kwargsdict

Keyword arguments passed to _gmres_householder or _gmres_mgs.

Returns:

array

Updated guess after k iterations to the solution of Ax = b.

int

Halting status

0	successful exit
>0	convergence to tolerance not achieved, return iteration count instead.
<0	numerical breakdown, or illegal input

Notes

The LinearOperator class is in scipy.sparse.linalg. Use this class if you prefer to define A or M as a mat-vec routine as opposed to explicitly constructing the matrix.

The orthogonalization method, orthog='householder', is more robust than orthog='mgs', however for the majority of problems your problem will converge before ‘mgs’ loses orthogonality in your basis.

orthog='householder' has been more rigorously tested, and is therefore currently the default.

The residual is the preconditioned residual.

References

[1]

Yousef Saad, “Iterative Methods for Sparse Linear Systems, Second Edition”, SIAM, pp. 151-172, pp. 272-275, 2003 http://www-users.cs.umn.edu/~saad/books.html

Examples

>>> from pyamg.krylov import gmres
>>> from pyamg.util.linalg import norm
>>> import numpy as np
>>> from pyamg.gallery import poisson
>>> A = poisson((10,10))
>>> b = np.ones((A.shape[0],))
>>> (x,flag) = gmres(A,b, maxiter=2, tol=1e-8)
>>> print(f'{norm(b - A@x):.6}')
6.54282

pyamg.krylov.gmres_householder(A, b, x0=None, tol=1e-05, restart=None, maxiter=None, M=None, callback=None, residuals=None, restrt=None)

Generalized Minimum Residual Method (GMRES) based on Housholder.

GMRES iteratively refines the initial solution guess to the system Ax = b. Householder reflections are used for orthogonalization. Left preconditioning, leading to preconditioned residuals.

Parameters:

Aarray, matrix, sparse matrix, LinearOperator

Linear system of size (n,n) to solve.

barray, matrix

Right hand side of size (n,) or (n,1).

x0array, matrix

Initial guess, default is a vector of zeros.

tolfloat

Tolerance for stopping criteria, let r=r_k

||M r|| < tol ||M b||

if ||b||=0, then set ||M b||=1 for these tests.

restartNone, int

• If int, restart is max number of inner iterations and maxiter is the max number of outer iterations.

• If None, do not restart GMRES, and max number of inner iterations is maxiter.

maxiterNone, int

• If restart is None, maxiter is the max number of inner iterations and GMRES does not restart.

• If restart is int, maxiter is the max number of outer iterations, and restart is the max number of inner iterations.

• Defaults to min(n,40) if restart=None.

Marray, matrix, sparse matrix, LinearOperator

Inverted preconditioner of size (n,n), i.e. solve M A x = M b.

callbackfunction

User-supplied function is called after each iteration as callback(xk), where xk is the current solution vector.

residualslist

Preconditioned residual history in the 2-norm, including the initial preconditioned residual.

restrtNone, int

Deprecated. See restart.

Returns:

array

Updated guess after k iterations to the solution of Ax = b.

int

Halting status

0	successful exit
>0	convergence to tolerance not achieved, return iteration count instead.
<0	numerical breakdown, or illegal input

Notes

The LinearOperator class is in scipy.sparse.linalg. Use this class if you prefer to define A or M as a mat-vec routine as opposed to explicitly constructing the matrix.

For robustness, modified Gram-Schmidt is used to orthogonalize the Krylov Space Givens Rotations are used to provide the residual norm each iteration

The residual is the preconditioned residual.

References

[1]

Yousef Saad, “Iterative Methods for Sparse Linear Systems, Second Edition”, SIAM, pp. 151-172, pp. 272-275, 2003 http://www-users.cs.umn.edu/~saad/books.html

Examples

>>> from pyamg.krylov import gmres
>>> from pyamg.util.linalg import norm
>>> import numpy as np
>>> from pyamg.gallery import poisson
>>> A = poisson((10, 10))
>>> b = np.ones((A.shape[0],))
>>> (x, flag) = gmres(A, b, maxiter=2, tol=1e-8, orthog='householder')
>>> print(f'{norm(b - A@x):.6}')
6.54282

pyamg.krylov.gmres_mgs(A, b, x0=None, tol=1e-05, restart=None, maxiter=None, M=None, callback=None, residuals=None, reorth=False, restrt=None)

Generalized Minimum Residual Method (GMRES) based on MGS.

GMRES iteratively refines the initial solution guess to the system Ax = b. Modified Gram-Schmidt version. Left preconditioning, leading to preconditioned residuals.

Parameters:

Aarray, matrix, sparse matrix, LinearOperator

Linear system of size (n,n) to solve.

barray, matrix

Right-hand side of size (n,) or (n,1).

x0array, matrix

Initial guess, default is a vector of zeros.

tolfloat

Tolerance for stopping criteria, let r=r_k

||M r|| < tol ||M b||

if ||b||=0, then set ||M b||=1 for these tests.

restartNone, int

• If int, restart is max number of inner iterations and maxiter is the max number of outer iterations.

• If None, do not restart GMRES, and max number of inner iterations is maxiter.

maxiterNone, int

• If restart is None, maxiter is the max number of inner iterations and GMRES does not restart.

• If restart is int, maxiter is the max number of outer iterations, and restart is the max number of inner iterations.

• Defaults to min(n,40) if restart=None.

Marray, matrix, sparse matrix, LinearOperator

Inverted preconditioner of size (n,n), i.e. solve M A x = M b.

callbackfunction

User-supplied function is called after each iteration as callback(xk), where xk is the current solution vector.

residualslist

Preconditioned residual history in the 2-norm, including the initial preconditioned residual.

reorthbool

If True, then a check is made whether to re-orthogonalize the Krylov space each GMRES iteration.

restrtNone, int

Deprecated. See restart.

Returns:

array

Updated guess after k iterations to the solution of Ax = b.

int

Halting status

0	successful exit
>0	convergence to tolerance not achieved, return iteration count instead.
<0	numerical breakdown, or illegal input

Notes

The LinearOperator class is in scipy.sparse.linalg. Use this class if you prefer to define A or M as a mat-vec routine as opposed to explicitly constructing the matrix.

For robustness, modified Gram-Schmidt is used to orthogonalize the Krylov Space Givens Rotations are used to provide the residual norm each iteration

The residual is the preconditioned residual.

References

[1]

Yousef Saad, “Iterative Methods for Sparse Linear Systems, Second Edition”, SIAM, pp. 151-172, pp. 272-275, 2003 http://www-users.cs.umn.edu/~saad/books.html

[2]

3. 20. Kelley, http://www4.ncsu.edu/~ctk/matlab_roots.html

Examples

>>> from pyamg.krylov import gmres
>>> from pyamg.util.linalg import norm
>>> import numpy as np
>>> from pyamg.gallery import poisson
>>> A = poisson((10,10))
>>> b = np.ones((A.shape[0],))
>>> (x,flag) = gmres(A,b, maxiter=2, tol=1e-8, orthog='mgs')
>>> print(f'{norm(b - A@x):.6}')
6.54282

pyamg.krylov.minimal_residual(A, b, x0=None, tol=1e-05, maxiter=None, M=None, callback=None, residuals=None)

Minimal residual (MR) algorithm. 1D projection method.

Solves the linear system Ax = b. Left preconditioning is supported.

Parameters:

Aarray, matrix, sparse matrix, LinearOperator

Linear system of size (n,n) to solve.

barray, matrix

Right-hand side of size (n,) or (n,1).

x0array, matrix

Initial guess, default is a vector of zeros.

tolfloat

Tolerance for stopping criteria, let r=r_k

||M r|| < tol ||M b||

if ||b||=0, then set ||M b||=1 for these tests.

maxiterint

Maximum number of iterations allowed.

Marray, matrix, sparse matrix, LinearOperator

Inverted preconditioner of size (n,n), i.e. solve M A x = M b.

callbackfunction

User-supplied function is called after each iteration as callback(xk), where xk is the current solution vector.

residualslist

Preconditioned residual history in the 2-norm, including the initial preconditioned residual.

Returns:

array

Updated guess after k iterations to the solution of Ax = b.

int

Halting status

0	successful exit
>0	convergence to tolerance not achieved, return iteration count instead.
<0	numerical breakdown, or illegal input

See also

_steepest_descent

Notes

The LinearOperator class is in scipy.sparse.linalg. Use this class if you prefer to define A or M as a mat-vec routine as opposed to explicitly constructing the matrix.

The algorithm:

minimal residual algorithm:      Preconditioned version:
r = b - A x                      r = b - A x, z = M r
while not converged:             while not converged:
    p = A r                          p = M A z
    alpha = (p,r) / (p,p)            alpha = (p, z) / (p, p)
    x = x + alpha r                  x = x + alpha z
    r = r - alpha p                  z = z - alpha p

References

[1]

Yousef Saad, “Iterative Methods for Sparse Linear Systems, Second Edition”, SIAM, pp. 137–142, 2003 http://www-users.cs.umn.edu/~saad/books.html

Examples

>>> from pyamg.krylov import minimal_residual
>>> from pyamg.util.linalg import norm
>>> import numpy as np
>>> from pyamg.gallery import poisson
>>> A = poisson((10,10))
>>> b = np.ones((A.shape[0],))
>>> (x,flag) = minimal_residual(A,b, maxiter=2, tol=1e-8)
>>> print(f'{norm(b - A@x):.6}')
7.26369

pyamg.krylov.steepest_descent(A, b, x0=None, tol=1e-05, criteria='rr', maxiter=None, M=None, callback=None, residuals=None)

Steepest descent algorithm.

Solves the linear system Ax = b. Left preconditioning is supported.

Parameters:

Aarray, matrix, sparse matrix, LinearOperator

Linear system of size (n,n) to solve.

barray, matrix

Right hand side of size (n,) or (n,1).

x0array, matrix

Initial guess, default is a vector of zeros.

tolfloat

Tolerance for stopping criteria.

criteriastr

Stopping criteria, let r=r_k, x=x_k

‘rr’	||r|| < tol ||b||
‘rr+’	||r|| < tol (||b|| + ||A||_F ||x||)
‘MrMr’	||M r|| < tol ||M b||
‘rMr’	<r, Mr>^1/2 < tol

if ||b||=0, then set ||b||=1 for these tests.

maxiterint

Maximum number of iterations allowed.

Marray, matrix, sparse matrix, LinearOperator

Inverted preconditioner of size (n,n), i.e. solve M A x = M b.

callbackfunction

User-supplied function is called after each iteration as callback(xk), where xk is the current solution vector.

residualslist

Residual history in the 2-norm, including the initial residual.

Returns:

array

Updated guess after k iterations to the solution of Ax = b.

int

Halting status of cg

0	successful exit
>0	convergence to tolerance not achieved, return iteration count instead.
<0	numerical breakdown, or illegal input

See also

_minimal_residual

Notes

The LinearOperator class is in scipy.sparse.linalg. Use this class if you prefer to define A or M as a mat-vec routine as opposed to explicitly constructing the matrix.

References

[1]

Yousef Saad, “Iterative Methods for Sparse Linear Systems, Second Edition”, SIAM, pp. 137–142, 2003 http://www-users.cs.umn.edu/~saad/books.html

Examples

>>> from pyamg.krylov import steepest_descent
>>> from pyamg.util.linalg import norm
>>> import numpy as np
>>> from pyamg.gallery import poisson
>>> A = poisson((10,10))
>>> b = np.ones((A.shape[0],))
>>> (x,flag) = steepest_descent(A,b, maxiter=2, tol=1e-8)
>>> print(f'{norm(b - A@x):.6}')
7.89436