![]() |
MAGMA
2.0.0
Matrix Algebra for GPU and Multicore Architectures
|
Functions | |
magma_int_t | magma_dposv (magma_uplo_t uplo, magma_int_t n, magma_int_t nrhs, double *A, magma_int_t lda, double *B, magma_int_t ldb, magma_int_t *info) |
DPOSV computes the solution to a real system of linear equations A * X = B, where A is an N-by-N symmetric positive definite matrix and X and B are N-by-NRHS matrices. More... | |
magma_int_t | magma_dposv_batched (magma_uplo_t uplo, magma_int_t n, magma_int_t nrhs, double **dA_array, magma_int_t ldda, double **dB_array, magma_int_t lddb, magma_int_t *dinfo_array, magma_int_t batchCount, magma_queue_t queue) |
DPOSV computes the solution to a real system of linear equations A * X = B, where A is an N-by-N symmetric positive definite matrix and X and B are N-by-NRHS matrices. More... | |
magma_int_t | magma_dposv_gpu (magma_uplo_t uplo, magma_int_t n, magma_int_t nrhs, magmaDouble_ptr dA, magma_int_t ldda, magmaDouble_ptr dB, magma_int_t lddb, magma_int_t *info) |
DPOSV computes the solution to a real system of linear equations A * X = B, where A is an N-by-N symmetric positive definite matrix and X and B are N-by-NRHS matrices. More... | |
magma_int_t | magma_dsposv_gpu (magma_uplo_t uplo, magma_int_t n, magma_int_t nrhs, magmaDouble_ptr dA, magma_int_t ldda, magmaDouble_ptr dB, magma_int_t lddb, magmaDouble_ptr dX, magma_int_t lddx, magmaDouble_ptr dworkd, magmaFloat_ptr dworks, magma_int_t *iter, magma_int_t *info) |
DSPOSV computes the solution to a real system of linear equations A * X = B, where A is an N-by-N symmetric positive definite matrix and X and B are N-by-NRHS matrices. More... | |
magma_int_t magma_dposv | ( | magma_uplo_t | uplo, |
magma_int_t | n, | ||
magma_int_t | nrhs, | ||
double * | A, | ||
magma_int_t | lda, | ||
double * | B, | ||
magma_int_t | ldb, | ||
magma_int_t * | info | ||
) |
DPOSV computes the solution to a real system of linear equations A * X = B, where A is an N-by-N symmetric positive definite matrix and X and B are N-by-NRHS matrices.
The Cholesky decomposition is used to factor A as A = U**H * U, if UPLO = MagmaUpper, or A = L * L**H, if UPLO = MagmaLower, where U is an upper triangular matrix and L is a lower triangular matrix. The factored form of A is then used to solve the system of equations A * X = B.
[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,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**H*U or A = L*L**H. |
[in] | lda | INTEGER The leading dimension of the array A. LDA >= max(1,N). |
[in,out] | B | DOUBLE_PRECISION array, dimension (LDB,NRHS) On entry, the right hand side matrix B. On exit, the solution matrix X. |
[in] | ldb | INTEGER The leading dimension of the array B. LDB >= max(1,N). |
[out] | info | INTEGER
|
magma_int_t magma_dposv_batched | ( | magma_uplo_t | uplo, |
magma_int_t | n, | ||
magma_int_t | nrhs, | ||
double ** | dA_array, | ||
magma_int_t | ldda, | ||
double ** | dB_array, | ||
magma_int_t | lddb, | ||
magma_int_t * | dinfo_array, | ||
magma_int_t | batchCount, | ||
magma_queue_t | queue | ||
) |
DPOSV computes the solution to a real system of linear equations A * X = B, where A is an N-by-N symmetric positive definite matrix and X and B are N-by-NRHS matrices.
The Cholesky decomposition is used to factor A as A = U**H * U, if UPLO = MagmaUpper, or A = L * L**H, if UPLO = MagmaLower, where U is an upper triangular matrix and L is a lower triangular matrix. The factored form of A is then used to solve the system of equations A * X = B.
[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,out] | dA_array | Array of pointers, dimension (batchCount). Each is a DOUBLE_PRECISION array on the GPU, dimension (LDDA,N) On entry, each pointer is a 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 dA 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 corresponding entry in dinfo_array = 0, each pointer is the factor U or L from the Cholesky factorization A = U**H*U or A = L*L**H. |
[in] | ldda | INTEGER The leading dimension of each array A. LDA >= max(1,N). |
[in,out] | dB_array | Array of pointers, dimension (batchCount). Each is a DOUBLE_PRECISION array on the GPU, dimension (LDB,NRHS) On entry, each pointer is a right hand side matrix B. On exit, each pointer is the corresponding solution matrix X. |
[in] | lddb | INTEGER The leading dimension of each array B. LDB >= max(1,N). |
[out] | dinfo_array | Array of INTEGERs, dimension (batchCount), for corresponding matrices.
|
[in] | batchCount | INTEGER The number of matrices to operate on. |
[in] | queue | magma_queue_t Queue to execute in. |
magma_int_t magma_dposv_gpu | ( | magma_uplo_t | uplo, |
magma_int_t | n, | ||
magma_int_t | nrhs, | ||
magmaDouble_ptr | dA, | ||
magma_int_t | ldda, | ||
magmaDouble_ptr | dB, | ||
magma_int_t | lddb, | ||
magma_int_t * | info | ||
) |
DPOSV computes the solution to a real system of linear equations A * X = B, where A is an N-by-N symmetric positive definite matrix and X and B are N-by-NRHS matrices.
The Cholesky decomposition is used to factor A as A = U**H * U, if UPLO = MagmaUpper, or A = L * L**H, if UPLO = MagmaLower, where U is an upper triangular matrix and L is a lower triangular matrix. The factored form of A is then used to solve the system of equations A * X = B.
[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,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**H*U or dA = L*L**H. |
[in] | ldda | INTEGER The leading dimension of the array A. LDA >= max(1,N). |
[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_dsposv_gpu | ( | magma_uplo_t | uplo, |
magma_int_t | n, | ||
magma_int_t | nrhs, | ||
magmaDouble_ptr | dA, | ||
magma_int_t | ldda, | ||
magmaDouble_ptr | dB, | ||
magma_int_t | lddb, | ||
magmaDouble_ptr | dX, | ||
magma_int_t | lddx, | ||
magmaDouble_ptr | dworkd, | ||
magmaFloat_ptr | dworks, | ||
magma_int_t * | iter, | ||
magma_int_t * | info | ||
) |
DSPOSV computes the solution to a real system of linear equations A * X = B, where A is an N-by-N symmetric positive definite matrix and X and B are N-by-NRHS matrices.
DSPOSV first attempts to factorize the matrix in real SINGLE PRECISION and use this factorization within an iterative refinement procedure to produce a solution with real DOUBLE PRECISION norm-wise backward error quality (see below). If the approach fails the method switches to a real DOUBLE PRECISION factorization and solve.
The iterative refinement is not going to be a winning strategy if the ratio real SINGLE PRECISION performance over real DOUBLE PRECISION performance is too small. A reasonable strategy should take the number of right-hand sides and the size of the matrix into account. This might be done with a call to ILAENV in the future. Up to now, we always try iterative refinement.
The iterative refinement process is stopped if ITER > ITERMAX or for all the RHS we have: RNRM < SQRT(N)*XNRM*ANRM*EPS*BWDMAX where o ITER is the number of the current iteration in the iterative refinement process o RNRM is the infinity-norm of the residual o XNRM is the infinity-norm of the solution o ANRM is the infinity-operator-norm of the matrix A o EPS is the machine epsilon returned by DLAMCH('Epsilon') The value ITERMAX and BWDMAX are fixed to 30 and 1.0D+00 respectively.
[in] | uplo | magma_uplo_t
|
[in] | n | INTEGER The number of linear equations, i.e., 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,out] | dA | DOUBLE PRECISION array on the GPU, dimension (LDDA,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 iterative refinement has been successfully used (INFO.EQ.0 and ITER.GE.0, see description below), then A is unchanged, if double factorization has been used (INFO.EQ.0 and ITER.LT.0, see description below), then the array dA contains the factor U or L from the Cholesky factorization A = U**T*U or A = L*L**T. |
[in] | ldda | INTEGER The leading dimension of the array dA. LDDA >= max(1,N). |
[in] | dB | DOUBLE PRECISION array on the GPU, dimension (LDDB,NRHS) The N-by-NRHS right hand side matrix B. |
[in] | lddb | INTEGER The leading dimension of the array dB. LDDB >= max(1,N). |
[out] | dX | DOUBLE PRECISION array on the GPU, dimension (LDDX,NRHS) If INFO = 0, the N-by-NRHS solution matrix X. |
[in] | lddx | INTEGER The leading dimension of the array dX. LDDX >= max(1,N). |
dworkd | (workspace) DOUBLE PRECISION array on the GPU, dimension (N*NRHS) This array is used to hold the residual vectors. | |
dworks | (workspace) SINGLE PRECISION array on the GPU, dimension (N*(N+NRHS)) This array is used to store the real single precision matrix and the right-hand sides or solutions in single precision. | |
[out] | iter | INTEGER
|
[out] | info | INTEGER
|