SLATE 2024.05.31
Software for Linear Algebra Targeting Exascale
Loading...
Searching...
No Matches
tbsm: Triangular solve band matrix

\(C = A^{-1} B\) or \(C = B A^{-1}\) where \(A\) is band triangular More...

Functions

template<typename scalar_t >
void slate::tbsm (Side side, scalar_t alpha, TriangularBandMatrix< scalar_t > &A, Matrix< scalar_t > &B, Options const &opts)
 Distributed parallel triangular band matrix-matrix solve.
 
template<typename scalar_t >
void slate::tbsm (blas::Side side, scalar_t alpha, TriangularBandMatrix< scalar_t > &A, Pivots &pivots, Matrix< scalar_t > &B, Options const &opts)
 Distributed parallel triangular band matrix-matrix solve.
 

Detailed Description

\(C = A^{-1} B\) or \(C = B A^{-1}\) where \(A\) is band triangular

Function Documentation

◆ tbsm() [1/2]

template<typename scalar_t >
void slate::tbsm ( blas::Side  side,
scalar_t  alpha,
TriangularBandMatrix< scalar_t > &  A,
Pivots &  pivots,
Matrix< scalar_t > &  B,
Options const &  opts 
)

Distributed parallel triangular band matrix-matrix solve.

Solves one of the triangular matrix equations

\[ A X = \alpha B, \]

or

\[ X A = \alpha B, \]

where alpha is a scalar, B is an m-by-n matrix and A is a unit or non-unit, upper or lower triangular band matrix. The matrix X overwrites B. Pivoting from tbtrf is applied during the solve. The matrices can be transposed or conjugate-transposed beforehand, e.g.,

auto AT = slate::transpose( A );
slate::tbsm( Side::Left, alpha, AT, pivots, B );
Template Parameters
scalar_tOne of float, double, std::complex<float>, std::complex<double>.
Parameters
[in]sideWhether A appears on the left or on the right of X:
  • Side::Left: solve \(A X = \alpha B\)
  • Side::Right: solve \(X A = \alpha B\)
[in]alphaThe scalar alpha.
[in]A
  • If side = left, the m-by-m triangular matrix A;
  • if side = right, the n-by-n triangular matrix A.
[in]pivotsPivot information from gbtrf. If pivots is an empty vector, no pivoting is applied.
[in,out]BOn entry, the m-by-n matrix B. On exit, overwritten by the result X.
[in]optsAdditional options, as map of name = value pairs. Possible options:
  • Option::Lookahead: Number of panels to overlap with matrix updates. lookahead >= 0. Default 1.
  • Option::Target: Implementation to target. Possible values:
    • HostTask: OpenMP tasks on CPU host [default].
    • HostNest: nested OpenMP parallel for loop on CPU host.
    • HostBatch: batched BLAS on CPU host.
    • Devices: batched BLAS on GPU device.

◆ tbsm() [2/2]

template<typename scalar_t >
void slate::tbsm ( Side  side,
scalar_t  alpha,
TriangularBandMatrix< scalar_t > &  A,
Matrix< scalar_t > &  B,
Options const &  opts 
)

Distributed parallel triangular band matrix-matrix solve.

Solves one of the triangular matrix equations

\[ A X = \alpha B, \]

or

\[ X A = \alpha B, \]

where alpha is a scalar, B is an m-by-n matrix and A is a unit or non-unit, upper or lower triangular band matrix. The matrix X overwrites B. Pivoting from tbtrf is applied during the solve. The matrices can be transposed or conjugate-transposed beforehand, e.g.,

auto AT = slate::transpose( A );
slate::tbsm( Side::Left, alpha, AT, B );
Template Parameters
scalar_tOne of float, double, std::complex<float>, std::complex<double>.
Parameters
[in]sideWhether A appears on the left or on the right of X:
  • Side::Left: solve \(A X = \alpha B\)
  • Side::Right: solve \(X A = \alpha B\)
[in]alphaThe scalar alpha.
[in]A
  • If side = left, the m-by-m triangular matrix A;
  • if side = right, the n-by-n triangular matrix A.
[in,out]BOn entry, the m-by-n matrix B. On exit, overwritten by the result X.
[in]optsAdditional options, as map of name = value pairs. Possible options:
  • Option::Lookahead: Number of panels to overlap with matrix updates. lookahead >= 0. Default 1.
  • Option::Target: Implementation to target. Possible values:
    • HostTask: OpenMP tasks on CPU host [default].
    • HostNest: nested OpenMP parallel for loop on CPU host.
    • HostBatch: batched BLAS on CPU host.
    • Devices: batched BLAS on GPU device.