LAPACK++
LAPACK C++ API
|
Functions | |
int64_t | lapack::hetrd (lapack::Uplo uplo, int64_t n, std::complex< double > *A, int64_t lda, double *D, double *E, std::complex< double > *tau) |
Reduces a Hermitian matrix A to real symmetric tridiagonal form T by a unitary similarity transformation: \(Q^H A Q = T\). More... | |
int64_t | lapack::hetrd (lapack::Uplo uplo, int64_t n, std::complex< float > *A, int64_t lda, float *D, float *E, std::complex< float > *tau) |
int64_t | lapack::hetrd_2stage (lapack::Job jobz, lapack::Uplo uplo, int64_t n, std::complex< double > *A, int64_t lda, double *D, double *E, std::complex< double > *tau, std::complex< double > *hous2, int64_t lhous2) |
Reduces a Hermitian matrix A to real symmetric tridiagonal form T by a unitary similarity transformation: \(Q1^H Q2^H A Q2 Q1 = T\). More... | |
int64_t | lapack::hetrd_2stage (lapack::Job jobz, lapack::Uplo uplo, int64_t n, std::complex< float > *A, int64_t lda, float *D, float *E, std::complex< float > *tau, std::complex< float > *hous2, int64_t lhous2) |
int64_t | lapack::orgtr (lapack::Uplo uplo, int64_t n, double *A, int64_t lda, double const *tau) |
int64_t | lapack::orgtr (lapack::Uplo uplo, int64_t n, float *A, int64_t lda, float const *tau) |
int64_t | lapack::ormtr (lapack::Side side, lapack::Uplo uplo, lapack::Op trans, int64_t m, int64_t n, double const *A, int64_t lda, double const *tau, double *C, int64_t ldc) |
int64_t | lapack::ormtr (lapack::Side side, lapack::Uplo uplo, lapack::Op trans, int64_t m, int64_t n, float const *A, int64_t lda, float const *tau, float *C, int64_t ldc) |
int64_t | lapack::sytrd (lapack::Uplo uplo, int64_t n, double *A, int64_t lda, double *D, double *E, double *tau) |
int64_t | lapack::sytrd (lapack::Uplo uplo, int64_t n, float *A, int64_t lda, float *D, float *E, float *tau) |
int64_t | lapack::sytrd_2stage (lapack::Job jobz, lapack::Uplo uplo, int64_t n, double *A, int64_t lda, double *D, double *E, double *tau, double *hous2, int64_t lhous2) |
int64_t | lapack::sytrd_2stage (lapack::Job jobz, lapack::Uplo uplo, int64_t n, float *A, int64_t lda, float *D, float *E, float *tau, float *hous2, int64_t lhous2) |
int64_t | lapack::ungtr (lapack::Uplo uplo, int64_t n, std::complex< double > *A, int64_t lda, std::complex< double > const *tau) |
Generates an n-by-n unitary matrix Q which is defined as the product of n-1 elementary reflectors of order n, as returned by lapack::hetrd : More... | |
int64_t | lapack::ungtr (lapack::Uplo uplo, int64_t n, std::complex< float > *A, int64_t lda, std::complex< float > const *tau) |
int64_t | lapack::unmtr (lapack::Side side, lapack::Uplo uplo, lapack::Op trans, int64_t m, int64_t n, std::complex< double > const *A, int64_t lda, std::complex< double > const *tau, std::complex< double > *C, int64_t ldc) |
Multiplies the general m-by-n matrix C by Q from lapack::hetrd as follows: More... | |
int64_t | lapack::unmtr (lapack::Side side, lapack::Uplo uplo, lapack::Op trans, int64_t m, int64_t n, std::complex< float > const *A, int64_t lda, std::complex< float > const *tau, std::complex< float > *C, int64_t ldc) |
int64_t | lapack::upmtr (lapack::Side side, lapack::Uplo uplo, lapack::Op trans, int64_t m, int64_t n, std::complex< double > const *AP, std::complex< double > const *tau, std::complex< double > *C, int64_t ldc) |
Multiplies the general m-by-n matrix C by Q from lapack::hptrd as follows: More... | |
int64_t | lapack::upmtr (lapack::Side side, lapack::Uplo uplo, lapack::Op trans, int64_t m, int64_t n, std::complex< float > const *AP, std::complex< float > const *tau, std::complex< float > *C, int64_t ldc) |
int64_t lapack::hetrd | ( | lapack::Uplo | uplo, |
int64_t | n, | ||
std::complex< double > * | A, | ||
int64_t | lda, | ||
double * | D, | ||
double * | E, | ||
std::complex< double > * | tau | ||
) |
Reduces a Hermitian matrix A to real symmetric tridiagonal form T by a unitary similarity transformation: \(Q^H A Q = T\).
Overloaded versions are available for float
, double
, std::complex<float>
, and std::complex<double>
. For real matrices, this is an alias for lapack::sytrd
.
[in] | uplo |
|
[in] | n | The order of the matrix A. n >= 0. |
[in,out] | A | The n-by-n matrix A, stored in an lda-by-n array. On entry, the Hermitian matrix A.
|
[in] | lda | The leading dimension of the array A. lda >= max(1,n). |
[out] | D | The vector D of length n. The diagonal elements of the tridiagonal matrix T: D(i) = A(i,i). |
[out] | E | The vector E of length n-1. The off-diagonal elements of the tridiagonal matrix T: E(i) = A(i,i+1) if uplo = Upper, E(i) = A(i+1,i) if uplo = Lower. |
[out] | tau | The vector tau of length n-1. The scalar factors of the elementary reflectors (see Further Details). |
If uplo = Upper, the matrix Q is represented as a product of elementary reflectors
\[ Q = H(n-1) . . . H(2) H(1). \]
Each H(i) has the form
\[ H(i) = I - \tau v v^H \]
where \(\tau\) is a scalar, and v is a vector with v(i+1:n) = 0 and v(i) = 1; v(1:i-1) is stored on exit in A(1:i-1,i+1), and \(\tau\) in tau(i).
If uplo = Lower, the matrix Q is represented as a product of elementary reflectors
\[ Q = H(1) H(2) . . . H(n-1). \]
Each H(i) has the form
\[ H(i) = I - \tau v v^H \]
where \(\tau\) is a scalar, and v is a vector with v(1:i) = 0 and v(i+1) = 1; v(i+2:n) is stored on exit in A(i+2:n,i), and \(\tau\) in tau(i).
The contents of A on exit are illustrated by the following examples with n = 5:
if uplo = Upper: if uplo = Lower: ( d e v2 v3 v4 ) ( d ) ( d e v3 v4 ) ( e d ) ( d e v4 ) ( v1 e d ) ( d e ) ( v1 v2 e d ) ( d ) ( v1 v2 v3 e d )
where d and e denote diagonal and off-diagonal elements of T, and vi denotes an element of the vector defining H(i).
int64_t lapack::hetrd_2stage | ( | lapack::Job | jobz, |
lapack::Uplo | uplo, | ||
int64_t | n, | ||
std::complex< double > * | A, | ||
int64_t | lda, | ||
double * | D, | ||
double * | E, | ||
std::complex< double > * | tau, | ||
std::complex< double > * | hous2, | ||
int64_t | lhous2 | ||
) |
Reduces a Hermitian matrix A to real symmetric tridiagonal form T by a unitary similarity transformation: \(Q1^H Q2^H A Q2 Q1 = T\).
Overloaded versions are available for float
, double
, std::complex<float>
, and std::complex<double>
. For real matrices, this is an alias for lapack::sytrd_2stage
.
[in] | jobz |
|
[in] | uplo |
|
[in] | n | The order of the matrix A. n >= 0. |
[in,out] | A | The n-by-n matrix A, stored in an lda-by-n array. On entry, the Hermitian matrix A.
|
[in] | lda | The leading dimension of the array A. lda >= max(1,n). |
[out] | D | The vector D of length n. The diagonal elements of the tridiagonal matrix T. |
[out] | E | The vector E of length n-1. The off-diagonal elements of the tridiagonal matrix T. |
[out] | tau | The vector tau of length n-kd. The scalar factors of the elementary reflectors of the first stage (see Further Details). |
[out] | hous2 | The vector hous2 of length lhous2. Stores the Householder representation of the stage2 band to tridiagonal. |
[in] | lhous2 | The dimension of the array hous2.
|
Implemented by Azzam Haidar.
All details are available on technical report, SC11, SC13 papers.
Azzam Haidar, Hatem Ltaief, and Jack Dongarra. Parallel reduction to condensed forms for symmetric eigenvalue problems using aggregated fine-grained and memory-aware kernels. In Proceedings of 2011 International Conference for High Performance Computing, Networking, Storage and Analysis (SC '11), New York, NY, USA, Article 8, 11 pages. http://doi.acm.org/10.1145/2063384.2063394
A. Haidar, J. Kurzak, P. Luszczek, 2013. An improved parallel singular value algorithm and its implementation for multicore hardware, In Proceedings of 2013 International Conference for High Performance Computing, Networking, Storage and Analysis (SC '13). Denver, Colorado, USA, 2013. Article 90, 12 pages. http://doi.acm.org/10.1145/2503210.2503292
A. Haidar, R. Solca, S. Tomov, T. Schulthess and J. Dongarra. A novel hybrid CPU-GPU generalized eigensolver for electronic structure calculations based on fine-grained memory aware tasks. International Journal of High Performance Computing Applications. Volume 28 Issue 2, Pages 196-209, May 2014. http://hpc.sagepub.com/content/28/2/196
int64_t lapack::orgtr | ( | lapack::Uplo | uplo, |
int64_t | n, | ||
double * | A, | ||
int64_t | lda, | ||
double const * | tau | ||
) |
int64_t lapack::ormtr | ( | lapack::Side | side, |
lapack::Uplo | uplo, | ||
lapack::Op | trans, | ||
int64_t | m, | ||
int64_t | n, | ||
double const * | A, | ||
int64_t | lda, | ||
double const * | tau, | ||
double * | C, | ||
int64_t | ldc | ||
) |
int64_t lapack::sytrd | ( | lapack::Uplo | uplo, |
int64_t | n, | ||
double * | A, | ||
int64_t | lda, | ||
double * | D, | ||
double * | E, | ||
double * | tau | ||
) |
int64_t lapack::sytrd_2stage | ( | lapack::Job | jobz, |
lapack::Uplo | uplo, | ||
int64_t | n, | ||
double * | A, | ||
int64_t | lda, | ||
double * | D, | ||
double * | E, | ||
double * | tau, | ||
double * | hous2, | ||
int64_t | lhous2 | ||
) |
int64_t lapack::ungtr | ( | lapack::Uplo | uplo, |
int64_t | n, | ||
std::complex< double > * | A, | ||
int64_t | lda, | ||
std::complex< double > const * | tau | ||
) |
Generates an n-by-n unitary matrix Q which is defined as the product of n-1 elementary reflectors of order n, as returned by lapack::hetrd
:
Overloaded versions are available for float
, double
, std::complex<float>
, and std::complex<double>
. For real matrices, this is an alias for lapack::orgtr
.
[in] | uplo |
|
[in] | n | The order of the matrix Q. n >= 0. |
[in,out] | A | The n-by-n matrix A, stored in an lda-by-n array. On entry, the vectors which define the elementary reflectors, as returned by lapack::hetrd . On exit, the n-by-n unitary matrix Q. |
[in] | lda | The leading dimension of the array A. lda >= n. |
[in] | tau | The vector tau of length n-1. tau(i) must contain the scalar factor of the elementary reflector H(i), as returned by lapack::hetrd . |
int64_t lapack::unmtr | ( | lapack::Side | side, |
lapack::Uplo | uplo, | ||
lapack::Op | trans, | ||
int64_t | m, | ||
int64_t | n, | ||
std::complex< double > const * | A, | ||
int64_t | lda, | ||
std::complex< double > const * | tau, | ||
std::complex< double > * | C, | ||
int64_t | ldc | ||
) |
Multiplies the general m-by-n matrix C by Q from lapack::hetrd
as follows:
where Q is a unitary matrix of order nq, with nq = m if side = Left and nq = n if side = Right. Q is defined as the product of nq-1 elementary reflectors, as returned by lapack::hetrd
:
Overloaded versions are available for float
, double
, std::complex<float>
, and std::complex<double>
. For real matrices, this is an alias for lapack::ormtr
.
[in] | side |
|
[in] | uplo |
|
[in] | trans |
|
[in] | m | The number of rows of the matrix C. m >= 0. |
[in] | n | The number of columns of the matrix C. n >= 0. |
[in] | A | The vectors which define the elementary reflectors, as returned by lapack::hetrd .
|
[in] | lda | The leading dimension of the array A.
|
[in] | tau | tau(i) must contain the scalar factor of the elementary reflector H(i), as returned by lapack::hetrd .
|
[in,out] | C | The m-by-n matrix C, stored in an ldc-by-n array. On entry, the m-by-n matrix C. On exit, C is overwritten by \(Q C\) or \(Q^H C\) or \(C Q^H\) or \(C Q\). |
[in] | ldc | The leading dimension of the array C. ldc >= max(1,m). |
int64_t lapack::upmtr | ( | lapack::Side | side, |
lapack::Uplo | uplo, | ||
lapack::Op | trans, | ||
int64_t | m, | ||
int64_t | n, | ||
std::complex< double > const * | AP, | ||
std::complex< double > const * | tau, | ||
std::complex< double > * | C, | ||
int64_t | ldc | ||
) |
Multiplies the general m-by-n matrix C by Q from lapack::hptrd
as follows:
where Q is a unitary matrix of order nq, with nq = m if side = Left and nq = n if side = Right. Q is defined as the product of nq-1 elementary reflectors, as returned by lapack::hptrd
using packed storage:
Overloaded versions are available for float
, double
, std::complex<float>
, and std::complex<double>
. For real matrices, this is an alias for lapack::opmtr
.
[in] | side |
|
[in] | uplo |
|
[in] | trans |
|
[in] | m | The number of rows of the matrix C. m >= 0. |
[in] | n | The number of columns of the matrix C. n >= 0. |
[in] | AP | The vectors which define the elementary reflectors, as returned by lapack::hptrd . AP is modified by the routine but restored on exit.
|
[in] | tau | tau(i) must contain the scalar factor of the elementary reflector H(i), as returned by lapack::hptrd .
|
[in,out] | C | The m-by-n matrix C, stored in an ldc-by-n array. On entry, the m-by-n matrix C. On exit, C is overwritten by \(Q C\) or \(Q^H C\) or \(C Q^H\) or \(C Q\). |
[in] | ldc | The leading dimension of the array C. ldc >= max(1,m). |