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

\(Ax = \lambda B x\), etc. More...

Functions

template<typename scalar_t >
void slate::hegv (int64_t itype, HermitianMatrix< scalar_t > &A, HermitianMatrix< scalar_t > &B, std::vector< blas::real_type< scalar_t > > &Lambda, Matrix< scalar_t > &Z, Options const &opts)
 Distributed parallel computation of all the eigenvalues, and optionally, the eigenvectors of a complex generalized Hermitian-definite eigenproblem, of the form:
 

Detailed Description

\(Ax = \lambda B x\), etc.

Function Documentation

◆ hegv()

template<typename scalar_t >
void slate::hegv ( int64_t  itype,
HermitianMatrix< scalar_t > &  A,
HermitianMatrix< scalar_t > &  B,
std::vector< blas::real_type< scalar_t > > &  Lambda,
Matrix< scalar_t > &  Z,
Options const &  opts 
)

Distributed parallel computation of all the eigenvalues, and optionally, the eigenvectors of a complex generalized Hermitian-definite eigenproblem, of the form:

itype Problem
itype = 1 \(A z = \lambda B z\)
itype = 2 \(A B z = \lambda z\)
itype = 3 \(B A z = \lambda z\)

Here A and B are assumed to be Hermitian and B is also positive definite.

Template Parameters
scalar_tOne of float, double, std::complex<float>, std::complex<double>.
Parameters
[in]itype
  • itype = 1: Compute \(A z = \lambda B z\);
  • itype = 2: Compute \(A B z = \lambda z\);
  • itype = 3: Compute \(B A z = \lambda z\).
[in,out]AOn entry, the n-by-n Hermitian matrix \(A\). On exit, contents are destroyed.
[in,out]BOn entry, the n-by-n Hermitian positive definite matrix \(B\). On exit, B is overwritten by the triangular factor U or L from the Cholesky factorization \(B = U^H U\) or \(B = L L^H\).
[out]LambdaThe vector Lambda of length n. If successful, the eigenvalues in ascending order.
[out]ZOn entry, if Z is empty, does not compute eigenvectors. On exit, orthonormal eigenvectors of the 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::InnerBlocking: Inner blocking to use for panel. Default 16.
  • Option::MaxPanelThreads: Number of threads to use for panel. Default omp_get_max_threads()/2.
  • 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