SLATE 2024.05.31
Software for Linear Algebra Targeting Exascale
Loading...
Searching...
No Matches
Computational

Inverse, multiply. More...

Functions

template<typename scalar_t >
void slate::trtri (TriangularMatrix< scalar_t > &A, Options const &opts)
 Distributed parallel inverse of a triangular matrix.
 
template<typename scalar_t >
void slate::trtrm (TriangularMatrix< scalar_t > &A, Options const &opts)
 Distributed parallel inverse of a triangular matrix.
 

Detailed Description

Inverse, multiply.

Function Documentation

◆ trtri()

template<typename scalar_t >
void slate::trtri ( TriangularMatrix< scalar_t > &  A,
Options const &  opts 
)

Distributed parallel inverse of a triangular matrix.

Computes the inverse of an upper or lower triangular matrix \(A\).

Complexity (in real): \(\approx \frac{1}{3} n^{3}\) flops.

Template Parameters
scalar_tOne of float, double, std::complex<float>, std::complex<double>.
Parameters
[in,out]AOn entry, the n-by-n triangular matrix \(A\). On exit, if return value = 0, the (triangular) inverse of the original matrix \(A\).
[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.

TODO: return value

Return values
0successful exit
>0for return value = \(i\), \(A(i,i)\) is exactly zero. The triangular matrix is singular and its inverse can not be computed.

◆ trtrm()

template<typename scalar_t >
void slate::trtrm ( TriangularMatrix< scalar_t > &  A,
Options const &  opts 
)

Distributed parallel inverse of a triangular matrix.

Computes the inverse of an upper or lower triangular matrix \(A\).

Template Parameters
scalar_tOne of float, double, std::complex<float>, std::complex<double>.
Parameters
[in,out]AOn entry, the n-by-n triangular matrix \(A\). On exit, if return value = 0, the (triangular) inverse of the original matrix \(A\).
[in]optsAdditional options, as map of name = value pairs. Possible options:
  • 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.

TODO: return value

Return values
0successful exit