LAPACK++
LAPACK C++ API
|
Functions | |
template<typename scalar_t > | |
void | lapack::generate_correlation_factor (Matrix< scalar_t > &A) |
Given matrix A with singular values such that sum(sigma_i^2) = n, returns A with columns of unit norm, with the same condition number. More... | |
template<typename scalar_t > | |
void | lapack::generate_geev (MatrixParams ¶ms, Dist dist, blas::real_type< scalar_t > cond, blas::real_type< scalar_t > sigma_max, Matrix< scalar_t > &A, Vector< blas::real_type< scalar_t > > &sigma) |
Generates matrix using general eigenvalue decomposition, \(A = V T V^H\), with orthogonal eigenvectors. More... | |
template<typename scalar_t > | |
void | lapack::generate_geevx (MatrixParams ¶ms, Dist dist, blas::real_type< scalar_t > cond, blas::real_type< scalar_t > sigma_max, Matrix< scalar_t > &A, Vector< blas::real_type< scalar_t > > &sigma) |
Generates matrix using general eigenvalue decomposition, \(A = X T X^{-1}\), with random eigenvectors. More... | |
template<typename scalar_t > | |
void | lapack::generate_heev (MatrixParams ¶ms, Dist dist, bool rand_sign, blas::real_type< scalar_t > cond, blas::real_type< scalar_t > condD, blas::real_type< scalar_t > sigma_max, Matrix< scalar_t > &A, Vector< blas::real_type< scalar_t > > &sigma) |
Generates matrix using Hermitian eigenvalue decomposition, \(A = U Sigma U^H\). More... | |
template<typename scalar_t > | |
void | lapack::generate_matrix (MatrixParams ¶ms, int64_t m, int64_t n, scalar_t *A_ptr, int64_t lda, blas::real_type< scalar_t > *sigma_ptr) |
Generates an m-by-n test matrix. More... | |
template<typename scalar_t > | |
void | lapack::generate_matrix (MatrixParams ¶ms, Matrix< scalar_t > &A, Vector< blas::real_type< scalar_t > > &sigma) |
Generates an m-by-n test matrix. More... | |
template<typename scalar_t > | |
void | lapack::generate_sigma (MatrixParams ¶ms, Dist dist, bool rand_sign, blas::real_type< scalar_t > cond, blas::real_type< scalar_t > sigma_max, Matrix< scalar_t > &A, Vector< blas::real_type< scalar_t > > &sigma) |
Generates sigma vector of singular or eigenvalues, according to distribution. More... | |
template<typename scalar_t > | |
void | lapack::generate_svd (MatrixParams ¶ms, Dist dist, blas::real_type< scalar_t > cond, blas::real_type< scalar_t > condD, blas::real_type< scalar_t > sigma_max, Matrix< scalar_t > &A, Vector< blas::real_type< scalar_t > > &sigma) |
Generates matrix using SVD, \(A = U Sigma V^H\). More... | |
void lapack::generate_correlation_factor | ( | Matrix< scalar_t > & | A | ) |
Given matrix A with singular values such that sum(sigma_i^2) = n, returns A with columns of unit norm, with the same condition number.
see: Davies and Higham, 2000, Numerically stable generation of correlation matrices and their factors.
Internal function, called from generate_matrix().
void lapack::generate_geev | ( | MatrixParams & | params, |
Dist | dist, | ||
blas::real_type< scalar_t > | cond, | ||
blas::real_type< scalar_t > | sigma_max, | ||
Matrix< scalar_t > & | A, | ||
Vector< blas::real_type< scalar_t > > & | sigma | ||
) |
Generates matrix using general eigenvalue decomposition, \(A = V T V^H\), with orthogonal eigenvectors.
Not yet implemented.
Internal function, called from generate_matrix().
void lapack::generate_geevx | ( | MatrixParams & | params, |
Dist | dist, | ||
blas::real_type< scalar_t > | cond, | ||
blas::real_type< scalar_t > | sigma_max, | ||
Matrix< scalar_t > & | A, | ||
Vector< blas::real_type< scalar_t > > & | sigma | ||
) |
Generates matrix using general eigenvalue decomposition, \(A = X T X^{-1}\), with random eigenvectors.
Not yet implemented.
Internal function, called from generate_matrix().
void lapack::generate_heev | ( | MatrixParams & | params, |
Dist | dist, | ||
bool | rand_sign, | ||
blas::real_type< scalar_t > | cond, | ||
blas::real_type< scalar_t > | condD, | ||
blas::real_type< scalar_t > | sigma_max, | ||
Matrix< scalar_t > & | A, | ||
Vector< blas::real_type< scalar_t > > & | sigma | ||
) |
Generates matrix using Hermitian eigenvalue decomposition, \(A = U Sigma U^H\).
Internal function, called from generate_matrix().
void lapack::generate_matrix | ( | MatrixParams & | params, |
int64_t | m, | ||
int64_t | n, | ||
scalar_t * | A_ptr, | ||
int64_t | lda, | ||
blas::real_type< scalar_t > * | sigma_ptr | ||
) |
Generates an m-by-n test matrix.
Traditional interface with m, n, lda instead of Matrix object.
void lapack::generate_matrix | ( | MatrixParams & | params, |
Matrix< scalar_t > & | A, | ||
Vector< blas::real_type< scalar_t > > & | sigma | ||
) |
Generates an m-by-n test matrix.
Similar to LAPACK's libtmg functionality, but a level 3 BLAS implementation.
[in] | params | Test matrix parameters. Uses matrix, cond, condD parameters; see further details. |
[out] | A | Complex array, dimension (lda, n). On output, the m-by-n test matrix A in an lda-by-n array. |
[in,out] | sigma | Real array, dimension (min(m,n))
|
The matrix parameter specifies the matrix kind according to the tables below. As indicated, kinds take an optional distribution suffix (^) and an optional scaling and modifier suffix (@). The default distribution is logrand. Examples: rand, rand_small, svd_arith, heev_geo_small.
The cond parameter specifies the condition number \(cond(S)\), where \(S\) is either the singular values \(\Sigma\) or the eigenvalues \(\Lambda\), as described by the distributions below. It does not apply to some matrices and distributions. For geev and geevx, cond(A) is generally much worse than cond(S). If _dominant is applied, cond(A) generally improves. By default, cond(S) = sqrt( 1/eps ) = 6.7e7 for double, 2.9e3 for single.
The condD parameter specifies the condition number cond(D), where D is a diagonal scaling matrix [1]. By default, condD = 1. If condD != 1, then:
Note using condD changes the singular or eigenvalues of \(A\); on output, sigma contains the singular or eigenvalues of \(A_0\), not of \(A\).
Notation used below: \(\Sigma\) is a diagonal matrix with entries \(\sigma_i\) for \(i = 1, \dots, n\); \(\Lambda\) is a diagonal matrix with entries \(\lambda_i = \pm \sigma_i\), with random sign; \(U\) and \(V\) are random orthogonal matrices from the Haar distribution [2], \(X\) is a random matrix.
See LAPACK Working Note (LAWN) 41:
Table 5 (Test matrices for the nonsymmetric eigenvalue problem)
Table 10 (Test matrices for the symmetric eigenvalue problem)
Table 11 (Test matrices for the singular value decomposition)
Matrix | Description |
---|---|
zero | all zero |
identity | ones on diagonal, rest zero |
jordan | ones on diagonal and first subdiagonal, rest zero |
– | – |
rand</td> | matrix entries random uniform on (0, 1) |
rands</td> | matrix entries random uniform on (-1, 1) |
randn</td> | matrix entries random normal with mean 0, std 1 |
– | – |
diag^</td> | \(A = \Sigma\) |
svd^</td> | \(A = U \Sigma V^H\) |
poev^</td> | \(A = V \Sigma V^H\) (eigenvalues positive, i.e., matrix SPD) |
spd^</td> | alias for poev |
heev^</td> | \(A = V \Lambda V^H\) (eigenvalues mixed signs) |
syev^</td> | alias for heev |
geev^</td> | \(A = V T V^H\), Schur-form \(T\) [not yet implemented] |
geevx^</td> | \(A = X T X^{-1}\), Schur-form \(T\), \(X\) ill-conditioned [not yet implemented] |
Note for geev that \(cond(\Lambda)\) is specified, where \(\Lambda = diag(T)\); while \(cond(T)\) and \(cond(A)\) are usually much worse.
^ and @ denote optional suffixes described below.
^ Distribution | Description |
---|---|
_logrand | \(\log(\sigma_i)\) random uniform on \([ \log(1/cond), \log(1) ]\); default |
_arith | \(\sigma_i = 1 - \frac{i - 1}{n - 1} (1 - 1/cond)\); arithmetic: \(\sigma_{i+1} - \sigma_i\) is constant |
_geo | \(\sigma_i = (cond)^{ -(i-1)/(n-1) }\); geometric: \(\sigma_{i+1} / \sigma_i\) is constant |
_cluster0 | \(\Sigma = [ 1, 1/cond, ..., 1/cond ]\); 1 unit value, \(n-1\) small values |
_cluster1 | \(\Sigma = [ 1, ..., 1, 1/cond ]\); \(n-1\) unit values, 1 small value |
_rarith | _arith, reversed order |
_rgeo | _geo, reversed order |
_rcluster0 | _cluster0, reversed order |
_rcluster1 | _cluster1, reversed order |
_specified | user specified sigma on input |
– | – |
_rand | \(\sigma_i\) random uniform on (0, 1) |
_rands | \(\sigma_i\) random uniform on (-1, 1) |
_randn | \(\sigma_i\) random normal with mean 0, std 1 |
Note _rand, _rands, _randn do not use cond; the condition number is random. For randn, Expected( log( cond ) ) = log( 4.65 n ) [Edelman, 1988].
Note for _rands and _randn, \(\Sigma\) contains negative values. This means poev_rands and poev_randn will not generate an SPD matrix.
@ Scaling | Description |
---|---|
_ufl | scale near underflow = 1e-308 for double |
_ofl | scale near overflow = 2e+308 for double |
_small | scale near sqrt( underflow ) = 1e-154 for double |
_large | scale near sqrt( overflow ) = 6e+153 for double |
Note scaling changes the singular or eigenvalues, but not the condition number.
@ Modifier | Description -------------—|----------— _dominant | diagonally dominant: set \(A_{i,i} = \pm \max_i( \sum_j |A_{i,j}|, \sum_j |A_{j,i}| )\).
Note _dominant changes the singular or eigenvalues, and the condition number.
[1] Demmel and Veselic, Jacobi's method is more accurate than QR, 1992.
[2] Stewart, The efficient generation of random orthogonal matrices with an application to condition estimators, 1980.
void lapack::generate_sigma | ( | MatrixParams & | params, |
Dist | dist, | ||
bool | rand_sign, | ||
blas::real_type< scalar_t > | cond, | ||
blas::real_type< scalar_t > | sigma_max, | ||
Matrix< scalar_t > & | A, | ||
Vector< blas::real_type< scalar_t > > & | sigma | ||
) |
Generates sigma vector of singular or eigenvalues, according to distribution.
Internal function, called from generate_matrix().
void lapack::generate_svd | ( | MatrixParams & | params, |
Dist | dist, | ||
blas::real_type< scalar_t > | cond, | ||
blas::real_type< scalar_t > | condD, | ||
blas::real_type< scalar_t > | sigma_max, | ||
Matrix< scalar_t > & | A, | ||
Vector< blas::real_type< scalar_t > > & | sigma | ||
) |
Generates matrix using SVD, \(A = U Sigma V^H\).
Internal function, called from generate_matrix().