LAPACK++
LAPACK C++ API
|
Functions | |
int64_t | lapack::gebak (lapack::Balance balance, lapack::Side side, int64_t n, int64_t ilo, int64_t ihi, double const *scale, int64_t m, std::complex< double > *V, int64_t ldv) |
Forms the right or left eigenvectors of a complex general matrix by backward transformation on the computed eigenvectors of the balanced matrix output by lapack::gebal . More... | |
int64_t | lapack::gebak (lapack::Balance balance, lapack::Side side, int64_t n, int64_t ilo, int64_t ihi, double const *scale, int64_t m, double *V, int64_t ldv) |
int64_t | lapack::gebak (lapack::Balance balance, lapack::Side side, int64_t n, int64_t ilo, int64_t ihi, float const *scale, int64_t m, float *V, int64_t ldv) |
int64_t | lapack::gebak (lapack::Balance balance, lapack::Side side, int64_t n, int64_t ilo, int64_t ihi, float const *scale, int64_t m, std::complex< float > *V, int64_t ldv) |
int64_t | lapack::gebal (lapack::Balance balance, int64_t n, double *A, int64_t lda, int64_t *ilo, int64_t *ihi, double *scale) |
int64_t | lapack::gebal (lapack::Balance balance, int64_t n, float *A, int64_t lda, int64_t *ilo, int64_t *ihi, float *scale) |
int64_t | lapack::gebal (lapack::Balance balance, int64_t n, std::complex< double > *A, int64_t lda, int64_t *ilo, int64_t *ihi, double *scale) |
Balances a general complex matrix A. More... | |
int64_t | lapack::gebal (lapack::Balance balance, int64_t n, std::complex< float > *A, int64_t lda, int64_t *ilo, int64_t *ihi, float *scale) |
int64_t | lapack::gehrd (int64_t n, int64_t ilo, int64_t ihi, double *A, int64_t lda, double *tau) |
int64_t | lapack::gehrd (int64_t n, int64_t ilo, int64_t ihi, float *A, int64_t lda, float *tau) |
int64_t | lapack::gehrd (int64_t n, int64_t ilo, int64_t ihi, std::complex< double > *A, int64_t lda, std::complex< double > *tau) |
Reduces a general matrix A to upper Hessenberg form H by an unitary similarity transformation: \(Q^H A Q = H\). More... | |
int64_t | lapack::gehrd (int64_t n, int64_t ilo, int64_t ihi, std::complex< float > *A, int64_t lda, std::complex< float > *tau) |
int64_t | lapack::hseqr (lapack::JobSchur jobschur, lapack::Job compz, int64_t n, int64_t ilo, int64_t ihi, double *H, int64_t ldh, std::complex< double > *W, double *Z, int64_t ldz) |
int64_t | lapack::hseqr (lapack::JobSchur jobschur, lapack::Job compz, int64_t n, int64_t ilo, int64_t ihi, float *H, int64_t ldh, std::complex< float > *W, float *Z, int64_t ldz) |
int64_t | lapack::hseqr (lapack::JobSchur jobschur, lapack::Job compz, int64_t n, int64_t ilo, int64_t ihi, std::complex< double > *H, int64_t ldh, std::complex< double > *W, std::complex< double > *Z, int64_t ldz) |
Computes the eigenvalues of a Hessenberg matrix H and, optionally, the matrices T and Z from the Schur decomposition. More... | |
int64_t | lapack::hseqr (lapack::JobSchur jobschur, lapack::Job compz, int64_t n, int64_t ilo, int64_t ihi, std::complex< float > *H, int64_t ldh, std::complex< float > *W, std::complex< float > *Z, int64_t ldz) |
int64_t | lapack::trevc (lapack::Sides side, lapack::HowMany howmany, bool *select, int64_t n, double const *T, int64_t ldt, double *VL, int64_t ldvl, double *VR, int64_t ldvr, int64_t mm, int64_t *m) |
int64_t | lapack::trevc (lapack::Sides side, lapack::HowMany howmany, bool *select, int64_t n, float const *T, int64_t ldt, float *VL, int64_t ldvl, float *VR, int64_t ldvr, int64_t mm, int64_t *m) |
int64_t | lapack::trevc (lapack::Sides side, lapack::HowMany howmany, bool const *select, int64_t n, std::complex< float > *T, int64_t ldt, std::complex< float > *VL, int64_t ldvl, std::complex< float > *VR, int64_t ldvr, int64_t mm, int64_t *m) |
int64_t | lapack::trevc (lapack::Sides sides, lapack::HowMany howmany, bool const *select, int64_t n, std::complex< double > *T, int64_t ldt, std::complex< double > *VL, int64_t ldvl, std::complex< double > *VR, int64_t ldvr, int64_t mm, int64_t *m) |
Computes some or all of the right and/or left eigenvectors of a complex upper triangular matrix T. More... | |
int64_t | lapack::trevc3 (lapack::Sides side, lapack::HowMany howmany, bool *select, int64_t n, double const *T, int64_t ldt, double *VL, int64_t ldvl, double *VR, int64_t ldvr, int64_t mm, int64_t *m) |
int64_t | lapack::trevc3 (lapack::Sides side, lapack::HowMany howmany, bool *select, int64_t n, float const *T, int64_t ldt, float *VL, int64_t ldvl, float *VR, int64_t ldvr, int64_t mm, int64_t *m) |
int64_t | lapack::trevc3 (lapack::Sides side, lapack::HowMany howmany, bool const *select, int64_t n, std::complex< float > *T, int64_t ldt, std::complex< float > *VL, int64_t ldvl, std::complex< float > *VR, int64_t ldvr, int64_t mm, int64_t *m) |
int64_t | lapack::trevc3 (lapack::Sides sides, lapack::HowMany howmany, bool const *select, int64_t n, std::complex< double > *T, int64_t ldt, std::complex< double > *VL, int64_t ldvl, std::complex< double > *VR, int64_t ldvr, int64_t mm, int64_t *m) |
Computes some or all of the right and/or left eigenvectors of a complex upper triangular matrix T. More... | |
int64_t | lapack::unghr (int64_t n, int64_t ilo, int64_t ihi, 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 ihi-ilo elementary reflectors of order n, as returned by lapack::gehrd : More... | |
int64_t | lapack::unghr (int64_t n, int64_t ilo, int64_t ihi, std::complex< float > *A, int64_t lda, std::complex< float > const *tau) |
int64_t | lapack::unmhr (lapack::Side side, lapack::Op trans, int64_t m, int64_t n, int64_t ilo, int64_t ihi, 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::gehrd as follows: More... | |
int64_t | lapack::unmhr (lapack::Side side, lapack::Op trans, int64_t m, int64_t n, int64_t ilo, int64_t ihi, std::complex< float > const *A, int64_t lda, std::complex< float > const *tau, std::complex< float > *C, int64_t ldc) |
int64_t lapack::gebak | ( | lapack::Balance | balance, |
lapack::Side | side, | ||
int64_t | n, | ||
int64_t | ilo, | ||
int64_t | ihi, | ||
double const * | scale, | ||
int64_t | m, | ||
std::complex< double > * | V, | ||
int64_t | ldv | ||
) |
Forms the right or left eigenvectors of a complex general matrix by backward transformation on the computed eigenvectors of the balanced matrix output by lapack::gebal
.
Overloaded versions are available for float
, double
, std::complex<float>
, and std::complex<double>
.
[in] | balance | Specifies the type of backward transformation required:
|
[in] | side |
|
[in] | n | The number of rows of the matrix V. n >= 0. |
[in] | ilo | |
[in] | ihi | The integers ilo and ihi determined by lapack::gebal .
|
[in] | scale | The vector scale of length n. Details of the permutation and scaling factors, as returned by lapack::gebal . |
[in] | m | The number of columns of the matrix V. m >= 0. |
[in,out] | V | The n-by-m matrix V, stored in an ldv-by-m array. On entry, the matrix of right or left eigenvectors to be transformed, as returned by lapack::hsein or lapack::trevc . On exit, V is overwritten by the transformed eigenvectors. |
[in] | ldv | The leading dimension of the array V. ldv >= max(1,n). |
int64_t lapack::gebal | ( | lapack::Balance | balance, |
int64_t | n, | ||
std::complex< double > * | A, | ||
int64_t | lda, | ||
int64_t * | ilo, | ||
int64_t * | ihi, | ||
double * | scale | ||
) |
Balances a general complex matrix A.
This involves, first, permuting A by a similarity transformation to isolate eigenvalues in the first 1 to ilo-1 and last ihi+1 to n elements on the diagonal; and second, applying a diagonal similarity transformation to rows and columns ilo to ihi to make the rows and columns as close in norm as possible. Both steps are optional.
Balancing may reduce the 1-norm of the matrix, and improve the accuracy of the computed eigenvalues and/or eigenvectors.
Overloaded versions are available for float
, double
, std::complex<float>
, and std::complex<double>
.
[in] | balance | Specifies the operations to be performed on A:
|
[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 input matrix A. On exit, A is overwritten by the balanced matrix. If balance = None, A is not referenced. See Further Details. |
[in] | lda | The leading dimension of the array A. lda >= max(1,n). |
[out] | ilo | |
[out] | ihi | ilo and ihi are set to integers such that on exit A(i,j) = 0 if i > j and j = 1, ..., ilo-1 or i = ihi+1, ..., n. If balance = None or Scale, ilo = 1 and ihi = n. |
[out] | scale | The vector scale of length n. Details of the permutations and scaling factors applied to A. If P(j) is the index of the row and column interchanged with row and column j and D(j) is the scaling factor applied to row and column j, then: scale(j) = P(j) for j = 1, ..., ilo-1; scale(j) = D(j) for j = ilo, ..., ihi; scale(j) = P(j) for j = ihi+1, ..., n. The order in which the interchanges are made is n to ihi+1, then 1 to ilo-1. |
The permutations consist of row and column interchanges which put the matrix in the form
\[ P A P = \begin{bmatrix} T1 & X & Y \\ 0 & B & Z \\ 0 & 0 & T2 \end{bmatrix}, \]
where T1 and T2 are upper triangular matrices whose eigenvalues lie along the diagonal. The column indices ilo and ihi mark the starting and ending columns of the submatrix B. Balancing consists of applying a diagonal similarity transformation \(D^{-1} B D\) to make the 1-norms of each row of B and its corresponding column nearly equal. The output matrix is
\[ \begin{bmatrix} T1 & X D & Y \\ 0 & D^{-1} B D & D^{-1} Z \\ 0 & 0 & T2 \end{bmatrix}. \]
Information about the permutations P and the diagonal matrix D is returned in the vector scale.
This subroutine is based on the EISPACK routine CBAL.
Modified by Tzu-Yi Chen, Computer Science Division, University of California at Berkeley, USA
int64_t lapack::gehrd | ( | int64_t | n, |
int64_t | ilo, | ||
int64_t | ihi, | ||
std::complex< double > * | A, | ||
int64_t | lda, | ||
std::complex< double > * | tau | ||
) |
Reduces a general matrix A to upper Hessenberg form H by an unitary similarity transformation: \(Q^H A Q = H\).
Overloaded versions are available for float
, double
, std::complex<float>
, and std::complex<double>
.
[in] | n | The order of the matrix A. n >= 0. |
[in] | ilo | |
[in] | ihi | It is assumed that A is already upper triangular in rows and columns 1:ilo-1 and ihi+1:n. ilo and ihi are normally set by a previous call to lapack::gebal ; otherwise they should be set to 1 and n respectively. See Further Details.
|
[in,out] | A | The n-by-n matrix A, stored in an lda-by-n array. On entry, the n-by-n general matrix to be reduced. On exit, the upper triangle and the first subdiagonal of A are overwritten with the upper Hessenberg matrix H, and the elements below the first subdiagonal, with the array tau, represent the unitary matrix Q as a product of elementary reflectors. See Further Details. |
[in] | lda | The leading dimension of the array A. lda >= max(1,n). |
[out] | tau | The vector tau of length n-1. The scalar factors of the elementary reflectors (see Further Details). Elements 1:ilo-1 and ihi:n-1 of tau are set to zero. |
The matrix Q is represented as a product of (ihi-ilo) elementary reflectors
\[ Q = H(ilo) H(ilo+1) . . . H(ihi-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, v(i+1) = 1 and v(ihi+1:n) = 0; v(i+2:ihi) is stored on exit in A(i+2:ihi,i), and \(\tau\) in tau(i).
The contents of A are illustrated by the following example, with n = 7, ilo = 2 and ihi = 6:
on entry, on exit, ( a a a a a a a ) ( a a h h h h a ) ( a a a a a a ) ( a h h h h a ) ( a a a a a a ) ( h h h h h h ) ( a a a a a a ) ( v2 h h h h h ) ( a a a a a a ) ( v2 v3 h h h h ) ( a a a a a a ) ( v2 v3 v4 h h h ) ( a ) ( a )
where a denotes an element of the original matrix A, h denotes a modified element of the upper Hessenberg matrix H, and vi denotes an element of the vector defining H(i).
This routine is a slight modification of the LAPACK 3.0 gehrd
subroutine incorporating improvements proposed by Quintana-Orti and Van de Geijn (2006). (See lapack::lahr2
.)
int64_t lapack::hseqr | ( | lapack::JobSchur | jobschur, |
lapack::Job | compz, | ||
int64_t | n, | ||
int64_t | ilo, | ||
int64_t | ihi, | ||
std::complex< double > * | H, | ||
int64_t | ldh, | ||
std::complex< double > * | W, | ||
std::complex< double > * | Z, | ||
int64_t | ldz | ||
) |
Computes the eigenvalues of a Hessenberg matrix H and, optionally, the matrices T and Z from the Schur decomposition.
\[ H = Z T Z^H, \]
where T is an upper triangular matrix (the Schur form), and Z is the unitary matrix of Schur vectors.
Optionally Z may be postmultiplied into an input unitary matrix Q so that this routine can give the Schur factorization of a matrix A which has been reduced to the Hessenberg form H by the unitary matrix Q: \(A = Q H Q^H = (QZ) T (QZ)^H\).
Overloaded versions are available for float
, double
, std::complex<float>
, and std::complex<double>
.
[in] | jobschur |
|
[in] | compz |
|
[in] | n | The order of the matrix H. n >= 0. |
[in] | ilo | |
[in] | ihi | It is assumed that H is already upper triangular in rows and columns 1:ilo-1 and ihi+1:n. ilo and ihi are normally set by a previous call to lapack::gebal , and then passed to lapack::gehrd when the matrix output by lapack::gebal is reduced to Hessenberg form. Otherwise ilo and ihi should be set to 1 and n respectively.
|
[in,out] | H | The n-by-n matrix H, stored in an ldh-by-n array. On entry, the upper Hessenberg matrix H. On exit, if successful and job = Schur, H contains the upper triangular matrix T from the Schur decomposition (the Schur form). If successful and job = None, the contents of H are unspecified on exit. (The output value of H when return value > 0 is given under the description of info below.) Unlike earlier versions of hseqr , this subroutine may explicitly H(i,j) = 0 for i > j and j = 1, 2, ... ilo-1 or j = ihi+1, ihi+2, ... n. |
[in] | ldh | The leading dimension of the array H. ldh >= max(1,n). |
[out] | W | The vector W of length n. The computed eigenvalues. If job = Schur, the eigenvalues are stored in the same order as on the diagonal of the Schur form returned in H, with W(i) = H(i,i). Note: In LAPACK++, W is always complex, whereas LAPACK with a real matrix H uses a split-complex representation (WR, WI) for W. |
[in,out] | Z | The n-by-n matrix Z, stored in an ldz-by-n array.
|
[in] | ldz | The leading dimension of the array Z. if compz = Vec or compz = Update, then ldz >= max(1,n). Otherwize, ldz >= 1. |
hseqr
failed to compute all of the eigenvalues. Elements 1:ilo-1 and i+1:n of W contain those eigenvalues which have been successfully computed. (Failures are rare.) int64_t lapack::trevc | ( | lapack::Sides | sides, |
lapack::HowMany | howmany, | ||
bool const * | select, | ||
int64_t | n, | ||
std::complex< double > * | T, | ||
int64_t | ldt, | ||
std::complex< double > * | VL, | ||
int64_t | ldvl, | ||
std::complex< double > * | VR, | ||
int64_t | ldvr, | ||
int64_t | mm, | ||
int64_t * | m | ||
) |
Computes some or all of the right and/or left eigenvectors of a complex upper triangular matrix T.
Matrices of this type are produced by the Schur factorization of a complex general matrix: \(A = Q T Q^H\), as computed by lapack::hseqr
.
The right eigenvector x and the left eigenvector y of T corresponding to an eigenvalue \(\lambda\) are defined by:
\[ T x = \lambda x, \]
\[ y^H T = \lambda y^H, \]
where \(y^H\) denotes the conjugate transpose of the vector y. The eigenvalues are not input to this routine, but are read directly from the diagonal of T.
This routine returns the matrices X and/or Y of right and left eigenvectors of T, or the products \(Q X\) and/or \(Q Y\), where Q is an input matrix. If Q is the unitary factor that reduces a matrix A to Schur form T, then \(Q X\) and \(Q Y\) are the matrices of right and left eigenvectors of A.
Overloaded versions are available for float
, double
, std::complex<float>
, and std::complex<double>
.
[in] | sides |
|
[in] | howmany |
|
[in] | select | The vector select of length n. If howmany = Select, select specifies the eigenvectors to be computed. The eigenvector corresponding to the j-th eigenvalue is computed if select(j) = true. Not referenced if howmany = All or Backtransform. TODO: updated in real case. See [sd]trevc. |
[in] | n | The order of the matrix T. n >= 0. |
[in,out] | T | The n-by-n matrix T, stored in an ldt-by-n array. The upper triangular matrix T. T is modified, but restored on exit. |
[in] | ldt | The leading dimension of the array T. ldt >= max(1,n). |
[in,out] | VL | The n-by-mm matrix VL, stored in an ldvl-by-mm array. On entry, if side = Left or Both and howmany = Backtransform, VL must contain an n-by-n matrix Q (usually the unitary matrix Q of Schur vectors returned by lapack::hseqr ). On exit, if side = Left or Both, VL contains:
|
[in] | ldvl | The leading dimension of the array VL. ldvl >= 1, and if side = Left or Both, ldvl >= n. |
[in,out] | VR | The n-by-mm matrix VR, stored in an ldvr-by-mm array. On entry, if side = Right or Both and howmany = Backtransform, VR must contain an n-by-n matrix Q (usually the unitary matrix Q of Schur vectors returned by lapack::hseqr ). On exit, if side = Right or Both, VR contains:
|
[in] | ldvr | The leading dimension of the array VR. ldvr >= 1, and if side = Right or Both, ldvr >= n. |
[in] | mm | The number of columns in the arrays VL and/or VR. mm >= m. |
[out] | m | The number of columns in the arrays VL and/or VR actually used to store the eigenvectors. If howmany = All or Backtransform, m is set to n. Each selected eigenvector occupies one column. |
The algorithm used in this program is basically backward (forward) substitution, with scaling to make the the code robust against possible overflow.
Each eigenvector is normalized so that the element of largest magnitude has magnitude 1; here the magnitude of a complex number (x,y) is taken to be |x| + |y|.
int64_t lapack::trevc3 | ( | lapack::Sides | sides, |
lapack::HowMany | howmany, | ||
bool const * | select, | ||
int64_t | n, | ||
std::complex< double > * | T, | ||
int64_t | ldt, | ||
std::complex< double > * | VL, | ||
int64_t | ldvl, | ||
std::complex< double > * | VR, | ||
int64_t | ldvr, | ||
int64_t | mm, | ||
int64_t * | m | ||
) |
Computes some or all of the right and/or left eigenvectors of a complex upper triangular matrix T.
Matrices of this type are produced by the Schur factorization of a complex general matrix: \(A = Q T Q^H\), as computed by lapack::hseqr
.
The right eigenvector x and the left eigenvector y of T corresponding to an eigenvalue \(\lambda\) are defined by:
\[ T x = \lambda x, \]
\[ y^H T = \lambda y^H, \]
where \(y^H\) denotes the conjugate transpose of the vector y. The eigenvalues are not input to this routine, but are read directly from the diagonal of T.
This routine returns the matrices X and/or Y of right and left eigenvectors of T, or the products \(Q X\) and/or \(Q Y\), where Q is an input matrix. If Q is the unitary factor that reduces a matrix A to Schur form T, then \(Q X\) and \(Q Y\) are the matrices of right and left eigenvectors of A.
This uses a Level 3 BLAS version of the back transformation.
Overloaded versions are available for float
, double
, std::complex<float>
, and std::complex<double>
.
[in] | sides |
|
[in] | howmany |
|
[in] | select | The vector select of length n. If howmany = Select, select specifies the eigenvectors to be computed. The eigenvector corresponding to the j-th eigenvalue is computed if select(j) = true. Not referenced if howmany = All or Backtransform. TODO: updated in real case. See [sd]trevc. |
[in] | n | The order of the matrix T. n >= 0. |
[in,out] | T | The n-by-n matrix T, stored in an ldt-by-n array. The upper triangular matrix T. T is modified, but restored on exit. |
[in] | ldt | The leading dimension of the array T. ldt >= max(1,n). |
[in,out] | VL | The n-by-mm matrix VL, stored in an ldvl-by-mm array. On entry, if side = Left or Both and howmany = Backtransform, VL must contain an n-by-n matrix Q (usually the unitary matrix Q of Schur vectors returned by lapack::hseqr ). On exit, if side = Left or Both, VL contains:
|
[in] | ldvl | The leading dimension of the array VL. ldvl >= 1, and if side = Left or Both, ldvl >= n. |
[in,out] | VR | The n-by-mm matrix VR, stored in an ldvr-by-mm array. On entry, if side = Right or Both and howmany = Backtransform, VR must contain an n-by-n matrix Q (usually the unitary matrix Q of Schur vectors returned by lapack::hseqr ). On exit, if side = Right or Both, VR contains:
|
[in] | ldvr | The leading dimension of the array VR. ldvr >= 1, and if side = Right or Both, ldvr >= n. |
[in] | mm | The number of columns in the arrays VL and/or VR. mm >= m. |
[out] | m | The number of columns in the arrays VL and/or VR actually used to store the eigenvectors. If howmany = All or Backtransform, m is set to n. Each selected eigenvector occupies one column. |
The algorithm used in this program is basically backward (forward) substitution, with scaling to make the the code robust against possible overflow.
Each eigenvector is normalized so that the element of largest magnitude has magnitude 1; here the magnitude of a complex number (x,y) is taken to be |x| + |y|.
int64_t lapack::unghr | ( | int64_t | n, |
int64_t | ilo, | ||
int64_t | ihi, | ||
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 ihi-ilo elementary reflectors of order n, as returned by lapack::gehrd
:
\[ Q = H(ilo) H(ilo+1) \dots H(ihi-1). \]
Overloaded versions are available for float
, double
, std::complex<float>
, and std::complex<double>
. For real matrices, this is an alias for lapack::orghr
.
[in] | n | The order of the matrix Q. n >= 0. |
[in] | ilo | |
[in] | ihi | ilo and ihi must have the same values as in the previous call of lapack::gehrd . Q is equal to the unit matrix except in the submatrix Q(ilo+1:ihi,ilo+1:ihi).
|
[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::gehrd . On exit, the n-by-n unitary matrix Q. |
[in] | lda | The leading dimension of the array A. lda >= max(1,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::gehrd . |
int64_t lapack::unmhr | ( | lapack::Side | side, |
lapack::Op | trans, | ||
int64_t | m, | ||
int64_t | n, | ||
int64_t | ilo, | ||
int64_t | ihi, | ||
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::gehrd
as follows:
where Q is a unitary matrix of order m if side = Left and order n if side = Right. Q is defined as the product of ihi-ilo elementary reflectors, as returned by lapack::gehrd
:
\[ Q = H(ilo) H(ilo+1) \dots H(ihi-1). \]
Overloaded versions are available for float
, double
, std::complex<float>
, and std::complex<double>
. For real matrices, this is an alias for lapack::ormhr
.
[in] | side |
|
[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] | ilo | |
[in] | ihi | ilo and ihi must have the same values as in the previous call of lapack::gehrd . Q is equal to the unit matrix except in the submatrix Q(ilo+1:ihi,ilo+1:ihi).
|
[in] | A | The vectors which define the elementary reflectors, as returned by lapack::gehrd .
|
[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::gehrd .
|
[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). |