SLATE 2024.05.31
Software for Linear Algebra Targeting Exascale
|
Factor, forward and back solve, invert. More...
Functions | |
template<typename scalar_t > | |
int64_t | slate::potrf (HermitianMatrix< scalar_t > &A, Options const &opts) |
Distributed parallel Cholesky factorization. | |
template<typename scalar_t > | |
void | slate::potri (HermitianMatrix< scalar_t > &A, Options const &opts) |
Distributed parallel Cholesky inversion. | |
template<typename scalar_t > | |
void | slate::potrs (HermitianMatrix< scalar_t > &A, Matrix< scalar_t > &B, Options const &opts) |
Distributed parallel Cholesky solve. | |
Factor, forward and back solve, invert.
int64_t slate::potrf | ( | HermitianMatrix< scalar_t > & | A, |
Options const & | opts | ||
) |
Distributed parallel Cholesky factorization.
Performs the Cholesky factorization of a Hermitian positive definite matrix \(A\).
The factorization has the form
\[ A = L L^H, \]
if \(A\) is stored lower, where \(L\) is a lower triangular matrix, or
\[ A = U^H U, \]
if \(A\) is stored upper, where \(U\) is an upper triangular matrix.
Complexity (in real): \(\approx \frac{1}{3} n^{3}\) flops.
scalar_t | One of float, double, std::complex<float>, std::complex<double>. |
[in,out] | A | On entry, the n-by-n Hermitian positive definite matrix \(A\). On exit, if return value = 0, the factor \(U\) or \(L\) from the Cholesky factorization \(A = U^H U\) or \(A = L L^H\). If scalar_t is real, \(A\) can be a SymmetricMatrix object. |
[in] | opts | Additional options, as map of name = value pairs. Possible options:
|
void slate::potri | ( | HermitianMatrix< scalar_t > & | A, |
Options const & | opts | ||
) |
Distributed parallel Cholesky inversion.
Computes the inverse of a complex Hermitian positive definite matrix \(A\) using the Cholesky factorization \(A = U^H U\) or \(A = L L^H\) computed by potrf
.
Complexity (in real): \(\approx \frac{2}{3} n^{3}\) flops.
scalar_t | One of float, double, std::complex<float>, std::complex<double>. |
[in,out] | A | On entry, the triangular factor \(U\) or \(L\) from the Cholesky factorization \(A = U^H U\) or \(A = L L^H\), as computed by potrf . On exit, the upper or lower triangle of the (Hermitian) inverse of \(A\), overwriting the input factor \(U\) or \(L\). If scalar_t is real, \(A\) can be a SymmetricMatrix object. |
[in] | opts | Additional options, as map of name = value pairs. Possible options:
|
TODO: return value
0 | successful exit |
>0 | for return value = \(i\), \(A(i,i)\) is exactly zero. The triangular matrix is singular and its inverse can not be computed. |
void slate::potrs | ( | HermitianMatrix< scalar_t > & | A, |
Matrix< scalar_t > & | B, | ||
Options const & | opts | ||
) |
Distributed parallel Cholesky solve.
Solves a system of linear equations
\[ A X = B \]
with a Hermitian positive definite matrix \(A\) using the Cholesky factorization \(A = U^H U\) or \(A = L L^H\) computed by potrf.
Complexity (in real): \(2 n^{2} r\) flops.
scalar_t | One of float, double, std::complex<float>, std::complex<double>. |
[in] | A | The n-by-n triangular factor \(U\) or \(L\) from the Cholesky factorization \(A = U^H U\) or \(A = L L^H\), computed by potrf. If scalar_t is real, \(A\) can be a SymmetricMatrix object. |
[in,out] | B | On entry, the n-by-nrhs right hand side matrix \(B\). On exit, the n-by-nrhs solution matrix \(X\). |
[in] | opts | Additional options, as map of name = value pairs. Possible options:
|