dsposv#

Functions

void dsposv(
    const char*           uplo,
    const INT             n,
    const INT             nrhs,
          f64*   restrict A,
    const INT             lda,
    const f64*   restrict B,
    const INT             ldb,
          f64*   restrict X,
    const INT             ldx,
          f64*   restrict work,
          float* restrict swork,
          INT*            iter,
          INT*            info
);
void dsposv(const char *uplo, const INT n, const INT nrhs, f64 *restrict A, const INT lda, const f64 *restrict B, const INT ldb, f64 *restrict X, const INT ldx, f64 *restrict work, float *restrict swork, INT *iter, INT *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 SINGLE PRECISION and use this factorization within an iterative refinement procedure to produce a solution with DOUBLE PRECISION normwise backward error quality (see below). If the approach fails the method switches to a DOUBLE PRECISION factorization and solve.

The iterative refinement is not going to be a winning strategy if the ratio SINGLE PRECISION performance over 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.

Parameters

in
uplo

  • 'U': Upper triangle of A is stored

  • 'L': Lower triangle of A is stored

in
n

The number of linear equations, i.e., the order of the matrix A. n>=0.

in
nrhs

The number of right hand sides, i.e., the number of columns of the matrix B. nrhs>=0.

inout
A

Array of dimension (lda, n). On entry, the symmetric matrix A. If uplo='U', the leading n-by-n upper triangular part of A contains the upper triangular part of the matrix A. If uplo='L', the leading n-by-n lower triangular part of A contains the lower triangular part of the matrix A. On exit, if iterative refinement has been successfully used (info=0 and iter>=0) then A is unchanged. If double precision factorization has been used (info=0 and iter<0) then the array A contains the factor U or L from the Cholesky factorization A = U**T*U or A = L*L**T.

in
lda

The leading dimension of the array A. lda>=max(1,n).

in
B

Array of dimension (ldb, nrhs). The n-by-nrhs right hand side matrix B.

in
ldb

The leading dimension of the array B. ldb>=max(1,n).

out
X

Array of dimension (ldx, nrhs). If info=0, the n-by-nrhs solution matrix X.

in
ldx

The leading dimension of the array X. ldx>=max(1,n).

out
work

Double precision workspace for residual vectors. Array of dimension (n, nrhs).

out
swork

Single precision workspace for the matrix and the right-hand sides or solutions in single precision. Array of dimension n*(n+nrhs).

out
iter

Iteration count:

  • iter<0: iterative refinement has failed, double precision factorization has been performed

    • -1 : the routine fell back to full precision for implementation- or machine-specific reasons

    • -2 : narrowing the precision induced an overflow, the routine fell back to full precision

    • -3 : failure of SPOTRF

    • -31: stop the iterative refinement after the 30th iterations

  • iter>0: iterative refinement has been successfully used. Returns the number of iterations

out
info

  • info=0: successful exit

  • info<0: if info=-i, the i-th argument had an illegal value

  • info>0: if info=i, the leading principal minor of order i of (DOUBLE PRECISION) A is not positive, so the factorization could not be completed, and the solution has not been computed.

On this page
Show Source