![]() |
MAGMA
1.5.0
Matrix Algebra for GPU and Multicore Architectures
|
Functions | |
magma_int_t | magma_dgetrf (magma_int_t m, magma_int_t n, double *A, magma_int_t lda, magma_int_t *ipiv, magma_int_t *info) |
DGETRF computes an LU factorization of a general M-by-N matrix A using partial pivoting with row interchanges. More... | |
magma_int_t | magma_dgetrf2_mgpu (magma_int_t num_gpus, magma_int_t m, magma_int_t n, magma_int_t nb, magma_int_t offset, double *d_lAT[], magma_int_t lddat, magma_int_t *ipiv, double *d_lAP[], double *w, magma_int_t ldw, magma_queue_t streaml[][2], magma_int_t *info) |
DGETRF computes an LU factorization of a general M-by-N matrix A using partial pivoting with row interchanges. More... | |
magma_int_t | magma_dgetrf_gpu (magma_int_t m, magma_int_t n, double *dA, magma_int_t ldda, magma_int_t *ipiv, magma_int_t *info) |
DGETRF computes an LU factorization of a general M-by-N matrix A using partial pivoting with row interchanges. More... | |
magma_int_t | magma_dgetrf_m (magma_int_t num_gpus0, magma_int_t m, magma_int_t n, double *A, magma_int_t lda, magma_int_t *ipiv, magma_int_t *info) |
DGETRF_m computes an LU factorization of a general M-by-N matrix A using partial pivoting with row interchanges. More... | |
magma_int_t | magma_dgetrf_mgpu (magma_int_t num_gpus, magma_int_t m, magma_int_t n, double **d_lA, magma_int_t ldda, magma_int_t *ipiv, magma_int_t *info) |
DGETRF computes an LU factorization of a general M-by-N matrix A using partial pivoting with row interchanges. More... | |
magma_int_t | magma_dgetrf_nopiv (magma_int_t *m, magma_int_t *n, double *A, magma_int_t *lda, magma_int_t *info) |
DGETRF_NOPIV computes an LU factorization of a general M-by-N matrix A without pivoting. More... | |
magma_int_t | magma_dgetrf_nopiv_gpu (magma_int_t m, magma_int_t n, double *dA, magma_int_t ldda, magma_int_t *info) |
DGETRF_NOPIV_GPU computes an LU factorization of a general M-by-N matrix A without any pivoting. More... | |
magma_int_t | magma_dgetri_gpu (magma_int_t n, double *dA, magma_int_t ldda, magma_int_t *ipiv, double *dwork, magma_int_t lwork, magma_int_t *info) |
DGETRI computes the inverse of a matrix using the LU factorization computed by DGETRF. More... | |
magma_int_t | magma_dgetrs_gpu (magma_trans_t trans, magma_int_t n, magma_int_t nrhs, double *dA, magma_int_t ldda, magma_int_t *ipiv, double *dB, magma_int_t lddb, magma_int_t *info) |
Solves a system of linear equations A * X = B, A**T * X = B, or A**T * X = B with a general N-by-N matrix A using the LU factorization computed by DGETRF_GPU. More... | |
magma_int_t | magma_dsgetrs_gpu (magma_trans_t trans, magma_int_t n, magma_int_t nrhs, float *dA, magma_int_t ldda, magma_int_t *dipiv, double *dB, magma_int_t lddb, double *dX, magma_int_t lddx, float *dSX, magma_int_t *info) |
DSGETRS solves a system of linear equations A * X = B, A**T * X = B, or A**T * X = B with a general N-by-N matrix A using the LU factorization computed by MAGMA_SGETRF_GPU. More... | |
magma_int_t magma_dgetrf | ( | magma_int_t | m, |
magma_int_t | n, | ||
double * | A, | ||
magma_int_t | lda, | ||
magma_int_t * | ipiv, | ||
magma_int_t * | info | ||
) |
DGETRF computes an LU factorization of a general M-by-N matrix A using partial pivoting with row interchanges.
This version does not require work space on the GPU passed as input. GPU memory is allocated in the routine.
The factorization has the form A = P * L * U where P is a permutation matrix, L is lower triangular with unit diagonal elements (lower trapezoidal if m > n), and U is upper triangular (upper trapezoidal if m < n).
This is the right-looking Level 3 BLAS version of the algorithm. If the current stream is NULL, this version replaces it with user defined stream to overlap computation with communication.
[in] | m | INTEGER The number of rows of the matrix A. M >= 0. |
[in] | n | INTEGER The number of columns of the matrix A. N >= 0. |
[in,out] | A | DOUBLE_PRECISION array, dimension (LDA,N) On entry, the M-by-N matrix to be factored. On exit, the factors L and U from the factorization A = P*L*U; the unit diagonal elements of L are not stored. Higher performance is achieved if A is in pinned memory, e.g. allocated using magma_malloc_pinned. |
[in] | lda | INTEGER The leading dimension of the array A. LDA >= max(1,M). |
[out] | ipiv | INTEGER array, dimension (min(M,N)) The pivot indices; for 1 <= i <= min(M,N), row i of the matrix was interchanged with row IPIV(i). |
[out] | info | INTEGER
|
magma_int_t magma_dgetrf2_mgpu | ( | magma_int_t | num_gpus, |
magma_int_t | m, | ||
magma_int_t | n, | ||
magma_int_t | nb, | ||
magma_int_t | offset, | ||
double * | d_lAT[], | ||
magma_int_t | lddat, | ||
magma_int_t * | ipiv, | ||
double * | d_lAP[], | ||
double * | w, | ||
magma_int_t | ldw, | ||
magma_queue_t | streaml[][2], | ||
magma_int_t * | info | ||
) |
DGETRF computes an LU factorization of a general M-by-N matrix A using partial pivoting with row interchanges.
The factorization has the form A = P * L * U where P is a permutation matrix, L is lower triangular with unit diagonal elements (lower trapezoidal if m > n), and U is upper triangular (upper trapezoidal if m < n).
This is the right-looking Level 3 BLAS version of the algorithm. Use two buffer to send panels.
[in] | num_gpus | INTEGER The number of GPUs to be used for the factorization. |
[in] | m | INTEGER The number of rows of the matrix A. M >= 0. |
[in] | n | INTEGER The number of columns of the matrix A. N >= 0. |
[in,out] | A | DOUBLE_PRECISION array on the GPU, dimension (LDDA,N). On entry, the M-by-N matrix to be factored. On exit, the factors L and U from the factorization A = P*L*U; the unit diagonal elements of L are not stored. |
[in] | ldda | INTEGER The leading dimension of the array A. LDDA >= max(1,M). |
[out] | ipiv | INTEGER array, dimension (min(M,N)) The pivot indices; for 1 <= i <= min(M,N), row i of the matrix was interchanged with row IPIV(i). |
[out] | info | INTEGER
|
magma_int_t magma_dgetrf_gpu | ( | magma_int_t | m, |
magma_int_t | n, | ||
double * | dA, | ||
magma_int_t | ldda, | ||
magma_int_t * | ipiv, | ||
magma_int_t * | info | ||
) |
DGETRF computes an LU factorization of a general M-by-N matrix A using partial pivoting with row interchanges.
The factorization has the form A = P * L * U where P is a permutation matrix, L is lower triangular with unit diagonal elements (lower trapezoidal if m > n), and U is upper triangular (upper trapezoidal if m < n).
This is the right-looking Level 3 BLAS version of the algorithm. If the current stream is NULL, this version replaces it with user defined stream to overlap computation with communication.
[in] | m | INTEGER The number of rows of the matrix A. M >= 0. |
[in] | n | INTEGER The number of columns of the matrix A. N >= 0. |
[in,out] | dA | DOUBLE_PRECISION array on the GPU, dimension (LDDA,N). On entry, the M-by-N matrix to be factored. On exit, the factors L and U from the factorization A = P*L*U; the unit diagonal elements of L are not stored. |
[in] | ldda | INTEGER The leading dimension of the array A. LDDA >= max(1,M). |
[out] | ipiv | INTEGER array, dimension (min(M,N)) The pivot indices; for 1 <= i <= min(M,N), row i of the matrix was interchanged with row IPIV(i). |
[out] | info | INTEGER
|
magma_int_t magma_dgetrf_m | ( | magma_int_t | num_gpus0, |
magma_int_t | m, | ||
magma_int_t | n, | ||
double * | A, | ||
magma_int_t | lda, | ||
magma_int_t * | ipiv, | ||
magma_int_t * | info | ||
) |
DGETRF_m computes an LU factorization of a general M-by-N matrix A using partial pivoting with row interchanges.
This version does not require work space on the GPU passed as input. GPU memory is allocated in the routine. The matrix may not fit entirely in the GPU memory.
The factorization has the form A = P * L * U where P is a permutation matrix, L is lower triangular with unit diagonal elements (lower trapezoidal if m > n), and U is upper triangular (upper trapezoidal if m < n).
This is the right-looking Level 3 BLAS version of the algorithm.
Note: The factorization of big panel is done calling multiple-gpu-interface. Pivots are applied on GPU within the big panel.
[in] | m | INTEGER The number of rows of the matrix A. M >= 0. |
[in] | n | INTEGER The number of columns of the matrix A. N >= 0. |
[in,out] | A | DOUBLE_PRECISION array, dimension (LDA,N) On entry, the M-by-N matrix to be factored. On exit, the factors L and U from the factorization A = P*L*U; the unit diagonal elements of L are not stored. Higher performance is achieved if A is in pinned memory, e.g. allocated using magma_malloc_pinned. |
[in] | lda | INTEGER The leading dimension of the array A. LDA >= max(1,M). |
[out] | ipiv | INTEGER array, dimension (min(M,N)) The pivot indices; for 1 <= i <= min(M,N), row i of the matrix was interchanged with row IPIV(i). |
[out] | info | INTEGER
|
magma_int_t magma_dgetrf_mgpu | ( | magma_int_t | num_gpus, |
magma_int_t | m, | ||
magma_int_t | n, | ||
double ** | d_lA, | ||
magma_int_t | ldda, | ||
magma_int_t * | ipiv, | ||
magma_int_t * | info | ||
) |
DGETRF computes an LU factorization of a general M-by-N matrix A using partial pivoting with row interchanges.
The factorization has the form A = P * L * U where P is a permutation matrix, L is lower triangular with unit diagonal elements (lower trapezoidal if m > n), and U is upper triangular (upper trapezoidal if m < n).
This is the right-looking Level 3 BLAS version of the algorithm.
[in] | num_gpus | INTEGER The number of GPUs to be used for the factorization. |
[in] | m | INTEGER The number of rows of the matrix A. M >= 0. |
[in] | n | INTEGER The number of columns of the matrix A. N >= 0. |
[in,out] | A | DOUBLE_PRECISION array on the GPU, dimension (LDDA,N). On entry, the M-by-N matrix to be factored. On exit, the factors L and U from the factorization A = P*L*U; the unit diagonal elements of L are not stored. |
[in] | ldda | INTEGER The leading dimension of the array A. LDDA >= max(1,M). |
[out] | ipiv | INTEGER array, dimension (min(M,N)) The pivot indices; for 1 <= i <= min(M,N), row i of the matrix was interchanged with row IPIV(i). |
[out] | info | INTEGER
|
magma_int_t magma_dgetrf_nopiv | ( | magma_int_t * | m, |
magma_int_t * | n, | ||
double * | A, | ||
magma_int_t * | lda, | ||
magma_int_t * | info | ||
) |
DGETRF_NOPIV computes an LU factorization of a general M-by-N matrix A without pivoting.
The factorization has the form A = L * U where L is lower triangular with unit diagonal elements (lower trapezoidal if m > n), and U is upper triangular (upper trapezoidal if m < n).
This is the right-looking Level 3 BLAS version of the algorithm.
[in] | m | INTEGER The number of rows of the matrix A. M >= 0. |
[in] | n | INTEGER The number of columns of the matrix A. N >= 0. |
[in,out] | A | DOUBLE_PRECISION array, dimension (LDA,N) On entry, the M-by-N matrix to be factored. On exit, the factors L and U from the factorization A = P*L*U; the unit diagonal elements of L are not stored. |
[in] | lda | INTEGER The leading dimension of the array A. LDA >= max(1,M). |
[out] | info | INTEGER
|
magma_int_t magma_dgetrf_nopiv_gpu | ( | magma_int_t | m, |
magma_int_t | n, | ||
double * | dA, | ||
magma_int_t | ldda, | ||
magma_int_t * | info | ||
) |
DGETRF_NOPIV_GPU computes an LU factorization of a general M-by-N matrix A without any pivoting.
The factorization has the form A = P * L * U where P is a permutation matrix, L is lower triangular with unit diagonal elements (lower trapezoidal if m > n), and U is upper triangular (upper trapezoidal if m < n).
This is the right-looking Level 3 BLAS version of the algorithm.
[in] | m | INTEGER The number of rows of the matrix A. M >= 0. |
[in] | n | INTEGER The number of columns of the matrix A. N >= 0. |
[in,out] | dA | DOUBLE_PRECISION array on the GPU, dimension (LDDA,N). On entry, the M-by-N matrix to be factored. On exit, the factors L and U from the factorization A = P*L*U; the unit diagonal elements of L are not stored. |
[in] | ldda | INTEGER The leading dimension of the array A. LDDA >= max(1,M). |
[out] | info | INTEGER
|
magma_int_t magma_dgetri_gpu | ( | magma_int_t | n, |
double * | dA, | ||
magma_int_t | ldda, | ||
magma_int_t * | ipiv, | ||
double * | dwork, | ||
magma_int_t | lwork, | ||
magma_int_t * | info | ||
) |
DGETRI computes the inverse of a matrix using the LU factorization computed by DGETRF.
This method inverts U and then computes inv(A) by solving the system inv(A)*L = inv(U) for inv(A).
Note that it is generally both faster and more accurate to use DGESV, or DGETRF and DGETRS, to solve the system AX = B, rather than inverting the matrix and multiplying to form X = inv(A)*B. Only in special instances should an explicit inverse be computed with this routine.
[in] | n | INTEGER The order of the matrix A. N >= 0. |
[in,out] | dA | DOUBLE_PRECISION array on the GPU, dimension (LDDA,N) On entry, the factors L and U from the factorization A = P*L*U as computed by DGETRF_GPU. On exit, if INFO = 0, the inverse of the original matrix A. |
[in] | ldda | INTEGER The leading dimension of the array A. LDDA >= max(1,N). |
[in] | ipiv | INTEGER array, dimension (N) The pivot indices from DGETRF; for 1 <= i <= N, row i of the matrix was interchanged with row IPIV(i). |
[out] | dwork | (workspace) DOUBLE_PRECISION array on the GPU, dimension (MAX(1,LWORK)) |
[in] | lwork | INTEGER The dimension of the array DWORK. LWORK >= N*NB, where NB is the optimal blocksize returned by magma_get_dgetri_nb(n). Unlike LAPACK, this version does not currently support a workspace query, because the workspace is on the GPU. |
[out] | info | INTEGER
|
magma_int_t magma_dgetrs_gpu | ( | magma_trans_t | trans, |
magma_int_t | n, | ||
magma_int_t | nrhs, | ||
double * | dA, | ||
magma_int_t | ldda, | ||
magma_int_t * | ipiv, | ||
double * | dB, | ||
magma_int_t | lddb, | ||
magma_int_t * | info | ||
) |
Solves a system of linear equations A * X = B, A**T * X = B, or A**T * X = B with a general N-by-N matrix A using the LU factorization computed by DGETRF_GPU.
[in] | trans | magma_trans_t Specifies the form of the system of equations:
|
[in] | n | INTEGER The order of the matrix A. N >= 0. |
[in] | nrhs | INTEGER The number of right hand sides, i.e., the number of columns of the matrix B. NRHS >= 0. |
[in] | dA | DOUBLE_PRECISION array on the GPU, dimension (LDA,N) The factors L and U from the factorization A = P*L*U as computed by DGETRF_GPU. |
[in] | ldda | INTEGER The leading dimension of the array A. LDA >= max(1,N). |
[in] | ipiv | INTEGER array, dimension (N) The pivot indices from DGETRF; for 1 <= i <= N, row i of the matrix was interchanged with row IPIV(i). |
[in,out] | dB | DOUBLE_PRECISION array on the GPU, dimension (LDB,NRHS) On entry, the right hand side matrix B. On exit, the solution matrix X. |
[in] | lddb | INTEGER The leading dimension of the array B. LDB >= max(1,N). |
[out] | info | INTEGER
|
magma_int_t magma_dsgetrs_gpu | ( | magma_trans_t | trans, |
magma_int_t | n, | ||
magma_int_t | nrhs, | ||
float * | dA, | ||
magma_int_t | ldda, | ||
magma_int_t * | dipiv, | ||
double * | dB, | ||
magma_int_t | lddb, | ||
double * | dX, | ||
magma_int_t | lddx, | ||
float * | dSX, | ||
magma_int_t * | info | ||
) |
DSGETRS solves a system of linear equations A * X = B, A**T * X = B, or A**T * X = B with a general N-by-N matrix A using the LU factorization computed by MAGMA_SGETRF_GPU.
B and X are in DOUBLE PRECISION, and A is in SINGLE PRECISION. This routine is used in the mixed precision iterative solver magma_dsgesv.
[in] | trans | magma_trans_t Specifies the form of the system of equations:
|
[in] | n | INTEGER The order of the matrix A. N >= 0. |
[in] | nrhs | INTEGER The number of right hand sides, i.e., the number of columns of the matrix B. NRHS >= 0. |
[in] | dA | SINGLE PRECISION array on the GPU, dimension (LDDA,N) The factors L and U from the factorization A = P*L*U as computed by CGETRF_GPU. |
[in] | ldda | INTEGER The leading dimension of the array dA. LDDA >= max(1,N). |
[in] | dipiv | INTEGER array on the GPU, dimension (N) The pivot indices; for 1 <= i <= N, after permuting, row i of the matrix was moved to row dIPIV(i). Note this is different than IPIV from DGETRF, where interchanges are applied one-after-another. |
[in] | dB | DOUBLE PRECISION array on the GPU, dimension (LDDB,NRHS) On entry, the right hand side matrix B. |
[in] | lddb | INTEGER The leading dimension of the arrays X and B. LDDB >= max(1,N). |
[out] | dX | DOUBLE PRECISION array on the GPU, dimension (LDDX, NRHS) On exit, the solution matrix dX. |
[in] | lddx | INTEGER The leading dimension of the array dX, LDDX >= max(1,N). |
dSX | (workspace) SINGLE PRECISION array on the GPU used as workspace, dimension (N, NRHS) | |
[out] | info | INTEGER
|