SRC/sspsvx.f(3)

📦
version 3.12.1 · 4 versions
, in arch rolling

💅 SRC/sspsvx.f — the girlies' guide to SSPSVX

okok so this is a LAPACK man page, bestie — it lives in SRC/sspsvx.f in the library functions manual, and honestly? it's kind of a baddie. it's a Fortran subroutine that solves symmetric packed linear systems and then checks its own work. main character energy fr fr. 🎀

✨ SYNOPSIS ✨

Functions/Subroutines

subroutine sspsvx (fact, uplo, n, nrhs, ap, afp, ipiv, b, ldb, x, ldx,
    rcond, ferr, berr, work, iwork, info)

SSPSVX computes the solution to system of linear equations A * X = B for OTHER matrices — that's the one-line summary, and yes it's giving "I contain multitudes." 💜

🦋 Function/Subroutine Documentation 🦋

subroutine sspsvx (character fact, character uplo, integer n, integer nrhs,
    real, dimension( * ) ap, real, dimension( * ) afp, integer, dimension( *
    ) ipiv, real, dimension( ldb, * ) b, integer ldb, real, dimension( ldx,
    * ) x, integer ldx, real rcond, real, dimension( * ) ferr, real,
    dimension( * ) berr, real, dimension( * ) work, integer, dimension( * )
    iwork, integer info)

SSPSVX computes the solution to system of linear equations A * X = B for OTHER matrices. iconic, no notes. 🌟

🌸 Purpose 🌸

SSPSVX uses the diagonal pivoting factorization A = U*D*UT or A = L*D*LT to compute the solution to a real system of linear equations A * X = B, where A is an N-by-N symmetric matrix stored in packed format and X and B are N-by-NRHS matrices. Error bounds on the solution and a condition estimate are also provided.

so basically: solve the thing, then tell you how much you should trust the answer. she's thorough. she's that girl. 💅

🧚 Description 🧚

The following steps are performed: 1. If FACT = 'N', the diagonal pivoting method is used to factor A as A = U * D * UT, if UPLO = 'U', or A = L * D * LT, if UPLO = 'L', where U (or L) is a product of permutation and unit upper (lower) triangular matrices and D is symmetric and block diagonal with 1-by-1 and 2-by-2 diagonal blocks. 2. If some D(i,i)=0, so that D is exactly singular, 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. 3. The system of equations is solved for X using the factored form of A. 4. Iterative refinement is applied to improve the computed solution matrix and calculate error bounds and backward error estimates for it.

four steps, zero skips, she does NOT cut corners. 🫶

👑 Parameters 👑

FACT

FACT is CHARACTER*1 Specifies whether or not the factored form of A has been supplied on entry. = 'F': On entry, AFP and IPIV contain the factored form of A. AP, AFP and IPIV will not be modified. = 'N': The matrix A will be copied to AFP and factored.

tl;dr — is the factorization already done or do we have to do it ourselves? 🤔

UPLO

UPLO is CHARACTER*1 = 'U': Upper triangle of A is stored; = 'L': Lower triangle of A is stored.

pick your triangle, babe. only one gets stored, she's not storing both. ✨

N

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

how big is the matrix, basically. N >= 0, bestie.

NRHS

NRHS is INTEGER The number of right hand sides, i.e., the number of columns of the matrices B and X. NRHS >= 0.

how many right-hand-side columns we're doing at once. NRHS >= 0.

AP

AP is REAL array, dimension (N*(N+1)/2) The upper or lower triangle of the symmetric matrix A, packed columnwise in a linear array. The j-th column of A is stored in the array AP as follows: if UPLO = 'U', AP(i + (j-1)j/2) = A(i,j) for 1<=i<=j; if UPLO = 'L', AP(i + (j-1)(2*n-j)/2) = A(i,j) for j<=i<=n. See below for further details.

packed storage my beloved 😭 the indexing looks scary but there's an example further down, I gotchu.

AFP

AFP is REAL array, dimension (N*(N+1)/2) If FACT = 'F', then AFP is an input argument and on entry contains the block diagonal matrix D and the multipliers used to obtain the factor U or L from the factorization A = U*D*UT or A = L*D*LT as computed by SSPTRF, stored as a packed triangular matrix in the same storage format as A. If FACT = 'N', then AFP is an output argument and on exit contains the block diagonal matrix D and the multipliers used to obtain the factor U or L from the factorization A = U*D*UT or A = L*D*LT as computed by SSPTRF, stored as a packed triangular matrix in the same storage format as A.

input if FACT = 'F', output if FACT = 'N' — she switches roles depending on the vibe check. 💀

IPIV

IPIV is INTEGER array, dimension (N) If FACT = 'F', then IPIV is an input argument and on entry contains details of the interchanges and the block structure of D, as determined by SSPTRF. If IPIV(k) > 0, then rows and columns k and IPIV(k) were interchanged and D(k,k) is a 1-by-1 diagonal block. If UPLO = 'U' and IPIV(k) = IPIV(k-1) < 0, then rows and columns k-1 and -IPIV(k) were interchanged and D(k-1:k,k-1:k) is a 2-by-2 diagonal block. If UPLO = 'L' and IPIV(k) = IPIV(k+1) < 0, then rows and columns k+1 and -IPIV(k) were interchanged and D(k:k+1,k:k+1) is a 2-by-2 diagonal block. If FACT = 'N', then IPIV is an output argument and on exit contains details of the interchanges and the block structure of D, as determined by SSPTRF.

the swap receipts 🧾 plus the block structure of D. negative values = 2-by-2 block situation, and which neighbor it pairs with depends on UPLO (k-1 for 'U', k+1 for 'L'). iykyk.

B

B is REAL array, dimension (LDB,NRHS) The N-by-NRHS right hand side matrix B.

the right hand side. the problem itself. 🍓

LDB

LDB is INTEGER The leading dimension of the array B. LDB >= max(1,N).

leading dimension of B, and yes she'll check: LDB >= max(1,N).

X

X is REAL array, dimension (LDX,NRHS) If INFO = 0 or INFO = N+1, the N-by-NRHS solution matrix X.

the answer matrix!! only filled in if INFO = 0 or INFO = N+1.

LDX

LDX is INTEGER The leading dimension of the array X. LDX >= max(1,N).

leading dimension of X: LDX >= max(1,N). same rules, different array. 💖

RCOND

RCOND is REAL The estimate of the reciprocal condition number of the matrix A. If RCOND is less than the machine precision (in particular, if RCOND = 0), the matrix is singular to working precision. This condition is indicated by a return code of INFO > 0.

the trustworthiness score. if it dips below machine precision (especially RCOND = 0), the matrix is singular to working precision and you get INFO > 0. not me warning you in advance 😭

FERR

FERR is REAL array, dimension (NRHS) The estimated forward error bound for each solution vector X(j) (the j-th column of the solution matrix X). If XTRUE is the true solution corresponding to X(j), FERR(j) is an estimated upper bound for the magnitude of the largest element in (X(j) - XTRUE) divided by the magnitude of the largest element in X(j). The estimate is as reliable as the estimate for RCOND, and is almost always a slight overestimate of the true error.

one forward error bound per solution column, and it's almost always a slight overestimate — she's cautious like that. 🎀

BERR

BERR is REAL array, dimension (NRHS) The componentwise relative backward error of each solution vector X(j) (i.e., the smallest relative change in any element of A or B that makes X(j) an exact solution).

how much you'd have to nudge A or B for X(j) to be exactly right. componentwise, relative, per column. 🌈

WORK

WORK is REAL array, dimension (3*N)

workspace, 3*N of it. do not ask questions, just pass the array. 💜

IWORK

IWORK is INTEGER array, dimension (N)

integer workspace, size N. 🙈

INFO

INFO is INTEGER = 0: successful exit < 0: if INFO = -i, the i-th argument had an illegal value > 0: if INFO = i, and i is <= N: D(i,i) is exactly zero. The factorization has been completed but the factor D is exactly singular, so the solution and error bounds could not be computed. RCOND = 0 is returned. = N+1: D 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.

the report card:

  • = 0 — slayed, successful exit ✨
  • < 0 — if INFO = -i, the i-th argument had an illegal value (you passed something cursed)
  • > 0 — if INFO = i:
    • <= N: D(i,i) is exactly zero. factorization finished but D is exactly singular, so solution and error bounds could NOT be computed, and RCOND = 0 is returned. 💀
    • = N+1: D is nonsingular but RCOND is less than machine precision — matrix is singular to working precision. still computes solution and error bounds anyway, because sometimes the computed solution is more accurate than RCOND would suggest. she's generous like that. 🫶

⭐ Author ⭐

  • Univ. of Tennessee
  • Univ. of California Berkeley
  • Univ. of Colorado Denver
  • NAG Ltd.

the whole squad contributed, we love a collab. 👑

🦋 Further Details: 🦋

The packed storage scheme is illustrated by the following example when N = 4, UPLO = 'U': Two-dimensional storage of the symmetric matrix A: a11 a12 a13 a14 a22 a23 a24 a33 a34 (aij = aji) a44 Packed storage of the upper triangle of A: AP = [ a11, a12, a22, a13, a23, a33, a14, a24, a34, a44 ]

see!! not scary at all once you see it laid out. column by column, upper triangle only, (aij = aji) because she's symmetric obviously. 🌸

📍 Where she lives

Definition at line 274 of file sspsvx.f.

💖 Author (of this page, not the subroutine)

Generated automatically by Doxygen for LAPACK from the source code.

ngl if you got this far you're lowkey a LAPACK girlie now. take the win. 🏆

LAPACK · Version 3.12.0 · SRC/sspsvx.f(3) — that's the version and date stamp at the bottom of the original page, mwah 💋 bye bestie!! 🎀✨

📦 where this page comes from