posvx#

Functions

void sposvx(
    const char*          fact,
    const char*          uplo,
    const INT            n,
    const INT            nrhs,
          f32*  restrict A,
    const INT            lda,
          f32*  restrict AF,
    const INT            ldaf,
          char*          equed,
          f32*  restrict S,
          f32*  restrict B,
    const INT            ldb,
          f32*  restrict X,
    const INT            ldx,
          f32*           rcond,
          f32*  restrict ferr,
          f32*  restrict berr,
          f32*  restrict work,
          INT*  restrict iwork,
          INT*           info
);
void sposvx(const char *fact, const char *uplo, const INT n, const INT nrhs, f32 *restrict A, const INT lda, f32 *restrict AF, const INT ldaf, char *equed, f32 *restrict S, f32 *restrict B, const INT ldb, f32 *restrict X, const INT ldx, f32 *rcond, f32 *restrict ferr, f32 *restrict berr, f32 *restrict work, INT *restrict iwork, INT *info)#

SPOSVX uses the Cholesky factorization A = U**T*U or A = L*L**T to compute 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.

Error bounds on the solution and a condition estimate are also provided.

The following steps are performed:

  1. If fact='E', real scaling factors are computed to equilibrate the system:

    diag(S)*A*diag(S) * inv(diag(S))*X = diag(S)*B
    

    Whether or not the system will be equilibrated depends on the scaling of the matrix A, but if equilibration is used, A is overwritten by diag(S)*A*diag(S) and B by diag(S)*B.

  2. If fact='N' or 'E', the Cholesky decomposition is used to factor the matrix A (after equilibration if fact='E') as:

    A = U**T * U,  if uplo = 'U', or
    A = L * L**T,  if uplo = 'L',
    

    where U is an upper triangular matrix and L is a lower triangular matrix.

  3. If the leading principal minor of order i is not positive, then the routine returns with info=i. Otherwise, the factored form of A is used to estimate the condition number of the matrix A. If the reciprocal of the condition number is less than machine precision, info=n+1 is returned as a warning, but the routine still goes on to solve for X and compute error bounds as described below.

  4. The system of equations is solved for X using the factored form of A.

  5. Iterative refinement is applied to improve the computed solution matrix and calculate error bounds and backward error estimates for it.

  6. If equilibration was used, the matrix X is premultiplied by diag(S) so that it solves the original system before equilibration.

equed is an input argument if fact='F'; otherwise it is an output argument.

Parameters

in
fact

  • 'F': AF contains the factored form of A. If equed='Y', A has been equilibrated with scaling factors given by S; A and AF will not be modified.

  • 'N': The matrix A will be copied to AF and factored.

  • 'E': The matrix A will be equilibrated if necessary, then copied to AF and factored.

in
uplo

  • 'U': Upper triangle of A is stored

  • 'L': Lower triangle of A is stored

in
n

The number of linear equations. n>=0.

in
nrhs

The number of right hand sides. nrhs>=0.

inout
A

Array of dimension (lda, n). On entry, the symmetric matrix A, except if fact='F' and equed='Y', then A must contain the equilibrated matrix diag(S)*A*diag(S). A is not modified if fact='F' or 'N', or if fact='E' and equed='N' on exit. On exit, if fact='E' and equed='Y', A is overwritten by diag(S)*A*diag(S).

in
lda

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

inout
AF

Array of dimension (ldaf, n). If fact='F', an input argument containing the triangular factor U or L from the Cholesky factorization, in the same storage format as A. If fact='N' or 'E', an output argument returning the triangular factor of the (possibly equilibrated) matrix A.

in
ldaf

The leading dimension of the array AF. ldaf>=max(1,n).

inout
equed

  • 'N': No equilibration (always true if fact='N')

  • 'Y': Equilibration was done, i.e., A has been replaced by diag(S) * A * diag(S)

inout
S

Array of dimension (n). The scale factors for A; not accessed if equed='N'. S is an input argument if fact='F'; otherwise it is an output argument. If fact='F' and equed='Y', each element of S must be positive.

inout
B

Array of dimension (ldb, nrhs). On entry, the n-by-nrhs right hand side matrix B. On exit, if equed='N', B is not modified; if equed='Y', B is overwritten by diag(S) * B.

in
ldb

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

out
X

Array of dimension (ldx, nrhs). If info=0 or info=n+1, the n-by-nrhs solution matrix X to the original system of equations. Note that if equed='Y', A and B are modified on exit, and the solution to the equilibrated system is inv(diag(S))*X.

in
ldx

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

out
rcond

The estimate of the reciprocal condition number of the matrix A after equilibration (if done).

out
ferr

Array of dimension (nrhs). The estimated forward error bound for each solution vector.

out
berr

Array of dimension (nrhs). The componentwise relative backward error of each solution vector.

out
work

Array of dimension 3*n.

out
iwork

Array of dimension n.

out
info

  • info=0: successful exit

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

  • info>0: if info=k, and k<=n, the leading principal minor of order k of A is not positive, so the factorization could not be completed, and the solution has not been computed. rcond=0 is returned.

  • info=n+1: U is nonsingular, but rcond is less than machine precision, meaning that the matrix is singular to working precision. Nevertheless, the solution and error bounds are computed because there are a number of situations where the computed solution can be more accurate than the value of rcond would suggest.

Functions

void dposvx(
    const char*          fact,
    const char*          uplo,
    const INT            n,
    const INT            nrhs,
          f64*  restrict A,
    const INT            lda,
          f64*  restrict AF,
    const INT            ldaf,
          char*          equed,
          f64*  restrict S,
          f64*  restrict B,
    const INT            ldb,
          f64*  restrict X,
    const INT            ldx,
          f64*           rcond,
          f64*  restrict ferr,
          f64*  restrict berr,
          f64*  restrict work,
          INT*  restrict iwork,
          INT*           info
);
void dposvx(const char *fact, const char *uplo, const INT n, const INT nrhs, f64 *restrict A, const INT lda, f64 *restrict AF, const INT ldaf, char *equed, f64 *restrict S, f64 *restrict B, const INT ldb, f64 *restrict X, const INT ldx, f64 *rcond, f64 *restrict ferr, f64 *restrict berr, f64 *restrict work, INT *restrict iwork, INT *info)#

DPOSVX uses the Cholesky factorization A = U**T*U or A = L*L**T to compute 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.

Error bounds on the solution and a condition estimate are also provided.

The following steps are performed:

  1. If fact='E', real scaling factors are computed to equilibrate the system:

    diag(S)*A*diag(S) * inv(diag(S))*X = diag(S)*B
    

    Whether or not the system will be equilibrated depends on the scaling of the matrix A, but if equilibration is used, A is overwritten by diag(S)*A*diag(S) and B by diag(S)*B.

  2. If fact='N' or 'E', the Cholesky decomposition is used to factor the matrix A (after equilibration if fact='E') as:

    A = U**T * U,  if uplo = 'U', or
    A = L * L**T,  if uplo = 'L',
    

    where U is an upper triangular matrix and L is a lower triangular matrix.

  3. If the leading principal minor of order i is not positive, then the routine returns with info=i. Otherwise, the factored form of A is used to estimate the condition number of the matrix A. If the reciprocal of the condition number is less than machine precision, info=n+1 is returned as a warning, but the routine still goes on to solve for X and compute error bounds as described below.

  4. The system of equations is solved for X using the factored form of A.

  5. Iterative refinement is applied to improve the computed solution matrix and calculate error bounds and backward error estimates for it.

  6. If equilibration was used, the matrix X is premultiplied by diag(S) so that it solves the original system before equilibration.

equed is an input argument if fact='F'; otherwise it is an output argument.

Parameters

in
fact

  • 'F': AF contains the factored form of A. If equed='Y', A has been equilibrated with scaling factors given by S; A and AF will not be modified.

  • 'N': The matrix A will be copied to AF and factored.

  • 'E': The matrix A will be equilibrated if necessary, then copied to AF and factored.

in
uplo

  • 'U': Upper triangle of A is stored

  • 'L': Lower triangle of A is stored

in
n

The number of linear equations. n>=0.

in
nrhs

The number of right hand sides. nrhs>=0.

inout
A

Array of dimension (lda, n). On entry, the symmetric matrix A, except if fact='F' and equed='Y', then A must contain the equilibrated matrix diag(S)*A*diag(S). A is not modified if fact='F' or 'N', or if fact='E' and equed='N' on exit. On exit, if fact='E' and equed='Y', A is overwritten by diag(S)*A*diag(S).

in
lda

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

inout
AF

Array of dimension (ldaf, n). If fact='F', an input argument containing the triangular factor U or L from the Cholesky factorization, in the same storage format as A. If fact='N' or 'E', an output argument returning the triangular factor of the (possibly equilibrated) matrix A.

in
ldaf

The leading dimension of the array AF. ldaf>=max(1,n).

inout
equed

  • 'N': No equilibration (always true if fact='N')

  • 'Y': Equilibration was done, i.e., A has been replaced by diag(S) * A * diag(S)

inout
S

Array of dimension (n). The scale factors for A; not accessed if equed='N'. S is an input argument if fact='F'; otherwise it is an output argument. If fact='F' and equed='Y', each element of S must be positive.

inout
B

Array of dimension (ldb, nrhs). On entry, the n-by-nrhs right hand side matrix B. On exit, if equed='N', B is not modified; if equed='Y', B is overwritten by diag(S) * B.

in
ldb

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

out
X

Array of dimension (ldx, nrhs). If info=0 or info=n+1, the n-by-nrhs solution matrix X to the original system of equations. Note that if equed='Y', A and B are modified on exit, and the solution to the equilibrated system is inv(diag(S))*X.

in
ldx

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

out
rcond

The estimate of the reciprocal condition number of the matrix A after equilibration (if done).

out
ferr

Array of dimension (nrhs). The estimated forward error bound for each solution vector.

out
berr

Array of dimension (nrhs). The componentwise relative backward error of each solution vector.

out
work

Array of dimension 3*n.

out
iwork

Array of dimension n.

out
info

  • info=0: successful exit

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

  • info>0: if info=k, and k<=n, the leading principal minor of order k of A is not positive, so the factorization could not be completed, and the solution has not been computed. rcond=0 is returned.

  • info=n+1: U is nonsingular, but rcond is less than machine precision, meaning that the matrix is singular to working precision. Nevertheless, the solution and error bounds are computed because there are a number of situations where the computed solution can be more accurate than the value of rcond would suggest.

Functions

void cposvx(
    const char*          fact,
    const char*          uplo,
    const INT            n,
    const INT            nrhs,
          c64*  restrict A,
    const INT            lda,
          c64*  restrict AF,
    const INT            ldaf,
          char*          equed,
          f32*  restrict S,
          c64*  restrict B,
    const INT            ldb,
          c64*  restrict X,
    const INT            ldx,
          f32*           rcond,
          f32*  restrict ferr,
          f32*  restrict berr,
          c64*  restrict work,
          f32*  restrict rwork,
          INT*           info
);
void cposvx(const char *fact, const char *uplo, const INT n, const INT nrhs, c64 *restrict A, const INT lda, c64 *restrict AF, const INT ldaf, char *equed, f32 *restrict S, c64 *restrict B, const INT ldb, c64 *restrict X, const INT ldx, f32 *rcond, f32 *restrict ferr, f32 *restrict berr, c64 *restrict work, f32 *restrict rwork, INT *info)#

CPOSVX uses the Cholesky factorization A = U**H*U or A = L*L**H to compute the solution to a complex system of linear equations.

A * X = B

where A is an n-by-n Hermitian positive definite matrix and X and B are n-by-nrhs matrices.

Error bounds on the solution and a condition estimate are also provided.

The following steps are performed:

  1. If fact='E', real scaling factors are computed to equilibrate the system:

    diag(S)*A*diag(S) * inv(diag(S))*X = diag(S)*B
    

    Whether or not the system will be equilibrated depends on the scaling of the matrix A, but if equilibration is used, A is overwritten by diag(S)*A*diag(S) and B by diag(S)*B.

  2. If fact='N' or 'E', the Cholesky decomposition is used to factor the matrix A (after equilibration if fact='E') as:

    A = U**H * U,  if uplo = 'U', or
    A = L * L**H,  if uplo = 'L',
    

    where U is an upper triangular matrix and L is a lower triangular matrix.

  3. If the leading principal minor of order i is not positive, then the routine returns with info=i. Otherwise, the factored form of A is used to estimate the condition number of the matrix A. If the reciprocal of the condition number is less than machine precision, info=n+1 is returned as a warning, but the routine still goes on to solve for X and compute error bounds as described below.

  4. The system of equations is solved for X using the factored form of A.

  5. Iterative refinement is applied to improve the computed solution matrix and calculate error bounds and backward error estimates for it.

  6. If equilibration was used, the matrix X is premultiplied by diag(S) so that it solves the original system before equilibration.

equed is an input argument if fact='F'; otherwise it is an output argument.

Parameters

in
fact

  • 'F': AF contains the factored form of A. If equed='Y', A has been equilibrated with scaling factors given by S; A and AF will not be modified.

  • 'N': The matrix A will be copied to AF and factored.

  • 'E': The matrix A will be equilibrated if necessary, then copied to AF and factored.

in
uplo

  • 'U': Upper triangle of A is stored

  • 'L': Lower triangle of A is stored

in
n

The number of linear equations. n>=0.

in
nrhs

The number of right hand sides. nrhs>=0.

inout
A

Array of dimension (lda, n). On entry, the Hermitian matrix A, except if fact='F' and equed='Y', then A must contain the equilibrated matrix diag(S)*A*diag(S). A is not modified if fact='F' or 'N', or if fact='E' and equed='N' on exit. On exit, if fact='E' and equed='Y', A is overwritten by diag(S)*A*diag(S).

in
lda

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

inout
AF

Array of dimension (ldaf, n). If fact='F', an input argument containing the triangular factor U or L from the Cholesky factorization, in the same storage format as A. If fact='N' or 'E', an output argument returning the triangular factor of the (possibly equilibrated) matrix A.

in
ldaf

The leading dimension of the array AF. ldaf>=max(1,n).

inout
equed

  • 'N': No equilibration (always true if fact='N')

  • 'Y': Equilibration was done, i.e., A has been replaced by diag(S) * A * diag(S)

inout
S

Array of dimension (n). The scale factors for A; not accessed if equed='N'. S is an input argument if fact='F'; otherwise it is an output argument. If fact='F' and equed='Y', each element of S must be positive.

inout
B

Array of dimension (ldb, nrhs). On entry, the n-by-nrhs right hand side matrix B. On exit, if equed='N', B is not modified; if equed='Y', B is overwritten by diag(S) * B.

in
ldb

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

out
X

Array of dimension (ldx, nrhs). If info=0 or info=n+1, the n-by-nrhs solution matrix X to the original system of equations. Note that if equed='Y', A and B are modified on exit, and the solution to the equilibrated system is inv(diag(S))*X.

in
ldx

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

out
rcond

The estimate of the reciprocal condition number of the matrix A after equilibration (if done).

out
ferr

Array of dimension (nrhs). The estimated forward error bound for each solution vector.

out
berr

Array of dimension (nrhs). The componentwise relative backward error of each solution vector.

out
work

Complex array of dimension 2*n.

out
rwork

Array of dimension n.

out
info

  • info=0: successful exit

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

  • info>0: if info=k, and k<=n, the leading principal minor of order k of A is not positive, so the factorization could not be completed, and the solution has not been computed. rcond=0 is returned.

  • info=n+1: U is nonsingular, but rcond is less than machine precision, meaning that the matrix is singular to working precision. Nevertheless, the solution and error bounds are computed because there are a number of situations where the computed solution can be more accurate than the value of rcond would suggest.

Functions

void zposvx(
    const char*          fact,
    const char*          uplo,
    const INT            n,
    const INT            nrhs,
          c128* restrict A,
    const INT            lda,
          c128* restrict AF,
    const INT            ldaf,
          char*          equed,
          f64*  restrict S,
          c128* restrict B,
    const INT            ldb,
          c128* restrict X,
    const INT            ldx,
          f64*           rcond,
          f64*  restrict ferr,
          f64*  restrict berr,
          c128* restrict work,
          f64*  restrict rwork,
          INT*           info
);
void zposvx(const char *fact, const char *uplo, const INT n, const INT nrhs, c128 *restrict A, const INT lda, c128 *restrict AF, const INT ldaf, char *equed, f64 *restrict S, c128 *restrict B, const INT ldb, c128 *restrict X, const INT ldx, f64 *rcond, f64 *restrict ferr, f64 *restrict berr, c128 *restrict work, f64 *restrict rwork, INT *info)#

ZPOSVX uses the Cholesky factorization A = U**H*U or A = L*L**H to compute the solution to a complex system of linear equations.

A * X = B

where A is an n-by-n Hermitian positive definite matrix and X and B are n-by-nrhs matrices.

Error bounds on the solution and a condition estimate are also provided.

The following steps are performed:

  1. If fact='E', real scaling factors are computed to equilibrate the system:

    diag(S)*A*diag(S) * inv(diag(S))*X = diag(S)*B
    

    Whether or not the system will be equilibrated depends on the scaling of the matrix A, but if equilibration is used, A is overwritten by diag(S)*A*diag(S) and B by diag(S)*B.

  2. If fact='N' or 'E', the Cholesky decomposition is used to factor the matrix A (after equilibration if fact='E') as:

    A = U**H * U,  if uplo = 'U', or
    A = L * L**H,  if uplo = 'L',
    

    where U is an upper triangular matrix and L is a lower triangular matrix.

  3. If the leading principal minor of order i is not positive, then the routine returns with info=i. Otherwise, the factored form of A is used to estimate the condition number of the matrix A. If the reciprocal of the condition number is less than machine precision, info=n+1 is returned as a warning, but the routine still goes on to solve for X and compute error bounds as described below.

  4. The system of equations is solved for X using the factored form of A.

  5. Iterative refinement is applied to improve the computed solution matrix and calculate error bounds and backward error estimates for it.

  6. If equilibration was used, the matrix X is premultiplied by diag(S) so that it solves the original system before equilibration.

equed is an input argument if fact='F'; otherwise it is an output argument.

Parameters

in
fact

  • 'F': AF contains the factored form of A. If equed='Y', A has been equilibrated with scaling factors given by S; A and AF will not be modified.

  • 'N': The matrix A will be copied to AF and factored.

  • 'E': The matrix A will be equilibrated if necessary, then copied to AF and factored.

in
uplo

  • 'U': Upper triangle of A is stored

  • 'L': Lower triangle of A is stored

in
n

The number of linear equations. n>=0.

in
nrhs

The number of right hand sides. nrhs>=0.

inout
A

Array of dimension (lda, n). On entry, the Hermitian matrix A, except if fact='F' and equed='Y', then A must contain the equilibrated matrix diag(S)*A*diag(S). A is not modified if fact='F' or 'N', or if fact='E' and equed='N' on exit. On exit, if fact='E' and equed='Y', A is overwritten by diag(S)*A*diag(S).

in
lda

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

inout
AF

Array of dimension (ldaf, n). If fact='F', an input argument containing the triangular factor U or L from the Cholesky factorization, in the same storage format as A. If fact='N' or 'E', an output argument returning the triangular factor of the (possibly equilibrated) matrix A.

in
ldaf

The leading dimension of the array AF. ldaf>=max(1,n).

inout
equed

  • 'N': No equilibration (always true if fact='N')

  • 'Y': Equilibration was done, i.e., A has been replaced by diag(S) * A * diag(S)

inout
S

Array of dimension (n). The scale factors for A; not accessed if equed='N'. S is an input argument if fact='F'; otherwise it is an output argument. If fact='F' and equed='Y', each element of S must be positive.

inout
B

Array of dimension (ldb, nrhs). On entry, the n-by-nrhs right hand side matrix B. On exit, if equed='N', B is not modified; if equed='Y', B is overwritten by diag(S) * B.

in
ldb

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

out
X

Array of dimension (ldx, nrhs). If info=0 or info=n+1, the n-by-nrhs solution matrix X to the original system of equations. Note that if equed='Y', A and B are modified on exit, and the solution to the equilibrated system is inv(diag(S))*X.

in
ldx

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

out
rcond

The estimate of the reciprocal condition number of the matrix A after equilibration (if done).

out
ferr

Array of dimension (nrhs). The estimated forward error bound for each solution vector.

out
berr

Array of dimension (nrhs). The componentwise relative backward error of each solution vector.

out
work

Complex array of dimension 2*n.

out
rwork

Array of dimension n.

out
info

  • info=0: successful exit

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

  • info>0: if info=k, and k<=n, the leading principal minor of order k of A is not positive, so the factorization could not be completed, and the solution has not been computed. rcond=0 is returned.

  • info=n+1: U is nonsingular, but rcond is less than machine precision, meaning that the matrix is singular to working precision. Nevertheless, the solution and error bounds are computed because there are a number of situations where the computed solution can be more accurate than the value of rcond would suggest.