SLATE 2024.05.31
Software for Linear Algebra Targeting Exascale
|
Factor, forward and back solve, invert. More...
Functions | |
template<typename scalar_t > | |
void | slate::gerbt (Matrix< scalar_t > &U_in, Matrix< scalar_t > &A, Matrix< scalar_t > &V) |
Applies a 2-sided RBT to the given matrix. | |
template<typename scalar_t > | |
void | slate::gerbt (Matrix< scalar_t > &Uin, Matrix< scalar_t > &B) |
Applies a 1-sided RBT to the given matrix on the left. | |
template<typename scalar_t > | |
int64_t | slate::getrf (Matrix< scalar_t > &A, Pivots &pivots, Options const &opts) |
Distributed parallel LU factorization. | |
template<typename scalar_t > | |
int64_t | slate::getrf_nopiv (Matrix< scalar_t > &A, Options const &opts) |
Distributed parallel LU factorization without pivoting. | |
template<typename scalar_t > | |
int64_t | slate::getrf_tntpiv (Matrix< scalar_t > &A, Pivots &pivots, Options const &opts) |
Distributed parallel LU factorization. | |
template<typename scalar_t > | |
void | slate::getri (Matrix< scalar_t > &A, Pivots &pivots, Options const &opts) |
Distributed parallel LU inversion. | |
template<typename scalar_t > | |
void | slate::getri (Matrix< scalar_t > &A, Pivots &pivots, Matrix< scalar_t > &B, Options const &opts) |
Distributed parallel LU inversion (out-of-place version). | |
template<typename scalar_t > | |
void | slate::getrs (Matrix< scalar_t > &A, Pivots &pivots, Matrix< scalar_t > &B, Options const &opts) |
Distributed parallel LU solve. | |
template<typename scalar_t > | |
void | slate::getrs_nopiv (Matrix< scalar_t > &A, Matrix< scalar_t > &B, Options const &opts=Options()) |
Factor, forward and back solve, invert.
void slate::gerbt | ( | Matrix< scalar_t > & | U_in, |
Matrix< scalar_t > & | A, | ||
Matrix< scalar_t > & | V | ||
) |
Applies a 2-sided RBT to the given matrix.
[in] | U_in | The left transform in packed storage. Should be transposed. |
[in,out] | A | The matrix to transform |
[in] | V | The right transform in packed storage. Should not be transposed. |
Applies a 1-sided RBT to the given matrix on the left.
[in] | U_in | The transform in packed storage |
[in,out] | A | The matrix to transform |
int64_t slate::getrf | ( | Matrix< scalar_t > & | A, |
Pivots & | pivots, | ||
Options const & | opts | ||
) |
Distributed parallel LU factorization.
Computes an LU factorization of a general m-by-n matrix \(A\) using partial pivoting with row interchanges.
The factorization has the form
\[ A = P L U \]
where \(P\) is a permutation matrix, \(L\) is lower triangular with unit diagonal elements (lower trapezoidal if m > n), and \(U\) is upper triangular (upper trapezoidal if m < n).
This is the right-looking Level 3 BLAS version of the algorithm.
Complexity (in real):
scalar_t | One of float, double, std::complex<float>, std::complex<double>. |
[in,out] | A | On entry, the matrix \(A\) to be factored. On exit, the factors \(L\) and \(U\) from the factorization \(A = P L U\); the unit diagonal elements of \(L\) are not stored. |
[out] | pivots | The pivot indices that define the permutation matrix \(P\). |
[in] | opts | Additional options, as map of name = value pairs. Possible options: |
int64_t slate::getrf_nopiv | ( | Matrix< scalar_t > & | A, |
Options const & | opts | ||
) |
Distributed parallel LU factorization without pivoting.
Computes an LU factorization without pivoting of a general m-by-n matrix \(A\)
The factorization has the form
\[ A = L U \]
where \(L\) is lower triangular with unit diagonal elements (lower trapezoidal if m > n), and \(U\) is upper triangular (upper trapezoidal if m < n).
This is the right-looking Level 3 BLAS version of the algorithm.
scalar_t | One of float, double, std::complex<float>, std::complex<double>. |
[in,out] | A | On entry, the matrix \(A\) to be factored. On exit, the factors \(L\) and \(U\) from the factorization \(A = P L U\); the unit diagonal elements of \(L\) are not stored. |
[in] | opts | Additional options, as map of name = value pairs. Possible options:
|
int64_t slate::getrf_tntpiv | ( | Matrix< scalar_t > & | A, |
Pivots & | pivots, | ||
Options const & | opts | ||
) |
Distributed parallel LU factorization.
Computes an LU factorization of a general m-by-n matrix \(A\) using partial pivoting with row interchanges.
The factorization has the form
\[ A = P L U \]
where \(P\) is a permutation matrix, \(L\) is lower triangular with unit diagonal elements (lower trapezoidal if m > n), and \(U\) is upper triangular (upper trapezoidal if m < n).
This is the right-looking Level 3 BLAS version of the algorithm.
scalar_t | One of float, double, std::complex<float>, std::complex<double>. |
[in,out] | A | On entry, the matrix \(A\) to be factored. On exit, the factors \(L\) and \(U\) from the factorization \(A = P L U\); the unit diagonal elements of \(L\) are not stored. |
[out] | pivots | The pivot indices that define the permutation matrix \(P\). |
[in] | opts | Additional options, as map of name = value pairs. Possible options:
|
void slate::getri | ( | Matrix< scalar_t > & | A, |
Pivots & | pivots, | ||
Matrix< scalar_t > & | B, | ||
Options const & | opts | ||
) |
Distributed parallel LU inversion (out-of-place version).
Computes the inverse of a matrix \(A\) using the LU factorization \(A = L*U\) computed by getrf
. Stores the result in \(B\). Does not change \(A\).
scalar_t | One of float, double, std::complex<float>, std::complex<double>. |
[in] | A | On entry, the factors \(L\) and \(U\) from the factorization \(A = P L U\) as computed by getrf. On exit, the inverse of the original matrix \(A\). |
[in] | pivots | The pivot indices that define the permutation matrix \(P\) as computed by getrf. |
[out] | B | On exit, if return value = 0, the n-by-n inverse of matrix \(A\). |
[in] | opts | Additional options, as map of name = value pairs. Possible options:
|
void slate::getri | ( | Matrix< scalar_t > & | A, |
Pivots & | pivots, | ||
Options const & | opts | ||
) |
Distributed parallel LU inversion.
Computes the inverse of a matrix \(A\) using the LU factorization \(A = L*U\) computed by getrf
.
Complexity (in real): \(\approx \frac{4}{3} n^{3}\) flops.
scalar_t | One of float, double, std::complex<float>, std::complex<double>. |
[in] | A | On entry, the factors \(L\) and \(U\) from the factorization \(A = P L U\) as computed by getrf. On exit, the inverse of the original matrix \(A\). |
[in] | pivots | The pivot indices that define the permutation matrix \(P\) as computed by getrf. |
[in] | opts | Additional options, as map of name = value pairs. Possible options:
|
void slate::getrs | ( | Matrix< scalar_t > & | A, |
Pivots & | pivots, | ||
Matrix< scalar_t > & | B, | ||
Options const & | opts | ||
) |
Distributed parallel LU solve.
Solves a system of linear equations
\[ A X = B \]
with a general n-by-n matrix \(A\) using the LU factorization computed by getrf. \(A\) can be transposed or conjugate-transposed.
Complexity (in real): \(\approx 2 n^{2} r\) flops.
scalar_t | One of float, double, std::complex<float>, std::complex<double>. |
[in] | A | The factors \(L\) and \(U\) from the factorization \(A = P L U\) as computed by getrf. |
[in] | pivots | The pivot indices that define the permutation matrix \(P\) as computed by gbtrf. |
[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: |
Option::MethodLU: Algorithm for LU factorization.
void slate::getrs_nopiv | ( | Matrix< scalar_t > & | A, |
Matrix< scalar_t > & | B, | ||
Options const & | opts | ||
) |
Distributed parallel LU solve.
Solves a system of linear equations
\[ A X = B \]
with a general n-by-n matrix \(A\) using the LU factorization computed by getrf. \(A\) can be transposed or conjugate-transposed.
scalar_t | One of float, double, std::complex<float>, std::complex<double>. |
[in] | A | The factors \(L\) and \(U\) from the factorization \(A = P L U\) as computed by getrf. |
[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:
|