Classes | |
class | Diagonal_Preconditioner |
class | Identity_Preconditioner |
class | Cholesky |
class | Eigenvalue |
class | LU |
class | QR |
class | SVD |
Functions | |
template<class MATRIX , class VECTOR , class PRECONDITIONER , class REAL > | |
int | PCG (const MATRIX &A, VECTOR &x, const VECTOR &b, const PRECONDITIONER &M, int &max_iter, REAL &tol) |
template<class const_MATRIX , class VECTOR , class const_VECTOR , class const_PRECONDITIONER , class const_TRANSPOSE_PRECONDITIONER , class REAL > | |
int | PBiCG (const_MATRIX &A, VECTOR &x, const_VECTOR &b, const_PRECONDITIONER &M, const_TRANSPOSE_PRECONDITIONER &Mt, int &max_iter, REAL &tol) |
template<class const_MATRIX , class VECTOR , class const_VECTOR , class const_PRECONDITIONER , class REAL > | |
int | BiCGSTAB (const_MATRIX &A, VECTOR &x, const_VECTOR &b, const_PRECONDITIONER &M, int &max_iter, REAL &tol) |
template<class const_MATRIX , class VECTOR , class const_VECTOR , class const_PRECONDITIONER , class REAL > | |
int | CGS (const_MATRIX &A, VECTOR &x, const_VECTOR &b, const_PRECONDITIONER &M, int &max_iter, REAL &tol) |
template<class const_MATRIX , class const_VECTOR , class VECTOR , class const_PRECONDITIONER , class REAL , class SCALAR > | |
int | CHEBY (const_MATRIX &A, VECTOR &x, const_VECTOR &b, const_PRECONDITIONER &M, int &max_iter, REAL &tol, SCALAR eigmin, SCALAR eigmax) |
template<class MATRIX , class VECTOR > | |
void | Update (VECTOR &x, int k, MATRIX &h, VECTOR &s, VECTOR v[]) |
template<class REAL > | |
REAL | abs (REAL x) |
template<class Operator , class VECTOR , class PRECONDITIONER , class MATRIX , class REAL > | |
int | GMRES (const Operator &A, VECTOR &x, const VECTOR &b, const PRECONDITIONER &M, MATRIX &H, int &m, int &max_iter, REAL &tol) |
template<class REAL > | |
void | GeneratePlaneRotation (REAL &dx, REAL &dy, REAL &cs, REAL &sn) |
template<class REAL > | |
void | ApplyPlaneRotation (REAL &dx, REAL &dy, REAL &cs, REAL &sn) |
template<class MATRIX , class VECTOR , class PRECONDITIONER , class REAL > | |
int | IR (const MATRIX &A, VECTOR &x, const VECTOR &b, const PRECONDITIONER &M, int &max_iter, REAL &tol) |
template<class const_MATRIX , class const_VECTOR , class VECTOR , class const_PRECONDITIONER1 , class const_PRECONDITIONER1t , class const_PRECONDITIONER2 , class const_PRECONDITIONER2t , class REAL > | |
int | QMR (const_MATRIX &A, VECTOR &x, const_VECTOR &b, const_PRECONDITIONER1 &M1, const_PRECONDITIONER1t &M1t, const_PRECONDITIONER2 &M2, const_PRECONDITIONER2t &M2t, int &max_iter, REAL &tol) |
REAL TNT::Linear_Algebra::abs | ( | REAL | x | ) | [inline] |
void TNT::Linear_Algebra::ApplyPlaneRotation | ( | REAL & | dx, | |
REAL & | dy, | |||
REAL & | cs, | |||
REAL & | sn | |||
) | [inline] |
int TNT::Linear_Algebra::BiCGSTAB | ( | const_MATRIX & | A, | |
VECTOR & | x, | |||
const_VECTOR & | b, | |||
const_PRECONDITIONER & | M, | |||
int & | max_iter, | |||
REAL & | tol | |||
) | [inline] |
BiCGSTAB : Iterative Method for solving system of linear equations.
PBiCG solves the unsymmetric linear system Ax = b using the Preconditioned BiConjugate Gradient Stabilized method. It is a stabilized version of BiCG. This implementation follows the algorithm described on p. 27 of the SIAM Templates book.
This algorithm requires a method for multiplying with A and the transpose of A.
A | (in) the matrix for which Ax=b is being solved. | |
x | (in/out) on input, the initial solution for Ax=b; upon return, the calculated value for the solution. | |
b | (in) the right hand side of the equation Ax=b. | |
M | (in) a function object that approximates A. Its use appears as M(y), to solve the system M*x = y. | |
max_iter | (in/out) on input, the maximum allowed number of iterations to perform, on return, the number of iterations actually performed. | |
tol | tolerace (in/out) in input, stop iteration when the scaled residual |(A*x-b)|/|b| is less than this number; on return, actual value of this scaled residual after the final iteration.. |
int TNT::Linear_Algebra::CGS | ( | const_MATRIX & | A, | |
VECTOR & | x, | |||
const_VECTOR & | b, | |||
const_PRECONDITIONER & | M, | |||
int & | max_iter, | |||
REAL & | tol | |||
) | [inline] |
CGS solves the unsymmetric linear system Ax = b using the preconditioned Conjugate Gradient Squared method.
CGS follows the algorithm described on p. 26 of the SIAM Templates book.
A | (in) the matrix for which Ax=b is being solved. | |
x | (in/out) on input, the initial solution for Ax=b; upon return, the calculated value for the solution. | |
b | (in) the right hand side of the equation Ax=b. | |
M | (in) a function object that approximates A. Its use appears as M(y), to solve the system M*x = y. | |
max_iter | (in/out) on input, the maximum allowed number of iterations to perform, on return, the number of iterations actually performed. | |
tol | tolerace (in/out) in input, stop iteration when the scaled residual |(A*x-b)|/|b| is less than this number; on return, actual value of this scaled residual after the final iteration.. |
int TNT::Linear_Algebra::CHEBY | ( | const_MATRIX & | A, | |
VECTOR & | x, | |||
const_VECTOR & | b, | |||
const_PRECONDITIONER & | M, | |||
int & | max_iter, | |||
REAL & | tol, | |||
SCALAR | eigmin, | |||
SCALAR | eigmax | |||
) | [inline] |
CHEBY solves the symmetric positive definite linear system Ax = b using the Preconditioned Chebyshev Method
CHEBY follows the algorithm described on p. 30 of the SIAM Templates book.
A | (in) the matrix for which Ax=b is being solved. | |
x | (in/out) on input, the initial solution for Ax=b; upon return, the calculated value for the solution. | |
b | (in) the right hand side of the equation Ax=b. | |
M | (in) a function object that approximates A. Its use appears as M(y), to solve the system M*x = y. | |
max_iter | (in/out) on input, the maximum allowed number of iterations to perform, on return, the number of iterations actually performed. | |
tol | tolerance (in/out) in input, stop iteration when the scaled residual |(A*x-b)|/|b| is less than this number; on return, actual value of this scaled residual after the final iteration.. | |
eigmin | minimum eigenvalue of A | |
eigmax | maximum eigenvalue of A |
void TNT::Linear_Algebra::GeneratePlaneRotation | ( | REAL & | dx, | |
REAL & | dy, | |||
REAL & | cs, | |||
REAL & | sn | |||
) | [inline] |
int TNT::Linear_Algebra::GMRES | ( | const Operator & | A, | |
VECTOR & | x, | |||
const VECTOR & | b, | |||
const PRECONDITIONER & | M, | |||
MATRIX & | H, | |||
int & | m, | |||
int & | max_iter, | |||
REAL & | tol | |||
) | [inline] |
Iterative template routine -- GMRES
GMRES solves the unsymmetric linear system Ax = b using the Generalized Minimum Residual method
GMRES follows the algorithm described on p. 20 of the SIAM Templates book.
A | (in) the matrix for which Ax=b is being solved. | |
x | (in/out) on input, the initial solution for Ax=b; upon return, the calculated value for the solution. | |
b | (in) the right hand side of the equation Ax=b. | |
M | (in) a function object that approximates A. Its use appears as M(y), to solve the system M*x = y. | |
H | ||
m | ||
max_iter | (in/out) on input, the maximum allowed number of iterations to perform, on return, the number of iterations actually performed. | |
tol | tolerance (in/out) in input, stop iteration when the scaled residual |(A*x-b)|/|b| is less than this number; on return, actual value of this scaled residual after the final iteration.. |
int TNT::Linear_Algebra::IR | ( | const MATRIX & | A, | |
VECTOR & | x, | |||
const VECTOR & | b, | |||
const PRECONDITIONER & | M, | |||
int & | max_iter, | |||
REAL & | tol | |||
) | [inline] |
Preconditioned Richardson Iteration (IR). IR solves the unsymmetric linear system Ax = b using Iterative Refinement (preconditioned Richardson iteration).
A | (in) the matrix for which Ax=b is being solved. | |
x | (in/out) on input, the initial solution for Ax=b; upon return, the calculated value for the solution. | |
b | (in) the right hand side of the equation Ax=b. | |
M | (in) a function object that approximates A. Its use appears as M(y), to solve the system M*x = y. | |
max_iter | (in/out) on input, the maximum allowed number of iterations to perform, on return, the number of iterations actually performed. | |
tol | tolerance (in/out) in input, stop iteration when the scaled residual |(A*x-b)|/|b| is less than this number; on return, actual value of this scaled residual after the final iteration.. |
int TNT::Linear_Algebra::PBiCG | ( | const_MATRIX & | A, | |
VECTOR & | x, | |||
const_VECTOR & | b, | |||
const_PRECONDITIONER & | M, | |||
const_TRANSPOSE_PRECONDITIONER & | Mt, | |||
int & | max_iter, | |||
REAL & | tol | |||
) | [inline] |
PBiCG solves the unsymmetric linear system Ax = b using the Preconditioned BiConjugate Gradient method BiCG follows the algorithm described on p. 22 of the SIAM Templates book.
This algorithm requires a method for multiplying with A and the transpose of A.
A | (in) the matrix for which Ax=b is being solved. | |
x | (in/out) on input, the initial solution for Ax=b; upon return, the calculated value for the solution. | |
b | (in) the right hand side of the equation Ax=b. | |
M | (in) a function object that approximates A. Its use appears as M(y), to solve the system M*x = y. | |
Mt | (in) a function object that approximates the transpose of A. Its use appears as Mt(y), to solve the system Mt*x = y. | |
max_iter | (in/out) on input, the maximum allowed number of iterations to perform, on return, the number of iterations actually performed. | |
tol | tolerance (in/out) in input, stop iteration when the scaled residual |(A*x-b)|/|b| is less than this number; on return, actual value of this scaled residual after the final iteration.. |
int TNT::Linear_Algebra::PCG | ( | const MATRIX & | A, | |
VECTOR & | x, | |||
const VECTOR & | b, | |||
const PRECONDITIONER & | M, | |||
int & | max_iter, | |||
REAL & | tol | |||
) | [inline] |
Preconditioned Conjuate Gradient (CG) solves the symmetric positive definite linear system Ax=b using the Conjugate Gradient method. (See p. 15 of the SIAM Template book.)
A | (in) the matrix for which Ax=b is being solved. | |
x | (in/out) on input, the initial solution for Ax=b; upon return, the calculated value for the solution. | |
b | (in) the right hand side of the equation Ax=b. | |
M | (in) a function object that approximates the inverse of A. Its use appears as M(y), to solve the system M*x = y. | |
max_iter | (in/out) on input, the maximum allowed number of iterations to perform, on return, the number of iterations actually performed. | |
tol | tolerance (in/out) in input, stop iteration when the scaled residual |(A*x-b)|/|b| is less than this number; on return, actual value of this scaled residual after the final iteration.. |
int TNT::Linear_Algebra::QMR | ( | const_MATRIX & | A, | |
VECTOR & | x, | |||
const_VECTOR & | b, | |||
const_PRECONDITIONER1 & | M1, | |||
const_PRECONDITIONER1t & | M1t, | |||
const_PRECONDITIONER2 & | M2, | |||
const_PRECONDITIONER2t & | M2t, | |||
int & | max_iter, | |||
REAL & | tol | |||
) | [inline] |
QMR solves the unsymmetric linear system Ax = b using the Quasi-Minimal Residual method following the algorithm as described on p. 24 in the SIAM Templates book.
A | (in) the matrix for which Ax=b is being solved. | |
x | (in/out) on input, the initial solution for Ax=b; upon return, the calculated value for the solution. | |
b | (in) the right hand side of the equation Ax=b. | |
M1 | (in) a function object that approximates A. Its use appears as M(y), to solve the system M*x = y. | |
M1t | (in) a function object that approximates A. Its use appears as M(y), to solve the system M*x = y. | |
M2 | (in) a function object that approximates A. Its use appears as M(y), to solve the system M*x = y. | |
M2t | (in) a function object that approximates A. Its use appears as M(y), to solve the system M*x = y. | |
max_iter | (in/out) on input, the maximum allowed number of iterations to perform, on return, the number of iterations actually performed. | |
tol | tolerance (in/out) in input, stop iteration when the scaled residual |(A*x-b)|/|b| is less than this number; on return, actual value of this scaled residual after the final iteration.. |
0 : | successful return | |
1 : | no convergence | |
2 : | rho | |
3 : | beta | |
4 : | gamma | |
5 : | delta | |
6 : | ep | |
7 : | xi |
void TNT::Linear_Algebra::Update | ( | VECTOR & | x, | |
int | k, | |||
MATRIX & | h, | |||
VECTOR & | s, | |||
VECTOR | v[] | |||
) | [inline] |