![]() |
fml
0.1-0
Fused Matrix Library
|
Linear algebra functions. More...
Functions | |
template<typename REAL > | |
void | rsvd (const int k, const int q, cpumat< REAL > &x, cpuvec< REAL > &s, cpumat< REAL > &u, cpumat< REAL > &vt) |
template<typename REAL > | |
void | add (const bool transx, const bool transy, const REAL alpha, const REAL beta, const cpumat< REAL > &x, const cpumat< REAL > &y, cpumat< REAL > &ret) |
Returns alpha*op(x) + beta*op(y) where op(A) is A or A^T. More... | |
template<typename REAL > | |
cpumat< REAL > | add (const bool transx, const bool transy, const REAL alpha, const REAL beta, const cpumat< REAL > &x, const cpumat< REAL > &y) |
This is an overloaded member function, provided for convenience. It differs from the above function only in what argument(s) it accepts. | |
template<typename REAL > | |
cpumat< REAL > | matmult (const bool transx, const bool transy, const REAL alpha, const cpumat< REAL > &x, const cpumat< REAL > &y) |
Returns alpha*op(x)*op(y) where op(A) is A or A^T. More... | |
template<typename REAL > | |
void | matmult (const bool transx, const bool transy, const REAL alpha, const cpumat< REAL > &x, const cpumat< REAL > &y, cpumat< REAL > &ret) |
Computes ret = alpha*op(x)*op(y) where op(A) is A or A^T. More... | |
template<typename REAL > | |
void | crossprod (const REAL alpha, const cpumat< REAL > &x, cpumat< REAL > &ret) |
Computes lower triangle of alpha*x^T*x. More... | |
template<typename REAL > | |
cpumat< REAL > | crossprod (const REAL alpha, const cpumat< REAL > &x) |
This is an overloaded member function, provided for convenience. It differs from the above function only in what argument(s) it accepts. | |
template<typename REAL > | |
void | tcrossprod (const REAL alpha, const cpumat< REAL > &x, cpumat< REAL > &ret) |
Computes lower triangle of alpha*x*x^T. More... | |
template<typename REAL > | |
cpumat< REAL > | tcrossprod (const REAL alpha, const cpumat< REAL > &x) |
template<typename REAL > | |
void | xpose (const cpumat< REAL > &x, cpumat< REAL > &tx) |
Computes the transpose out-of-place (i.e. in a copy). More... | |
template<typename REAL > | |
cpumat< REAL > | xpose (const cpumat< REAL > &x) |
This is an overloaded member function, provided for convenience. It differs from the above function only in what argument(s) it accepts. | |
template<typename REAL > | |
void | lu (cpumat< REAL > &x, cpuvec< int > &p, int &info) |
Computes the PLU factorization with partial pivoting. More... | |
template<typename REAL > | |
void | lu (cpumat< REAL > &x) |
This is an overloaded member function, provided for convenience. It differs from the above function only in what argument(s) it accepts. | |
template<typename REAL > | |
void | det (cpumat< REAL > &x, int &sign, REAL &modulus) |
Computes the determinant in logarithmic form. More... | |
template<typename REAL > | |
REAL | trace (const cpumat< REAL > &x) |
Computes the trace, i.e. the sum of the diagonal. More... | |
template<typename REAL > | |
void | svd (cpumat< REAL > &x, cpuvec< REAL > &s) |
Computes the singular value decomposition. More... | |
template<typename REAL > | |
void | svd (cpumat< REAL > &x, cpuvec< REAL > &s, cpumat< REAL > &u, cpumat< REAL > &vt) |
This is an overloaded member function, provided for convenience. It differs from the above function only in what argument(s) it accepts. | |
template<typename REAL > | |
void | eigen_sym (cpumat< REAL > &x, cpuvec< REAL > &values) |
Compute the eigenvalues and optionally the eigenvectors for a symmetric matrix. More... | |
template<typename REAL > | |
void | eigen_sym (cpumat< REAL > &x, cpuvec< REAL > &values, cpumat< REAL > &vectors) |
This is an overloaded member function, provided for convenience. It differs from the above function only in what argument(s) it accepts. | |
template<typename REAL > | |
void | invert (cpumat< REAL > &x) |
Compute the matrix inverse. More... | |
template<typename REAL > | |
void | solve (cpumat< REAL > &x, cpuvec< REAL > &y) |
Solve a system of equations. More... | |
template<typename REAL > | |
void | solve (cpumat< REAL > &x, cpumat< REAL > &y) |
This is an overloaded member function, provided for convenience. It differs from the above function only in what argument(s) it accepts. | |
template<typename REAL > | |
void | qr (const bool pivot, cpumat< REAL > &x, cpuvec< REAL > &qraux) |
Computes the QR decomposition. More... | |
template<typename REAL > | |
void | qr_Q (const cpumat< REAL > &QR, const cpuvec< REAL > &qraux, cpumat< REAL > &Q, cpuvec< REAL > &work) |
Recover the Q matrix from a QR decomposition. More... | |
template<typename REAL > | |
void | qr_R (const cpumat< REAL > &QR, cpumat< REAL > &R) |
Recover the R matrix from a QR decomposition. More... | |
template<typename REAL > | |
void | lq (cpumat< REAL > &x, cpuvec< REAL > &lqaux) |
Computes the LQ decomposition. More... | |
template<typename REAL > | |
void | lq_L (const cpumat< REAL > &LQ, cpumat< REAL > &L) |
Recover the L matrix from an LQ decomposition. More... | |
template<typename REAL > | |
void | lq_Q (const cpumat< REAL > &LQ, const cpuvec< REAL > &lqaux, cpumat< REAL > &Q, cpuvec< REAL > &work) |
Recover the Q matrix from an LQ decomposition. More... | |
template<typename REAL > | |
void | tssvd (cpumat< REAL > &x, cpuvec< REAL > &s, cpumat< REAL > &u, cpumat< REAL > &vt) |
Computes the singular value decomposition for tall/skinny data. The number of rows must be greater than the number of columns. If the number of rows is not significantly larger than the number of columns, this may not be more efficient than simply calling linalg::svd() . More... | |
template<typename REAL > | |
void | tssvd (cpumat< REAL > &x, cpuvec< REAL > &s) |
This is an overloaded member function, provided for convenience. It differs from the above function only in what argument(s) it accepts. | |
template<typename REAL > | |
void | cpsvd (const cpumat< REAL > &x, cpuvec< REAL > &s, cpumat< REAL > &u, cpumat< REAL > &vt) |
Computes the singular value decomposition using the "crossproducts SVD". This method is not numerically stable. More... | |
template<typename REAL > | |
void | cpsvd (const cpumat< REAL > &x, cpuvec< REAL > &s) |
This is an overloaded member function, provided for convenience. It differs from the above function only in what argument(s) it accepts. | |
template<typename REAL > | |
void | chol (cpumat< REAL > &x) |
Compute the Choleski factorization. More... | |
template<typename REAL > | |
void | add (const bool transx, const bool transy, const REAL alpha, const REAL beta, const gpumat< REAL > &x, const gpumat< REAL > &y, gpumat< REAL > &ret) |
Returns alpha*op(x) + beta*op(y) where op(A) is A or A^T. More... | |
template<typename REAL > | |
gpumat< REAL > | add (const bool transx, const bool transy, const REAL alpha, const REAL beta, const gpumat< REAL > &x, const gpumat< REAL > &y) |
This is an overloaded member function, provided for convenience. It differs from the above function only in what argument(s) it accepts. | |
template<typename REAL > | |
gpumat< REAL > | matmult (const bool transx, const bool transy, const REAL alpha, const gpumat< REAL > &x, const gpumat< REAL > &y) |
Returns alpha*op(x)*op(y) where op(A) is A or A^T. More... | |
template<typename REAL > | |
void | matmult (const bool transx, const bool transy, const REAL alpha, const gpumat< REAL > &x, const gpumat< REAL > &y, gpumat< REAL > &ret) |
Computes ret = alpha*op(x)*op(y) where op(A) is A or A^T. More... | |
template<typename REAL > | |
void | crossprod (const REAL alpha, const gpumat< REAL > &x, gpumat< REAL > &ret) |
Computes lower triangle of alpha*x^T*x. More... | |
template<typename REAL > | |
gpumat< REAL > | crossprod (const REAL alpha, const gpumat< REAL > &x) |
This is an overloaded member function, provided for convenience. It differs from the above function only in what argument(s) it accepts. | |
template<typename REAL > | |
void | tcrossprod (const REAL alpha, const gpumat< REAL > &x, gpumat< REAL > &ret) |
Computes lower triangle of alpha*x*x^T. More... | |
template<typename REAL > | |
gpumat< REAL > | tcrossprod (const REAL alpha, const gpumat< REAL > &x) |
This is an overloaded member function, provided for convenience. It differs from the above function only in what argument(s) it accepts. | |
template<typename REAL > | |
void | xpose (const gpumat< REAL > &x, gpumat< REAL > &tx) |
Computes the transpose out-of-place (i.e. in a copy). More... | |
template<typename REAL > | |
gpumat< REAL > | xpose (const gpumat< REAL > &x) |
This is an overloaded member function, provided for convenience. It differs from the above function only in what argument(s) it accepts. | |
template<typename REAL > | |
void | lu (gpumat< REAL > &x, gpuvec< int > &p, int &info) |
Computes the PLU factorization with partial pivoting. More... | |
template<typename REAL > | |
void | lu (gpumat< REAL > &x) |
This is an overloaded member function, provided for convenience. It differs from the above function only in what argument(s) it accepts. | |
template<typename REAL > | |
void | det (gpumat< REAL > &x, int &sign, REAL &modulus) |
Computes the determinant in logarithmic form. More... | |
template<typename REAL > | |
REAL | trace (const gpumat< REAL > &x) |
Computes the trace, i.e. the sum of the diagonal. More... | |
template<typename REAL > | |
void | svd (gpumat< REAL > &x, gpuvec< REAL > &s) |
Computes the singular value decomposition. More... | |
template<typename REAL > | |
void | svd (gpumat< REAL > &x, gpuvec< REAL > &s, gpumat< REAL > &u, gpumat< REAL > &vt) |
This is an overloaded member function, provided for convenience. It differs from the above function only in what argument(s) it accepts. | |
template<typename REAL > | |
void | eigen_sym (gpumat< REAL > &x, gpuvec< REAL > &values) |
Compute the eigenvalues and optionally the eigenvectors for a symmetric matrix. More... | |
template<typename REAL > | |
void | eigen_sym (gpumat< REAL > &x, gpuvec< REAL > &values, gpumat< REAL > &vectors) |
This is an overloaded member function, provided for convenience. It differs from the above function only in what argument(s) it accepts. | |
template<typename REAL > | |
void | invert (gpumat< REAL > &x) |
Compute the matrix inverse. More... | |
template<typename REAL > | |
void | solve (gpumat< REAL > &x, gpuvec< REAL > &y) |
Solve a system of equations. More... | |
template<typename REAL > | |
void | solve (gpumat< REAL > &x, gpumat< REAL > &y) |
This is an overloaded member function, provided for convenience. It differs from the above function only in what argument(s) it accepts. | |
template<typename REAL > | |
void | qr (const bool pivot, gpumat< REAL > &x, gpuvec< REAL > &qraux) |
Computes the QR decomposition. More... | |
template<typename REAL > | |
void | qr_Q (const gpumat< REAL > &QR, const gpuvec< REAL > &qraux, gpumat< REAL > &Q, gpuvec< REAL > &work) |
Recover the Q matrix from a QR decomposition. More... | |
template<typename REAL > | |
void | qr_R (const gpumat< REAL > &QR, gpumat< REAL > &R) |
Recover the R matrix from a QR decomposition. More... | |
template<typename REAL > | |
void | lq (gpumat< REAL > &x, gpuvec< REAL > &lqaux) |
Computes the LQ decomposition. More... | |
template<typename REAL > | |
void | lq_L (const gpumat< REAL > &LQ, gpumat< REAL > &L) |
Recover the L matrix from an LQ decomposition. More... | |
template<typename REAL > | |
void | lq_Q (const gpumat< REAL > &LQ, const gpuvec< REAL > &lqaux, gpumat< REAL > &Q, gpuvec< REAL > &work) |
Recover the Q matrix from an LQ decomposition. More... | |
template<typename REAL > | |
void | tssvd (gpumat< REAL > &x, gpuvec< REAL > &s, gpumat< REAL > &u, gpumat< REAL > &vt) |
Computes the singular value decomposition for tall/skinny data. The number of rows must be greater than the number of columns. If the number of rows is not significantly larger than the number of columns, this may not be more efficient than simply calling linalg::svd() . More... | |
template<typename REAL > | |
void | tssvd (gpumat< REAL > &x, gpuvec< REAL > &s) |
This is an overloaded member function, provided for convenience. It differs from the above function only in what argument(s) it accepts. | |
template<typename REAL > | |
void | cpsvd (const gpumat< REAL > &x, gpuvec< REAL > &s, gpumat< REAL > &u, gpumat< REAL > &vt) |
Computes the singular value decomposition using the "crossproducts SVD". This method is not numerically stable. More... | |
template<typename REAL > | |
void | cpsvd (const gpumat< REAL > &x, gpuvec< REAL > &s) |
This is an overloaded member function, provided for convenience. It differs from the above function only in what argument(s) it accepts. | |
template<typename REAL > | |
void | chol (gpumat< REAL > &x) |
Compute the Choleski factorization. More... | |
template<typename REAL > | |
void | add (const bool transx, const bool transy, const REAL alpha, const REAL beta, const mpimat< REAL > &x, const mpimat< REAL > &y, mpimat< REAL > &ret) |
Returns alpha*op(x) + beta*op(y) where op(A) is A or A^T. More... | |
template<typename REAL > | |
mpimat< REAL > | add (const bool transx, const bool transy, const REAL alpha, const REAL beta, const mpimat< REAL > &x, const mpimat< REAL > &y) |
This is an overloaded member function, provided for convenience. It differs from the above function only in what argument(s) it accepts. | |
template<typename REAL > | |
mpimat< REAL > | matmult (const bool transx, const bool transy, const REAL alpha, const mpimat< REAL > &x, const mpimat< REAL > &y) |
Returns alpha*op(x)*op(y) where op(A) is A or A^T. More... | |
template<typename REAL > | |
void | matmult (const bool transx, const bool transy, const REAL alpha, const mpimat< REAL > &x, const mpimat< REAL > &y, mpimat< REAL > &ret) |
Computes ret = alpha*op(x)*op(y) where op(A) is A or A^T. More... | |
template<typename REAL > | |
void | crossprod (const REAL alpha, const mpimat< REAL > &x, mpimat< REAL > &ret) |
Computes lower triangle of alpha*x^T*x. More... | |
template<typename REAL > | |
mpimat< REAL > | crossprod (const REAL alpha, const mpimat< REAL > &x) |
This is an overloaded member function, provided for convenience. It differs from the above function only in what argument(s) it accepts. | |
template<typename REAL > | |
void | tcrossprod (const REAL alpha, const mpimat< REAL > &x, mpimat< REAL > &ret) |
Computes lower triangle of alpha*x*x^T. More... | |
template<typename REAL > | |
mpimat< REAL > | tcrossprod (const REAL alpha, const mpimat< REAL > &x) |
This is an overloaded member function, provided for convenience. It differs from the above function only in what argument(s) it accepts. | |
template<typename REAL > | |
void | xpose (const mpimat< REAL > &x, mpimat< REAL > &tx) |
Computes the transpose out-of-place (i.e. in a copy). More... | |
template<typename REAL > | |
mpimat< REAL > | xpose (const mpimat< REAL > &x) |
This is an overloaded member function, provided for convenience. It differs from the above function only in what argument(s) it accepts. | |
template<typename REAL > | |
void | lu (mpimat< REAL > &x, cpuvec< int > &p, int &info) |
Computes the PLU factorization with partial pivoting. More... | |
template<typename REAL > | |
void | lu (mpimat< REAL > &x) |
This is an overloaded member function, provided for convenience. It differs from the above function only in what argument(s) it accepts. | |
template<typename REAL > | |
void | det (mpimat< REAL > &x, int &sign, REAL &modulus) |
Computes the determinant in logarithmic form. More... | |
template<typename REAL > | |
REAL | trace (const mpimat< REAL > &x) |
Computes the trace, i.e. the sum of the diagonal. More... | |
template<typename REAL > | |
void | svd (mpimat< REAL > &x, cpuvec< REAL > &s) |
Computes the singular value decomposition. More... | |
template<typename REAL > | |
void | svd (mpimat< REAL > &x, cpuvec< REAL > &s, mpimat< REAL > &u, mpimat< REAL > &vt) |
This is an overloaded member function, provided for convenience. It differs from the above function only in what argument(s) it accepts. | |
template<typename REAL > | |
void | eigen_sym (mpimat< REAL > &x, cpuvec< REAL > &values) |
Compute the eigenvalues and optionally the eigenvectors for a symmetric matrix. More... | |
template<typename REAL > | |
void | eigen_sym (mpimat< REAL > &x, cpuvec< REAL > &values, mpimat< REAL > &vectors) |
This is an overloaded member function, provided for convenience. It differs from the above function only in what argument(s) it accepts. | |
template<typename REAL > | |
void | invert (mpimat< REAL > &x) |
Compute the matrix inverse. More... | |
template<typename REAL > | |
void | solve (mpimat< REAL > &x, mpimat< REAL > &y) |
Solve a system of equations. More... | |
template<typename REAL > | |
void | qr (const bool pivot, mpimat< REAL > &x, cpuvec< REAL > &qraux) |
Computes the QR decomposition. More... | |
template<typename REAL > | |
void | qr_Q (const mpimat< REAL > &QR, const cpuvec< REAL > &qraux, mpimat< REAL > &Q, cpuvec< REAL > &work) |
Recover the Q matrix from a QR decomposition. More... | |
template<typename REAL > | |
void | qr_R (const mpimat< REAL > &QR, mpimat< REAL > &R) |
Recover the R matrix from a QR decomposition. More... | |
template<typename REAL > | |
void | lq (mpimat< REAL > &x, cpuvec< REAL > &lqaux) |
Computes the LQ decomposition. More... | |
template<typename REAL > | |
void | lq_L (const mpimat< REAL > &LQ, mpimat< REAL > &L) |
Recover the L matrix from a LQ decomposition. More... | |
template<typename REAL > | |
void | lq_Q (const mpimat< REAL > &LQ, const cpuvec< REAL > &lqaux, mpimat< REAL > &Q, cpuvec< REAL > &work) |
Recover the Q matrix from a LQ decomposition. More... | |
template<typename REAL > | |
void | tssvd (mpimat< REAL > &x, cpuvec< REAL > &s, mpimat< REAL > &u, mpimat< REAL > &vt) |
Computes the singular value decomposition for tall/skinny data. The number of rows must be greater than the number of columns. If the number of rows is not significantly larger than the number of columns, this may not be more efficient than simply calling linalg::svd() . More... | |
template<typename REAL > | |
void | tssvd (mpimat< REAL > &x, cpuvec< REAL > &s) |
This is an overloaded member function, provided for convenience. It differs from the above function only in what argument(s) it accepts. | |
template<typename REAL > | |
void | cpsvd (const mpimat< REAL > &x, cpuvec< REAL > &s, mpimat< REAL > &u, mpimat< REAL > &vt) |
Computes the singular value decomposition using the "crossproducts SVD". This method is not numerically stable. More... | |
template<typename REAL > | |
void | cpsvd (const mpimat< REAL > &x, cpuvec< REAL > &s) |
This is an overloaded member function, provided for convenience. It differs from the above function only in what argument(s) it accepts. | |
template<typename REAL > | |
void | chol (mpimat< REAL > &x) |
Compute the Choleski factorization. More... | |
template<typename REAL > | |
void | crossprod (const REAL alpha, const parmat_cpu< REAL > &x, cpumat< REAL > &ret) |
Computes lower triangle of alpha*x^T*x. | |
template<typename REAL > | |
cpumat< REAL > | crossprod (const REAL alpha, const parmat_cpu< REAL > &x) |
This is an overloaded member function, provided for convenience. It differs from the above function only in what argument(s) it accepts. | |
Linear algebra functions.
void linalg::add | ( | const bool | transx, |
const bool | transy, | ||
const REAL | alpha, | ||
const REAL | beta, | ||
const cpumat< REAL > & | x, | ||
const cpumat< REAL > & | y, | ||
cpumat< REAL > & | ret | ||
) |
Returns alpha*op(x) + beta*op(y) where op(A) is A or A^T.
[in] | transx | Should x^T be used? |
[in] | transy | Should y^T be used? |
[in] | alpha,beta | Scalars. |
[in] | x,y | The inputs to the sum. |
[out] | ret | The sum. |
REAL | should be 'float' or 'double'. |
void linalg::add | ( | const bool | transx, |
const bool | transy, | ||
const REAL | alpha, | ||
const REAL | beta, | ||
const gpumat< REAL > & | x, | ||
const gpumat< REAL > & | y, | ||
gpumat< REAL > & | ret | ||
) |
Returns alpha*op(x) + beta*op(y) where op(A) is A or A^T.
[in] | transx | Should x^T be used? |
[in] | transy | Should y^T be used? |
[in] | alpha,beta | Scalars. |
[in] | x,y | The inputs to the sum. |
[out] | ret | The sum. |
REAL | should be '__half', 'float', or 'double'. |
void linalg::add | ( | const bool | transx, |
const bool | transy, | ||
const REAL | alpha, | ||
const REAL | beta, | ||
const mpimat< REAL > & | x, | ||
const mpimat< REAL > & | y, | ||
mpimat< REAL > & | ret | ||
) |
Returns alpha*op(x) + beta*op(y) where op(A) is A or A^T.
[in] | transx | Should x^T be used? |
[in] | transy | Should y^T be used? |
[in] | alpha,beta | Scalars. |
[in] | x,y | The inputs to the sum. |
[out] | ret | The sum. |
REAL | should be 'float' or 'double'. |
void linalg::chol | ( | cpumat< REAL > & | x | ) |
Compute the Choleski factorization.
The matrix should be 1. square, 2. symmetric, 3. positive-definite. Failure of any of these conditions can lead to a runtime exception. The input is replaced by its lower-triangular Choleski factor.
[in,out] | x | Input data matrix, replaced by its lower-triangular Choleski factor. |
REAL | should be 'float' or 'double'. |
void linalg::chol | ( | gpumat< REAL > & | x | ) |
Compute the Choleski factorization.
The matrix should be 1. square, 2. symmetric, 3. positive-definite. Failure of any of these conditions can lead to a runtime exception. The input is replaced by its lower-triangular Choleski factor.
[in,out] | x | Input data matrix, replaced by its lower-triangular Choleski factor. |
Uses the cuSOLVER function cusolverDnXpotrf()
.
REAL | should be 'float' or 'double'. |
void linalg::chol | ( | mpimat< REAL > & | x | ) |
Compute the Choleski factorization.
The matrix should be 1. square, 2. symmetric, 3. positive-definite. Failure of any of these conditions can lead to a runtime exception. The input is replaced by its lower-triangular Choleski factor.
[in,out] | x | Input data matrix, replaced by its lower-triangular Choleski factor. |
REAL | should be 'float' or 'double'. |
void linalg::cpsvd | ( | const cpumat< REAL > & | x, |
cpuvec< REAL > & | s, | ||
cpumat< REAL > & | u, | ||
cpumat< REAL > & | vt | ||
) |
Computes the singular value decomposition using the "crossproducts SVD". This method is not numerically stable.
The operation works by computing the crossproducts matrix X^T * X or X * X^T (whichever is smaller) and then computing the eigenvalue decomposition.
[in,out] | x | Input data matrix. |
[out] | s | Vector of singular values. |
[out] | u | Matrix of left singular vectors. |
[out] | vt | Matrix of (transposed) right singular vectors. |
eigen_sym()
.REAL | should be 'float' or 'double'. |
void linalg::cpsvd | ( | const gpumat< REAL > & | x, |
gpuvec< REAL > & | s, | ||
gpumat< REAL > & | u, | ||
gpumat< REAL > & | vt | ||
) |
Computes the singular value decomposition using the "crossproducts SVD". This method is not numerically stable.
The operation works by computing the crossproducts matrix X^T * X or X * X^T (whichever is smaller) and then computing the eigenvalue decomposition.
[in,out] | x | Input data matrix. |
[out] | s | Vector of singular values. |
[out] | u | Matrix of left singular vectors. |
[out] | vt | Matrix of (transposed) right singular vectors. |
eigen_sym()
.REAL | should be 'float' or 'double'. |
void linalg::cpsvd | ( | const mpimat< REAL > & | x, |
cpuvec< REAL > & | s, | ||
mpimat< REAL > & | u, | ||
mpimat< REAL > & | vt | ||
) |
Computes the singular value decomposition using the "crossproducts SVD". This method is not numerically stable.
The operation works by computing the crossproducts matrix X^T * X or X * X^T (whichever is smaller) and then computing the eigenvalue decomposition.
[in,out] | x | Input data matrix. |
[out] | s | Vector of singular values. |
[out] | u | Matrix of left singular vectors. |
[out] | vt | Matrix of (transposed) right singular vectors. |
eigen_sym()
.REAL | should be 'float' or 'double'. |
void linalg::crossprod | ( | const REAL | alpha, |
const cpumat< REAL > & | x, | ||
cpumat< REAL > & | ret | ||
) |
Computes lower triangle of alpha*x^T*x.
[in] | alpha | Scalar. |
[in] | x | Input data matrix. |
[out] | ret | The product. |
REAL | should be 'float' or 'double'. |
void linalg::crossprod | ( | const REAL | alpha, |
const gpumat< REAL > & | x, | ||
gpumat< REAL > & | ret | ||
) |
Computes lower triangle of alpha*x^T*x.
[in] | alpha | Scalar. |
[in] | x | Input data matrix. |
[out] | ret | The product. |
REAL | should be '__half', 'float', or 'double'. |
void linalg::crossprod | ( | const REAL | alpha, |
const mpimat< REAL > & | x, | ||
mpimat< REAL > & | ret | ||
) |
Computes lower triangle of alpha*x^T*x.
[in] | alpha | Scalar. |
[in] | x | Input data matrix. |
[out] | ret | The product. |
REAL | should be 'float' or 'double'. |
void linalg::det | ( | cpumat< REAL > & | x, |
int & | sign, | ||
REAL & | modulus | ||
) |
Computes the determinant in logarithmic form.
The input is replaced by its LU factorization.
[in,out] | x | Input data matrix, replaced by its LU factorization. |
[out] | sign | The sign of the determinant. |
[out] | modulus | Log of the modulus. |
REAL | should be 'float' or 'double'. |
void linalg::det | ( | gpumat< REAL > & | x, |
int & | sign, | ||
REAL & | modulus | ||
) |
Computes the determinant in logarithmic form.
The input is replaced by its LU factorization.
[in,out] | x | Input data matrix, replaced by its LU factorization. |
[out] | sign | The sign of the determinant. |
[out] | modulus | Log of the modulus. |
REAL | should be 'float' or 'double'. |
void linalg::det | ( | mpimat< REAL > & | x, |
int & | sign, | ||
REAL & | modulus | ||
) |
Computes the determinant in logarithmic form.
The input is replaced by its LU factorization.
[in,out] | x | Input data matrix, replaced by its LU factorization. |
[out] | sign | The sign of the determinant. |
[out] | modulus | Log of the modulus. |
REAL | should be 'float' or 'double'. |
Compute the eigenvalues and optionally the eigenvectors for a symmetric matrix.
The input data is overwritten.
[in,out] | x | Input data matrix. Should be square. |
[out] | values | Eigenvalues. |
[out] | vectors | Eigenvectors. |
bad_alloc
exception will be thrown.REAL | should be 'float' or 'double'. |
Compute the eigenvalues and optionally the eigenvectors for a symmetric matrix.
The input data is overwritten.
[in,out] | x | Input data matrix. Should be square. |
[out] | values | Eigenvalues. |
[out] | vectors | Eigenvectors. |
bad_alloc
exception will be thrown.REAL | should be '__half', 'float', or 'double'. |
Compute the eigenvalues and optionally the eigenvectors for a symmetric matrix.
The input data is overwritten.
[in,out] | x | Input data matrix. Should be square. |
[out] | values | Eigenvalues. |
[out] | vectors | Eigenvectors. |
bad_alloc
exception will be thrown.REAL | should be 'float' or 'double'. |
void linalg::invert | ( | cpumat< REAL > & | x | ) |
Compute the matrix inverse.
The input is replaced by its inverse, computed via a PLU.
[in,out] | x | Input data matrix. Should be square. |
bad_alloc
exception will be thrown.REAL | should be 'float' or 'double'. |
void linalg::invert | ( | gpumat< REAL > & | x | ) |
Compute the matrix inverse.
The input is replaced by its inverse, computed via a PLU.
[in,out] | x | Input data matrix. Should be square. |
cusolverDnXgetrs()
(solve).bad_alloc
exception will be thrown.REAL | should be '__half', 'float', or 'double'. |
void linalg::invert | ( | mpimat< REAL > & | x | ) |
Compute the matrix inverse.
The input is replaced by its inverse, computed via a PLU.
[in,out] | x | Input data matrix. Should be square. |
bad_alloc
exception will be thrown.REAL | should be 'float' or 'double'. |
Computes the LQ decomposition.
The factorization works mostly in-place by modifying the input data. After execution, the matrix will be the LAPACK-like compact LQ representation.
[in,out] | x | Input data matrix. Values are overwritten. |
[out] | lqaux | Auxiliary data for compact LQ. |
REAL | should be 'float' or 'double'. |
Computes the LQ decomposition.
The factorization works mostly in-place by modifying the input data. After execution, the matrix will be the LAPACK-like compact LQ representation.
[in,out] | x | Input data matrix. Values are overwritten. |
[out] | lqaux | Auxiliary data for compact LQ. |
REAL | should be 'float' or 'double'. |
Computes the LQ decomposition.
The factorization works mostly in-place by modifying the input data. After execution, the matrix will be the LAPACK-like compact LQ representation.
[in,out] | x | Input data matrix. Values are overwritten. |
[out] | lqaux | Auxiliary data for compact LQ. |
REAL | should be 'float' or 'double'. |
Recover the L matrix from an LQ decomposition.
[in] | LQ | The compact LQ factorization, as computed via lq() . |
[out] | L | The L matrix. |
REAL | should be 'float' or 'double'. |
Recover the L matrix from an LQ decomposition.
[in] | LQ | The compact LQ factorization, as computed via lq() . |
[out] | L | The L matrix. |
REAL | should be 'float' or 'double'. |
Recover the L matrix from a LQ decomposition.
[in] | LQ | The compact LQ factorization, as computed via lq() . |
[out] | L | The L matrix. |
REAL | should be 'float' or 'double'. |
void linalg::lq_Q | ( | const cpumat< REAL > & | LQ, |
const cpuvec< REAL > & | lqaux, | ||
cpumat< REAL > & | Q, | ||
cpuvec< REAL > & | work | ||
) |
Recover the Q matrix from an LQ decomposition.
[in] | LQ | The compact LQ factorization, as computed via lq() . |
[in] | lqaux | Auxiliary data for compact LQ. |
[out] | Q | The Q matrix. |
[out] | work | Workspace array. Will be resized as necessary. |
REAL | should be 'float' or 'double'. |
void linalg::lq_Q | ( | const gpumat< REAL > & | LQ, |
const gpuvec< REAL > & | lqaux, | ||
gpumat< REAL > & | Q, | ||
gpuvec< REAL > & | work | ||
) |
Recover the Q matrix from an LQ decomposition.
[in] | LQ | The compact LQ factorization, as computed via lq() . |
[in] | lqaux | Auxiliary data for compact LQ. |
[out] | Q | The Q matrix. |
[out] | work | Workspace array. Will be resized as necessary. |
REAL | should be 'float' or 'double'. |
void linalg::lq_Q | ( | const mpimat< REAL > & | LQ, |
const cpuvec< REAL > & | lqaux, | ||
mpimat< REAL > & | Q, | ||
cpuvec< REAL > & | work | ||
) |
Recover the Q matrix from a LQ decomposition.
[in] | LQ | The compact LQ factorization, as computed via lq() . |
[in] | lqaux | Auxiliary data for compact LQ. |
[out] | Q | The Q matrix. |
[out] | work | Workspace array. Will be resized as necessary. |
REAL | should be 'float' or 'double'. |
Computes the PLU factorization with partial pivoting.
The input is replaced by its LU factorization, with L unit-diagonal.
[in,out] | x | Input data matrix, replaced by its LU factorization. |
[out] | p | Vector of pivots, representing the diagonal matrix P in the PLU. |
[out] | info | The LAPACK return number. |
REAL | should be 'float' or 'double'. |
Computes the PLU factorization with partial pivoting.
The input is replaced by its LU factorization, with L unit-diagonal.
[in,out] | x | Input data matrix, replaced by its LU factorization. |
[out] | p | Vector of pivots, representing the diagonal matrix P in the PLU. |
[out] | info | The LAPACK return number. |
REAL | should be '__half', 'float', or 'double'. |
Computes the PLU factorization with partial pivoting.
The input is replaced by its LU factorization, with L unit-diagonal.
[in,out] | x | Input data matrix, replaced by its LU factorization. |
[out] | p | Vector of pivots, representing the diagonal matrix P in the PLU. |
[out] | info | The ScaLAPACK return number. |
REAL | should be 'float' or 'double'. |
cpumat<REAL> linalg::matmult | ( | const bool | transx, |
const bool | transy, | ||
const REAL | alpha, | ||
const cpumat< REAL > & | x, | ||
const cpumat< REAL > & | y | ||
) |
Returns alpha*op(x)*op(y) where op(A) is A or A^T.
[in] | transx | Should x^T be used? |
[in] | transy | Should y^T be used? |
[in] | alpha | Scalar. |
[in] | x | Left multiplicand. |
[in] | y | Right multiplicand. |
REAL | should be 'float' or 'double'. |
void linalg::matmult | ( | const bool | transx, |
const bool | transy, | ||
const REAL | alpha, | ||
const cpumat< REAL > & | x, | ||
const cpumat< REAL > & | y, | ||
cpumat< REAL > & | ret | ||
) |
Computes ret = alpha*op(x)*op(y) where op(A) is A or A^T.
[in] | transx | Should x^T be used? |
[in] | transy | Should y^T be used? |
[in] | alpha | Scalar. |
[in] | x | Left multiplicand. |
[in] | y | Right multiplicand. |
[out] | ret | The product. |
REAL | should be 'float' or 'double'. |
gpumat<REAL> linalg::matmult | ( | const bool | transx, |
const bool | transy, | ||
const REAL | alpha, | ||
const gpumat< REAL > & | x, | ||
const gpumat< REAL > & | y | ||
) |
Returns alpha*op(x)*op(y) where op(A) is A or A^T.
[in] | transx | Should x^T be used? |
[in] | transy | Should y^T be used? |
[in] | alpha | Scalar. |
[in] | x | Left multiplicand. |
[in] | y | Right multiplicand. |
REAL | should be '__half', 'float', or 'double'. |
void linalg::matmult | ( | const bool | transx, |
const bool | transy, | ||
const REAL | alpha, | ||
const gpumat< REAL > & | x, | ||
const gpumat< REAL > & | y, | ||
gpumat< REAL > & | ret | ||
) |
Computes ret = alpha*op(x)*op(y) where op(A) is A or A^T.
[in] | transx | Should x^T be used? |
[in] | transy | Should y^T be used? |
[in] | alpha | Scalar. |
[in] | x | Left multiplicand. |
[in] | y | Right multiplicand. |
[out] | ret | The product. |
REAL | should be '__half', 'float', or 'double'. |
mpimat<REAL> linalg::matmult | ( | const bool | transx, |
const bool | transy, | ||
const REAL | alpha, | ||
const mpimat< REAL > & | x, | ||
const mpimat< REAL > & | y | ||
) |
Returns alpha*op(x)*op(y) where op(A) is A or A^T.
[in] | transx | Should x^T be used? |
[in] | transy | Should y^T be used? |
[in] | alpha | Scalar. |
[in] | x | Left multiplicand. |
[in] | y | Right multiplicand. |
runtime_exception
is thrown.REAL | should be 'float' or 'double'. |
void linalg::matmult | ( | const bool | transx, |
const bool | transy, | ||
const REAL | alpha, | ||
const mpimat< REAL > & | x, | ||
const mpimat< REAL > & | y, | ||
mpimat< REAL > & | ret | ||
) |
Computes ret = alpha*op(x)*op(y) where op(A) is A or A^T.
[in] | transx | Should x^T be used? |
[in] | transy | Should y^T be used? |
[in] | alpha | Scalar. |
[in] | x | Left multiplicand. |
[in] | y | Right multiplicand. |
[out] | ret | The product. |
runtime_exception
is thrown.REAL | should be 'float' or 'double'. |
void linalg::qr | ( | const bool | pivot, |
cpumat< REAL > & | x, | ||
cpuvec< REAL > & | qraux | ||
) |
Computes the QR decomposition.
The factorization works mostly in-place by modifying the input data. After execution, the matrix will be the LAPACK-like compact QR representation.
[in] | pivot | Should the factorization use column pivoting? |
[in,out] | x | Input data matrix. Values are overwritten. |
[out] | qraux | Auxiliary data for compact QR. |
REAL | should be 'float' or 'double'. |
void linalg::qr | ( | const bool | pivot, |
gpumat< REAL > & | x, | ||
gpuvec< REAL > & | qraux | ||
) |
Computes the QR decomposition.
The factorization works mostly in-place by modifying the input data. After execution, the matrix will be the LAPACK-like compact QR representation.
[in] | pivot | NOTE Pivoting does not yet work on GPU. Should the factorization use column pivoting? |
[in,out] | x | Input data matrix. Values are overwritten. |
[out] | qraux | Auxiliary data for compact QR. |
REAL | should be 'float' or 'double'. |
void linalg::qr | ( | const bool | pivot, |
mpimat< REAL > & | x, | ||
cpuvec< REAL > & | qraux | ||
) |
Computes the QR decomposition.
The factorization works mostly in-place by modifying the input data. After execution, the matrix will be the LAPACK-like compact QR representation.
[in] | pivot | Should the factorization use column pivoting? |
[in,out] | x | Input data matrix. Values are overwritten. |
[out] | qraux | Auxiliary data for compact QR. |
REAL | should be 'float' or 'double'. |
void linalg::qr_Q | ( | const cpumat< REAL > & | QR, |
const cpuvec< REAL > & | qraux, | ||
cpumat< REAL > & | Q, | ||
cpuvec< REAL > & | work | ||
) |
Recover the Q matrix from a QR decomposition.
[in] | QR | The compact QR factorization, as computed via qr() . |
[in] | qraux | Auxiliary data for compact QR. |
[out] | Q | The Q matrix. |
[out] | work | Workspace array. Will be resized as necessary. |
REAL | should be 'float' or 'double'. |
void linalg::qr_Q | ( | const gpumat< REAL > & | QR, |
const gpuvec< REAL > & | qraux, | ||
gpumat< REAL > & | Q, | ||
gpuvec< REAL > & | work | ||
) |
Recover the Q matrix from a QR decomposition.
[in] | QR | The compact QR factorization, as computed via qr() . |
[in] | qraux | Auxiliary data for compact QR. |
[out] | Q | The Q matrix. |
[out] | work | Workspace array. Will be resized as necessary. |
REAL | should be 'float' or 'double'. |
void linalg::qr_Q | ( | const mpimat< REAL > & | QR, |
const cpuvec< REAL > & | qraux, | ||
mpimat< REAL > & | Q, | ||
cpuvec< REAL > & | work | ||
) |
Recover the Q matrix from a QR decomposition.
[in] | QR | The compact QR factorization, as computed via qr() . |
[in] | qraux | Auxiliary data for compact QR. |
[out] | Q | The Q matrix. |
[out] | work | Workspace array. Will be resized as necessary. |
REAL | should be 'float' or 'double'. |
Recover the R matrix from a QR decomposition.
[in] | QR | The compact QR factorization, as computed via qr() . |
[out] | R | The R matrix. |
REAL | should be 'float' or 'double'. |
Recover the R matrix from a QR decomposition.
[in] | QR | The compact QR factorization, as computed via qr() . |
[out] | R | The R matrix. |
REAL | should be 'float' or 'double'. |
Recover the R matrix from a QR decomposition.
[in] | QR | The compact QR factorization, as computed via qr() . |
[out] | R | The R matrix. |
REAL | should be 'float' or 'double'. |
Solve a system of equations.
The input is replaced by its LU factorization.
[in,out] | x | Input LHS. Should be square. Overwritten by LU. |
[in,out] | y | Input RHS. Overwritten by solution. |
runtime_error
exception is thrown. If an allocation fails, a bad_alloc
exception will be thrown.REAL | should be 'float' or 'double'. |
Solve a system of equations.
The input is replaced by its PLU factorization.
[in,out] | x | Input LHS. Should be square. Overwritten by LU. |
[in,out] | y | Input RHS. Overwritten by solution. |
cusolverDnXgetrs()
(solve).runtime_error
exception is thrown. If an allocation fails, a bad_alloc
exception will be thrown.REAL | should be '__half', 'float', or 'double'. |
Solve a system of equations.
The input is replaced by its LU factorization.
[in,out] | x | Input LHS. Should be square. Overwritten by LU. |
[in,out] | y | Input RHS. Overwritten by solution. |
runtime_error
exception is thrown. If the inputs are distributed on different grids, a runtime_exception
is thrown. If an allocation fails, a bad_alloc
exception will be thrown.REAL | should be 'float' or 'double'. |
Computes the singular value decomposition.
[in,out] | x | Input data matrix. Values are overwritten. |
[out] | s | Vector of singular values. |
[out] | u | Matrix of left singular vectors. |
[out] | vt | Matrix of (transposed) right singular vectors. |
REAL | should be 'float' or 'double'. |
Computes the singular value decomposition.
[in,out] | x | Input data matrix. Values are overwritten. |
[out] | s | Vector of singular values. |
[out] | u | Matrix of left singular vectors. |
[out] | vt | Matrix of (transposed) right singnular vectors. |
REAL | should be '__half', 'float', or 'double'. |
Computes the singular value decomposition.
[in,out] | x | Input data matrix. Values are overwritten. |
[out] | s | Vector of singular values. |
[out] | u | Matrix of left singular vectors. |
[out] | vt | Matrix of (transposed) right singnular vectors. |
REAL | should be 'float' or 'double'. |
void linalg::tcrossprod | ( | const REAL | alpha, |
const cpumat< REAL > & | x, | ||
cpumat< REAL > & | ret | ||
) |
Computes lower triangle of alpha*x*x^T.
[in] | alpha | Scalar. |
[in] | x | Input data matrix. |
[out] | ret | The product. |
REAL | should be 'float' or 'double'. |
void linalg::tcrossprod | ( | const REAL | alpha, |
const gpumat< REAL > & | x, | ||
gpumat< REAL > & | ret | ||
) |
Computes lower triangle of alpha*x*x^T.
[in] | alpha | Scalar. |
[in] | x | Input data matrix. |
[out] | ret | The product. |
REAL | should be '__half', 'float', or 'double'. |
void linalg::tcrossprod | ( | const REAL | alpha, |
const mpimat< REAL > & | x, | ||
mpimat< REAL > & | ret | ||
) |
Computes lower triangle of alpha*x*x^T.
[in] | alpha | Scalar. |
[in] | x | Input data matrix. |
[out] | ret | The product. |
REAL | should be 'float' or 'double'. |
REAL linalg::trace | ( | const cpumat< REAL > & | x | ) |
Computes the trace, i.e. the sum of the diagonal.
[in] | x | Input data matrix. |
REAL | should be 'float' or 'double'. |
REAL linalg::trace | ( | const gpumat< REAL > & | x | ) |
Computes the trace, i.e. the sum of the diagonal.
[in] | x | Input data matrix. |
REAL | should be '__half', 'float', or 'double'. |
REAL linalg::trace | ( | const mpimat< REAL > & | x | ) |
Computes the trace, i.e. the sum of the diagonal.
[in] | x | Input data matrix. |
REAL | should be 'float' or 'double'. |
void linalg::tssvd | ( | cpumat< REAL > & | x, |
cpuvec< REAL > & | s, | ||
cpumat< REAL > & | u, | ||
cpumat< REAL > & | vt | ||
) |
Computes the singular value decomposition for tall/skinny data. The number of rows must be greater than the number of columns. If the number of rows is not significantly larger than the number of columns, this may not be more efficient than simply calling linalg::svd()
.
The operation works by computing a QR and then taking the SVD of the R matrix. The left singular vectors are Q times the left singular vectors from R's SVD, and the singular value and the right singular vectors are those from R's SVD.
[in,out] | x | Input data matrix. Values are overwritten. |
[out] | s | Vector of singular values. |
[out] | u | Matrix of left singular vectors. |
[out] | vt | Matrix of (transposed) right singular vectors. |
linalg::qr_R()
and linalg::qr_Q()
.REAL | should be 'float' or 'double'. |
void linalg::tssvd | ( | gpumat< REAL > & | x, |
gpuvec< REAL > & | s, | ||
gpumat< REAL > & | u, | ||
gpumat< REAL > & | vt | ||
) |
Computes the singular value decomposition for tall/skinny data. The number of rows must be greater than the number of columns. If the number of rows is not significantly larger than the number of columns, this may not be more efficient than simply calling linalg::svd()
.
The operation works by computing a QR and then taking the SVD of the R matrix. The left singular vectors are Q times the left singular vectors from R's SVD, and the singular value and the right singular vectors are those from R's SVD.
[in,out] | x | Input data matrix. Values are overwritten. |
[out] | s | Vector of singular values. |
[out] | u | Matrix of left singular vectors. |
[out] | vt | Matrix of (transposed) right singular vectors. |
linalg::qr_R()
and linalg::qr_Q()
.REAL | should be 'float' or 'double'. |
void linalg::tssvd | ( | mpimat< REAL > & | x, |
cpuvec< REAL > & | s, | ||
mpimat< REAL > & | u, | ||
mpimat< REAL > & | vt | ||
) |
Computes the singular value decomposition for tall/skinny data. The number of rows must be greater than the number of columns. If the number of rows is not significantly larger than the number of columns, this may not be more efficient than simply calling linalg::svd()
.
The operation works by computing a QR and then taking the SVD of the R matrix. The left singular vectors are Q times the left singular vectors from R's SVD, and the singular value and the right singular vectors are those from R's SVD.
[in,out] | x | Input data matrix. Values are overwritten. |
[out] | s | Vector of singular values. |
[out] | u | Matrix of left singular vectors. |
[out] | vt | Matrix of (transposed) right singular vectors. |
linalg::qr_R()
and linalg::qr_Q()
.REAL | should be 'float' or 'double'. |
Computes the transpose out-of-place (i.e. in a copy).
[in] | x | Input data matrix. |
[out] | tx | The transpose. |
REAL | should be 'float' or 'double'. |
Computes the transpose out-of-place (i.e. in a copy).
[in] | x | Input data matrix. |
[out] | tx | The transpose. |
REAL | should be '__half', 'float', or 'double'. |
Computes the transpose out-of-place (i.e. in a copy).
[in] | x | Input data matrix. |
[out] | tx | The transpose. |
REAL | should be 'float' or 'double'. |