PLASMA
Parallel Linear Algebra Software for Multicore Architectures
sy/hetrs: Symmetric indefinite forward and back solves

Functions

int plasma_chetrs (plasma_enum_t uplo, int n, int nrhs, plasma_complex32_t *pA, int lda, int *ipiv, plasma_complex32_t *pT, int ldt, int *ipiv2, plasma_complex32_t *pB, int ldb)
 
void plasma_omp_chetrs (plasma_enum_t uplo, plasma_desc_t A, int *ipiv, plasma_desc_t T, int *ipiv2, plasma_desc_t B, plasma_sequence_t *sequence, plasma_request_t *request)
 
int plasma_dsytrs (plasma_enum_t uplo, int n, int nrhs, double *pA, int lda, int *ipiv, double *pT, int ldt, int *ipiv2, double *pB, int ldb)
 
void plasma_omp_dsytrs (plasma_enum_t uplo, plasma_desc_t A, int *ipiv, plasma_desc_t T, int *ipiv2, plasma_desc_t B, plasma_sequence_t *sequence, plasma_request_t *request)
 
int plasma_ssytrs (plasma_enum_t uplo, int n, int nrhs, float *pA, int lda, int *ipiv, float *pT, int ldt, int *ipiv2, float *pB, int ldb)
 
void plasma_omp_ssytrs (plasma_enum_t uplo, plasma_desc_t A, int *ipiv, plasma_desc_t T, int *ipiv2, plasma_desc_t B, plasma_sequence_t *sequence, plasma_request_t *request)
 
int plasma_zhetrs (plasma_enum_t uplo, int n, int nrhs, plasma_complex64_t *pA, int lda, int *ipiv, plasma_complex64_t *pT, int ldt, int *ipiv2, plasma_complex64_t *pB, int ldb)
 
void plasma_omp_zhetrs (plasma_enum_t uplo, plasma_desc_t A, int *ipiv, plasma_desc_t T, int *ipiv2, plasma_desc_t B, plasma_sequence_t *sequence, plasma_request_t *request)
 

Detailed Description

Function Documentation

int plasma_chetrs ( plasma_enum_t  uplo,
int  n,
int  nrhs,
plasma_complex32_t *  pA,
int  lda,
int *  ipiv,
plasma_complex32_t *  pT,
int  ldt,
int *  ipiv2,
plasma_complex32_t *  pB,
int  ldb 
)

Solves a system of linear equations A * X = B with LTLt factorization computed by plasma_chetrf.

Parameters
[in]uplo
  • PlasmaUpper: Upper triangle of A is stored;
  • PlasmaLower: Lower triangle of A is stored. TODO: only support Lower for now
[in]nThe order of the matrix A. n >= 0.
[in]nrhsThe number of right hand sides, i.e., the number of columns of the matrix B. nrhs >= 0.
[in,out]ADetails of the LTL factorization of the Hermitian matrix A, as computed by plasma_chetrf.
[in]ldaThe leading dimension of the array A.
[in,out]TDetails of the LU factorization of the band matrix A, as computed by plasma_cgbtrf.
[in]ldtThe leading dimension of the array T.
[in]ipivThe pivot indices used for chetrf; for 1 <= i <= min(m,n), row i of the matrix was interchanged with row ipiv(i).
[in]ipiv2The pivot indices used for cgbtrf; for 1 <= i <= min(m,n), row i of the matrix was interchanged with row ipiv(i).
[in,out]BOn entry, the n-by-nrhs right hand side matrix B. On exit, if return value = 0, the n-by-nrhs solution matrix X.
[in]ldbThe leading dimension of the array B. ldb >= max(1,n).
Return values
PlasmaSuccesssuccessful exit
<0 if -i, the i-th argument had an illegal value
See also
plasma_omp_chetrs
plasma_chetrs
plasma_dsytrs
plasma_ssytrs
plasma_chetrf
void plasma_omp_chetrs ( plasma_enum_t  uplo,
plasma_desc_t  A,
int *  ipiv,
plasma_desc_t  T,
int *  ipiv2,
plasma_desc_t  B,
plasma_sequence_t *  sequence,
plasma_request_t *  request 
)

Solves a system of linear equations using previously computed factorization. Non-blocking tile version of plasma_chetrs(). May return before the computation is finished. Operates on matrices stored by tiles. All matrices are passed through descriptors. All dimensions are taken from the descriptors. Allows for pipelining of operations at runtime.

Parameters
[in]uplo
  • PlasmaUpper: Upper triangle of A is stored;
  • PlasmaLower: Lower triangle of A is stored.
[in]AThe triangular factor U or L from the Cholesky factorization A = U^H*U or A = L*L^H, computed by plasma_cpotrf.
[in,out]BOn entry, the n-by-nrhs right hand side matrix B. On exit, if return value = 0, the n-by-nrhs solution matrix X.
[in]sequenceIdentifies the sequence of function calls that this call belongs to (for completion checks and exception handling purposes). Check the sequence->status for errors.
[out]requestIdentifies this function call (for exception handling purposes).
Return values
voidErrors are returned by setting sequence->status and request->status to error values. The sequence->status and request->status should never be set to PlasmaSuccess (the initial values) since another async call may be setting a failure value at the same time.
See also
plasma_chetrs
plasma_omp_chetrs
plasma_omp_chetrs
plasma_omp_dsytrs
plasma_omp_ssytrs
plasma_omp_chetrf
int plasma_dsytrs ( plasma_enum_t  uplo,
int  n,
int  nrhs,
double *  pA,
int  lda,
int *  ipiv,
double *  pT,
int  ldt,
int *  ipiv2,
double *  pB,
int  ldb 
)

Solves a system of linear equations A * X = B with LTLt factorization computed by plasma_dsytrf.

Parameters
[in]uplo
  • PlasmaUpper: Upper triangle of A is stored;
  • PlasmaLower: Lower triangle of A is stored. TODO: only support Lower for now
[in]nThe order of the matrix A. n >= 0.
[in]nrhsThe number of right hand sides, i.e., the number of columns of the matrix B. nrhs >= 0.
[in,out]ADetails of the LTL factorization of the symmetric matrix A, as computed by plasma_dsytrf.
[in]ldaThe leading dimension of the array A.
[in,out]TDetails of the LU factorization of the band matrix A, as computed by plasma_dgbtrf.
[in]ldtThe leading dimension of the array T.
[in]ipivThe pivot indices used for dsytrf; for 1 <= i <= min(m,n), row i of the matrix was interchanged with row ipiv(i).
[in]ipiv2The pivot indices used for dgbtrf; for 1 <= i <= min(m,n), row i of the matrix was interchanged with row ipiv(i).
[in,out]BOn entry, the n-by-nrhs right hand side matrix B. On exit, if return value = 0, the n-by-nrhs solution matrix X.
[in]ldbThe leading dimension of the array B. ldb >= max(1,n).
Return values
PlasmaSuccesssuccessful exit
<0 if -i, the i-th argument had an illegal value
See also
plasma_omp_dsytrs
plasma_chetrs
plasma_dsytrs
plasma_ssytrs
plasma_dsytrf
void plasma_omp_dsytrs ( plasma_enum_t  uplo,
plasma_desc_t  A,
int *  ipiv,
plasma_desc_t  T,
int *  ipiv2,
plasma_desc_t  B,
plasma_sequence_t *  sequence,
plasma_request_t *  request 
)

Solves a system of linear equations using previously computed factorization. Non-blocking tile version of plasma_dsytrs(). May return before the computation is finished. Operates on matrices stored by tiles. All matrices are passed through descriptors. All dimensions are taken from the descriptors. Allows for pipelining of operations at runtime.

Parameters
[in]uplo
  • PlasmaUpper: Upper triangle of A is stored;
  • PlasmaLower: Lower triangle of A is stored.
[in]AThe triangular factor U or L from the Cholesky factorization A = U^T*U or A = L*L^T, computed by plasma_dpotrf.
[in,out]BOn entry, the n-by-nrhs right hand side matrix B. On exit, if return value = 0, the n-by-nrhs solution matrix X.
[in]sequenceIdentifies the sequence of function calls that this call belongs to (for completion checks and exception handling purposes). Check the sequence->status for errors.
[out]requestIdentifies this function call (for exception handling purposes).
Return values
voidErrors are returned by setting sequence->status and request->status to error values. The sequence->status and request->status should never be set to PlasmaSuccess (the initial values) since another async call may be setting a failure value at the same time.
See also
plasma_dsytrs
plasma_omp_dsytrs
plasma_omp_chetrs
plasma_omp_dsytrs
plasma_omp_ssytrs
plasma_omp_dsytrf
int plasma_ssytrs ( plasma_enum_t  uplo,
int  n,
int  nrhs,
float *  pA,
int  lda,
int *  ipiv,
float *  pT,
int  ldt,
int *  ipiv2,
float *  pB,
int  ldb 
)

Solves a system of linear equations A * X = B with LTLt factorization computed by plasma_ssytrf.

Parameters
[in]uplo
  • PlasmaUpper: Upper triangle of A is stored;
  • PlasmaLower: Lower triangle of A is stored. TODO: only support Lower for now
[in]nThe order of the matrix A. n >= 0.
[in]nrhsThe number of right hand sides, i.e., the number of columns of the matrix B. nrhs >= 0.
[in,out]ADetails of the LTL factorization of the symmetric matrix A, as computed by plasma_ssytrf.
[in]ldaThe leading dimension of the array A.
[in,out]TDetails of the LU factorization of the band matrix A, as computed by plasma_sgbtrf.
[in]ldtThe leading dimension of the array T.
[in]ipivThe pivot indices used for ssytrf; for 1 <= i <= min(m,n), row i of the matrix was interchanged with row ipiv(i).
[in]ipiv2The pivot indices used for sgbtrf; for 1 <= i <= min(m,n), row i of the matrix was interchanged with row ipiv(i).
[in,out]BOn entry, the n-by-nrhs right hand side matrix B. On exit, if return value = 0, the n-by-nrhs solution matrix X.
[in]ldbThe leading dimension of the array B. ldb >= max(1,n).
Return values
PlasmaSuccesssuccessful exit
<0 if -i, the i-th argument had an illegal value
See also
plasma_omp_ssytrs
plasma_chetrs
plasma_dsytrs
plasma_ssytrs
plasma_ssytrf
void plasma_omp_ssytrs ( plasma_enum_t  uplo,
plasma_desc_t  A,
int *  ipiv,
plasma_desc_t  T,
int *  ipiv2,
plasma_desc_t  B,
plasma_sequence_t *  sequence,
plasma_request_t *  request 
)

Solves a system of linear equations using previously computed factorization. Non-blocking tile version of plasma_ssytrs(). May return before the computation is finished. Operates on matrices stored by tiles. All matrices are passed through descriptors. All dimensions are taken from the descriptors. Allows for pipelining of operations at runtime.

Parameters
[in]uplo
  • PlasmaUpper: Upper triangle of A is stored;
  • PlasmaLower: Lower triangle of A is stored.
[in]AThe triangular factor U or L from the Cholesky factorization A = U^T*U or A = L*L^T, computed by plasma_spotrf.
[in,out]BOn entry, the n-by-nrhs right hand side matrix B. On exit, if return value = 0, the n-by-nrhs solution matrix X.
[in]sequenceIdentifies the sequence of function calls that this call belongs to (for completion checks and exception handling purposes). Check the sequence->status for errors.
[out]requestIdentifies this function call (for exception handling purposes).
Return values
voidErrors are returned by setting sequence->status and request->status to error values. The sequence->status and request->status should never be set to PlasmaSuccess (the initial values) since another async call may be setting a failure value at the same time.
See also
plasma_ssytrs
plasma_omp_ssytrs
plasma_omp_chetrs
plasma_omp_dsytrs
plasma_omp_ssytrs
plasma_omp_ssytrf
int plasma_zhetrs ( plasma_enum_t  uplo,
int  n,
int  nrhs,
plasma_complex64_t *  pA,
int  lda,
int *  ipiv,
plasma_complex64_t *  pT,
int  ldt,
int *  ipiv2,
plasma_complex64_t *  pB,
int  ldb 
)

Solves a system of linear equations A * X = B with LTLt factorization computed by plasma_zhetrf.

Parameters
[in]uplo
  • PlasmaUpper: Upper triangle of A is stored;
  • PlasmaLower: Lower triangle of A is stored. TODO: only support Lower for now
[in]nThe order of the matrix A. n >= 0.
[in]nrhsThe number of right hand sides, i.e., the number of columns of the matrix B. nrhs >= 0.
[in,out]ADetails of the LTL factorization of the Hermitian matrix A, as computed by plasma_zhetrf.
[in]ldaThe leading dimension of the array A.
[in,out]TDetails of the LU factorization of the band matrix A, as computed by plasma_zgbtrf.
[in]ldtThe leading dimension of the array T.
[in]ipivThe pivot indices used for zhetrf; for 1 <= i <= min(m,n), row i of the matrix was interchanged with row ipiv(i).
[in]ipiv2The pivot indices used for zgbtrf; for 1 <= i <= min(m,n), row i of the matrix was interchanged with row ipiv(i).
[in,out]BOn entry, the n-by-nrhs right hand side matrix B. On exit, if return value = 0, the n-by-nrhs solution matrix X.
[in]ldbThe leading dimension of the array B. ldb >= max(1,n).
Return values
PlasmaSuccesssuccessful exit
<0 if -i, the i-th argument had an illegal value
See also
plasma_omp_zhetrs
plasma_chetrs
plasma_dsytrs
plasma_ssytrs
plasma_zhetrf
void plasma_omp_zhetrs ( plasma_enum_t  uplo,
plasma_desc_t  A,
int *  ipiv,
plasma_desc_t  T,
int *  ipiv2,
plasma_desc_t  B,
plasma_sequence_t *  sequence,
plasma_request_t *  request 
)

Solves a system of linear equations using previously computed factorization. Non-blocking tile version of plasma_zhetrs(). May return before the computation is finished. Operates on matrices stored by tiles. All matrices are passed through descriptors. All dimensions are taken from the descriptors. Allows for pipelining of operations at runtime.

Parameters
[in]uplo
  • PlasmaUpper: Upper triangle of A is stored;
  • PlasmaLower: Lower triangle of A is stored.
[in]AThe triangular factor U or L from the Cholesky factorization A = U^H*U or A = L*L^H, computed by plasma_zpotrf.
[in,out]BOn entry, the n-by-nrhs right hand side matrix B. On exit, if return value = 0, the n-by-nrhs solution matrix X.
[in]sequenceIdentifies the sequence of function calls that this call belongs to (for completion checks and exception handling purposes). Check the sequence->status for errors.
[out]requestIdentifies this function call (for exception handling purposes).
Return values
voidErrors are returned by setting sequence->status and request->status to error values. The sequence->status and request->status should never be set to PlasmaSuccess (the initial values) since another async call may be setting a failure value at the same time.
See also
plasma_zhetrs
plasma_omp_zhetrs
plasma_omp_chetrs
plasma_omp_dsytrs
plasma_omp_ssytrs
plasma_omp_zhetrf