MAGMA  2.7.1
Matrix Algebra for GPU and Multicore Architectures
 All Classes Files Functions Friends Groups Pages
getri: LU inverse

Functions

magma_int_t magma_cgetri_outofplace_batched (magma_int_t n, magmaFloatComplex **dA_array, magma_int_t ldda, magma_int_t **dipiv_array, magmaFloatComplex **dinvA_array, magma_int_t lddia, magma_int_t *info_array, magma_int_t batchCount, magma_queue_t queue)
 CGETRI computes the inverse of a matrix using the LU factorization computed by CGETRF. More...
 
magma_int_t magma_dgetri_outofplace_batched (magma_int_t n, double **dA_array, magma_int_t ldda, magma_int_t **dipiv_array, double **dinvA_array, magma_int_t lddia, magma_int_t *info_array, magma_int_t batchCount, magma_queue_t queue)
 DGETRI computes the inverse of a matrix using the LU factorization computed by DGETRF. More...
 
magma_int_t magma_sgetri_outofplace_batched (magma_int_t n, float **dA_array, magma_int_t ldda, magma_int_t **dipiv_array, float **dinvA_array, magma_int_t lddia, magma_int_t *info_array, magma_int_t batchCount, magma_queue_t queue)
 SGETRI computes the inverse of a matrix using the LU factorization computed by SGETRF. More...
 
magma_int_t magma_zgetri_outofplace_batched (magma_int_t n, magmaDoubleComplex **dA_array, magma_int_t ldda, magma_int_t **dipiv_array, magmaDoubleComplex **dinvA_array, magma_int_t lddia, magma_int_t *info_array, magma_int_t batchCount, magma_queue_t queue)
 ZGETRI computes the inverse of a matrix using the LU factorization computed by ZGETRF. More...
 

Detailed Description

Function Documentation

magma_int_t magma_cgetri_outofplace_batched ( magma_int_t  n,
magmaFloatComplex **  dA_array,
magma_int_t  ldda,
magma_int_t **  dipiv_array,
magmaFloatComplex **  dinvA_array,
magma_int_t  lddia,
magma_int_t *  info_array,
magma_int_t  batchCount,
magma_queue_t  queue 
)

CGETRI computes the inverse of a matrix using the LU factorization computed by CGETRF.

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 CGESV, or CGETRF and CGETRS, 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.

Parameters
[in]nINTEGER The order of the matrix A. N >= 0.
[in,out]dA_arrayArray of pointers, dimension (batchCount). Each is a COMPLEX array on the GPU, dimension (LDDA,N) On entry, the factors L and U from the factorization A = P*L*U as computed by CGETRF_GPU. On exit, if INFO = 0, the inverse of the original matrix A.
[in]lddaINTEGER The leading dimension of the array A. LDDA >= max(1,N).
[in]dipiv_arrayArray of pointers, dimension (batchCount), for corresponding matrices. Each is an INTEGER array, dimension (N) The pivot indices; for 1 <= i <= min(M,N), row i of the matrix was interchanged with row IPIV(i).
[out]dinvA_arrayArray of pointers, dimension (batchCount). Each is a COMPLEX array on the GPU, dimension (LDDIA,N) It contains the inverse of the matrix
[in]lddiaINTEGER The leading dimension of the array invA_array. LDDIA >= max(1,N).
[out]info_arrayArray of INTEGERs, dimension (batchCount), for corresponding matrices.
  • = 0: successful exit
  • < 0: if INFO = -i, the i-th argument had an illegal value or another error occured, such as memory allocation failed.
  • > 0: if INFO = i, U(i,i) is exactly zero. The factorization has been completed, but the factor U is exactly singular, and division by zero will occur if it is used to solve a system of equations.
[in]batchCountINTEGER The number of matrices to operate on.
[in]queuemagma_queue_t Queue to execute in.
magma_int_t magma_dgetri_outofplace_batched ( magma_int_t  n,
double **  dA_array,
magma_int_t  ldda,
magma_int_t **  dipiv_array,
double **  dinvA_array,
magma_int_t  lddia,
magma_int_t *  info_array,
magma_int_t  batchCount,
magma_queue_t  queue 
)

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.

Parameters
[in]nINTEGER The order of the matrix A. N >= 0.
[in,out]dA_arrayArray of pointers, dimension (batchCount). Each is a 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]lddaINTEGER The leading dimension of the array A. LDDA >= max(1,N).
[in]dipiv_arrayArray of pointers, dimension (batchCount), for corresponding matrices. Each is an INTEGER array, dimension (N) The pivot indices; for 1 <= i <= min(M,N), row i of the matrix was interchanged with row IPIV(i).
[out]dinvA_arrayArray of pointers, dimension (batchCount). Each is a DOUBLE PRECISION array on the GPU, dimension (LDDIA,N) It contains the inverse of the matrix
[in]lddiaINTEGER The leading dimension of the array invA_array. LDDIA >= max(1,N).
[out]info_arrayArray of INTEGERs, dimension (batchCount), for corresponding matrices.
  • = 0: successful exit
  • < 0: if INFO = -i, the i-th argument had an illegal value or another error occured, such as memory allocation failed.
  • > 0: if INFO = i, U(i,i) is exactly zero. The factorization has been completed, but the factor U is exactly singular, and division by zero will occur if it is used to solve a system of equations.
[in]batchCountINTEGER The number of matrices to operate on.
[in]queuemagma_queue_t Queue to execute in.
magma_int_t magma_sgetri_outofplace_batched ( magma_int_t  n,
float **  dA_array,
magma_int_t  ldda,
magma_int_t **  dipiv_array,
float **  dinvA_array,
magma_int_t  lddia,
magma_int_t *  info_array,
magma_int_t  batchCount,
magma_queue_t  queue 
)

SGETRI computes the inverse of a matrix using the LU factorization computed by SGETRF.

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 SGESV, or SGETRF and SGETRS, 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.

Parameters
[in]nINTEGER The order of the matrix A. N >= 0.
[in,out]dA_arrayArray of pointers, dimension (batchCount). Each is a REAL array on the GPU, dimension (LDDA,N) On entry, the factors L and U from the factorization A = P*L*U as computed by SGETRF_GPU. On exit, if INFO = 0, the inverse of the original matrix A.
[in]lddaINTEGER The leading dimension of the array A. LDDA >= max(1,N).
[in]dipiv_arrayArray of pointers, dimension (batchCount), for corresponding matrices. Each is an INTEGER array, dimension (N) The pivot indices; for 1 <= i <= min(M,N), row i of the matrix was interchanged with row IPIV(i).
[out]dinvA_arrayArray of pointers, dimension (batchCount). Each is a REAL array on the GPU, dimension (LDDIA,N) It contains the inverse of the matrix
[in]lddiaINTEGER The leading dimension of the array invA_array. LDDIA >= max(1,N).
[out]info_arrayArray of INTEGERs, dimension (batchCount), for corresponding matrices.
  • = 0: successful exit
  • < 0: if INFO = -i, the i-th argument had an illegal value or another error occured, such as memory allocation failed.
  • > 0: if INFO = i, U(i,i) is exactly zero. The factorization has been completed, but the factor U is exactly singular, and division by zero will occur if it is used to solve a system of equations.
[in]batchCountINTEGER The number of matrices to operate on.
[in]queuemagma_queue_t Queue to execute in.
magma_int_t magma_zgetri_outofplace_batched ( magma_int_t  n,
magmaDoubleComplex **  dA_array,
magma_int_t  ldda,
magma_int_t **  dipiv_array,
magmaDoubleComplex **  dinvA_array,
magma_int_t  lddia,
magma_int_t *  info_array,
magma_int_t  batchCount,
magma_queue_t  queue 
)

ZGETRI computes the inverse of a matrix using the LU factorization computed by ZGETRF.

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 ZGESV, or ZGETRF and ZGETRS, 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.

Parameters
[in]nINTEGER The order of the matrix A. N >= 0.
[in,out]dA_arrayArray of pointers, dimension (batchCount). Each is a COMPLEX_16 array on the GPU, dimension (LDDA,N) On entry, the factors L and U from the factorization A = P*L*U as computed by ZGETRF_GPU. On exit, if INFO = 0, the inverse of the original matrix A.
[in]lddaINTEGER The leading dimension of the array A. LDDA >= max(1,N).
[in]dipiv_arrayArray of pointers, dimension (batchCount), for corresponding matrices. Each is an INTEGER array, dimension (N) The pivot indices; for 1 <= i <= min(M,N), row i of the matrix was interchanged with row IPIV(i).
[out]dinvA_arrayArray of pointers, dimension (batchCount). Each is a COMPLEX_16 array on the GPU, dimension (LDDIA,N) It contains the inverse of the matrix
[in]lddiaINTEGER The leading dimension of the array invA_array. LDDIA >= max(1,N).
[out]info_arrayArray of INTEGERs, dimension (batchCount), for corresponding matrices.
  • = 0: successful exit
  • < 0: if INFO = -i, the i-th argument had an illegal value or another error occured, such as memory allocation failed.
  • > 0: if INFO = i, U(i,i) is exactly zero. The factorization has been completed, but the factor U is exactly singular, and division by zero will occur if it is used to solve a system of equations.
[in]batchCountINTEGER The number of matrices to operate on.
[in]queuemagma_queue_t Queue to execute in.