This page contains the Krylov Package documentation.
Biconjugate Gradient Algorithm with Stabilization
Solves the linear system Ax = b. Left preconditioning is supported.
Parameters : | A : {array, matrix, sparse matrix, LinearOperator}
b : {array, matrix}
x0 : {array, matrix}
tol : float
maxiter : int
xtype : type
M : {array, matrix, sparse matrix, LinearOperator}
callback : function
residuals : list
|
||||||
---|---|---|---|---|---|---|---|
Returns : | (xNew, info) : xNew : an updated guess to the solution of Ax = b info : halting status of bicgstab
|
Notes
The LinearOperator class is in scipy.sparse.linalg.interface. Use this class if you prefer to define A or M as a mat-vec routine as opposed to explicitly constructing the matrix. A.psolve(..) is still supported as a legacy.
References
[R19] | 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.bicgstab import bicgstab
>>> from pyamg.util.linalg import norm
>>> import numpy
>>> from pyamg.gallery import poisson
>>> A = poisson((10,10))
>>> b = numpy.ones((A.shape[0],))
>>> (x,flag) = bicgstab(A,b, maxiter=2, tol=1e-8)
>>> print norm(b - A*x)
4.68163045309
Conjugate Gradient algorithm
Solves the linear system Ax = b. Left preconditioning is supported.
Parameters : | A : {array, matrix, sparse matrix, LinearOperator}
b : {array, matrix}
x0 : {array, matrix}
tol : float
maxiter : int
xtype : type
M : {array, matrix, sparse matrix, LinearOperator}
callback : function
residuals : list
|
||||||
---|---|---|---|---|---|---|---|
Returns : | (xNew, info) : xNew : an updated guess to the solution of Ax = b info : halting status of cg
|
Notes
The LinearOperator class is in scipy.sparse.linalg.interface. Use this class if you prefer to define A or M as a mat-vec routine as opposed to explicitly constructing the matrix. A.psolve(..) is still supported as a legacy.
The residual in the preconditioner norm is both used for halting and returned in the residuals list.
References
[R20] | 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.cg import cg
>>> from pyamg.util.linalg import norm
>>> import numpy
>>> from pyamg.gallery import poisson
>>> A = poisson((10,10))
>>> b = numpy.ones((A.shape[0],))
>>> (x,flag) = cg(A,b, maxiter=2, tol=1e-8)
>>> print norm(b - A*x)
10.9370700187
Conjugate Gradient, Normal Error algorithm
Applies CG to the normal equations, A.H A x = b. Left preconditioning is supported. Note that unless A is well-conditioned, the use of CGNE is inadvisable
Parameters : | A : {array, matrix, sparse matrix, LinearOperator}
b : {array, matrix}
x0 : {array, matrix}
tol : float
maxiter : int
xtype : type
M : {array, matrix, sparse matrix, LinearOperator}
callback : function
residuals : list
|
||||||
---|---|---|---|---|---|---|---|
Returns : | (xNew, info) : xNew : an updated guess to the solution of Ax = b info : halting status of cgne
|
Notes
References
[R21] | 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.cgne import cgne
>>> from pyamg.util.linalg import norm
>>> import numpy
>>> from pyamg.gallery import poisson
>>> A = poisson((10,10))
>>> b = numpy.ones((A.shape[0],))
>>> (x,flag) = cgne(A,b, maxiter=2, tol=1e-8)
>>> print norm(b - A*x)
46.1547104367
Conjugate Gradient, Normal Residual algorithm
Applies CG to the normal equations, A.H A x = b. Left preconditioning is supported. Note that unless A is well-conditioned, the use of CGNR is inadvisable
Parameters : | A : {array, matrix, sparse matrix, LinearOperator}
b : {array, matrix}
x0 : {array, matrix}
tol : float
maxiter : int
xtype : type
M : {array, matrix, sparse matrix, LinearOperator}
callback : function
residuals : list
|
||||||
---|---|---|---|---|---|---|---|
Returns : | (xNew, info) : xNew : an updated guess to the solution of Ax = b info : halting status of cgnr
|
Notes
The LinearOperator class is in scipy.sparse.linalg.interface. Use this class if you prefer to define A or M as a mat-vec routine as opposed to explicitly constructing the matrix. A.psolve(..) is still supported as a legacy.
References
[R22] | 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.cgnr import cgnr
>>> from pyamg.util.linalg import norm
>>> import numpy
>>> from pyamg.gallery import poisson
>>> A = poisson((10,10))
>>> b = numpy.ones((A.shape[0],))
>>> (x,flag) = cgnr(A,b, maxiter=2, tol=1e-8)
>>> print norm(b - A*x)
9.3910201849
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 : | A : {array, matrix, sparse matrix, LinearOperator}
b : {array, matrix}
x0 : {array, matrix}
tol : float
restrt : {None, int}
maxiter : {None, int}
xtype : type
M : {array, matrix, sparse matrix, LinearOperator}
callback : function
residuals : list
|
||||||
---|---|---|---|---|---|---|---|
Returns : | (xNew, info) : xNew : an updated guess to the solution of Ax = b info : halting status of gmres
|
Notes
References
[R23] | 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.fgmres import fgmres
>>> from pyamg.util.linalg import norm
>>> import numpy
>>> from pyamg.gallery import poisson
>>> A = poisson((10,10))
>>> b = numpy.ones((A.shape[0],))
>>> (x,flag) = fgmres(A,b, maxiter=2, tol=1e-8)
>>> print norm(b - A*x)
6.5428213057
Parameters : | A : {array, matrix, sparse matrix, LinearOperator}
b : {array, matrix}
x0 : {array, matrix}
tol : float
restrt : {None, int}
maxiter : {None, int}
xtype : type
M : {array, matrix, sparse matrix, LinearOperator}
callback : function
residuals : list
orthog : string
|
||||||
---|---|---|---|---|---|---|---|
Returns : | (xNew, info) : xNew : an updated guess to the solution of Ax = b info : halting status of gmres
|
Notes
References
[R24] | 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
>>> from pyamg.gallery import poisson
>>> A = poisson((10,10))
>>> b = numpy.ones((A.shape[0],))
>>> (x,flag) = gmres(A,b, maxiter=2, tol=1e-8)
>>> print norm(b - A*x)
6.5428213057
Parameters : | A : {array, matrix, sparse matrix, LinearOperator}
b : {array, matrix}
x0 : {array, matrix}
tol : float
restrt : {None, int}
maxiter : {None, int}
xtype : type
M : {array, matrix, sparse matrix, LinearOperator}
callback : function
residuals : list
|
||||||
---|---|---|---|---|---|---|---|
Returns : | (xNew, info) : xNew : an updated guess to the solution of Ax = b info : halting status of gmres
|
Notes
References
[R25] | 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
>>> from pyamg.gallery import poisson
>>> A = poisson((10,10))
>>> b = numpy.ones((A.shape[0],))
>>> (x,flag) = gmres(A,b, maxiter=2, tol=1e-8, orthog='householder')
>>> print norm(b - A*x)
6.5428213057
Parameters : | A : {array, matrix, sparse matrix, LinearOperator}
b : {array, matrix}
x0 : {array, matrix}
tol : float
restrt : {None, int}
maxiter : {None, int}
xtype : type
M : {array, matrix, sparse matrix, LinearOperator}
callback : function
residuals : list
reorth : boolean
|
||||||
---|---|---|---|---|---|---|---|
Returns : | (xNew, info) : xNew : an updated guess to the solution of Ax = b info : halting status of gmres
|
Notes
References
[R26] | 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 |
[R27] |
Examples
>>> from pyamg.krylov import gmres
>>> from pyamg.util.linalg import norm
>>> import numpy
>>> from pyamg.gallery import poisson
>>> A = poisson((10,10))
>>> b = numpy.ones((A.shape[0],))
>>> (x,flag) = gmres(A,b, maxiter=2, tol=1e-8, orthog='mgs')
>>> print norm(b - A*x)
>>> 6.5428213057