SLATE 2024.05.31
Software for Linear Algebra Targeting Exascale
|
Solve \(AX = B\). More...
Functions | |
template<typename scalar_t > | |
int64_t | slate::gbsv (BandMatrix< scalar_t > &A, Pivots &pivots, Matrix< scalar_t > &B, Options const &opts) |
Distributed parallel band LU factorization and solve. | |
Solve \(AX = B\).
int64_t slate::gbsv | ( | BandMatrix< scalar_t > & | A, |
Pivots & | pivots, | ||
Matrix< scalar_t > & | B, | ||
Options const & | opts | ||
) |
Distributed parallel band LU factorization and solve.
Computes the solution to a system of linear equations
\[ A X = B, \]
where \(A\) is an n-by-n band matrix and \(X\) and \(B\) are n-by-nrhs matrices.
The LU decomposition with partial pivoting and row interchanges is used to factor \(A\) as
\[ A = L U, \]
where \(L\) is a product of permutation and unit lower triangular matrices, and \(U\) is upper triangular. The factored form of \(A\) is then used to solve the system of equations \(A X = B\).
scalar_t | One of float, double, std::complex<float>, std::complex<double>. |
[in,out] | A | On entry, the n-by-n band matrix \(A\) to be factored. Tiles outside the bandwidth do not need to exist. For tiles that are partially outside the bandwidth, data outside the bandwidth should be explicitly set to zero. On exit, the factors \(L\) and \(U\) from the factorization \(A = L U\); the unit diagonal elements of \(L\) are not stored. The upper bandwidth is increased to accommodate fill-in of \(U\). |
[out] | pivots | The pivot indices that define the permutation matrix \(P\). |
[in,out] | B | On entry, the n-by-nrhs right hand side matrix \(B\). On exit, if return value = 0, the n-by-nrhs solution matrix \(X\). |
[in] | opts | Additional options, as map of name = value pairs. Possible options:
|