LAPACK++
2022.05.00
LAPACK C++ API
|
Functions | |
int64_t | lapack::gbbrd (lapack::Vect vect, int64_t m, int64_t n, int64_t ncc, int64_t kl, int64_t ku, double *AB, int64_t ldab, double *D, double *E, double *Q, int64_t ldq, double *PT, int64_t ldpt, double *C, int64_t ldc) |
int64_t | lapack::gbbrd (lapack::Vect vect, int64_t m, int64_t n, int64_t ncc, int64_t kl, int64_t ku, float *AB, int64_t ldab, float *D, float *E, float *Q, int64_t ldq, float *PT, int64_t ldpt, float *C, int64_t ldc) |
int64_t | lapack::gbbrd (lapack::Vect vect, int64_t m, int64_t n, int64_t ncc, int64_t kl, int64_t ku, std::complex< double > *AB, int64_t ldab, double *D, double *E, std::complex< double > *Q, int64_t ldq, std::complex< double > *PT, int64_t ldpt, std::complex< double > *C, int64_t ldc) |
Reduces a general m-by-n band matrix A to real upper bidiagonal form B by a unitary transformation: \(Q^H A P = B\). More... | |
int64_t | lapack::gbbrd (lapack::Vect vect, int64_t m, int64_t n, int64_t ncc, int64_t kl, int64_t ku, std::complex< float > *AB, int64_t ldab, float *D, float *E, std::complex< float > *Q, int64_t ldq, std::complex< float > *PT, int64_t ldpt, std::complex< float > *C, int64_t ldc) |
int64_t | lapack::gebrd (int64_t m, int64_t n, double *A, int64_t lda, double *D, double *E, double *tauq, double *taup) |
int64_t | lapack::gebrd (int64_t m, int64_t n, float *A, int64_t lda, float *D, float *E, float *tauq, float *taup) |
int64_t | lapack::gebrd (int64_t m, int64_t n, std::complex< double > *A, int64_t lda, double *D, double *E, std::complex< double > *tauq, std::complex< double > *taup) |
Reduces a general m-by-n matrix A to upper or lower bidiagonal form B by a unitary transformation: \(Q^H A P = B\). More... | |
int64_t | lapack::gebrd (int64_t m, int64_t n, std::complex< float > *A, int64_t lda, float *D, float *E, std::complex< float > *tauq, std::complex< float > *taup) |
int64_t | lapack::orgbr (lapack::Vect vect, int64_t m, int64_t n, int64_t k, double *A, int64_t lda, double const *tau) |
int64_t | lapack::orgbr (lapack::Vect vect, int64_t m, int64_t n, int64_t k, float *A, int64_t lda, float const *tau) |
int64_t | lapack::ormbr (lapack::Vect vect, lapack::Side side, lapack::Op trans, int64_t m, int64_t n, int64_t k, double const *A, int64_t lda, double const *tau, double *C, int64_t ldc) |
int64_t | lapack::ormbr (lapack::Vect vect, lapack::Side side, lapack::Op trans, int64_t m, int64_t n, int64_t k, float const *A, int64_t lda, float const *tau, float *C, int64_t ldc) |
int64_t | lapack::ungbr (lapack::Vect vect, int64_t m, int64_t n, int64_t k, std::complex< double > *A, int64_t lda, std::complex< double > const *tau) |
Generates one of the complex unitary matrices \(Q\) or \(P^H\) determined by lapack::gebrd when reducing a complex matrix A to bidiagonal form: \(A = Q B P^H.\) \(Q\) and \(P^H\) are defined as products of elementary reflectors H(i) or G(i) respectively. More... | |
int64_t | lapack::ungbr (lapack::Vect vect, int64_t m, int64_t n, int64_t k, std::complex< float > *A, int64_t lda, std::complex< float > const *tau) |
int64_t | lapack::unmbr (lapack::Vect vect, lapack::Side side, lapack::Op trans, int64_t m, int64_t n, int64_t k, 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 P or Q from lapack::gebrd as follows: More... | |
int64_t | lapack::unmbr (lapack::Vect vect, lapack::Side side, lapack::Op trans, int64_t m, int64_t n, int64_t k, std::complex< float > const *A, int64_t lda, std::complex< float > const *tau, std::complex< float > *C, int64_t ldc) |
int64_t lapack::gbbrd | ( | lapack::Vect | vect, |
int64_t | m, | ||
int64_t | n, | ||
int64_t | ncc, | ||
int64_t | kl, | ||
int64_t | ku, | ||
std::complex< double > * | AB, | ||
int64_t | ldab, | ||
double * | D, | ||
double * | E, | ||
std::complex< double > * | Q, | ||
int64_t | ldq, | ||
std::complex< double > * | PT, | ||
int64_t | ldpt, | ||
std::complex< double > * | C, | ||
int64_t | ldc | ||
) |
Reduces a general m-by-n band matrix A to real upper bidiagonal form B by a unitary transformation: \(Q^H A P = B\).
The routine computes B, and optionally forms \(Q\) or \(P^H\), or computes \(Q^H C\) for a given matrix C.
Overloaded versions are available for float
, double
, std::complex<float>
, and std::complex<double>
.
[in] | vect | Whether or not the matrices Q and P^H are to be formed.
|
[in] | m | The number of rows of the matrix A. m >= 0. |
[in] | n | The number of columns of the matrix A. n >= 0. |
[in] | ncc | The number of columns of the matrix C. ncc >= 0. |
[in] | kl | The number of subdiagonals of the matrix A. kl >= 0. |
[in] | ku | The number of superdiagonals of the matrix A. ku >= 0. |
[in,out] | AB | The m-by-n band matrix AB, stored in an ldab-by-n array. On entry, the m-by-n band matrix A, stored in rows 1 to kl+ku+1. The j-th column of A is stored in the j-th column of the array AB as follows: AB(ku+1+i-j,j) = A(i,j) for max(1,j-ku) <= i <= min(m,j+kl). On exit, A is overwritten by values generated during the reduction. |
[in] | ldab | The leading dimension of the array A. ldab >= kl+ku+1. |
[out] | D | The vector D of length min(m,n). The diagonal elements of the bidiagonal matrix B. |
[out] | E | The vector E of length min(m,n)-1. The superdiagonal elements of the bidiagonal matrix B. |
[out] | Q | The m-by-m matrix Q, stored in an ldq-by-m array.
|
[in] | ldq | The leading dimension of the array Q.
|
[out] | PT | The n-by-n matrix PT, stored in an ldpt-by-n array.
|
[in] | ldpt | The leading dimension of the array PT.
|
[in,out] | C | The m-by-ncc matrix C, stored in an ldc-by-ncc array. On entry, an m-by-ncc matrix C. On exit, C is overwritten by \(Q^H C\). C is not referenced if ncc = 0. |
[in] | ldc | The leading dimension of the array C.
|
int64_t lapack::gebrd | ( | int64_t | m, |
int64_t | n, | ||
std::complex< double > * | A, | ||
int64_t | lda, | ||
double * | D, | ||
double * | E, | ||
std::complex< double > * | tauq, | ||
std::complex< double > * | taup | ||
) |
Reduces a general m-by-n matrix A to upper or lower bidiagonal form B by a unitary transformation: \(Q^H A P = B\).
If m >= n, B is upper bidiagonal; if m < n, B is lower bidiagonal.
Overloaded versions are available for float
, double
, std::complex<float>
, and std::complex<double>
.
[in] | m | The number of rows in the matrix A. m >= 0. |
[in] | n | The number of columns in the matrix A. n >= 0. |
[in,out] | A | The m-by-n matrix A, stored in an lda-by-n array. On entry, the m-by-n general matrix to be reduced. On exit:
|
[in] | lda | The leading dimension of the array A. lda >= max(1,m). |
[out] | D | The vector D of length min(m,n). The diagonal elements of the bidiagonal matrix B: D(i) = A(i,i). |
[out] | E | The vector E of length min(m,n)-1. The off-diagonal elements of the bidiagonal matrix B: if m >= n, E(i) = A(i,i+1) for i = 1, 2, ..., n-1; if m < n, E(i) = A(i+1,i) for i = 1, 2, ..., m-1. |
[out] | tauq | The vector tauq of length min(m,n). The scalar factors of the elementary reflectors which represent the unitary matrix Q. See Further Details. |
[out] | taup | The vector taup of length min(m,n). The scalar factors of the elementary reflectors which represent the unitary matrix P. See Further Details. |
The matrices Q and P are represented as products of elementary reflectors:
If m >= n,
\[ Q = H(1) H(2) . . . H(n) \]
\[ P = G(1) G(2) . . . G(n-1). \]
Each H(i) and G(i) has the form:
\[ H(i) = I - \tau_q v v^H \]
and
\[ G(i) = I - \tau_p u u^H \]
where \(\tau_q\) and \(\tau_p\) are scalars, and v and u are vectors; v(1:i-1) = 0, v(i) = 1, and v(i+1:m) is stored on exit in A(i+1:m,i); u(1:i) = 0, u(i+1) = 1, and u(i+2:n) is stored on exit in A(i,i+2:n); \(\tau_q\) is stored in tauq(i) and \(\tau_p\) in taup(i).
If m < n,
\[ Q = H(1) H(2) . . . H(m-1) \]
and
\[ P = G(1) G(2) . . . G(m) \]
Each H(i) and G(i) has the form:
\[ H(i) = I - \tau_q v v^H \]
and
\[ G(i) = I - \tau_p u u^H \]
where \(\tau_q\) and \(\tau_p\) are scalars, and v and u are vectors; v(1:i) = 0, v(i+1) = 1, and v(i+2:m) is stored on exit in A(i+2:m,i); u(1:i-1) = 0, u(i) = 1, and u(i+1:n) is stored on exit in A(i,i+1:n); \(\tau_q\) is stored in tauq(i) and \(\tau_p\) in taup(i).
The contents of A on exit are illustrated by the following examples:
m = 6 and n = 5 (m >= n): m = 5 and n = 6 (m < n): ( d e u1 u1 u1 ) ( d u1 u1 u1 u1 u1 ) ( v1 d e u2 u2 ) ( e d u2 u2 u2 u2 ) ( v1 v2 d e u3 ) ( v1 e d u3 u3 u3 ) ( v1 v2 v3 d e ) ( v1 v2 e d u4 u4 ) ( v1 v2 v3 v4 d ) ( v1 v2 v3 e d u5 ) ( v1 v2 v3 v4 v5 )
where d and e denote diagonal and off-diagonal elements of B, vi denotes an element of the vector defining H(i), and ui an element of the vector defining G(i).
int64_t lapack::orgbr | ( | lapack::Vect | vect, |
int64_t | m, | ||
int64_t | n, | ||
int64_t | k, | ||
double * | A, | ||
int64_t | lda, | ||
double const * | tau | ||
) |
int64_t lapack::ormbr | ( | lapack::Vect | vect, |
lapack::Side | side, | ||
lapack::Op | trans, | ||
int64_t | m, | ||
int64_t | n, | ||
int64_t | k, | ||
double const * | A, | ||
int64_t | lda, | ||
double const * | tau, | ||
double * | C, | ||
int64_t | ldc | ||
) |
int64_t lapack::ungbr | ( | lapack::Vect | vect, |
int64_t | m, | ||
int64_t | n, | ||
int64_t | k, | ||
std::complex< double > * | A, | ||
int64_t | lda, | ||
std::complex< double > const * | tau | ||
) |
Generates one of the complex unitary matrices \(Q\) or \(P^H\) determined by lapack::gebrd
when reducing a complex matrix A to bidiagonal form: \(A = Q B P^H.\) \(Q\) and \(P^H\) are defined as products of elementary reflectors H(i) or G(i) respectively.
ungbr
returns the first n columns of Q, where m >= n >= k;ungbr
returns Q as an m-by-m matrix.ungbr
returns the first m rows of \(P^H,\) where n >= m >= k;ungbr
returns \(P^H\) as an n-by-n matrix.Overloaded versions are available for float
, double
, std::complex<float>
, and std::complex<double>
. For real matrices, this is an alias for lapack::orgbr
.
[in] | vect | Specifies whether the matrix \(Q\) or the matrix \(P^H\) is required, as defined in the transformation applied by lapack::gebrd :
|
[in] | m | The number of rows of the matrix Q or \(P^H\) to be returned. m >= 0. |
[in] | n | The number of columns of the matrix Q or \(P^H\) to be returned. n >= 0.
|
[in] | k |
|
[in,out] | A | The m-by-n matrix A, stored in an lda-by-n array. On entry, the vectors which define the elementary reflectors, as returned by lapack::gebrd . On exit, the m-by-n matrix \(Q\) or \(P^H,\) |
[in] | lda | The leading dimension of the array A. lda >= m. |
[in] | tau | tau(i) must contain the scalar factor of the elementary reflector H(i) or G(i), which determines \(Q\) or \(P^H,\) as returned by lapack::gebrd in its array argument tauq or taup.
|
int64_t lapack::unmbr | ( | lapack::Vect | vect, |
lapack::Side | side, | ||
lapack::Op | trans, | ||
int64_t | m, | ||
int64_t | n, | ||
int64_t | k, | ||
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 P or Q from lapack::gebrd
as follows:
Here \(Q\) and \(P^H\) are the unitary matrices determined by lapack::gebrd
when reducing a complex matrix A to bidiagonal form: \(A = Q B P^H\). \(Q\) and \(P^H\) are defined as products of elementary reflectors H(i) and G(i) respectively.
Let nq = m if side = Left and nq = n if side = Right. Thus nq is the order of the unitary matrix \(Q\) or \(P^H\) that is applied.
Overloaded versions are available for float
, double
, std::complex<float>
, and std::complex<double>
. For real matrices, this is an alias for lapack::ormbr
.
[in] | vect |
|
[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] | k |
|
[in] | A | The vector A of length lda,min(nq,k) if vect = Q; lda,nq if vect = P. (lda,min(nq,k)) if vect = Q (lda,nq) if vect = P The vectors which define the elementary reflectors H(i) and G(i), whose products determine the matrices Q and P, as returned by lapack::gebrd .
|
[in] | lda | The leading dimension of the array A.
|
[in] | tau | The vector tau of length min(nq,k). tau(i) must contain the scalar factor of the elementary reflector H(i) or G(i) which determines Q or P, as returned by lapack::gebrd in the array argument tauq or taup. |
[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 one of \(Q C\), \(Q^H C\), \(C Q^H\), \(C Q\), \(P C\), \(P^H C\), \(C P^H\), or \(C P\). |
[in] | ldc | The leading dimension of the array C. ldc >= max(1,m). |