![]() |
MAGMA
1.5.0
Matrix Algebra for GPU and Multicore Architectures
|
Functions | |
magma_int_t | magma_dpotrf (magma_uplo_t uplo, magma_int_t n, double *A, magma_int_t lda, magma_int_t *info) |
DPOTRF computes the Cholesky factorization of a real symmetric positive definite matrix A. More... | |
magma_int_t | magma_dpotrf2_mgpu (int num_gpus, magma_uplo_t uplo, magma_int_t m, magma_int_t n, magma_int_t off_i, magma_int_t off_j, magma_int_t nb, double **d_lA, magma_int_t ldda, double **d_lP, magma_int_t lddp, double *A, magma_int_t lda, magma_int_t h, magma_queue_t stream[][3], magma_event_t event[][5], magma_int_t *info) |
DPOTRF computes the Cholesky factorization of a real symmetric positive definite matrix dA. More... | |
magma_int_t | magma_dpotrf3_mgpu (magma_int_t num_gpus, magma_uplo_t uplo, magma_int_t m, magma_int_t n, magma_int_t off_i, magma_int_t off_j, magma_int_t nb, double *d_lA[], magma_int_t ldda, double *d_lP[], magma_int_t lddp, double *A, magma_int_t lda, magma_int_t h, magma_queue_t stream[][3], magma_event_t event[][5], magma_int_t *info) |
DPOTRF computes the Cholesky factorization of a real symmetric positive definite matrix dA. More... | |
magma_int_t | magma_dpotrf_gpu (magma_uplo_t uplo, magma_int_t n, double *dA, magma_int_t ldda, magma_int_t *info) |
DPOTRF computes the Cholesky factorization of a real symmetric positive definite matrix dA. More... | |
magma_int_t | magma_dpotrf_m (magma_int_t num_gpus, magma_uplo_t uplo, magma_int_t n, double *A, magma_int_t lda, magma_int_t *info) |
DPOTRF_OOC computes the Cholesky factorization of a real symmetric positive definite matrix A. More... | |
magma_int_t | magma_dpotrf_mgpu (magma_int_t num_gpus, magma_uplo_t uplo, magma_int_t n, double **d_lA, magma_int_t ldda, magma_int_t *info) |
DPOTRF computes the Cholesky factorization of a real symmetric positive definite matrix dA. More... | |
magma_int_t | magma_dpotri (magma_uplo_t uplo, magma_int_t n, double *A, magma_int_t lda, magma_int_t *info) |
DPOTRI computes the inverse of a real symmetric positive definite matrix A using the Cholesky factorization A = U**T*U or A = L*L**T computed by DPOTRF. More... | |
magma_int_t | magma_dpotri_gpu (magma_uplo_t uplo, magma_int_t n, double *dA, magma_int_t ldda, magma_int_t *info) |
DPOTRI computes the inverse of a real symmetric positive definite matrix A using the Cholesky factorization A = U**T*U or A = L*L**T computed by DPOTRF. More... | |
magma_int_t | magma_dpotrs_gpu (magma_uplo_t uplo, magma_int_t n, magma_int_t nrhs, double *dA, magma_int_t ldda, double *dB, magma_int_t lddb, magma_int_t *info) |
DPOTRS solves a system of linear equations A*X = B with a symmetric positive definite matrix A using the Cholesky factorization A = U**T*U or A = L*L**T computed by DPOTRF. More... | |
magma_int_t magma_dpotrf | ( | magma_uplo_t | uplo, |
magma_int_t | n, | ||
double * | A, | ||
magma_int_t | lda, | ||
magma_int_t * | info | ||
) |
DPOTRF computes the Cholesky factorization of a real symmetric positive definite matrix A.
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 = U**T * U, if uplo = MagmaUpper, or A = L * L**T, if uplo = MagmaLower, where U is an upper triangular matrix and L is lower triangular.
This is the block version of the algorithm, calling Level 3 BLAS. If the current stream is NULL, this version replaces it with user defined stream to overlap computation with communication.
[in] | uplo | magma_uplo_t
|
[in] | n | INTEGER The order of the matrix A. N >= 0. |
[in,out] | A | DOUBLE_PRECISION array, dimension (LDA,N) On entry, the symmetric matrix A. If uplo = MagmaUpper, the leading N-by-N upper triangular part of A contains the upper triangular part of the matrix A, and the strictly lower triangular part of A is not referenced. If uplo = MagmaLower, the leading N-by-N lower triangular part of A contains the lower triangular part of the matrix A, and the strictly upper triangular part of A is not referenced. On exit, if INFO = 0, the factor U or L from the Cholesky factorization A = U**T * U or A = L * L**T. 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,N). |
[out] | info | INTEGER
|
magma_int_t magma_dpotrf2_mgpu | ( | int | num_gpus, |
magma_uplo_t | uplo, | ||
magma_int_t | m, | ||
magma_int_t | n, | ||
magma_int_t | off_i, | ||
magma_int_t | off_j, | ||
magma_int_t | nb, | ||
double ** | d_lA, | ||
magma_int_t | ldda, | ||
double ** | d_lP, | ||
magma_int_t | lddp, | ||
double * | A, | ||
magma_int_t | lda, | ||
magma_int_t | h, | ||
magma_queue_t | stream[][3], | ||
magma_event_t | event[][5], | ||
magma_int_t * | info | ||
) |
DPOTRF computes the Cholesky factorization of a real symmetric positive definite matrix dA.
The factorization has the form dA = U**T * U, if UPLO = MagmaUpper, or dA = L * L**T, if UPLO = MagmaLower, where U is an upper triangular matrix and L is lower triangular.
This is the block version of the algorithm, calling Level 3 BLAS.
[in] | uplo | magma_uplo_t
|
[in] | n | INTEGER The order of the matrix dA. N >= 0. |
[in,out] | dA | DOUBLE_PRECISION array on the GPU, dimension (LDDA,N) On entry, the symmetric matrix dA. If UPLO = MagmaUpper, the leading N-by-N upper triangular part of dA contains the upper triangular part of the matrix dA, and the strictly lower triangular part of dA is not referenced. If UPLO = MagmaLower, the leading N-by-N lower triangular part of dA contains the lower triangular part of the matrix dA, and the strictly upper triangular part of dA is not referenced. On exit, if INFO = 0, the factor U or L from the Cholesky factorization dA = U**T * U or dA = L * L**T. |
[in] | ldda | INTEGER The leading dimension of the array dA. LDDA >= max(1,N). To benefit from coalescent memory accesses LDDA must be divisible by 16. |
[out] | info | INTEGER
|
magma_int_t magma_dpotrf3_mgpu | ( | magma_int_t | num_gpus, |
magma_uplo_t | uplo, | ||
magma_int_t | m, | ||
magma_int_t | n, | ||
magma_int_t | off_i, | ||
magma_int_t | off_j, | ||
magma_int_t | nb, | ||
double * | d_lA[], | ||
magma_int_t | ldda, | ||
double * | d_lP[], | ||
magma_int_t | lddp, | ||
double * | A, | ||
magma_int_t | lda, | ||
magma_int_t | h, | ||
magma_queue_t | stream[][3], | ||
magma_event_t | event[][5], | ||
magma_int_t * | info | ||
) |
DPOTRF computes the Cholesky factorization of a real symmetric positive definite matrix dA.
Auxiliary subroutine for dpotrf2_ooc. It is multiple gpu interface to compute Cholesky of a "rectangular" matrix.
The factorization has the form dA = U**T * U, if UPLO = MagmaUpper, or dA = L * L**T, if UPLO = MagmaLower, where U is an upper triangular matrix and L is lower triangular.
This is the block version of the algorithm, calling Level 3 BLAS.
[in] | uplo | magma_uplo_t
|
[in] | n | INTEGER The order of the matrix dA. N >= 0. |
[in,out] | dA | DOUBLE_PRECISION array on the GPU, dimension (LDDA,N) On entry, the symmetric matrix dA. If UPLO = MagmaUpper, the leading N-by-N upper triangular part of dA contains the upper triangular part of the matrix dA, and the strictly lower triangular part of dA is not referenced. If UPLO = MagmaLower, the leading N-by-N lower triangular part of dA contains the lower triangular part of the matrix dA, and the strictly upper triangular part of dA is not referenced. On exit, if INFO = 0, the factor U or L from the Cholesky factorization dA = U**T * U or dA = L * L**T. |
[in] | ldda | INTEGER The leading dimension of the array dA. LDDA >= max(1,N). To benefit from coalescent memory accesses LDDA must be divisible by 16. |
[out] | info | INTEGER
|
magma_int_t magma_dpotrf_gpu | ( | magma_uplo_t | uplo, |
magma_int_t | n, | ||
double * | dA, | ||
magma_int_t | ldda, | ||
magma_int_t * | info | ||
) |
DPOTRF computes the Cholesky factorization of a real symmetric positive definite matrix dA.
The factorization has the form dA = U**T * U, if UPLO = MagmaUpper, or dA = L * L**T, if UPLO = MagmaLower, where U is an upper triangular matrix and L is lower triangular.
This is the block version of the algorithm, calling Level 3 BLAS. If the current stream is NULL, this version replaces it with user defined stream to overlap computation with communication.
[in] | uplo | magma_uplo_t
|
[in] | n | INTEGER The order of the matrix dA. N >= 0. |
[in,out] | dA | DOUBLE_PRECISION array on the GPU, dimension (LDDA,N) On entry, the symmetric matrix dA. If UPLO = MagmaUpper, the leading N-by-N upper triangular part of dA contains the upper triangular part of the matrix dA, and the strictly lower triangular part of dA is not referenced. If UPLO = MagmaLower, the leading N-by-N lower triangular part of dA contains the lower triangular part of the matrix dA, and the strictly upper triangular part of dA is not referenced. On exit, if INFO = 0, the factor U or L from the Cholesky factorization dA = U**T * U or dA = L * L**T. |
[in] | ldda | INTEGER The leading dimension of the array dA. LDDA >= max(1,N). To benefit from coalescent memory accesses LDDA must be divisible by 16. |
[out] | info | INTEGER
|
magma_int_t magma_dpotrf_m | ( | magma_int_t | num_gpus, |
magma_uplo_t | uplo, | ||
magma_int_t | n, | ||
double * | A, | ||
magma_int_t | lda, | ||
magma_int_t * | info | ||
) |
DPOTRF_OOC computes the Cholesky factorization of a real symmetric positive definite matrix A.
This version does not require work space on the GPU passed as input. GPU memory is allocated in the routine. The matrix A may not fit entirely in the GPU memory.
The factorization has the form A = U**T * U, if UPLO = MagmaUpper, or A = L * L**T, if UPLO = MagmaLower, where U is an upper triangular matrix and L is lower triangular.
This is the block version of the algorithm, calling Level 3 BLAS.
[in] | num_gpus | INTEGER The number of GPUs. num_gpus > 0. |
[in] | uplo | magma_uplo_t
|
[in] | n | INTEGER The order of the matrix A. N >= 0. |
[in,out] | A | DOUBLE_PRECISION array, dimension (LDA,N) On entry, the symmetric matrix A. If UPLO = MagmaUpper, the leading N-by-N upper triangular part of A contains the upper triangular part of the matrix A, and the strictly lower triangular part of A is not referenced. If UPLO = MagmaLower, the leading N-by-N lower triangular part of A contains the lower triangular part of the matrix A, and the strictly upper triangular part of A is not referenced. On exit, if INFO = 0, the factor U or L from the Cholesky factorization A = U**T * U or A = L * L**T. 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,N). |
[out] | info | INTEGER
|
magma_int_t magma_dpotrf_mgpu | ( | magma_int_t | num_gpus, |
magma_uplo_t | uplo, | ||
magma_int_t | n, | ||
double ** | d_lA, | ||
magma_int_t | ldda, | ||
magma_int_t * | info | ||
) |
DPOTRF computes the Cholesky factorization of a real symmetric positive definite matrix dA.
The factorization has the form dA = U**T * U, if UPLO = MagmaUpper, or dA = L * L**T, if UPLO = MagmaLower, where U is an upper triangular matrix and L is lower triangular.
This is the block version of the algorithm, calling Level 3 BLAS.
[in] | uplo | magma_uplo_t
|
[in] | n | INTEGER The order of the matrix dA. N >= 0. |
[in,out] | dA | DOUBLE_PRECISION array on the GPU, dimension (LDDA,N) On entry, the symmetric matrix dA. If UPLO = MagmaUpper, the leading N-by-N upper triangular part of dA contains the upper triangular part of the matrix dA, and the strictly lower triangular part of dA is not referenced. If UPLO = MagmaLower, the leading N-by-N lower triangular part of dA contains the lower triangular part of the matrix dA, and the strictly upper triangular part of dA is not referenced. On exit, if INFO = 0, the factor U or L from the Cholesky factorization dA = U**T * U or dA = L * L**T. |
[in] | ldda | INTEGER The leading dimension of the array dA. LDDA >= max(1,N). To benefit from coalescent memory accesses LDDA must be divisible by 16. |
[out] | info | INTEGER
|
magma_int_t magma_dpotri | ( | magma_uplo_t | uplo, |
magma_int_t | n, | ||
double * | A, | ||
magma_int_t | lda, | ||
magma_int_t * | info | ||
) |
DPOTRI computes the inverse of a real symmetric positive definite matrix A using the Cholesky factorization A = U**T*U or A = L*L**T computed by DPOTRF.
[in] | uplo | magma_uplo_t
|
[in] | n | INTEGER The order of the matrix A. N >= 0. |
[in,out] | A | DOUBLE_PRECISION array, dimension (LDA,N) On entry, the triangular factor U or L from the Cholesky factorization A = U**T*U or A = L*L**T, as computed by DPOTRF. On exit, the upper or lower triangle of the (symmetric) inverse of A, overwriting the input factor U or L. |
[in] | lda | INTEGER The leading dimension of the array A. LDA >= max(1,N). |
[out] | info | INTEGER
|
magma_int_t magma_dpotri_gpu | ( | magma_uplo_t | uplo, |
magma_int_t | n, | ||
double * | dA, | ||
magma_int_t | ldda, | ||
magma_int_t * | info | ||
) |
DPOTRI computes the inverse of a real symmetric positive definite matrix A using the Cholesky factorization A = U**T*U or A = L*L**T computed by DPOTRF.
[in] | uplo | magma_uplo_t
|
[in] | n | INTEGER The order of the matrix A. N >= 0. |
[in,out] | dA | DOUBLE_PRECISION array on the GPU, dimension (LDA,N) On entry, the triangular factor U or L from the Cholesky factorization A = U**T*U or A = L*L**T, as computed by DPOTRF. On exit, the upper or lower triangle of the (symmetric) inverse of A, overwriting the input factor U or L. |
[in] | ldda | INTEGER The leading dimension of the array dA. LDDA >= max(1,N). |
[out] | info | INTEGER
|
magma_int_t magma_dpotrs_gpu | ( | magma_uplo_t | uplo, |
magma_int_t | n, | ||
magma_int_t | nrhs, | ||
double * | dA, | ||
magma_int_t | ldda, | ||
double * | dB, | ||
magma_int_t | lddb, | ||
magma_int_t * | info | ||
) |
DPOTRS solves a system of linear equations A*X = B with a symmetric positive definite matrix A using the Cholesky factorization A = U**T*U or A = L*L**T computed by DPOTRF.
[in] | uplo | magma_uplo_t
|
[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 (LDDA,N) The triangular factor U or L from the Cholesky factorization A = U**T*U or A = L*L**T, as computed by DPOTRF. |
[in] | ldda | INTEGER The leading dimension of the array A. LDDA >= max(1,N). |
[in,out] | dB | DOUBLE_PRECISION array on the GPU, dimension (LDDB,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. LDDB >= max(1,N). |
[out] | info | INTEGER
|