complex16HEcomputational(3)

📦

version 3.11.0 · 2 versions

3.11.0 · debian bookworm (latest)
3.10.0 · ubuntu jammy

✨ complex16HEcomputational ✨

okay so bestie, welcome to complex16HEcomputational(3) — this is the LAPACK manual page for the entire squad of complex16 computational functions for HE matrices 💖 it's giving big group-chat energy, it's giving "let me introduce you to all my friends", and honestly? every single one of them ate. 👑

🌸 NAME 🌸

complex16HEcomputational — complex16 🦋

✨ SYNOPSIS ✨

🫶 Functions 🫶

subroutine zhecon (UPLO, N, A, LDA, IPIV, ANORM, RCOND, WORK, INFO) — ZHECON
subroutine zhecon_3 (UPLO, N, A, LDA, E, IPIV, ANORM, RCOND, WORK, INFO) — ZHECON_3
subroutine zhecon_rook (UPLO, N, A, LDA, IPIV, ANORM, RCOND, WORK, INFO) — ZHECON_ROOK estimates the reciprocal of the condition number for HE matrices using the factorization obtained with one of the bounded diagonal pivoting methods (max 2 interchanges).
subroutine zheequb (UPLO, N, A, LDA, S, SCOND, AMAX, WORK, INFO) — ZHEEQUB
subroutine zhegs2 (ITYPE, UPLO, N, A, LDA, B, LDB, INFO) — ZHEGS2 reduces a Hermitian definite generalized eigenproblem to standard form, using the factorization results obtained from cpotrf (unblocked algorithm).
subroutine zhegst (ITYPE, UPLO, N, A, LDA, B, LDB, INFO) — ZHEGST
subroutine zherfs (UPLO, N, NRHS, A, LDA, AF, LDAF, IPIV, B, LDB, X, LDX, FERR, BERR, WORK, RWORK, INFO) — ZHERFS
subroutine zherfsx (UPLO, EQUED, N, NRHS, A, LDA, AF, LDAF, IPIV, S, B, LDB, X, LDX, RCOND, BERR, N_ERR_BNDS, ERR_BNDS_NORM, ERR_BNDS_COMP, NPARAMS, PARAMS, WORK, RWORK, INFO) — ZHERFSX
subroutine zhetd2 (UPLO, N, A, LDA, D, E, TAU, INFO) — ZHETD2 reduces a Hermitian matrix to real symmetric tridiagonal form by an unitary similarity transformation (unblocked algorithm).
subroutine zhetf2 (UPLO, N, A, LDA, IPIV, INFO) — ZHETF2 computes the factorization of a complex Hermitian matrix, using the diagonal pivoting method (unblocked algorithm, calling Level 2 BLAS).
subroutine zhetf2_rk (UPLO, N, A, LDA, E, IPIV, INFO) — ZHETF2_RK computes the factorization of a complex Hermitian indefinite matrix using the bounded Bunch-Kaufman (rook) diagonal pivoting method (BLAS2 unblocked algorithm).
subroutine zhetf2_rook (UPLO, N, A, LDA, IPIV, INFO) — ZHETF2_ROOK computes the factorization of a complex Hermitian indefinite matrix using the bounded Bunch-Kaufman ('rook') diagonal pivoting method (unblocked algorithm).
subroutine zhetrd (UPLO, N, A, LDA, D, E, TAU, WORK, LWORK, INFO) — ZHETRD
subroutine zhetrd_2stage (VECT, UPLO, N, A, LDA, D, E, TAU, HOUS2, LHOUS2, WORK, LWORK, INFO) — ZHETRD_2STAGE
subroutine zhetrd_he2hb (UPLO, N, KD, A, LDA, AB, LDAB, TAU, WORK, LWORK, INFO) — ZHETRD_HE2HB
subroutine zhetrf (UPLO, N, A, LDA, IPIV, WORK, LWORK, INFO) — ZHETRF
subroutine zhetrf_aa (UPLO, N, A, LDA, IPIV, WORK, LWORK, INFO) — ZHETRF_AA
subroutine zhetrf_rk (UPLO, N, A, LDA, E, IPIV, WORK, LWORK, INFO) — ZHETRF_RK computes the factorization of a complex Hermitian indefinite matrix using the bounded Bunch-Kaufman (rook) diagonal pivoting method (BLAS3 blocked algorithm).
subroutine zhetrf_rook (UPLO, N, A, LDA, IPIV, WORK, LWORK, INFO) — ZHETRF_ROOK computes the factorization of a complex Hermitian indefinite matrix using the bounded Bunch-Kaufman ('rook') diagonal pivoting method (blocked algorithm, calling Level 3 BLAS).
subroutine zhetri (UPLO, N, A, LDA, IPIV, WORK, INFO) — ZHETRI
subroutine zhetri2 (UPLO, N, A, LDA, IPIV, WORK, LWORK, INFO) — ZHETRI2
subroutine zhetri2x (UPLO, N, A, LDA, IPIV, WORK, NB, INFO) — ZHETRI2X
subroutine zhetri_3 (UPLO, N, A, LDA, E, IPIV, WORK, LWORK, INFO) — ZHETRI_3
subroutine zhetri_3x (UPLO, N, A, LDA, E, IPIV, WORK, NB, INFO) — ZHETRI_3X
subroutine zhetri_rook (UPLO, N, A, LDA, IPIV, WORK, INFO) — ZHETRI_ROOK computes the inverse of HE matrix using the factorization obtained with the bounded Bunch-Kaufman ('rook') diagonal pivoting method.
subroutine zhetrs (UPLO, N, NRHS, A, LDA, IPIV, B, LDB, INFO) — ZHETRS
subroutine zhetrs2 (UPLO, N, NRHS, A, LDA, IPIV, B, LDB, WORK, INFO) — ZHETRS2
subroutine zhetrs_3 (UPLO, N, NRHS, A, LDA, E, IPIV, B, LDB, INFO) — ZHETRS_3
subroutine zhetrs_aa (UPLO, N, NRHS, A, LDA, IPIV, B, LDB, WORK, LWORK, INFO) — ZHETRS_AA
subroutine zhetrs_rook (UPLO, N, NRHS, A, LDA, IPIV, B, LDB, INFO) — ZHETRS_ROOK computes the solution to a system of linear equations A * X = B for HE matrices using the factorization obtained with one of the bounded diagonal pivoting methods (max 2 interchanges).
subroutine zla_heamv (UPLO, N, ALPHA, A, LDA, X, INCX, BETA, Y, INCY) — ZLA_HEAMV computes a matrix-vector product using a Hermitian indefinite matrix to calculate error bounds.
double precision function zla_hercond_c (UPLO, N, A, LDA, AF, LDAF, IPIV, C, CAPPLY, INFO, WORK, RWORK) — ZLA_HERCOND_C computes the infinity norm condition number of op(A)*inv(diag(c)) for Hermitian indefinite matrices.
double precision function zla_hercond_x (UPLO, N, A, LDA, AF, LDAF, IPIV, X, INFO, WORK, RWORK) — ZLA_HERCOND_X computes the infinity norm condition number of op(A)*diag(x) for Hermitian indefinite matrices.
subroutine zla_herfsx_extended (PREC_TYPE, UPLO, N, NRHS, A, LDA, AF, LDAF, IPIV, COLEQU, C, B, LDB, Y, LDY, BERR_OUT, N_NORMS, ERR_BNDS_NORM, ERR_BNDS_COMP, RES, AYB, DY, Y_TAIL, RCOND, ITHRESH, RTHRESH, DZ_UB, IGNORE_CWISE, INFO) — ZLA_HERFSX_EXTENDED improves the computed solution to a system of linear equations for Hermitian indefinite matrices by performing extra-precise iterative refinement and provides error bounds and backward error estimates for the solution.
double precision function zla_herpvgrw (UPLO, N, INFO, A, LDA, AF, LDAF, IPIV, WORK) — ZLA_HERPVGRW
subroutine zlahef (UPLO, N, NB, KB, A, LDA, IPIV, W, LDW, INFO) — ZLAHEF computes a partial factorization of a complex Hermitian indefinite matrix using the Bunch-Kaufman diagonal pivoting method (blocked algorithm, calling Level 3 BLAS).
subroutine zlahef_aa (UPLO, J1, M, NB, A, LDA, IPIV, H, LDH, WORK) — ZLAHEF_AA
subroutine zlahef_rk (UPLO, N, NB, KB, A, LDA, E, IPIV, W, LDW, INFO) — ZLAHEF_RK computes a partial factorization of a complex Hermitian indefinite matrix using bounded Bunch-Kaufman (rook) diagonal pivoting method.
subroutine zlahef_rook (UPLO, N, NB, KB, A, LDA, IPIV, W, LDW, INFO) — ZLAHEF_ROOK

💜 Detailed Description 💜

This is the group of complex16 computational functions for HE matrices. that's it. that's the whole description. short, sweet, and absolutely iconic of them 💅

👑 Function Documentation 👑

🌟 subroutine zhecon 🌟

subroutine zhecon (character UPLO, integer N, complex*16, dimension( lda, * ) A, integer LDA, integer, dimension( * ) IPIV, double precision ANORM, double precision RCOND, complex*16, dimension( * ) WORK, integer INFO)

ZHECON

🎀 Purpose 🎀

okay so ZHECON estimates the reciprocal of the condition number of a complex Hermitian matrix A using the factorization A = U*D*U**H or A = L*D*L**H computed by ZHETRF.

An estimate is obtained for norm(inv(A)), and the reciprocal of the condition number is computed as:

RCOND = 1 / (ANORM * norm(inv(A)))

🧚 Parameters 🧚

UPLO — UPLO is CHARACTER*1. Tells you whether the details of the factorization are stored as an upper or lower triangular matrix.
- = 'U': Upper triangular, form is A = U*D*U**H;
- = 'L': Lower triangular, form is A = L*D*L**H.
N — N is INTEGER. The order of the matrix A. N >= 0.
A — A is COMPLEX*16 array, dimension (LDA,N). The block diagonal matrix D and the multipliers used to obtain the factor U or L as computed by ZHETRF.
LDA — LDA is INTEGER. The leading dimension of the array A. LDA >= max(1,N).
IPIV — IPIV is INTEGER array, dimension (N). Details of the interchanges and the block structure of D as determined by ZHETRF.
ANORM — ANORM is DOUBLE PRECISION. The 1-norm of the original matrix A.
RCOND — RCOND is DOUBLE PRECISION. The reciprocal of the condition number of the matrix A, computed as RCOND = 1/(ANORM * AINVNM), where AINVNM is an estimate of the 1-norm of inv(A) computed in this routine.
WORK — WORK is COMPLEX*16 array, dimension (2*N).
INFO — INFO is INTEGER.
- = 0: successful exit
- < 0: if INFO = -i, the i-th argument had an illegal value

🫶 Author 🫶

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

🌟 subroutine zhecon_3 🌟

subroutine zhecon_3 (character UPLO, integer N, complex*16, dimension( lda, * ) A, integer LDA, complex*16, dimension( * ) E, integer, dimension( * ) IPIV, double precision ANORM, double precision RCOND, complex*16, dimension( * ) WORK, integer INFO)

ZHECON_3

🎀 Purpose 🎀

ZHECON_3 estimates the reciprocal of the condition number (in the 1-norm) of a complex Hermitian matrix A using the factorization computed by ZHETRF_RK or ZHETRF_BK:

A = P*U*D*(U**H)*(P**T) or A = P*L*D*(L**H)*(P**T),

where U (or L) is unit upper (or lower) triangular matrix, U**H (or L**H) is the conjugate of U (or L), P is a permutation matrix, P**T is the transpose of P, and D is Hermitian and block diagonal with 1-by-1 and 2-by-2 diagonal blocks.

An estimate is obtained for norm(inv(A)), and the reciprocal of the condition number is computed as RCOND = 1 / (ANORM * norm(inv(A))). This routine uses BLAS3 solver ZHETRS_3.

🧚 Parameters 🧚

UPLO — UPLO is CHARACTER*1. Specifies whether the details of the factorization are stored as an upper or lower triangular matrix:
- = 'U': Upper triangular, form is A = P*U*D*(U**H)*(P**T);
- = 'L': Lower triangular, form is A = P*L*D*(L**H)*(P**T).
N — N is INTEGER. The order of the matrix A. N >= 0.
A — A is COMPLEX*16 array, dimension (LDA,N). Diagonal of the block diagonal matrix D and factors U or L as computed by ZHETRF_RK and ZHETRF_BK:
1. ONLY diagonal elements of the Hermitian block diagonal matrix D on the diagonal of A, i.e. D(k,k) = A(k,k); (superdiagonal (or subdiagonal) elements of D should be provided on entry in array E), and
2. If UPLO = 'U': factor U in the superdiagonal part of A. If UPLO = 'L': factor L in the subdiagonal part of A.
LDA — LDA is INTEGER. The leading dimension of the array A. LDA >= max(1,N).
E — and bestie, the description of E is literally where this part of the page cuts off, so hold that thought 💖
E — COMPLEX*16 array, dimension (N). On entry it's holding the superdiagonal (or subdiagonal, depending on which side you're on) elements of the Hermitian block diagonal matrix D, which has 1-by-1 or 2-by-2 diagonal blocks. Specifically:
- If UPLO = 'U': E(i) = D(i-1,i), i=2:N, and E(1) is not referenced;
- If UPLO = 'L': E(i) = D(i+1,i), i=1:N-1, and E(N) is not referenced.

NOTE: for a 1-by-1 diagonal block D(k), where 1 <= k <= N, the element E(k) is not referenced — doesn't matter if UPLO = 'U' or UPLO = 'L', it just sits there unbothered. 🙈

IPIV — INTEGER array, dimension (N). All the deets on the interchanges and the block structure of D, as determined by ZHETRF_RK or ZHETRF_BK. Like, the receipts. 🧾
ANORM — DOUBLE PRECISION. The 1-norm of the original matrix A.
RCOND — DOUBLE PRECISION. The reciprocal of the condition number of the matrix A, computed as RCOND = 1/(ANORM * AINVNM), where AINVNM is an estimate of the 1-norm of inv(A) computed in this routine. It's giving "how invertible are we actually". 💅
WORK — COMPLEX*16 array, dimension (2*N).
INFO — INTEGER.
- = 0: successful exit ✨
- < 0: if INFO = -i, the i-th argument had an illegal value 😭

Author

Univ. of Tennessee

Univ. of California Berkeley

Univ. of Colorado Denver

NAG Ltd.

Contributors:

June 2017,  Igor Kozachenko,
                Computer Science Division,
                University of California, Berkeley

September 2007, Sven Hammarling, Nicholas J. Higham, Craig Lucas,
                School of Mathematics,
                University of Manchester

🦋 `zhecon_rook` 🦋

subroutine zhecon_rook (character UPLO, integer N, complex*16, dimension(
    lda, * ) A, integer LDA, integer, dimension( * ) IPIV, double precision
    ANORM, double precision RCOND, complex*16, dimension( * ) WORK, integer
    INFO)

ZHECON_ROOK estimates the reciprocal of the condition number fort HE matrices using factorization obtained with one of the bounded diagonal pivoting methods (max 2 interchanges)

Purpose:

ZHECON_ROOK estimates the reciprocal of the condition number of a complex
Hermitian matrix A using the factorization A = U*D*U**H or
A = L*D*L**H computed by CHETRF_ROOK.

An estimate is obtained for norm(inv(A)), and the reciprocal of the
condition number is computed as RCOND = 1 / (ANORM * norm(inv(A))).

Parameters

UPLO — CHARACTER*1. Tells you whether the factorization details are stored as an upper or lower triangular matrix.
- = 'U': Upper triangular, form is A = U*D*U**H;
- = 'L': Lower triangular, form is A = L*D*L**H.
N — INTEGER. The order of the matrix A. N >= 0.
A — COMPLEX*16 array, dimension (LDA,N). The block diagonal matrix D plus the multipliers used to obtain the factor U or L, as computed by CHETRF_ROOK.
LDA — INTEGER. The leading dimension of the array A. LDA >= max(1,N).
IPIV — INTEGER array, dimension (N). The details of the interchanges and the block structure of D, as determined by CHETRF_ROOK. 👑
ANORM — DOUBLE PRECISION. The 1-norm of the original matrix A.
RCOND — DOUBLE PRECISION. The reciprocal of the condition number of the matrix A, computed as RCOND = 1/(ANORM * AINVNM), where AINVNM is an estimate of the 1-norm of inv(A) computed in this routine.
WORK — COMPLEX*16 array, dimension (2*N).
INFO — INTEGER.
- = 0: successful exit ✨
- < 0: if INFO = -i, the i-th argument had an illegal value 💀

Author

Univ. of Tennessee

Univ. of California Berkeley

Univ. of Colorado Denver

NAG Ltd.

Contributors:

June 2017,  Igor Kozachenko,
                Computer Science Division,
                University of California, Berkeley

September 2007, Sven Hammarling, Nicholas J. Higham, Craig Lucas,
                School of Mathematics,
                University of Manchester

🌸 `zheequb` 🌸

subroutine zheequb (character UPLO, integer N, complex*16, dimension( lda, *
    ) A, integer LDA, double precision, dimension( * ) S, double precision
    SCOND, double precision AMAX, complex*16, dimension( * ) WORK, integer
    INFO)

ZHEEQUB

Purpose:

ZHEEQUB computes row and column scalings intended to equilibrate a
Hermitian matrix A (with respect to the Euclidean norm) and reduce
its condition number. The scale factors S are computed by the BIN
algorithm (see references) so that the scaled matrix B with elements
B(i,j) = S(i)*A(i,j)*S(j) has a condition number within a factor N of
the smallest possible condition number over all possible diagonal
scalings.

Parameters

UPLO — CHARACTER*1.
- = 'U': Upper triangle of A is stored;
- = 'L': Lower triangle of A is stored.
N — INTEGER. The order of the matrix A. N >= 0.
A — COMPLEX*16 array, dimension (LDA,N). The N-by-N Hermitian matrix whose scaling factors are about to be computed. Bestie, we're finding its glow-up numbers. ✨
LDA — INTEGER. The leading dimension of the array A. LDA >= max(1,N).
S — DOUBLE PRECISION array, dimension (N). If INFO = 0, S contains the scale factors for A.
SCOND — DOUBLE PRECISION. If INFO = 0, S contains the ratio of the smallest S(i) to the largest S(i). If SCOND >= 0.1 and AMAX is neither too large nor too small, it is not worth scaling by S. Don't do it, it's a waste of time fr fr.
AMAX — DOUBLE PRECISION. Largest absolute value of any matrix element. If AMAX is very close to overflow or very close to underflow, the matrix should be scaled. Like, genuinely. That's not a vibe, that's a requirement. 🚨
WORK — COMPLEX*16 array, dimension (2*N).
INFO — INTEGER.
- = 0: successful exit ✨
- < 0: if INFO = -i, the i-th argument had an illegal value 😭
- > 0: if INFO = i, the i-th diagonal element is nonpositive.

Author

Univ. of Tennessee

Univ. of California Berkeley

Univ. of Colorado Denver

NAG Ltd.

References:

Livne, O.E. and Golub, G.H., 'Scaling by Binormalization',
 Numerical Algorithms, vol. 35, no. 1, pp. 97-120, January 2004.
 DOI 10.1023/B:NUMA.0000016606.32820.69
 Tech report version:
http://citeseerx.ist.psu.edu/viewdoc/summary?doi=10.1.1.3.1679

🎀 `zhegs2` 🎀

subroutine zhegs2 (integer ITYPE, character UPLO, integer N, complex*16,
    dimension( lda, * ) A, integer LDA, complex*16, dimension( ldb, * ) B,
    integer LDB, integer INFO)

ZHEGS2 reduces a Hermitian definite generalized eigenproblem to standard form, using the factorization results obtained from cpotrf (unblocked algorithm).

Purpose:

ZHEGS2 reduces a complex Hermitian-definite generalized
eigenproblem to standard form.

If ITYPE = 1, the problem is A*x = lambda*B*x,
and A is overwritten by inv(U**H)*A*inv(U) or inv(L)*A*inv(L**H)

If ITYPE = 2 or 3, the problem is A*B*x = lambda*x or
B*A*x = lambda*x, and A is overwritten by U*A*U**H or L**H *A*L.

B must have been previously factorized as U**H *U or L*L**H by ZPOTRF.

Parameters

ITYPE — INTEGER.
- = 1: compute inv(U**H)*A*inv(U) or inv(L)*A*inv(L**H);
- = 2 or 3: compute U*A*U**H or L**H *A*L.
UPLO — CHARACTER*1. Specifies whether the upper or lower triangular part of the Hermitian matrix A is stored, and how B has been factorized.
- = 'U': Upper triangular
- = 'L': Lower triangular
N — INTEGER. The order of the matrices A and B. N >= 0.
A — COMPLEX*16 array, dimension (LDA,N). On entry, the Hermitian matrix A. If UPLO = 'U', the leading n by n upper triangular part of A contains the upper triangular part of the matrix A, and the strictly lower triangular part of A is not referenced. If UPLO = 'L', the leading n by n lower triangular part of A contains the lower triangular part of the matrix A, and the strictly upper triangular part of A is not referenced.

On exit, if INFO = 0, the transformed matrix, stored in the same format as A. She's different now. 💅

LDA — INTEGER. The leading dimension of the array A. LDA >= max(1,N).
B — COMPLEX*16 array, dimension (LDB,N). The triangular factor from the Cholesky factorization of B, as returned by ZPOTRF. B is modified by the routine but restored on exit — so don't panic, she comes back. 🫶
LDB — INTEGER. The leading dimension of the array B. LDB >= max(1,N).
INFO — INTEGER.
- = 0: successful exit.
- < 0: if INFO = -i, the i-th argument had an illegal value.

Author

Univ. of Tennessee

Univ. of California Berkeley

Univ. of Colorado Denver

NAG Ltd.

🌈 `zhegst` 🌈

subroutine zhegst (integer ITYPE, character UPLO, integer N, complex*16,
    dimension( lda, * ) A, integer LDA, complex*16, dimension( ldb, * ) B,
    integer LDB, integer INFO)

ZHEGST

Purpose:

ZHEGST reduces a complex Hermitian-definite generalized
eigenproblem to standard form.

If ITYPE = 1, the problem is A*x = lambda*B*x,
and A is overwritten by inv(U**H)*A*inv(U) or inv(L)*A*inv(L**H)

If ITYPE = 2 or 3, the problem is A*B*x = lambda*x or
B*A*x = lambda*x, and A is overwritten by U*A*U**H or L**H*A*L.

B must have been previously factorized as U**H*U or L*L**H by ZPOTRF.

Parameters

ITYPE
ITYPE — INTEGER, the vibe-setter for what exactly we're computing:
- = 1: compute inv(U**H)*A*inv(U) or inv(L)*A*inv(L**H);
- = 2 or 3: compute U*A*U**H or L**H*A*L.
UPLO — CHARACTER1*, which triangle of A we're even looking at:
- = 'U': upper triangle of A is stored and B is factored as U**H*U;
- = 'L': lower triangle of A is stored and B is factored as L*L**H.
N — INTEGER, the order of the matrices A and B. Constraint check bestie: N >= 0 💅
A — COMPLEX*16 array, dimension (LDA,N). On entry, this is the Hermitian matrix A. If UPLO = 'U', the leading N-by-N upper triangular part of A holds the upper triangular part of the matrix A, and the strictly lower triangular part of A is not referenced (aka don't bother, it's not invited to the party). If UPLO = 'L', the leading N-by-N lower triangular part of A holds the lower triangular part of the matrix A, and the strictly upper triangular part of A is not referenced.

On exit, if INFO = 0, this is the transformed matrix, stored in the same format as A. glow up complete ✨

LDA — INTEGER, the leading dimension of the array A. LDA >= max(1,N). No cap 🧢
B — COMPLEX*16 array, dimension (LDB,N). The triangular factor from the Cholesky factorization of B, as returned by ZPOTRF.

Heads up: B is modified by the routine but restored on exit. She's a temporary guest star, not a permanent resident 🎀

LDB — INTEGER, the leading dimension of the array B. LDB >= max(1,N).
INFO — INTEGER, the receipts 📜:
- = 0: successful exit 🎉
- < 0: if INFO = -i, the i-th argument had an illegal value. not slay

Author

Univ. of Tennessee

Univ. of California Berkeley

Univ. of Colorado Denver

NAG Ltd.

`zherfs`

subroutine zherfs (character UPLO, integer N, integer NRHS, complex*16,
    dimension( lda, * ) A, integer LDA, complex*16, dimension( ldaf, * ) AF,
    integer LDAF, integer, dimension( * ) IPIV, complex*16, dimension( ldb,
    * ) B, integer LDB, complex*16, dimension( ldx, * ) X, integer LDX,
    double precision, dimension( * ) FERR, double precision, dimension( * )
    BERR, complex*16, dimension( * ) WORK, double precision, dimension( * )
    RWORK, integer INFO)

Purpose

ZHERFS improves the computed solution to a system of linear equations when the coefficient matrix is Hermitian indefinite, and provides error bounds and backward error estimates for the solution. She's not just giving you an answer, she's giving you an answer with receipts 😌✨

Parameters

UPLO — CHARACTER1*, which triangle situation we're in:
- = 'U': Upper triangle of A is stored;
- = 'L': Lower triangle of A is stored.
N — INTEGER, the order of the matrix A. N >= 0.
NRHS — INTEGER, the number of right hand sides, i.e., the number of columns of the matrices B and X. NRHS >= 0. Like, how many problems are we solving at once bestie 💖
A — COMPLEX*16 array, dimension (LDA,N). The Hermitian matrix A. If UPLO = 'U', the leading N-by-N upper triangular part of A contains the upper triangular part of the matrix A, and the strictly lower triangular part of A is not referenced. If UPLO = 'L', the leading N-by-N lower triangular part of A contains the lower triangular part of the matrix A, and the strictly upper triangular part of A is not referenced. Same energy as before: only one triangle gets a seat at the table 🪑
LDA — INTEGER, the leading dimension of the array A. LDA >= max(1,N).
AF — COMPLEX*16 array, dimension (LDAF,N). The factored form of the matrix A. AF contains the block diagonal matrix D and the multipliers used to obtain the factor U or L from the factorization A = U*D*U**H or A = L*D*L**H as computed by ZHETRF.
LDAF — INTEGER, the leading dimension of the array AF. LDAF >= max(1,N).
IPIV — INTEGER array, dimension (N). Details of the interchanges and the block structure of D as determined by ZHETRF. It's the lore 📖
B — COMPLEX*16 array, dimension (LDB,NRHS). The right hand side matrix B.
LDB — INTEGER, the leading dimension of the array B. LDB >= max(1,N).
X — COMPLEX*16 array, dimension (LDX,NRHS). On entry, the solution matrix X, as computed by ZHETRS. On exit, the improved solution matrix X. glow-up arc, but for numbers 💅
LDX — INTEGER, the leading dimension of the array X. LDX >= max(1,N).
FERR — DOUBLE PRECISION 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. She's a little dramatic, she overestimates, but she's trustworthy fr fr 🫶
BERR — DOUBLE PRECISION 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). Tiny lil nudge and suddenly it's perfect ✨
WORK — COMPLEX*16 array, dimension (2*N)
RWORK — DOUBLE PRECISION array, dimension (N)
INFO — INTEGER, the receipts:
- = 0: successful exit 🎉
- < 0: if INFO = -i, the i-th argument had an illegal value.

Internal Parameters

ITMAX is the maximum number of steps of iterative refinement. She keeps going until it's giving perfection, but there's a limit bestie 🎀

Author

Univ. of Tennessee

Univ. of California Berkeley

Univ. of Colorado Denver

NAG Ltd.

`zherfsx`

subroutine zherfsx (character UPLO, character EQUED, integer N, integer
    NRHS, complex*16, dimension( lda, * ) A, integer LDA, complex*16,
    dimension( ldaf, * ) AF, integer LDAF, integer, dimension( * ) IPIV,
    double precision, dimension( * ) S, complex*16, dimension( ldb, * ) B,
    integer LDB, complex*16, dimension( ldx, * ) X, integer LDX, double
    precision RCOND, double precision, dimension( * ) BERR, integer
    N_ERR_BNDS, double precision, dimension( nrhs, * ) ERR_BNDS_NORM, double
    precision, dimension( nrhs, * ) ERR_BNDS_COMP, integer NPARAMS, double
    precision, dimension( * ) PARAMS, complex*16, dimension( * ) WORK,
    double precision, dimension( * ) RWORK, integer INFO)

Purpose

ZHERFSX improves the computed solution to a system of linear equations when the coefficient matrix is Hermitian indefinite, and provides error bounds and backward error estimates for the solution. In addition to the normwise error bound, the code provides the maximum componentwise error bound if possible. See the comments for ERR_BNDS_NORM and ERR_BNDS_COMP for details of the error bounds. Like, she brought TWO sets of receipts, not just one 👑👑

The original system of linear equations may have been equilibrated before calling this routine, as described by arguments EQUED and S below. In this case, the solution and error bounds returned are for the original unequilibrated system.

Some optional parameters are bundled in the PARAMS array. These settings determine how refinement is performed, but often the defaults are acceptable. If the defaults are acceptable, users can pass NPARAMS = 0 which prevents the source code from accessing the PARAMS argument. Basically: defaults are the moment, so if you're chill with them, just say NPARAMS = 0 and skip the whole argument 😌

Parameters

UPLO — CHARACTER1*
- = 'U': Upper triangle of A is stored;
- = 'L': Lower triangle of A is stored.
EQUED — CHARACTER1*. Specifies the form of equilibration that was done to A before calling this routine. This is needed to compute the solution and error bounds correctly.
- = 'N': No equilibration
- = 'Y': Both row and column equilibration, i.e., A has been replaced by diag(S) * A * diag(S). The right hand side B has been changed accordingly.
N — INTEGER, the order of the matrix A. N >= 0.
NRHS — INTEGER, the number of right hand sides, i.e., the number of columns of the matrices B and X. NRHS >= 0.
A — COMPLEX*16 array, dimension (LDA,N). The Hermitian matrix A. If UPLO = 'U', the leading N-by-N upper triangular part of A contains the upper triangular part of the matrix A, and the strictly lower triangular part of A is not referenced. If UPLO = 'L', the leading N-by-N lower triangular part of A contains the lower triangular part of the matrix A, and the strictly upper triangular part of A is not referenced.
LDA — INTEGER, the leading dimension of the array A. LDA >= max(1,N).
AF — COMPLEX*16 array, dimension (LDAF,N). The factored form of the matrix A. AF contains the block diagonal matrix D and the multipliers used to obtain the factor U or L from the factorization A = U*D*U**H or A = L*D*L**H as computed by ZHETRF.
LDAF — INTEGER, the leading dimension of the array AF. LDAF >= max(1,N).
IPIV — INTEGER array, dimension (N). Details of the interchanges and the block structure of D as determined by ZHETRF.
S
S — S is a DOUBLE PRECISION array, dimension (N) 💖 the scale factors for A. If EQUED = 'Y', A gets multiplied on the left and right by diag(S). S is an input argument if FACT = 'F'; otherwise S is an output argument. If FACT = 'F' and EQUED = 'Y', every element of S must be positive. If S is output, each element of S is a power of the radix. If S is input, each element of S should be a power of the radix to ensure a reliable solution and error estimates. Scaling by powers of the radix does not cause rounding errors unless the result underflows or overflows. Rounding errors during scaling lead to refining with a matrix that is not equivalent to the input matrix, producing error estimates that may not be reliable. so real, bestie — powers of the radix or it gets messy 😭
B — B is COMPLEX*16 array, dimension (LDB,NRHS) 🦋 the right hand side matrix B.
LDB — LDB is INTEGER the leading dimension of the array B. LDB >= max(1,N).
X — X is COMPLEX*16 array, dimension (LDX,NRHS) ✨ on entry, the solution matrix X, as computed by ZHETRS. On exit, the improved solution matrix X. glow up 💅
LDX — LDX is INTEGER the leading dimension of the array X. LDX >= max(1,N).
RCOND — RCOND is DOUBLE PRECISION reciprocal scaled condition number. This is an estimate of the reciprocal Skeel condition number of the matrix A after equilibration (if done). If this is less than the machine precision (in particular, if it is zero), the matrix is singular to working precision. Note that the error may still be small even if this number is very small and the matrix appears ill-conditioned. don't panic if it's tiny — sometimes the answer still eats 💅
BERR — BERR is DOUBLE PRECISION array, dimension (NRHS) componentwise relative backward error. This is 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). the tiniest change that makes it perfect — no cap ✨
N_ERR_BNDS — N_ERR_BNDS is INTEGER number of error bounds to return for each right hand side and each type (normwise or componentwise). See ERR_BNDS_NORM and ERR_BNDS_COMP below.
ERR_BNDS_NORM — ERR_BNDS_NORM is DOUBLE PRECISION array, dimension (NRHS, N_ERR_BNDS) for each right-hand side, this array contains information about various error bounds and condition numbers corresponding to the normwise relative error, which is defined as follows:

Normwise relative error in the ith solution vector:
        max_j (abs(XTRUE(j,i) - X(j,i)))
       ------------------------------
             max_j abs(X(j,i))

The array is indexed by the type of error information as described below. There currently are up to three pieces of information returned. three little data gems 💎

The first index in ERR_BNDS_NORM(i,:) corresponds to the ith right-hand side.

The second index in ERR_BNDS_NORM(:,err) contains the following three fields:

err = 1 'Trust/don't trust' boolean. Trust the answer if the reciprocal condition number is less than the threshold sqrt(n) * dlamch('Epsilon'). vibe check: trust or don't 🫶
err = 2 'Guaranteed' error bound: The estimated forward error, almost certainly within a factor of 10 of the true error so long as the next entry is greater than the threshold sqrt(n) * dlamch('Epsilon'). This error bound should only be trusted if the previous boolean is true. only trust it if the boolean said trust — periodt 💀
err = 3 Reciprocal condition number: Estimated normwise reciprocal condition number. Compared with the threshold sqrt(n) * dlamch('Epsilon') to determine if the error estimate is 'guaranteed'. These reciprocal condition numbers are 1 / (norm(Z^{-1},inf) * norm(Z,inf)) for some appropriately scaled matrix Z. Let Z = S*A, where S scales each row by a power of the radix so all absolute row sums of Z are approximately 1. math but make it iconic 🧚

See Lapack Working Note 165 for further details and extra cautions. iykyk 📖

ERR_BNDS_COMP — ERR_BNDS_COMP is DOUBLE PRECISION array, dimension (NRHS, N_ERR_BNDS) for each right-hand side, this array contains information about various error bounds and condition numbers corresponding to the componentwise relative error, which is defined as follows:

Componentwise relative error in the ith solution vector:
               abs(XTRUE(j,i) - X(j,i))
        max_j ----------------------
                    abs(X(j,i))

The array is indexed by the right-hand side i (on which the componentwise relative error depends), and the type of error information as described below. There currently are up to three pieces of information returned for each right-hand side. If componentwise accuracy is not requested (PARAMS(3) = 0.0), then ERR_BNDS_COMP is not accessed. If N_ERR_BNDS < 3, then at most the first (:,N_ERR_BNDS) entries are returned. so if you don't ask for componentwise, she stays unbothered 🙈

The first index in ERR_BNDS_COMP(i,:) corresponds to the ith right-hand side.

The second index in ERR_BNDS_COMP(:,err) contains the following three fields:

err = 1 'Trust/don't trust' boolean. Trust the answer if the reciprocal condition number is less than the threshold sqrt(n) * dlamch('Epsilon'). same vibe check, componentwise edition 💅
err = 2 'Guaranteed' error bound: The estimated forward error, almost certainly within a factor of 10 of the true error so long as the next entry is greater than the threshold sqrt(n) * dlamch('Epsilon'). This error bound should only be trusted if the previous boolean is true. again — only trust if boolean says so, bestie 😭
err = 3 Reciprocal condition number: Estimated componentwise reciprocal condition number. Compared with the threshold sqrt(n) * dlamch('Epsilon') to determine if the error estimate is 'guaranteed'. These reciprocal condition numbers are 1 / (norm(Z^{-1},inf) * norm(Z,inf)) for some appropriately scaled matrix Z. Let Z = S*(A*diag(x)), where x is the solution for the current right-hand side and S scales each row of A*diag(x) by a power of the radix so all absolute row sums of Z are approximately 1. more math, more slay 🧚

See Lapack Working Note 165 for further details and extra cautions. iykyk again 📖

NPARAMS — NPARAMS is INTEGER specifies the number of parameters set in PARAMS. If <= 0, the PARAMS array is never referenced and default values are used. no params, no problem ✨
PARAMS — PARAMS is DOUBLE PRECISION array, dimension NPARAMS specifies algorithm parameters. If an entry is < 0.0, then that entry will be filled with default value used for that parameter. Only positions up to NPARAMS are accessed; defaults are used for higher-numbered parameters.
- PARAMS(LA_LINRX_ITREF_I = 1) : Whether to perform iterative refinement or not. Default: 1.0D+0
 - = 0.0: No refinement is performed, and no error bounds are computed. just raw dogging it, no refinement 😳
 - = 1.0: Use the double-precision refinement algorithm, possibly with doubled-single computations if the compilation environment does not support DOUBLE PRECISION. glow up with extra precision ✨
 - (other values are reserved for future use) mystery future vibes 👀
- PARAMS(LA_LINRX_ITHRESH_I = 2) : Maximum number of residual computations allowed for refinement. Default: 10 Aggressive: Set to 100 to permit convergence using approximate factorizations or factorizations other than LU. If the factorization uses a technique other than Gaussian elimination, the guarantees in err_bnds_norm and err_bnds_comp may no longer be trustworthy. so if you go aggressive, those guarantees might ghost you 💀
- PARAMS(LA_LINRX_CWISE_I = 3) : Flag determining if the code will attempt to find a solution with small componentwise relative error in the double-precision algorithm. Positive is true, 0.0 is false. Default: 1.0 (attempt componentwise convergence) yes queen, go componentwise 💅
WORK — WORK is COMPLEX*16 array, dimension (2*N)
RWORK — RWORK is DOUBLE PRECISION array, dimension (2*N)
INFO — INFO is INTEGER
- = 0: Successful exit. The solution to every right-hand side is guaranteed. slay, all guaranteed 🎀
- < 0: If INFO = -i, the i-th argument had an illegal value. whoops, bad arg 😭
- > 0 and <= N: U(INFO,INFO) is exactly zero. The factorization has been completed, but the factor U is exactly singular, so the solution and error bounds could not be computed. RCOND = 0 is returned. singular matrix, no solution for you 💀
- = N+J: The solution corresponding to the Jth right-hand side is not guaranteed. The solutions corresponding to other right-hand sides K with K > J may not be guaranteed as well, but only the first such right-hand side is reported. If a small componentwise error is not requested (PARAMS(3) = 0.0) then the Jth right-hand side is the first with a normwise error bound that is not guaranteed (the smallest J such that ERR_BNDS_NORM(J,1) = 0.0). By default (PARAMS(3) = 1.0) the Jth right-hand side is the first with either a normwise or componentwise error bound that is not guaranteed (the smallest J such that either ERR_BNDS_NORM(J,1) = 0.0 or ERR_BNDS_COMP(J,1) = 0.0). See the definition of ERR_BNDS_NORM(:,1) and ERR_BNDS_COMP(:,1). To get information about all of the right-hand sides check ERR_BNDS_NORM or ERR_BNDS_COMP. so basically the first messy right-hand side gets called out, and you gotta check the arrays yourself for the rest 👑

Author 👑

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

subroutine `zhetd2` 🦋✨

subroutine zhetd2 (character UPLO, integer N, complex*16, dimension( lda, * ) A, integer LDA, double precision, dimension( * ) D, double precision, dimension( * ) E, complex*16, dimension( * ) TAU, integer INFO)

ZHETD2 reduces a Hermitian matrix to real symmetric tridiagonal form by an unitary similarity transformation (unblocked algorithm). no cap, it's the glow-up pipeline for Hermitian matrices 💅

Purpose 🌸

ZHETD2 reduces a complex Hermitian matrix A to real symmetric tridiagonal form T by a unitary similarity transformation:

Q**H * A * Q = T.

Parameters 🎀

UPLO — CHARACTER*1. Specifies whether the upper or lower triangular part of the Hermitian matrix A is stored:
- = 'U': Upper triangular
- = 'L': Lower triangular
N — INTEGER. The order of the matrix A. N >= 0.
A — COMPLEX*16 array, dimension (LDA,N).

On entry, the Hermitian matrix A. If UPLO = 'U', the leading n-by-n upper triangular part of A contains the upper triangular part of the matrix A, and the strictly lower triangular part of A is not referenced. If UPLO = 'L', the leading n-by-n lower triangular part of A contains the lower triangular part of the matrix A, and the strictly upper triangular part of A is not referenced.

On exit, if UPLO = 'U', the diagonal and first superdiagonal of A are overwritten by the corresponding elements of the tridiagonal matrix T, and the elements above the first superdiagonal, with the array TAU, represent the unitary matrix Q as a product of elementary reflectors; if UPLO = 'L', the diagonal and first subdiagonal of A are over-written by the corresponding elements of the tridiagonal matrix T, and the elements below the first subdiagonal, with the array TAU, represent the unitary matrix Q as a product of elementary reflectors. See Further Details.

LDA — INTEGER. The leading dimension of the array A. LDA >= max(1,N).
D — DOUBLE PRECISION array, dimension (N). The diagonal elements of the tridiagonal matrix T: D(i) = A(i,i).
E — DOUBLE PRECISION array, dimension (N-1). The off-diagonal elements of the tridiagonal matrix T: E(i) = A(i,i+1) if UPLO = 'U', E(i) = A(i+1,i) if UPLO = 'L'.
TAU — COMPLEX*16 array, dimension (N-1). The scalar factors of the elementary reflectors (see Further Details).
INFO — INTEGER.
- = 0: successful exit
- < 0: if INFO = -i, the i-th argument had an illegal value.

Author 👑

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

Further Details 🧚

If UPLO = 'U', the matrix Q is represented as a product of elementary reflectors

Q = H(n-1) . . . H(2) H(1).

Each H(i) has the form

H(i) = I - tau * v * v**H

where tau is a complex scalar, and v is a complex vector with v(i+1:n) = 0 and v(i) = 1; v(1:i-1) is stored on exit in A(1:i-1,i+1), and tau in TAU(i).

If UPLO = 'L', the matrix Q is represented as a product of elementary reflectors

Q = H(1) H(2) . . . H(n-1).

Each H(i) has the form

H(i) = I - tau * v * v**H

where tau is a complex scalar, and v is a complex vector with v(1:i) = 0 and v(i+1) = 1; v(i+2:n) is stored on exit in A(i+2:n,i), and tau in TAU(i).

The contents of A on exit are illustrated by the following examples with n = 5:

if UPLO = 'U':                       if UPLO = 'L':

  (  d   e   v2  v3  v4 )              (  d                  )
  (      d   e   v3  v4 )              (  e   d              )
  (          d   e   v4 )              (  v1  e   d          )
  (              d   e  )              (  v1  v2  e   d      )
  (                  d  )              (  v1  v2  v3  e   d  )

where d and e denote diagonal and off-diagonal elements of T, and vi denotes an element of the vector defining H(i). obsessed with this little map honestly 🗺️

subroutine `zhetf2` 💖

subroutine zhetf2 (character UPLO, integer N, complex*16, dimension( lda, * ) A, integer LDA, integer, dimension( * ) IPIV, integer INFO)

ZHETF2 computes the factorization of a complex Hermitian matrix, using the diagonal pivoting method (unblocked algorithm, calling Level 2 BLAS).

Purpose ✨

ZHETF2 computes the factorization of a complex Hermitian matrix A using the Bunch-Kaufman diagonal pivoting method:

A = U*D*U**H  or  A = L*D*L**H

where U (or L) is a product of permutation and unit upper (lower) triangular matrices, U**H is the conjugate transpose of U, and D is Hermitian and block diagonal with 1-by-1 and 2-by-2 diagonal blocks.

This is the unblocked version of the algorithm, calling Level 2 BLAS.

Parameters 🎀

UPLO — CHARACTER*1. Specifies whether the upper or lower triangular part of the Hermitian matrix A is stored:
- = 'U': Upper triangular
- = 'L': Lower triangular
N — INTEGER. The order of the matrix A. N >= 0.
A — COMPLEX*16 array, dimension (LDA,N).

On exit, the block diagonal matrix D and the multipliers used to obtain the factor U or L (see below for further details).

LDA — INTEGER. The leading dimension of the array A. LDA >= max(1,N).
IPIV — INTEGER array, dimension (N). Details of the interchanges and the block structure of D.

If UPLO = 'U':

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 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':

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 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.

INFO — INTEGER.
- = 0: successful exit
- < 0: if INFO = -k, the k-th argument had an illegal value
- > 0: if INFO = k, D(k,k) is exactly zero. The factorization has been completed, but the block diagonal matrix D is exactly singular, and division by zero will occur if it is used to solve a system of equations.

⚠️ bestie that > 0 case is not a drill — D is exactly singular and if you try to solve with it you WILL divide by zero. do not pass go. 😭

Author 👑

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

Further Details 🧚

If UPLO = 'U', then A = U*D*U**H, where

U = P(n)*U(n)* ... *P(k)U(k)* ...,

i.e., U is a product of terms P(k)*U(k), where k decreases from n to 1 in steps of 1 or 2, and D is a block diagonal matrix with 1-by-1 and 2-by-2 diagonal blocks D(k). P(k) is a permutation matrix as defined by IPIV(k), and U(k) is a unit upper triangular matrix, such that if the diagonal block D(k) is of order s (s = 1 or 2), then

           (   I    v    0   )   k-s
    U(k) = (   0    I    0   )   s
           (   0    0    I   )   n-k
              k-s   s   n-k

If s = 1, D(k) overwrites A(k,k), and v overwrites A(1:k-1,k). If s = 2, the upper triangle of D(k) overwrites A(k-1,k-1), A(k-1,k), and A(k,k), and v overwrites A(1:k-2,k-1:k).

If UPLO = 'L', then A = L*D*L**H, where

L = P(1)*L(1)* ... *P(k)*L(k)* ...,

i.e., L is a product of terms P(k)*L(k), where k increases from 1 to n in steps of 1 or 2, and D is a block diagonal matrix with 1-by-1 and 2-by-2 diagonal blocks D(k). P(k) is a permutation matrix as defined by IPIV(k), and L(k) is a unit lower triangular matrix, such that if the diagonal block D(k) is of order s (s = 1 or 2), then

           (   I    0     0   )  k-1
    L(k) = (   0    I     0   )  s
           (   0    v     I   )  n-k-s+1
              k-1   s  n-k-s+1

If s = 1, D(k) overwrites A(k,k), and v overwrites A(k+1:n,k). If s = 2, the lower triangle of D(k) overwrites A(k,k), A(k+1,k), and A(k+1,k+1), and v overwrites A(k+2:n,k:k+1).

Contributors: 🫶

09-29-06 — patch from Bobby Cheng, MathWorks

Replace l.210 and l.393

     IF( MAX( ABSAKK, COLMAX ).EQ.ZERO ) THEN

     IF( (MAX( ABSAKK, COLMAX ).EQ.ZERO) .OR. DISNAN(ABSAKK) ) THEN

01-01-96 — Based on modifications by J. Lewis, Boeing Computer Services Company, A. Petitet, Computer Science Dept., Univ. of Tenn., Knoxville, USA

subroutine `zhetf2_rk` 🌈

subroutine zhetf2_rk (character UPLO, integer N, complex*16, dimension( lda, * ) A, integer LDA, complex*16, dimension( * ) E, integer, dimension( * ) IPIV, integer INFO)

ZHETF2_RK computes the factorization of a complex Hermitian indefinite matrix using the bounded Bunch-Kaufman (rook) diagonal pivoting method (BLAS2 unblocked algorithm). the rook variant!! ♜✨

Purpose 💜

ZHETF2_RK computes the factorization of a complex Hermitian matrix A using the bounded Bunch-Kaufman (rook) diagonal pivoting method:

A = P*U*D*(U**H)*(P**T) or A = P*L*D*(L**H)*(P**T),

This is the unblocked version of the algorithm, calling Level 2 BLAS. For more info, go peek at the Further Details section, bestie. 💅

✨ Parameters ✨

UPLO — UPLO is CHARACTER*1. It tells us which triangular half of the Hermitian matrix A is actually chillin' in memory:
- = 'U': Upper triangular
- = 'L': Lower triangular
N — N is INTEGER. The order of the matrix A. So like, N >= 0, obviously. 🎀
A — A is COMPLEX*16 array, dimension (LDA,N). On entry: it's the Hermitian matrix A, no cap.
- If UPLO = 'U': the leading N-by-N upper triangular part of A holds the upper triangular part of the matrix A, and the strictly lower triangular part of A? not referenced, we literally don't look at it. 🙈
- If UPLO = 'L': the leading N-by-N lower triangular part of A holds the lower triangular part of the matrix A, and the strictly upper triangular part of A is not referenced either.

On exit, A contains:

a) ONLY the diagonal elements of the Hermitian block diagonal matrix D on the diagonal of A, i.e. D(k,k) = A(k,k). (The superdiagonal — or subdiagonal — elements of D get stashed on exit in array E instead, iykyk.)
b) If UPLO = 'U': the factor U in the superdiagonal part of A. If UPLO = 'L': the factor L in the subdiagonal part of A.

LDA — LDA is INTEGER. The leading dimension of the array A. Gotta be LDA >= max(1,N), periodt.
E — E is COMPLEX*16 array, dimension (N). On exit, it holds the superdiagonal (or subdiagonal) elements of the Hermitian block diagonal matrix D with 1-by-1 or 2-by-2 diagonal blocks, where:
- If UPLO = 'U': E(i) = D(i-1,i), i=2:N, and E(1) is set to 0;
- If UPLO = 'L': E(i) = D(i+1,i), i=1:N-1, and E(N) is set to 0.

NOTE: for a 1-by-1 diagonal block D(k), where 1 <= k <= N, the element E(k) is set to 0 in both the UPLO = 'U' and UPLO = 'L' cases. 💖

IPIV — IPIV is INTEGER array, dimension (N). IPIV describes the permutation matrix P in the factorization of matrix A, like this: the absolute value of IPIV(k) represents the index of the row and column that got interchanged with the k-th row and column. The value of UPLO describes the order in which those interchanges were applied. Also — and this is the iconic part — the sign of IPIV represents the block structure of the Hermitian block diagonal matrix D with 1-by-1 or 2-by-2 diagonal blocks, which correspond to 1 or 2 interchanges at each factorization step. For more info see the Further Details section. 🦋

If UPLO = 'U' (in factorization order, k decreases from N to 1):

a) A single positive entry IPIV(k) > 0 means: D(k,k) is a 1-by-1 diagonal block.
- If IPIV(k) != k, rows and columns k and IPIV(k) were interchanged in the matrix A(1:N,1:N);
- If IPIV(k) = k, no interchange occurred.
b) A pair of consecutive negative entries IPIV(k) < 0 and IPIV(k-1) < 0 means: D(k-1:k,k-1:k) is a 2-by-2 diagonal block. (NOTE: negative entries in IPIV appear ONLY in pairs.)
1. If -IPIV(k) != k, rows and columns k and -IPIV(k) were interchanged in the matrix A(1:N,1:N). If -IPIV(k) = k, no interchange occurred.
2. If -IPIV(k-1) != k-1, rows and columns k-1 and -IPIV(k-1) were interchanged in the matrix A(1:N,1:N). If -IPIV(k-1) = k-1, no interchange occurred.
c) In both cases a) and b), always ABS( IPIV(k) ) <= k.
d) NOTE: any entry IPIV(k) is always NONZERO on output. no cap.

If UPLO = 'L' (in factorization order, k increases from 1 to N):

a) A single positive entry IPIV(k) > 0 means: D(k,k) is a 1-by-1 diagonal block.
- If IPIV(k) != k, rows and columns k and IPIV(k) were interchanged in the matrix A(1:N,1:N).
- If IPIV(k) = k, no interchange occurred.
b) A pair of consecutive negative entries IPIV(k) < 0 and IPIV(k+1) < 0 means: D(k:k+1,k:k+1) is a 2-by-2 diagonal block. (NOTE: negative entries in IPIV appear ONLY in pairs.)
1. If -IPIV(k) != k, rows and columns k and -IPIV(k) were interchanged in the matrix A(1:N,1:N). If -IPIV(k) = k, no interchange occurred.
2. If -IPIV(k+1) != k+1, rows and columns k-1 and -IPIV(k-1) were interchanged in the matrix A(1:N,1:N). If -IPIV(k+1) = k+1, no interchange occurred.
c) In both cases a) and b), always ABS( IPIV(k) ) >= k.
d) NOTE: any entry IPIV(k) is always NONZERO on output.

INFO — INFO is INTEGER.
- = 0: successful exit, we ate. 💅
- < 0: if INFO = -k, the k-th argument had an illegal value.
- > 0: if INFO = k, the matrix A is singular, because:
 - If UPLO = 'U': column k in the upper triangular part of A contains all zeros.
 - If UPLO = 'L': column k in the lower triangular part of A contains all zeros.

Therefore D(k,k) is exactly zero, and the superdiagonal elements of column k of U (or the subdiagonal elements of column k of L) are all zeros. The factorization has been completed, but the block diagonal matrix D is exactly singular, and division by zero will occur if it's used to solve a system of equations. 😭

NOTE: INFO only stores the first occurrence of a singularity — any subsequent occurrence of singularity is not stored in INFO, even though the factorization always completes. lowkey important.

👑 Author 👑

Univ. of Tennessee

Univ. of California Berkeley

Univ. of Colorado Denver

NAG Ltd.

🔍 Further Details: 🔍

TODO: put further details

🫶 Contributors: 🫶

December 2016, Igor Kozachenko, Computer Science Division, University of California, Berkeley

September 2007, Sven Hammarling, Nicholas J. Higham, Craig Lucas, School of Mathematics, University of Manchester

01-01-96 - Based on modifications by J. Lewis, Boeing Computer Services Company A. Petitet, Computer Science Dept., Univ. of Tenn., Knoxville abd , USA

🌈 subroutine zhetf2_rook 🌈

subroutine zhetf2_rook (character UPLO, integer N, complex*16, dimension(
    lda, * ) A, integer LDA, integer, dimension( * ) IPIV, integer INFO)

ZHETF2_ROOK computes the factorization of a complex Hermitian indefinite matrix using the bounded Bunch-Kaufman ('rook') diagonal pivoting method (unblocked algorithm).

✨ Purpose: ✨

ZHETF2_ROOK computes the factorization of a complex Hermitian matrix A using the bounded Bunch-Kaufman ('rook') diagonal pivoting method:

A = U*D*U**H  or  A = L*D*L**H

This is the unblocked version of the algorithm, calling Level 2 BLAS. 💅

✨ Parameters ✨

UPLO — UPLO is CHARACTER*1. Which triangular half of the Hermitian matrix A is stored, bestie:
- = 'U': Upper triangular
- = 'L': Lower triangular
N — N is INTEGER. The order of the matrix A. N >= 0.
A — A is COMPLEX*16 array, dimension (LDA,N). On entry, the Hermitian matrix A. If UPLO = 'U', the leading n-by-n upper triangular part of A contains the upper triangular part of the matrix A, and the strictly lower triangular part of A is not referenced. If UPLO = 'L', the leading n-by-n lower triangular part of A contains the lower triangular part of the matrix A, and the strictly upper triangular part of A is not referenced. 🙈

On exit, the block diagonal matrix D and the multipliers used to obtain the factor U or L (see below for further details).

LDA — LDA is INTEGER. The leading dimension of the array A. LDA >= max(1,N).
IPIV — IPIV is INTEGER array, dimension (N). All the deets on the interchanges and the block structure of D:
- If UPLO = 'U':
 - 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 IPIV(k) < 0 and IPIV(k-1) < 0, then rows and columns k and -IPIV(k) were interchanged and rows and columns k-1 and -IPIV(k-1) were inerchaged, D(k-1:k,k-1:k) is a 2-by-2 diagonal block.
- If UPLO = 'L':
 - 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 IPIV(k) < 0 and IPIV(k+1) < 0, then rows and columns k and -IPIV(k) were interchanged and rows and columns k+1 and -IPIV(k+1) were inerchaged, D(k:k+1,k:k+1) is a 2-by-2 diagonal block.
INFO — hold up bestie, INFO gets its whole glow-up in the next part 💅
INFO is INTEGER
- = 0: successful exit ✨
- < 0: if INFO = -k, the k-th argument had an illegal value 😬
- > 0: if INFO = k, D(k,k) is exactly zero. The factorization has been completed, but the block diagonal matrix D is exactly singular, and division by zero will occur if it is used to solve a system of equations. girl run, that's a singular block, do NOT divide by it 💀

👑 Author

Univ. of Tennessee

Univ. of California Berkeley

Univ. of Colorado Denver

NAG Ltd.

🧚 Further Details

If UPLO = 'U', then A = U*D*U**H, where

U = P(n)*U(n)* ... *P(k)U(k)* ...,

         (   I    v    0   )   k-s
U(k) =   (   0    I    0   )   s
         (   0    0    I   )   n-k
            k-s   s   n-k

If s = 1, D(k) overwrites A(k,k), and v overwrites A(1:k-1,k).

If s = 2, the upper triangle of D(k) overwrites A(k-1,k-1), A(k-1,k), and A(k,k), and v overwrites A(1:k-2,k-1:k).

If UPLO = 'L', then A = L*D*L**H, where

L = P(1)*L(1)* ... *P(k)*L(k)* ...,

         (   I    0     0   )  k-1
L(k) =   (   0    I     0   )  s
         (   0    v     I   )  n-k-s+1
            k-1   s  n-k-s+1

If s = 1, D(k) overwrites A(k,k), and v overwrites A(k+1:n,k).

If s = 2, the lower triangle of D(k) overwrites A(k,k), A(k+1,k), and A(k+1,k+1), and v overwrites A(k+2:n,k:k+1).

🌸 Contributors

November 2013, Igor Kozachenko, Computer Science Division, University of California, Berkeley

September 2007, Sven Hammarling, Nicholas J. Higham, Craig Lucas, School of Mathematics, University of Manchester

01-01-96 - Based on modifications by J. Lewis, Boeing Computer Services Company A. Petitet, Computer Science Dept., Univ. of Tenn., Knoxville, USA

🦋 subroutine zhetrd (character UPLO, integer N, complex16, dimension( lda, ) A, integer LDA, double precision, dimension( * ) D, double precision, dimension( * ) E, complex16, dimension( ) TAU, complex16, dimension( ) WORK, integer LWORK, integer INFO)

ZHETRD

✨ Purpose

ZHETRD reduces a complex Hermitian matrix A to real symmetric tridiagonal form T by a unitary similarity transformation:

Q**H * A * Q = T.

🎀 Parameters

UPLO

UPLO is CHARACTER*1

= 'U': Upper triangle of A is stored;
= 'L': Lower triangle of A is stored.

N

N is INTEGER The order of the matrix A. N >= 0.

A

A is COMPLEX*16 array, dimension (LDA,N) On entry, the Hermitian matrix A. If UPLO = 'U', the leading N-by-N upper triangular part of A contains the upper triangular part of the matrix A, and the strictly lower triangular part of A is not referenced. If UPLO = 'L', the leading N-by-N lower triangular part of A contains the lower triangular part of the matrix A, and the strictly upper triangular part of A is not referenced. On exit, if UPLO = 'U', the diagonal and first superdiagonal of A are overwritten by the corresponding elements of the tridiagonal matrix T, and the elements above the first superdiagonal, with the array TAU, represent the unitary matrix Q as a product of elementary reflectors; if UPLO = 'L', the diagonal and first subdiagonal of A are over- written by the corresponding elements of the tridiagonal matrix T, and the elements below the first subdiagonal, with the array TAU, represent the unitary matrix Q as a product of elementary reflectors. See Further Details. 👀

LDA

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

D

D is DOUBLE PRECISION array, dimension (N) The diagonal elements of the tridiagonal matrix T: D(i) = A(i,i).

E

E is DOUBLE PRECISION array, dimension (N-1) The off-diagonal elements of the tridiagonal matrix T: E(i) = A(i,i+1) if UPLO = 'U', E(i) = A(i+1,i) if UPLO = 'L'.

TAU

TAU is COMPLEX*16 array, dimension (N-1) The scalar factors of the elementary reflectors (see Further Details).

WORK

WORK is COMPLEX*16 array, dimension (MAX(1,LWORK)) On exit, if INFO = 0, WORK(1) returns the optimal LWORK. 💅

LWORK

LWORK is INTEGER The dimension of the array WORK. LWORK >= 1. For optimum performance LWORK >= N*NB, where NB is the optimal blocksize.

If LWORK = -1, then a workspace query is assumed; the routine only calculates the optimal size of the WORK array, returns this value as the first entry of the WORK array, and no error message related to LWORK is issued by XERBLA.

INFO

INFO is INTEGER

= 0: successful exit ✨
< 0: if INFO = -i, the i-th argument had an illegal value 😭

👑 Author

Univ. of Tennessee

Univ. of California Berkeley

Univ. of Colorado Denver

NAG Ltd.

🧚 Further Details

If UPLO = 'U', the matrix Q is represented as a product of elementary reflectors

Q = H(n-1) . . . H(2) H(1).

Each H(i) has the form

H(i) = I - tau * v * v**H

where tau is a complex scalar, and v is a complex vector with v(i+1:n) = 0 and v(i) = 1; v(1:i-1) is stored on exit in A(1:i-1,i+1), and tau in TAU(i).

If UPLO = 'L', the matrix Q is represented as a product of elementary reflectors

Q = H(1) H(2) . . . H(n-1).

Each H(i) has the form

H(i) = I - tau * v * v**H

where tau is a complex scalar, and v is a complex vector with v(1:i) = 0 and v(i+1) = 1; v(i+2:n) is stored on exit in A(i+2:n,i), and tau in TAU(i).

The contents of A on exit are illustrated by the following examples with n = 5:

if UPLO = 'U':                       if UPLO = 'L':

  (  d   e   v2  v3  v4 )              (  d                  )
  (      d   e   v3  v4 )              (  e   d              )
  (          d   e   v4 )              (  v1  e   d          )
  (              d   e  )              (  v1  v2  e   d      )
  (                  d  )              (  v1  v2  v3  e   d  )

where d and e denote diagonal and off-diagonal elements of T, and vi denotes an element of the vector defining H(i).

🌈 subroutine zhetrd_2stage (character VECT, character UPLO, integer N, complex16, dimension( lda, ) A, integer LDA, double precision, dimension( * ) D, double precision, dimension( * ) E, complex16, dimension( ) TAU, complex16, dimension( ) HOUS2, integer LHOUS2, complex16, dimension( ) WORK, integer LWORK, integer INFO)

ZHETRD_2STAGE

✨ Purpose

ZHETRD_2STAGE reduces a complex Hermitian matrix A to real symmetric tridiagonal form T by a unitary similarity transformation:

Q1**H Q2**H* A * Q2 * Q1 = T.

🎀 Parameters

VECT

VECT is CHARACTER*1

= 'N': No need for the Housholder representation, in particular for the second stage (Band to tridiagonal) and thus LHOUS2 is of size max(1, 4*N);
= 'V': the Householder representation is needed to either generate Q1 Q2 or to apply Q1 Q2, then LHOUS2 is to be queried and computed. (NOT AVAILABLE IN THIS RELEASE). 🫠

UPLO

UPLO is CHARACTER*1

= 'U': Upper triangle of A is stored;
= 'L': Lower triangle of A is stored.

N

N is INTEGER The order of the matrix A. N >= 0.

A

A is COMPLEX*16 array, dimension (LDA,N) On entry, the Hermitian matrix A. If UPLO = 'U', the leading N-by-N upper triangular part of A contains the upper triangular part of the matrix A, and the strictly lower triangular part of A is not referenced. If UPLO = 'L', the leading N-by-N lower triangular part of A contains the lower triangular part of the matrix A, and the strictly upper triangular part of A is not referenced. On exit, if UPLO = 'U', the band superdiagonal of A are overwritten by the corresponding elements of the internal band-diagonal matrix AB, and the elements above the KD superdiagonal, with the array TAU, represent the unitary matrix Q1 as a product of elementary reflectors; if UPLO = 'L', the diagonal and band subdiagonal of A are over- written by the corresponding elements of the internal band-diagonal matrix AB, and the elements below the KD subdiagonal, with the array TAU, represent the unitary matrix Q1 as a product of elementary reflectors. See Further Details. 👀

LDA

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

D

D is DOUBLE PRECISION array, dimension (N) The diagonal elements of the tridiagonal matrix T.

E

E is DOUBLE PRECISION array, dimension (N-1) The off-diagonal elements of the tridiagonal matrix T.

TAU
TAU — COMPLEX*16 array, dimension (N-KD). It's holding the scalar factors of the elementary reflectors from the first stage (go peek at Further Details, bestie 💅).
HOUS2 — COMPLEX*16 array, dimension (LHOUS2). This one stores the Householder representation of the stage2 band-to-tridiagonal glow-up ✨.
LHOUS2 — INTEGER. The dimension of the array HOUS2, no biggie.
- If LWORK = -1, or LHOUS2 = -1, then a query is assumed — the routine only calculates the optimal size of the HOUS2 array, returns this value as the first entry of the HOUS2 array, and no error message related to LHOUS2 is issued by XERBLA. Periodt.
- If VECT='N', LHOUS2 = max(1, 4*n);
- if VECT='V', option not yet available. 😭 (it's giving "coming soon", not us though)
WORK — COMPLEX*16 array, dimension (LWORK).
LWORK — INTEGER. The dimension of the array WORK. LWORK = MAX(1, dimension).
- If LWORK = -1, or LHOUS2=-1, then a workspace query is assumed — the routine only calculates the optimal size of the WORK array, returns this value as the first entry of the WORK array, and no error message related to LWORK is issued by XERBLA. 🫶
- LWORK = MAX(1, dimension) where

dimension   = max(stage1,stage2) + (KD+1)*N
            = N*KD + N*max(KD+1,FACTOPTNB)
              + max(2*KD*KD, KD*NTHREADS)
              + (KD+1)*N

where KD is the blocking size of the reduction, FACTOPTNB is the blocking used by the QR or LQ algorithm, usually FACTOPTNB=128 is a good choice, and NTHREADS is the number of threads used when openMP compilation is enabled, otherwise =1. Slay 💅

INFO — INTEGER.
- = 0: successful exit 🎉
- < 0: if INFO = -i, the i-th argument had an illegal value 💀

👑 Author

Univ. of Tennessee

Univ. of California Berkeley

Univ. of Colorado Denver

NAG Ltd.

🌈 Further Details:

Implemented by Azzam Haidar.

All details are available on technical report, SC11, SC13 papers. 📚

Azzam Haidar, Hatem Ltaief, and Jack Dongarra. Parallel reduction to condensed forms for symmetric eigenvalue problems using aggregated fine-grained and memory-aware kernels. In Proceedings of 2011 International Conference for High Performance Computing, Networking, Storage and Analysis (SC '11), New York, NY, USA, Article 8 , 11 pages. http://doi.acm.org/10.1145/2063384.2063394

A. Haidar, J. Kurzak, P. Luszczek, 2013. An improved parallel singular value algorithm and its implementation for multicore hardware, In Proceedings of 2013 International Conference for High Performance Computing, Networking, Storage and Analysis (SC '13). Denver, Colorado, USA, 2013. Article 90, 12 pages. http://doi.acm.org/10.1145/2503210.2503292

A. Haidar, R. Solca, S. Tomov, T. Schulthess and J. Dongarra. A novel hybrid CPU-GPU generalized eigensolver for electronic structure calculations based on fine-grained memory aware tasks. International Journal of High Performance Computing Applications. Volume 28 Issue 2, Pages 196-209, May 2014. http://hpc.sagepub.com/content/28/2/196

✨ `zhetrd_he2hb` — the next main character 🌟

subroutine zhetrd_he2hb (character UPLO, integer N, integer KD, complex*16,
    dimension( lda, * ) A, integer LDA, complex*16, dimension( ldab, * ) AB,
    integer LDAB, complex*16, dimension( * ) TAU, complex*16, dimension( * )
    WORK, integer LWORK, integer INFO)

ZHETRD_HE2HB 🦋

💖 Purpose:

ZHETRD_HE2HB reduces a complex Hermitian matrix A to complex Hermitian band-diagonal form AB by a unitary similarity transformation:

Q**H * A * Q = AB.

okay so it literally takes your Hermitian matrix and gives it the band-diagonal glow-up, no cap ✨

🎀 Parameters

UPLO — UPLO is CHARACTER*1.
- = 'U': Upper triangle of A is stored;
- = 'L': Lower triangle of A is stored.
N — N is INTEGER. The order of the matrix A. N >= 0.
KD — KD is INTEGER. The number of superdiagonals of the reduced matrix if UPLO = 'U', or the number of subdiagonals if UPLO = 'L'. KD >= 0. The reduced matrix is stored in the array AB.
A — A is COMPLEX*16 array, dimension (LDA,N).
- On entry, the Hermitian matrix A. If UPLO = 'U', the leading N-by-N upper triangular part of A contains the upper triangular part of the matrix A, and the strictly lower triangular part of A is not referenced. If UPLO = 'L', the leading N-by-N lower triangular part of A contains the lower triangular part of the matrix A, and the strictly upper triangular part of A is not referenced.
- On exit, if UPLO = 'U', the diagonal and first superdiagonal of A are overwritten by the corresponding elements of the tridiagonal matrix T, and the elements above the first superdiagonal, with the array TAU, represent the unitary matrix Q as a product of elementary reflectors; if UPLO = 'L', the diagonal and first subdiagonal of A are over-written by the corresponding elements of the tridiagonal matrix T, and the elements below the first subdiagonal, with the array TAU, represent the unitary matrix Q as a product of elementary reflectors. See Further Details. 👀
LDA — LDA is INTEGER. The leading dimension of the array A. LDA >= max(1,N).
AB — AB is COMPLEX*16 array, dimension (LDAB,N).
- On exit, the upper or lower triangle of the Hermitian band matrix A, stored in the first KD+1 rows of the array. The j-th column of A is stored in the j-th column of the array AB as follows:
- if UPLO = 'U', AB(kd+1+i-j,j) = A(i,j) for max(1,j-kd)<=i<=j;
- if UPLO = 'L', AB(1+i-j,j) = A(i,j) for j<=i<=min(n,j+kd).
LDAB — LDAB is INTEGER. The leading dimension of the array AB. LDAB >= KD+1.
TAU — TAU is COMPLEX*16 array, dimension (N-KD). The scalar factors of the elementary reflectors (see Further Details). ✨
WORK — WORK is COMPLEX*16 array, dimension (LWORK). On exit, if INFO = 0, or if LWORK=-1, WORK(1) returns the size of LWORK.
LWORK — LWORK is INTEGER. The dimension of the array WORK which should be calculated by a workspace query. LWORK = MAX(1, LWORK_QUERY).
- If LWORK = -1, then a workspace query is assumed; the routine only calculates the optimal size of the WORK array, returns this value as the first entry of the WORK array, and no error message related to LWORK is issued by XERBLA.
- LWORK_QUERY = N*KD + N*max(KD,FACTOPTNB) + 2*KD*KD
- where FACTOPTNB is the blocking used by the QR or LQ algorithm, usually FACTOPTNB=128 is a good choice otherwise putting LWORK=-1 will provide the size of WORK. 💅
INFO — INFO is INTEGER.
- = 0: successful exit 🎉
- < 0: if INFO = -i, the i-th argument had an illegal value 💀

👑 Author

Univ. of Tennessee

Univ. of California Berkeley

Univ. of Colorado Denver

NAG Ltd.

🌈 Further Details:

Implemented by Azzam Haidar.

All details are available on technical report, SC11, SC13 papers. 📚

If UPLO = 'U', the matrix Q is represented as a product of elementary reflectors

Q = H(k)**H . . . H(2)**H H(1)**H, where k = n-kd.

Each H(i) has the form

H(i) = I - tau * v * v**H

where tau is a complex scalar, and v is a complex vector with v(1:i+kd-1) = 0 and v(i+kd) = 1; conjg(v(i+kd+1:n)) is stored on exit in A(i,i+kd+1:n), and tau in TAU(i).

If UPLO = 'L', the matrix Q is represented as a product of elementary reflectors

Q = H(1) H(2) . . . H(k), where k = n-kd.

Each H(i) has the form

H(i) = I - tau * v * v**H

where tau is a complex scalar, and v is a complex vector with v(kd+1:i) = 0 and v(i+kd+1) = 1; v(i+kd+2:n) is stored on exit in A(i+kd+2:n,i), and tau in TAU(i).

The contents of A on exit are illustrated by the following examples with n = 5:

if UPLO = 'U':                       if UPLO = 'L':

  (  ab  ab/v1  v1      v1     v1    )              (  ab                            )
  (      ab     ab/v2   v2     v2    )              (  ab/v1  ab                     )
  (             ab      ab/v3  v3    )              (  v1     ab/v2  ab              )
  (                     ab     ab/v4 )              (  v1     v2     ab/v3  ab       )
  (                            ab    )              (  v1     v2     v3     ab/v4 ab )

where d and e denote diagonal and off-diagonal elements of T, and vi denotes an element of the vector defining H(i)..fi

🦋 ZHETRF 🦋

subroutine zhetrf (character UPLO, integer N, complex*16, dimension( lda, * ) A, integer LDA, integer, dimension( * ) IPIV, complex*16, dimension( * ) WORK, integer LWORK, integer INFO)

Purpose:

ZHETRF computes the factorization of a complex Hermitian matrix A using the Bunch-Kaufman diagonal pivoting method. The form of the factorization is

   A = U*D*U**H  or  A = L*D*L**H

where U (or L) is a product of permutation and unit upper (lower) triangular matrices, and D is Hermitian and block diagonal with 1-by-1 and 2-by-2 diagonal blocks. 💅

This is the blocked version of the algorithm, calling Level 3 BLAS. slay.

Parameters

UPLO — UPLO is CHARACTER*1
- = 'U': Upper triangle of A is stored;
- = 'L': Lower triangle of A is stored.
N — N is INTEGER. The order of the matrix A. N >= 0.
A — A is COMPLEX*16 array, dimension (LDA,N). On entry, the Hermitian matrix A. If UPLO = 'U', the leading N-by-N upper triangular part of A contains the upper triangular part of the matrix A, and the strictly lower triangular part of A is not referenced. If UPLO = 'L', the leading N-by-N lower triangular part of A contains the lower triangular part of the matrix A, and the strictly upper triangular part of A is not referenced.

On exit, the block diagonal matrix D and the multipliers used to obtain the factor U or L (see below for further details). 👑

LDA — LDA is INTEGER. The leading dimension of the array A. LDA >= max(1,N).
IPIV — IPIV is INTEGER array, dimension (N). Details of the interchanges and the block structure of D. 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. 🫶
WORK — WORK is COMPLEX*16 array, dimension (MAX(1,LWORK)). On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
LWORK — LWORK is INTEGER. The length of WORK. LWORK >=1. For best performance LWORK >= N*NB, where NB is the block size returned by ILAENV.
INFO — INFO is INTEGER
- = 0: successful exit ✨
- < 0: if INFO = -i, the i-th argument had an illegal value
- > 0: if INFO = i, D(i,i) is exactly zero. The factorization has been completed, but the block diagonal matrix D is exactly singular, and division by zero will occur if it is used to solve a system of equations. 💀 no cap

Author

Univ. of Tennessee

Univ. of California Berkeley

Univ. of Colorado Denver

NAG Ltd.

Further Details:

If UPLO = 'U', then A = U*D*U**H, where

      U = P(n)*U(n)* ... *P(k)U(k)* ...,

           (   I    v    0   )   k-s
   U(k) =  (   0    I    0   )   s
           (   0    0    I   )   n-k
              k-s   s   n-k

If s = 1, D(k) overwrites A(k,k), and v overwrites A(1:k-1,k).

If s = 2, the upper triangle of D(k) overwrites A(k-1,k-1), A(k-1,k), and A(k,k), and v overwrites A(1:k-2,k-1:k).

If UPLO = 'L', then A = L*D*L**H, where

      L = P(1)*L(1)* ... *P(k)*L(k)* ...,

           (   I    0     0   )  k-1
   L(k) =  (   0    I     0   )  s
           (   0    v     I   )  n-k-s+1
              k-1   s  n-k-s+1

If s = 1, D(k) overwrites A(k,k), and v overwrites A(k+1:n,k).

If s = 2, the lower triangle of D(k) overwrites A(k,k), A(k+1,k), and A(k+1,k+1), and v overwrites A(k+2:n,k:k+1).

🧚 ZHETRF_AA 🧚

subroutine zhetrf_aa (character UPLO, integer N, complex*16, dimension( lda, * ) A, integer LDA, integer, dimension( * ) IPIV, complex*16, dimension( * ) WORK, integer LWORK, integer INFO)

Purpose:

ZHETRF_AA computes the factorization of a complex hermitian matrix A using the Aasen's algorithm. The form of the factorization is

   A = U**H*T*U  or  A = L*T*L**H

where U (or L) is a product of permutation and unit upper (lower) triangular matrices, and T is a hermitian tridiagonal matrix. 🌸

This is the blocked version of the algorithm, calling Level 3 BLAS.

Parameters

UPLO — UPLO is CHARACTER*1
- = 'U': Upper triangle of A is stored;
- = 'L': Lower triangle of A is stored.
N — N is INTEGER. The order of the matrix A. N >= 0.
A — A is COMPLEX*16 array, dimension (LDA,N). On entry, the hermitian matrix A. If UPLO = 'U', the leading N-by-N upper triangular part of A contains the upper triangular part of the matrix A, and the strictly lower triangular part of A is not referenced. If UPLO = 'L', the leading N-by-N lower triangular part of A contains the lower triangular part of the matrix A, and the strictly upper triangular part of A is not referenced.

On exit, the tridiagonal matrix is stored in the diagonals and the subdiagonals of A just below (or above) the diagonals, and L is stored below (or above) the subdiagonals, when UPLO is 'L' (or 'U'). ✨

LDA — LDA is INTEGER. The leading dimension of the array A. LDA >= max(1,N).
IPIV — IPIV is INTEGER array, dimension (N). On exit, it contains the details of the interchanges, i.e., the row and column k of A were interchanged with the row and column IPIV(k).
WORK — WORK is COMPLEX*16 array, dimension (MAX(1,LWORK)). On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
LWORK — LWORK is INTEGER. The length of WORK. LWORK >= MAX(1,2*N). For optimum performance LWORK >= N*(1+NB), where NB is the optimal blocksize.

INFO — INFO is INTEGER
- = 0: successful exit 💖
- < 0: if INFO = -i, the i-th argument had an illegal value.

Author

Univ. of Tennessee

Univ. of California Berkeley

Univ. of Colorado Denver

NAG Ltd.

🔥 ZHETRF_RK 🔥

subroutine zhetrf_rk (character UPLO, integer N, complex*16, dimension( lda, * ) A, integer LDA, complex*16, dimension( * ) E, integer, dimension( * ) IPIV, complex*16, dimension( * ) WORK, integer LWORK, integer INFO)

ZHETRF_RK computes the factorization of a complex Hermitian indefinite matrix using the bounded Bunch-Kaufman (rook) diagonal pivoting method (BLAS3 blocked algorithm).

Purpose:

ZHETRF_RK computes the factorization of a complex Hermitian matrix A using the bounded Bunch-Kaufman (rook) diagonal pivoting method:

   A = P*U*D*(U**H)*(P**T) or A = P*L*D*(L**H)*(P**T),

This is the blocked version of the algorithm, calling Level 3 BLAS. For more information see Further Details section.

Parameters

UPLO — UPLO is CHARACTER*1. Specifies whether the upper or lower triangular part of the Hermitian matrix A is stored:
- = 'U': Upper triangular
- = 'L': Lower triangular
N — N is INTEGER. The order of the matrix A. N >= 0.
A — A is COMPLEX*16 array, dimension (LDA,N). On entry, the Hermitian matrix A.
- If UPLO = 'U': the leading N-by-N upper triangular part of A contains the upper triangular part of the matrix A, and the strictly lower triangular part of A is not referenced.
- If UPLO = 'L': the leading N-by-N lower triangular part of A contains the lower triangular part of the matrix A, and the strictly upper triangular part of A is not referenced.

On exit, contains:

ONLY diagonal elements of the Hermitian block diagonal matrix D on the diagonal of A, i.e. D(k,k) = A(k,k); (superdiagonal (or subdiagonal) elements of D are stored on exit in array E), and
If UPLO = 'U': factor U in the superdiagonal part of A. If UPLO = 'L': factor L in the subdiagonal part of A. 💅

LDA — LDA is INTEGER. The leading dimension of the array A. LDA >= max(1,N).
E — ⏳ to be continued…

E — a COMPLEX*16 array, dimension (N) 💜

On exit, bestie, it's holding the superdiagonal (or subdiagonal) elements of the Hermitian block diagonal matrix D, the one with all those 1-by-1 or 2-by-2 diagonal blocks, where:

If UPLO = 'U': E(i) = D(i-1,i), i=2:N, E(1) is set to 0;
If UPLO = 'L': E(i) = D(i+1,i), i=1:N-1, E(N) is set to 0.

NOTE: For a 1-by-1 diagonal block D(k), where 1 <= k <= N, the element E(k) is set to 0 in both the UPLO = 'U' and the UPLO = 'L' cases. No exceptions, no loopholes 💅

IPIV — INTEGER array, dimension (N) 🎀

IPIV is the whole story of the permutation matrix P in the factorization of matrix A, and here's how it goes. The absolute value of IPIV(k) represents the index of the row and column that were interchanged with the k-th row and column. The value of UPLO describes the order in which the interchanges were applied. And the sign of IPIV represents the block structure of the Hermitian block diagonal matrix D with 1-by-1 or 2-by-2 diagonal blocks, which correspond to 1 or 2 interchanges at each factorization step. Want the whole lore drop? See the Further Details section, it's all there ✨

If UPLO = 'U', and in factorization order k decreases from N to 1:

a) A single positive entry IPIV(k) > 0 means:
- D(k,k) is a 1-by-1 diagonal block.
- If IPIV(k) != k, rows and columns k and IPIV(k) were interchanged in the matrix A(1:N,1:N);
- If IPIV(k) = k, no interchange occurred.
b) A pair of consecutive negative entries IPIV(k) < 0 and IPIV(k-1) < 0 means:
- D(k-1:k,k-1:k) is a 2-by-2 diagonal block.
- (NOTE: negative entries in IPIV appear ONLY in pairs 👀)
- 1) If -IPIV(k) != k, rows and columns k and -IPIV(k) were interchanged in the matrix A(1:N,1:N). If -IPIV(k) = k, no interchange occurred.
- 2) If -IPIV(k-1) != k-1, rows and columns k-1 and -IPIV(k-1) were interchanged in the matrix A(1:N,1:N). If -IPIV(k-1) = k-1, no interchange occurred.
c) In both cases a) and b), always ABS( IPIV(k) ) <= k.
d) NOTE: Any entry IPIV(k) is always NONZERO on output. Every single one, periodt.

If UPLO = 'L', and in factorization order k increases from 1 to N:

a) A single positive entry IPIV(k) > 0 means:
- D(k,k) is a 1-by-1 diagonal block.
- If IPIV(k) != k, rows and columns k and IPIV(k) were interchanged in the matrix A(1:N,1:N).
- If IPIV(k) = k, no interchange occurred.
b) A pair of consecutive negative entries IPIV(k) < 0 and IPIV(k+1) < 0 means:
- D(k:k+1,k:k+1) is a 2-by-2 diagonal block.
- (NOTE: negative entries in IPIV appear ONLY in pairs 👀)
- 1) If -IPIV(k) != k, rows and columns k and -IPIV(k) were interchanged in the matrix A(1:N,1:N). If -IPIV(k) = k, no interchange occurred.
- 2) If -IPIV(k+1) != k+1, rows and columns k-1 and -IPIV(k-1) were interchanged in the matrix A(1:N,1:N). If -IPIV(k+1) = k+1, no interchange occurred.
c) In both cases a) and b), always ABS( IPIV(k) ) >= k.
d) NOTE: Any entry IPIV(k) is always NONZERO on output. Same energy as above, just flipped 🔄

WORK — COMPLEX*16 array, dimension ( MAX(1,LWORK) ) 🦋

On exit, if INFO = 0, WORK(1) returns the optimal LWORK. It's giving "here's the number you actually wanted" 💖

LWORK — INTEGER

The length of WORK. LWORK >=1. For best performance, LWORK >= N*NB, where NB is the block size returned by ILAENV. Go big or go home, honestly.

INFO — INTEGER

= 0: successful exit 🎉
< 0: If INFO = -k, the k-th argument had an illegal value 🚨
> 0: If INFO = k, the matrix A is singular, because:
- If UPLO = 'U': column k in the upper triangular part of A contains all zeros.
- If UPLO = 'L': column k in the lower triangular part of A contains all zeros.

Therefore D(k,k) is exactly zero, and superdiagonal elements of column k of U (or subdiagonal elements of column k of L) are all zeros. The factorization has been completed, but the block diagonal matrix D is exactly singular, and division by zero will occur if it is used to solve a system of equations. Like, she finished the job but the result is a landmine 💀

NOTE: INFO only stores the first occurrence of a singularity, any subsequent occurrence of singularity is not stored in INFO even though the factorization always completes. One warning, that's the whole alert system bestie.

Author 👑

Univ. of Tennessee

Univ. of California Berkeley

Univ. of Colorado Denver

NAG Ltd.

Further Details: 📝

TODO: put correct description

(not me leaving the homework blank, iconic honestly)

Contributors: 🌸

December 2016, Igor Kozachenko, Computer Science Division, University of California, Berkeley
September 2007, Sven Hammarling, Nicholas J. Higham, Craig Lucas, School of Mathematics, University of Manchester

`subroutine zhetrf_rook` 💅

subroutine zhetrf_rook (character UPLO, integer N, complex*16, dimension(
    lda, * ) A, integer LDA, integer, dimension( * ) IPIV, complex*16,
    dimension( * ) WORK, integer LWORK, integer INFO)

ZHETRF_ROOK computes the factorization of a complex Hermitian indefinite matrix using the bounded Bunch-Kaufman ('rook') diagonal pivoting method (blocked algorithm, calling Level 3 BLAS). Yes, she does all of that 🫶

Purpose: ✨

ZHETRF_ROOK computes the factorization of a complex Hermitian matrix A using the bounded Bunch-Kaufman ('rook') diagonal pivoting method. The form of the factorization is

A = U*D*U**T  or  A = L*D*L**T

where U (or L) is a product of permutation and unit upper (lower) triangular matrices, and D is Hermitian and block diagonal with 1-by-1 and 2-by-2 diagonal blocks.

This is the blocked version of the algorithm, calling Level 3 BLAS. Main character energy, obviously 🔥

Parameters

UPLO — CHARACTER*1 🌈

= 'U': Upper triangle of A is stored;
= 'L': Lower triangle of A is stored.

N — INTEGER

The order of the matrix A. N >= 0.

A — COMPLEX*16 array, dimension (LDA,N) 💜

On entry, the Hermitian matrix A. If UPLO = 'U', the leading N-by-N upper triangular part of A contains the upper triangular part of the matrix A, and the strictly lower triangular part of A is not referenced. If UPLO = 'L', the leading N-by-N lower triangular part of A contains the lower triangular part of the matrix A, and the strictly upper triangular part of A is not referenced. Mind ya business, other triangle 🙈

On exit, the block diagonal matrix D and the multipliers used to obtain the factor U or L (see below for further details).

LDA — INTEGER

The leading dimension of the array A. LDA >= max(1,N).

IPIV — INTEGER array, dimension (N) 🎀

Details of the interchanges and the block structure of D.

If UPLO = 'U':

Only the last KB elements of IPIV are set.
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 IPIV(k) < 0 and IPIV(k-1) < 0, then rows and columns k and -IPIV(k) were interchanged and rows and columns k-1 and -IPIV(k-1) were interchanged, D(k-1:k,k-1:k) is a 2-by-2 diagonal block.

If UPLO = 'L':

Only the first KB elements of IPIV are set.
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 IPIV(k) < 0 and IPIV(k+1) < 0, then rows and columns k and -IPIV(k) were interchanged and rows and columns k+1 and -IPIV(k+1) were interchanged, D(k:k+1,k:k+1) is a 2-by-2 diagonal block.

WORK — COMPLEX*16 array, dimension (MAX(1,LWORK)) 🦋

On exit, if INFO = 0, WORK(1) returns the optimal LWORK.

LWORK — INTEGER

The length of WORK. LWORK >=1. For best performance LWORK >= N*NB, where NB is the block size returned by ILAENV.

INFO — INTEGER

= 0: successful exit 🎉
< 0: if INFO = -i, the i-th argument had an illegal value 🚨
> 0: if INFO = i, D(i,i) is exactly zero. The factorization has been completed, but the block diagonal matrix D is exactly singular, and division by zero will occur if it is used to solve a system of equations. Yeah. Watch out for that one 💀

Author 👑

Univ. of Tennessee

Univ. of California Berkeley

Univ. of Colorado Denver

NAG Ltd.

Further Details: 📝

okay so if UPLO = 'U', then bestie, A = U*D*U**T, where

   U = P(n)*U(n)* ... *P(k)U(k)* ...,

i.e., U is literally a product of terms P(k)*U(k), where k decreases from n to 1 in steps of 1 or 2, and D is a block diagonal matrix with 1-by-1 and 2-by-2 diagonal blocks D(k). P(k) is a permutation matrix as defined by IPIV(k), and U(k) is a unit upper triangular matrix, such that if the diagonal block D(k) is of order s (s = 1 or 2), then

         (   I    v    0   )   k-s
U(k) =   (   0    I    0   )   s
         (   0    0    I   )   n-k
            k-s   s   n-k

If s = 1, D(k) overwrites A(k,k), and v overwrites A(1:k-1,k). If s = 2, the upper triangle of D(k) overwrites A(k-1,k-1), A(k-1,k), and A(k,k), and v overwrites A(1:k-2,k-1:k). periodt.

If UPLO = 'L', then A = L*D*L**T, where

   L = P(1)*L(1)* ... *P(k)*L(k)* ...,

i.e., L is a product of terms P(k)*L(k), where k increases from 1 to n in steps of 1 or 2, and D is a block diagonal matrix with 1-by-1 and 2-by-2 diagonal blocks D(k). P(k) is a permutation matrix as defined by IPIV(k), and L(k) is a unit lower triangular matrix, such that if the diagonal block D(k) is of order s (s = 1 or 2), then

         (   I    0     0   )  k-1
L(k) =   (   0    I     0   )  s
         (   0    v     I   )  n-k-s+1
            k-1   s  n-k-s+1

If s = 1, D(k) overwrites A(k,k), and v overwrites A(k+1:n,k). If s = 2, the lower triangle of D(k) overwrites A(k,k), A(k+1,k), and A(k+1,k+1), and v overwrites A(k+2:n,k:k+1). slay. 💅

Contributors:

June 2016, Igor Kozachenko, Computer Science Division, University of California, Berkeley

September 2007, Sven Hammarling, Nicholas J. Higham, Craig Lucas, School of Mathematics, University of Manchester

subroutine zhetri (character UPLO, integer N, complex16, dimension( lda, ) A, integer LDA, integer, dimension( * ) IPIV, complex16, dimension( ) WORK, integer INFO)

ZHETRI

okay so this one is the moment 💅. ZHETRI computes the inverse of a complex Hermitian indefinite matrix A using the factorization A = U*D*U**H or A = L*D*L**H computed by ZHETRF.

Parameters

UPLO — UPLO is CHARACTER*1. It's the vibe check for whether the details of the factorization are stored as an upper or lower triangular matrix. = 'U': Upper triangular, form is A = U*D*U**H; = 'L': Lower triangular, form is A = L*D*L**H.
N — N is INTEGER. The order of the matrix A. N >= 0. No cap.
A — A is COMPLEX*16 array, dimension (LDA,N). On entry, the block diagonal matrix D and the multipliers used to obtain the factor U or L as computed by ZHETRF.

On exit, if INFO = 0, the (Hermitian) inverse of the original matrix. If UPLO = 'U', the upper triangular part of the inverse is formed and the part of A below the diagonal is not referenced; if UPLO = 'L' the lower triangular part of the inverse is formed and the part of A above the diagonal is not referenced.

LDA — LDA is INTEGER. The leading dimension of the array A. LDA >= max(1,N).
IPIV — IPIV is INTEGER array, dimension (N). Details of the interchanges and the block structure of D as determined by ZHETRF.
WORK — WORK is COMPLEX*16 array, dimension (N).
INFO — INFO is INTEGER. = 0: successful exit < 0: if INFO = -i, the i-th argument had an illegal value

0: if INFO = i, D(i,i) = 0; the matrix is singular and its inverse could not be computed.

Author

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

subroutine zhetri2 (character UPLO, integer N, complex16, dimension( lda, ) A, integer LDA, integer, dimension( * ) IPIV, complex16, dimension( ) WORK, integer LWORK, integer INFO)

ZHETRI2

ZHETRI2 computes the inverse of a COMPLEX*16 hermitian indefinite matrix A using the factorization A = U*D*U**T or A = L*D*L**T computed by ZHETRF. ZHETRI2 set the LEADING DIMENSION of the workspace before calling ZHETRI2X that actually computes the inverse. It's giving very much "I prep the workspace, you do the math" ✨.

Parameters

UPLO — UPLO is CHARACTER*1. Specifies whether the details of the factorization are stored as an upper or lower triangular matrix. = 'U': Upper triangular, form is A = U*D*U**T; = 'L': Lower triangular, form is A = L*D*L**T.
N — N is INTEGER. The order of the matrix A. N >= 0.
A — A is COMPLEX*16 array, dimension (LDA,N). On entry, the block diagonal matrix D and the multipliers used to obtain the factor U or L as computed by ZHETRF.

On exit, if INFO = 0, the (symmetric) inverse of the original matrix. If UPLO = 'U', the upper triangular part of the inverse is formed and the part of A below the diagonal is not referenced; if UPLO = 'L' the lower triangular part of the inverse is formed and the part of A above the diagonal is not referenced.

LDA — LDA is INTEGER. The leading dimension of the array A. LDA >= max(1,N).
IPIV — IPIV is INTEGER array, dimension (N). Details of the interchanges and the block structure of D as determined by ZHETRF.
WORK — WORK is COMPLEX*16 array, dimension (N+NB+1)*(NB+3).
LWORK — LWORK is INTEGER. The dimension of the array WORK. WORK is size >= (N+NB+1)*(NB+3). If LWORK = -1, then a workspace query is assumed; the routine calculates:
- the optimal size of the WORK array, returns this value as the first entry of the WORK array,
- and no error message related to LWORK is issued by XERBLA.
INFO — INFO is INTEGER. = 0: successful exit < 0: if INFO = -i, the i-th argument had an illegal value

0: if INFO = i, D(i,i) = 0; the matrix is singular and its inverse could not be computed.

Author

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

subroutine zhetri2x (character UPLO, integer N, complex16, dimension( lda, ) A, integer LDA, integer, dimension( * ) IPIV, complex16, dimension( n+nb+1, ) WORK, integer NB, integer INFO)

ZHETRI2X

ZHETRI2X computes the inverse of a COMPLEX*16 Hermitian indefinite matrix A using the factorization A = U*D*U**H or A = L*D*L**H computed by ZHETRF. This is the actual computation, no cap 🧢.

Parameters

UPLO — UPLO is CHARACTER*1. Specifies whether the details of the factorization are stored as an upper or lower triangular matrix. = 'U': Upper triangular, form is A = U*D*U**H; = 'L': Lower triangular, form is A = L*D*L**H.
N — N is INTEGER. The order of the matrix A. N >= 0.
A — A is COMPLEX*16 array, dimension (LDA,N). On entry, the NNB diagonal matrix D and the multipliers used to obtain the factor U or L as computed by ZHETRF.

LDA — LDA is INTEGER. The leading dimension of the array A. LDA >= max(1,N).
IPIV — IPIV is INTEGER array, dimension (N). Details of the interchanges and the NNB structure of D as determined by ZHETRF.
WORK — WORK is COMPLEX*16 array, dimension (N+NB+1,NB+3).
NB — NB is INTEGER. Block size.
INFO — INFO is INTEGER. = 0: successful exit < 0: if INFO = -i, the i-th argument had an illegal value

0: if INFO = i, D(i,i) = 0; the matrix is singular and its inverse could not be computed.

Author

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

subroutine zhetri_3 (character UPLO, integer N, complex16, dimension( lda, ) A, integer LDA, complex16, dimension( ) E, integer, dimension( * ) IPIV, complex16, dimension( ) WORK, integer LWORK, integer INFO)

ZHETRI_3

ZHETRI_3 computes the inverse of a complex Hermitian indefinite matrix A using the factorization computed by ZHETRF_RK or ZHETRF_BK:

    A = P*U*D*(U**H)*(P**T) or A = P*L*D*(L**H)*(P**T),

where U (or L) is unit upper (or lower) triangular matrix, U**H (or L**H) is the conjugate of U (or L), P is a permutation matrix, P**T is the transpose of P, and D is Hermitian and block diagonal with 1-by-1 and 2-by-2 diagonal blocks.

ZHETRI_3 sets the leading dimension of the workspace before calling ZHETRI_3X that actually computes the inverse. This is the blocked version of the algorithm, calling Level 3 BLAS. It's the glow-up version 💅.

Parameters

UPLO — UPLO is CHARACTER*1. Specifies whether the details of the factorization are stored as an upper or lower triangular matrix. = 'U': Upper triangle of A is stored; = 'L': Lower triangle of A is stored.
N — N is INTEGER. The order of the matrix A. N >= 0.
A —
A — A is COMPLEX*16 array, dimension (LDA,N) 💖

On entry, bestie, you're handing over the diagonal of the block diagonal matrix D plus the factors U or L exactly as ZHETRF_RK and ZHETRF_BK cooked them up:

ONLY the diagonal elements of the Hermitian block diagonal matrix D live on the diagonal of A, i.e. D(k,k) = A(k,k); (the superdiagonal (or subdiagonal) elements of D should be provided on entry in array E), and
If UPLO = 'U': factor U sits in the superdiagonal part of A. If UPLO = 'L': factor L sits in the subdiagonal part of A.

On exit, if INFO = 0, you get the Hermitian inverse of the original matrix, slay:

If UPLO = 'U': the upper triangular part of the inverse is formed and the part of A below the diagonal is not referenced;
If UPLO = 'L': the lower triangular part of the inverse is formed and the part of A above the diagonal is not referenced.

LDA — LDA is INTEGER. The leading dimension of the array A. LDA >= max(1,N). 🫶
E — E is COMPLEX*16 array, dimension (N) 🦋

On entry, it holds the superdiagonal (or subdiagonal) elements of the Hermitian block diagonal matrix D with 1-by-1 or 2-by-2 diagonal blocks, where:

If UPLO = 'U': E(i) = D(i-1,i),i=2:N, E(1) not referenced;
If UPLO = 'L': E(i) = D(i+1,i),i=1:N-1, E(N) not referenced.

NOTE: For a 1-by-1 diagonal block D(k), where 1 <= k <= N, the element E(k) is not referenced in both the UPLO = 'U' or UPLO = 'L' cases.

IPIV — IPIV is INTEGER array, dimension (N). All the deets on the interchanges and the block structure of D, exactly as determined by ZHETRF_RK or ZHETRF_BK. 👑
WORK — WORK is COMPLEX*16 array, dimension (N+NB+1)*(NB+3). On exit, if INFO = 0, WORK(1) returns the optimal LWORK. ✨
LWORK — LWORK is INTEGER. The length of WORK. LWORK >= (N+NB+1)*(NB+3). 📏

If LDWORK = -1, then a workspace query is assumed; the routine only calculates the optimal size of the optimal size of the WORK array (yes bestie, the original really does say "the optimal size of the optimal size" twice, and honestly? iconic redundancy, we stan a queen who repeats herself 💅), returns this value as the first entry of the WORK array, and no error message related to LWORK is issued by XERBLA.

INFO — INFO is INTEGER
- = 0: successful exit 🌸
- < 0: if INFO = -i, the i-th argument had an illegal value 😭
- > 0: if INFO = i, D(i,i) = 0; the matrix is singular and its inverse could not be computed. 💀

Author 🌈

Univ. of Tennessee

Univ. of California Berkeley

Univ. of Colorado Denver

NAG Ltd.

Contributors: 🎀

November 2017,  Igor Kozachenko,
                Computer Science Division,
                University of California, Berkeley

👑 subroutine `zhetri_3x` 👑

subroutine zhetri_3x (character UPLO, integer N, complex*16, dimension( lda,
* ) A, integer LDA, complex*16, dimension( * ) E, integer, dimension( *
) IPIV, complex*16, dimension( n+nb+1, * ) WORK, integer NB, integer
INFO)

ZHETRI_3X ✨

Purpose 💜

ZHETRI_3X computes the inverse of a complex Hermitian indefinite matrix A using the factorization computed by ZHETRF_RK or ZHETRF_BK:

A = P*U*D*(U**H)*(P**T) or A = P*L*D*(L**H)*(P**T),

This is the blocked version of the algorithm, calling Level 3 BLAS. No cap, she's out here doing the most. 💅

Parameters 🦋

UPLO — UPLO is CHARACTER*1. Specifies whether the details of the factorization are stored as an upper or lower triangular matrix.
- = 'U': Upper triangle of A is stored;
- = 'L': Lower triangle of A is stored.
N — N is INTEGER. The order of the matrix A. N >= 0. ⭐
A — A is COMPLEX*16 array, dimension (LDA,N)

On entry, diagonal of the block diagonal matrix D and factors U or L as computed by ZHETRF_RK and ZHETRF_BK:

ONLY diagonal elements of the Hermitian block diagonal matrix D on the diagonal of A, i.e. D(k,k) = A(k,k); (superdiagonal (or subdiagonal) elements of D should be provided on entry in array E), and
If UPLO = 'U': factor U in the superdiagonal part of A. If UPLO = 'L': factor L in the subdiagonal part of A.

On exit, if INFO = 0, the Hermitian inverse of the original matrix. If UPLO = 'U': the upper triangular part of the inverse is formed and the part of A below the diagonal is not referenced; if UPLO = 'L': the lower triangular part of the inverse is formed and the part of A above the diagonal is not referenced.

LDA — LDA is INTEGER. The leading dimension of the array A. LDA >= max(1,N). 🫶
E — E is COMPLEX*16 array, dimension (N)

On entry, contains the superdiagonal (or subdiagonal) elements of the Hermitian block diagonal matrix D with 1-by-1 or 2-by-2 diagonal blocks, where:

If UPLO = 'U': E(i) = D(i-1,i), i=2:N, E(1) not referenced;
If UPLO = 'L': E(i) = D(i+1,i), i=1:N-1, E(N) not referenced.

NOTE: For 1-by-1 diagonal block D(k), where 1 <= k <= N, the element E(k) is not referenced in both UPLO = 'U' or UPLO = 'L' cases.

IPIV — IPIV is INTEGER array, dimension (N). Details of the interchanges and the block structure of D as determined by ZHETRF_RK or ZHETRF_BK. 👑
WORK — WORK is COMPLEX*16 array, dimension (N+NB+1,NB+3). ✨
NB — NB is INTEGER. Block size. 🧱
INFO — INFO is INTEGER
- = 0: successful exit 🌸
- < 0: if INFO = -i, the i-th argument had an illegal value 😭
- > 0: if INFO = i, D(i,i) = 0; the matrix is singular and its inverse could not be computed. 💀

Author 🌈

Univ. of Tennessee

Univ. of California Berkeley

Univ. of Colorado Denver

NAG Ltd.

Contributors: 🎀

June 2017,  Igor Kozachenko,
                Computer Science Division,
                University of California, Berkeley

💖 subroutine `zhetri_rook` 💖

subroutine zhetri_rook (character UPLO, integer N, complex*16, dimension(
lda, * ) A, integer LDA, integer, dimension( * ) IPIV, complex*16,
dimension( * ) WORK, integer INFO)

ZHETRI_ROOK computes the inverse of HE matrix using the factorization obtained with the bounded Bunch-Kaufman ('rook') diagonal pivoting method. 🏰 (rook as in the chess piece, not the bird, bestie — although honestly both are valid energy)

Purpose 💜

ZHETRI_ROOK computes the inverse of a complex Hermitian indefinite matrix A using the factorization A = U*D*U**H or A = L*D*L**H computed by ZHETRF_ROOK.

Parameters 🦋

UPLO — UPLO is CHARACTER*1. Specifies whether the details of the factorization are stored as an upper or lower triangular matrix.
- = 'U': Upper triangular, form is A = U*D*U**H;
- = 'L': Lower triangular, form is A = L*D*L**H.
N — N is INTEGER. The order of the matrix A. N >= 0. ⭐
A — A is COMPLEX*16 array, dimension (LDA,N)

On entry, the block diagonal matrix D and the multipliers used to obtain the factor U or L as computed by ZHETRF_ROOK.

On exit, if INFO = 0, the (Hermitian) inverse of the original matrix. If UPLO = 'U', the upper triangular part of the inverse is formed and the part of A below the diagonal is not referenced; if UPLO = 'L' the lower triangular part of the inverse is formed and the part of A above the diagonal is not referenced. 🙈

LDA — LDA is INTEGER. The leading dimension of the array A. LDA >= max(1,N). 🫶
IPIV — IPIV is INTEGER array, dimension (N). Details of the interchanges and the block structure of D as determined by ZHETRF_ROOK. 👑
WORK — WORK is COMPLEX*16 array, dimension (N). ✨
INFO — INFO is INTEGER
- = 0: successful exit 🌸
- < 0: if INFO = -i, the i-th argument had an illegal value 😭
- > 0: if INFO = i, D(i,i) = 0; the matrix is singular and its inverse could not be computed. 💀

Author 🌈

Univ. of Tennessee

Univ. of California Berkeley

Univ. of Colorado Denver

NAG Ltd.

Contributors: 🎀

November 2013,  Igor Kozachenko,
                Computer Science Division,
                University of California, Berkeley

September 2007, Sven Hammarling, Nicholas J. Higham, Craig Lucas,
                School of Mathematics,
                University of Manchester

🌟 subroutine `zhetrs` 🌟

subroutine zhetrs (character UPLO, integer N, integer NRHS, complex*16,
dimension( lda, * ) A, integer LDA, integer, dimension( * ) IPIV,
complex*16, dimension( ldb, * ) B, integer LDB, integer INFO)

ZHETRS

Purpose 💜

ZHETRS solves a system of linear equations A*X = B with a complex Hermitian matrix A using the factorization A = U*D*U**H or A = L*D*L**H computed by ZHETRF. Solving systems? She does that. Periodt. 💅

Parameters 🦋

UPLO — UPLO is CHARACTER*1. Specifies whether the details of the factorization are stored as an upper or lower triangular matrix.
- = 'U': Upper triangular, form is A = U*D*U**H;
- = 'L': Lower triangular, form is A = L*D*L**H.
N — N is INTEGER. The order of the matrix A. N >= 0. ⭐
NRHS — NRHS is INTEGER. The number of right hand sides, i.e., the number of columns of the matrix B. NRHS >= 0. 🍓
A —
A — A is a COMPLEX*16 array, dimension (LDA,N). It's holding the block diagonal matrix D plus all the multipliers that got you the factor U or L, exactly the way ZHETRF computed them. 💖
LDA — LDA is INTEGER. It's the leading dimension of the array A, and the rule is LDA >= max(1,N). bestie pls respect the leading dimension, she is not decorative. 🙈
IPIV — IPIV is INTEGER array, dimension (N). It's got the details of the interchanges and the block structure of D, as determined by ZHETRF. 🎀
B — B is COMPLEX*16 array, dimension (LDB,NRHS). On entry it's the right hand side matrix B. On exit? It's the solution matrix X. the way it just becomes the answer, main character energy fr. ✨
LDB — LDB is INTEGER. The leading dimension of the array B. LDB >= max(1,N). 💅
INFO — INFO is INTEGER.
- = 0: successful exit, we ate, no notes. 👑
- < 0: if INFO = -i, the i-th argument had an illegal value. something you passed in was simply not it. 💀

👑 Author

Univ. of Tennessee

Univ. of California Berkeley

Univ. of Colorado Denver

NAG Ltd.

🦋 subroutine zhetrs2

subroutine zhetrs2 (character UPLO, integer N, integer NRHS, complex*16,
    dimension( lda, * ) A, integer LDA, integer, dimension( * ) IPIV,
    complex*16, dimension( ldb, * ) B, integer LDB, complex*16, dimension( *
    ) WORK, integer INFO)

ZHETRS2

Purpose:

ZHETRS2 solves a system of linear equations A*X = B with a complex Hermitian matrix A using the factorization A = U*D*U**H or A = L*D*L**H computed by ZHETRF and converted by ZSYCONV. it's giving "we already did the hard part, now we just solve." 😌

Parameters

UPLO — UPLO is CHARACTER*1. It says whether the details of the factorization are stored as an upper or lower triangular matrix.
- = 'U': Upper triangular, form is A = U*D*U**H;
- = 'L': Lower triangular, form is A = L*D*L**H.
N — N is INTEGER. The order of the matrix A. N >= 0. 🦋
NRHS — NRHS is INTEGER. The number of right hand sides, i.e., the number of columns of the matrix B. NRHS >= 0. 🌸
A — A is COMPLEX*16 array, dimension (LDA,N). The block diagonal matrix D and the multipliers used to obtain the factor U or L as computed by ZHETRF. ✨
LDA — LDA is INTEGER. The leading dimension of the array A. LDA >= max(1,N). 💅
IPIV — IPIV is INTEGER array, dimension (N). Details of the interchanges and the block structure of D as determined by ZHETRF. 🎀
B — B is COMPLEX*16 array, dimension (LDB,NRHS). On entry, the right hand side matrix B. On exit, the solution matrix X. 🔥
LDB — LDB is INTEGER. The leading dimension of the array B. LDB >= max(1,N). 💖
WORK — WORK is COMPLEX*16 array, dimension (N). 🧚
INFO — INFO is INTEGER.
- = 0: successful exit
- < 0: if INFO = -i, the i-th argument had an illegal value

👑 Author

Univ. of Tennessee

Univ. of California Berkeley

Univ. of Colorado Denver

NAG Ltd.

🌈 subroutine zhetrs_3

subroutine zhetrs_3 (character UPLO, integer N, integer NRHS, complex*16,
    dimension( lda, * ) A, integer LDA, complex*16, dimension( * ) E,
    integer, dimension( * ) IPIV, complex*16, dimension( ldb, * ) B, integer
    LDB, integer INFO)

ZHETRS_3

Purpose:

ZHETRS_3 solves a system of linear equations A * X = B with a complex Hermitian matrix A using the factorization computed by ZHETRF_RK or ZHETRF_BK:

A = P*U*D*(U**H)*(P**T) or A = P*L*D*(L**H)*(P**T),

This algorithm is using Level 3 BLAS. the way it's out here doing Level 3 BLAS like it's nothing — iconic. 🔥

Parameters

UPLO — UPLO is CHARACTER*1. Specifies whether the details of the factorization are stored as an upper or lower triangular matrix:
- = 'U': Upper triangular, form is A = P*U*D*(U**H)*(P**T);
- = 'L': Lower triangular, form is A = P*L*D*(L**H)*(P**T).
N — N is INTEGER. The order of the matrix A. N >= 0. 🦋
NRHS — NRHS is INTEGER. The number of right hand sides, i.e., the number of columns of the matrix B. NRHS >= 0. 🌸
A — A is COMPLEX*16 array, dimension (LDA,N). Diagonal of the block diagonal matrix D and factors U or L as computed by ZHETRF_RK and ZHETRF_BK:
1. ONLY diagonal elements of the Hermitian block diagonal matrix D on the diagonal of A, i.e. D(k,k) = A(k,k); (superdiagonal (or subdiagonal) elements of D should be provided on entry in array E), and
2. If UPLO = 'U': factor U in the superdiagonal part of A. If UPLO = 'L': factor L in the subdiagonal part of A.
LDA — LDA is INTEGER. The leading dimension of the array A. LDA >= max(1,N). 💅
E — E is COMPLEX*16 array, dimension (N). On entry, contains the superdiagonal (or subdiagonal) elements of the Hermitian block diagonal matrix D with 1-by-1 or 2-by-2 diagonal blocks, where

If UPLO = 'U': E(i) = D(i-1,i),i=2:N, E(1) not referenced;
If UPLO = 'L': E(i) = D(i+1,i),i=1:N-1, E(N) not referenced.

NOTE: For 1-by-1 diagonal block D(k), where 1 <= k <= N, the element E(k) is not referenced in both UPLO = 'U' or UPLO = 'L' cases. so don't come for E(k), she's simply not in the conversation for that block. 🙈

IPIV — IPIV is INTEGER array, dimension (N). Details of the interchanges and the block structure of D as determined by ZHETRF_RK or ZHETRF_BK. 🎀
B — B is COMPLEX*16 array, dimension (LDB,NRHS). On entry, the right hand side matrix B. On exit, the solution matrix X. ✨
LDB — LDB is INTEGER. The leading dimension of the array B. LDB >= max(1,N). 💖
INFO — INFO is INTEGER.
- = 0: successful exit
- < 0: if INFO = -i, the i-th argument had an illegal value

👑 Author

Univ. of Tennessee

Univ. of California Berkeley

Univ. of Colorado Denver

NAG Ltd.

🌟 Contributors:

June 2017, Igor Kozachenko, Computer Science Division, University of California, Berkeley
September 2007, Sven Hammarling, Nicholas J. Higham, Craig Lucas, School of Mathematics, University of Manchester

🫶 subroutine zhetrs_aa

subroutine zhetrs_aa (character UPLO, integer N, integer NRHS, complex*16,
    dimension( lda, * ) A, integer LDA, integer, dimension( * ) IPIV,
    complex*16, dimension( ldb, * ) B, integer LDB, complex*16, dimension( *
    ) WORK, integer LWORK, integer INFO)

ZHETRS_AA

Purpose:

ZHETRS_AA solves a system of linear equations A*X = B with a complex hermitian matrix A using the factorization A = U**H*T*U or A = L*T*L**H computed by ZHETRF_AA. it really said "A double-A, no cap." 😌

Parameters

UPLO — UPLO is CHARACTER*1. Specifies whether the details of the factorization are stored as an upper or lower triangular matrix.
- = 'U': Upper triangular, form is A = U**H*T*U;
- = 'L': Lower triangular, form is A = L*T*L**H.
N — N is INTEGER. The order of the matrix A. N >= 0. 🦋
NRHS — NRHS is INTEGER. The number of right hand sides, i.e., the number of columns of the matrix B. NRHS >= 0. 🌸
A — A is COMPLEX*16 array, dimension (LDA,N). Details of factors computed by ZHETRF_AA. ✨
LDA — LDA is INTEGER. The leading dimension of the array A. LDA >= max(1,N). 💅
IPIV — IPIV is INTEGER array, dimension (N). Details of the interchanges as computed by ZHETRF_AA. 🎀
B — B is COMPLEX*16 array, dimension (LDB,NRHS). On entry, the right hand side matrix B. On exit, the solution matrix X. 🔥
LDB — LDB is INTEGER. The leading dimension of the array B. LDB >= max(1,N). 💖
WORK — WORK is COMPLEX*16 array, dimension (MAX(1,LWORK)). 🧚
LWORK — LWORK is INTEGER. The dimension of the array WORK. LWORK >= max(1,3*N-2). ngl if you undersize this one she will simply not participate. 👀
INFO — INFO is INTEGER.
- = 0: successful exit
- < 0: if INFO = -i, the i-th argument had an illegal value

👑 Author

Univ. of Tennessee

Univ. of California Berkeley

Univ. of Colorado Denver

NAG Ltd.

💜 subroutine zhetrs_rook

subroutine zhetrs_rook (character UPLO, integer N, integer NRHS, complex*16,
    dimension( lda, * ) A, integer LDA, integer, dimension( * ) IPIV,
    complex*16, dimension( ldb, * ) B, integer LDB, integer INFO)

ZHETRS_ROOK computes the solution to a system of linear equations A * X = B for HE matrices using factorization obtained with one of the bounded diagonal pivoting methods (max 2 interchanges). rook said two interchanges MAX and honestly, respect. ♟️

Purpose:

ZHETRS_ROOK solves a system of linear equations A*X = B with a complex Hermitian matrix A using the factorization A = U*D*U**H or A = L*D*L**H computed by ZHETRF_ROOK.

Parameters

UPLO — UPLO is CHARACTER*1. Specifies whether the details of the factorization are stored as an upper or lower triangular matrix.
- = 'U': Upper triangular, form is A = U*D*U**H;
- = 'L': Lower triangular, form is A = L*D*L**H.
N — N is INTEGER. The order of the matrix A. N >= 0. 🦋
NRHS — NRHS is INTEGER. The number of right hand sides, i.e., the number of columns of the matrix B. NRHS >= 0. 🌸
A — A is COMPLEX*16 array, dimension (LDA,N). The block diagonal matrix D and the multipliers used to obtain the factor U or L as computed by ZHETRF_ROOK. ✨
LDA — LDA is INTEGER. The leading dimension of the array A. LDA >= max(1,N). 💅
IPIV — it's an INTEGER array, dimension (N). It's holding all the deets on the interchanges and the block structure of D, exactly as ZHETRF_ROOK decided. 💅
B — COMPLEX*16 array, dimension (LDB,NRHS). On entry, it's the right hand side matrix B. On exit? Bestie, it becomes the solution matrix X. ✨
LDB — INTEGER. The leading dimension of array B. LDB >= max(1,N). No cap.
INFO — INTEGER. = 0: successful exit, we love to see it. < 0: if INFO = -i, then the i-th argument had an illegal value, so something was giving messy. 🚫

Author

Univ. of Tennessee

Univ. of California Berkeley

Univ. of Colorado Denver

NAG Ltd.

Contributors:

November 2013, Igor Kozachenko, Computer Science Division, University of California, Berkeley
September 2007, Sven Hammarling, Nicholas J. Higham, Craig Lucas, School of Mathematics, University of Manchester

✨ `subroutine zla_heamv (integer UPLO, integer N, double precision ALPHA, complex16, dimension( lda, ) A, integer LDA, complex16, dimension( ) X, integer INCX, double precision BETA, double precision, dimension( * ) Y, integer INCY)` ✨

ZLA_HEAMV computes a matrix-vector product using a Hermitian indefinite matrix to calculate error bounds.

Purpose:

ZLA_SYAMV performs the matrix-vector operation

y := alpha*abs(A)*abs(x) + beta*abs(y),

where alpha and beta are scalars, x and y are vectors and A is an n by n symmetric matrix.

This function is primarily used in calculating error bounds. To protect against underflow during evaluation, components in the resulting vector are perturbed away from zero by (N+1) times the underflow threshold. To prevent unnecessarily large errors for block-structure embedded in general matrices, 'symbolically' zero components are not perturbed. A zero entry is considered 'symbolic' if all multiplications involved in computing that entry have at least one zero multiplicand.

Parameters

UPLO — INTEGER. On entry, UPLO specifies whether the upper or lower triangular part of the array A is to be referenced as follows:
- UPLO = BLAS_UPPER — Only the upper triangular part of A is to be referenced.
- UPLO = BLAS_LOWER — Only the lower triangular part of A is to be referenced.

Unchanged on exit. 👑

N — INTEGER. On entry, N specifies the number of columns of the matrix A. N must be at least zero. Unchanged on exit.
ALPHA — DOUBLE PRECISION. On entry, ALPHA specifies the scalar alpha. Unchanged on exit.
A — COMPLEX*16 array, dimension ( LDA, n ). Before entry, the leading m by n part of the array A must contain the matrix of coefficients. Unchanged on exit.
LDA — INTEGER. On entry, LDA specifies the first dimension of A as declared in the calling (sub) program. LDA must be at least max( 1, n ). Unchanged on exit.
X — COMPLEX*16 array, dimension at least ( 1 + ( n - 1 )*abs( INCX ) ). Before entry, the incremented array X must contain the vector x. Unchanged on exit.
INCX — INTEGER. On entry, INCX specifies the increment for the elements of X. INCX must not be zero. Unchanged on exit.
BETA — DOUBLE PRECISION. On entry, BETA specifies the scalar beta. When BETA is supplied as zero then Y need not be set on input. Unchanged on exit.
Y — DOUBLE PRECISION array, dimension ( 1 + ( n - 1 )*abs( INCY ) ). Before entry with BETA non-zero, the incremented array Y must contain the vector y. On exit, Y is overwritten by the updated vector y.
INCY — INTEGER. On entry, INCY specifies the increment for the elements of Y. INCY must not be zero. Unchanged on exit. 🦋

Author

Univ. of Tennessee

Univ. of California Berkeley

Univ. of Colorado Denver

NAG Ltd.

Further Details:

Level 2 Blas routine.

Written on 22-October-1986.
- Jack Dongarra, Argonne National Lab.
- Jeremy Du Croz, Nag Central Office.
- Sven Hammarling, Nag Central Office.
- Richard Hanson, Sandia National Labs.
Modified for the absolute-value product, April 2006
- Jason Riedy, UC Berkeley

✨ `double precision function zla_hercond_c (character UPLO, integer N, complex16, dimension( lda, ) A, integer LDA, complex16, dimension( ldaf, ) AF, integer LDAF, integer, dimension( * ) IPIV, double precision, dimension ( * ) C, logical CAPPLY, integer INFO, complex16, dimension( ) WORK, double precision, dimension( * ) RWORK)` ✨

ZLA_HERCOND_C computes the infinity norm condition number of op(A)*inv(diag(c)) for Hermitian indefinite matrices.

Purpose:

ZLA_HERCOND_C computes the infinity norm condition number of op(A) * inv(diag(C)) where C is a DOUBLE PRECISION vector.

Parameters

UPLO — CHARACTER*1
- = 'U': Upper triangle of A is stored;
- = 'L': Lower triangle of A is stored.
N — INTEGER. The number of linear equations, i.e., the order of the matrix A. N >= 0.
A — COMPLEX*16 array, dimension (LDA,N). On entry, the N-by-N matrix A
LDA — INTEGER. The leading dimension of the array A. LDA >= max(1,N).
AF — COMPLEX*16 array, dimension (LDAF,N). The block diagonal matrix D and the multipliers used to obtain the factor U or L as computed by ZHETRF.
LDAF — INTEGER. The leading dimension of the array AF. LDAF >= max(1,N).
IPIV — INTEGER array, dimension (N). Details of the interchanges and the block structure of D as determined by CHETRF.
C — DOUBLE PRECISION array, dimension (N). The vector C in the formula op(A) * inv(diag(C)).
CAPPLY — LOGICAL. If .TRUE. then access the vector C in the formula above.
INFO — INTEGER
- = 0: Successful exit.
- i > 0: The ith argument is invalid.
WORK — COMPLEX*16 array, dimension (2*N). Workspace.
RWORK — DOUBLE PRECISION array, dimension (N). Workspace. 🌸

Author

Univ. of Tennessee

Univ. of California Berkeley

Univ. of Colorado Denver

NAG Ltd.

✨ `double precision function zla_hercond_x (character UPLO, integer N, complex16, dimension( lda, ) A, integer LDA, complex16, dimension( ldaf, ) AF, integer LDAF, integer, dimension( * ) IPIV, complex16, dimension( ) X, integer INFO, complex16, dimension( ) WORK, double precision, dimension( * ) RWORK)` ✨

ZLA_HERCOND_X computes the infinity norm condition number of op(A)*diag(x) for Hermitian indefinite matrices.

Purpose:

ZLA_HERCOND_X computes the infinity norm condition number of op(A) * diag(X) where X is a COMPLEX*16 vector.

Parameters

UPLO — CHARACTER*1
- = 'U': Upper triangle of A is stored;
- = 'L': Lower triangle of A is stored.
N — INTEGER. The number of linear equations, i.e., the order of the matrix A. N >= 0.
A — COMPLEX*16 array, dimension (LDA,N). On entry, the N-by-N matrix A.
LDA — INTEGER. The leading dimension of the array A. LDA >= max(1,N).
AF — COMPLEX*16 array, dimension (LDAF,N). The block diagonal matrix D and the multipliers used to obtain the factor U or L as computed by ZHETRF.
LDAF — INTEGER. The leading dimension of the array AF. LDAF >= max(1,N).
IPIV — INTEGER array, dimension (N). Details of the interchanges and the block structure of D as determined by CHETRF.
X — COMPLEX*16 array, dimension (N). The vector X in the formula op(A) * diag(X).
INFO — INTEGER
- = 0: Successful exit.
- i > 0: The ith argument is invalid.
WORK — COMPLEX*16 array, dimension (2*N). Workspace.
RWORK — DOUBLE PRECISION array, dimension (N). Workspace. 💖

Author

Univ. of Tennessee

Univ. of California Berkeley

Univ. of Colorado Denver

NAG Ltd.

💅 `subroutine zla_herfsx_extended` 💅

subroutine zla_herfsx_extended (integer PREC_TYPE, character UPLO, integer
    N, integer NRHS, complex*16, dimension( lda, * ) A, integer LDA,
    complex*16, dimension( ldaf, * ) AF, integer LDAF, integer, dimension( *
    ) IPIV, logical COLEQU, double precision, dimension( * ) C, complex*16,
    dimension( ldb, * ) B, integer LDB, complex*16, dimension( ldy, * ) Y,
    integer LDY, double precision, dimension( * ) BERR_OUT, integer N_NORMS,
    double precision, dimension( nrhs, * ) ERR_BNDS_NORM, double precision,
    dimension( nrhs, * ) ERR_BNDS_COMP, complex*16, dimension( * ) RES,
    double precision, dimension( * ) AYB, complex*16, dimension( * ) DY,
    complex*16, dimension( * ) Y_TAIL, double precision RCOND, integer
    ITHRESH, double precision RTHRESH, double precision DZ_UB, logical
    IGNORE_CWISE, integer INFO)

ZLA_HERFSX_EXTENDED improves the computed solution to a system of linear equations for Hermitian indefinite matrices by performing extra-precise iterative refinement and provides error bounds and backward error estimates for the solution.

✨ Purpose ✨

So basically, ZLA_HERFSX_EXTENDED glow-ups your computed solution to a system of linear equations by doing that extra-precise iterative refinement thing, and it hands you error bounds + backward error estimates for the solution. 💖

This subroutine is called by ZHERFSX to perform iterative refinement.

On top of the normwise error bound, the code also gives you the maximum componentwise error bound if possible — peek at the comments for ERR_BNDS_NORM and ERR_BNDS_COMP for the deets on those error bounds. 🫶 And heads up: this subroutine is only responsible for setting the second fields of ERR_BNDS_NORM and ERR_BNDS_COMP. That's it. That's the job. 💅

🎀 Parameters 🎀

PREC_TYPE — INTEGER. Specifies the intermediate precision to be used in refinement. The value is defined by ILAPREC(P) where P is a CHARACTER, and P:
- = 'S': Single
- = 'D': Double
- = 'I': Indigenous
- = 'X' or 'E': Extra
UPLO — CHARACTER*1. Which triangle we're vibing with:
- = 'U': Upper triangle of A is stored;
- = 'L': Lower triangle of A is stored.
N — INTEGER. The number of linear equations, i.e., the order of the matrix A. N >= 0.
NRHS — INTEGER. The number of right-hand-sides, i.e., the number of columns of the matrix B.
A — COMPLEX*16 array, dimension (LDA,N). On entry, the N-by-N matrix A.
LDA — INTEGER. The leading dimension of the array A. LDA >= max(1,N).
AF — COMPLEX*16 array, dimension (LDAF,N). The block diagonal matrix D and the multipliers used to obtain the factor U or L as computed by ZHETRF.
LDAF — INTEGER. The leading dimension of the array AF. LDAF >= max(1,N).
IPIV — INTEGER array, dimension (N). Details of the interchanges and the block structure of D as determined by ZHETRF.
COLEQU — LOGICAL. If .TRUE. then column equilibration was done to A before calling this routine. This is needed to compute the solution and error bounds correctly — so don't skip it, bestie.
C — DOUBLE PRECISION array, dimension (N). The column scale factors for A. If COLEQU = .FALSE., C is not accessed. If C is input, each element of C should be a power of the radix to ensure a reliable solution and error estimates. Scaling by powers of the radix does not cause rounding errors unless the result underflows or overflows. Rounding errors during scaling lead to refining with a matrix that is not equivalent to the input matrix, producing error estimates that may not be reliable. 😭

⚠️ That last bit is a real caveat, not a vibe — if your scaling isn't by powers of the radix, your error estimates can go sideways and stop being trustworthy.

B — COMPLEX*16 array, dimension (LDB,NRHS). The right-hand-side matrix B.
LDB — INTEGER. The leading dimension of the array B. LDB >= max(1,N).
Y — COMPLEX*16 array, dimension (LDY,NRHS). On entry, the solution matrix X, as computed by ZHETRS. On exit, the improved solution matrix Y. ✨
LDY — INTEGER. The leading dimension of the array Y. LDY >= max(1,N).
BERR_OUT — DOUBLE PRECISION array, dimension (NRHS). On exit, BERR_OUT(j) contains the componentwise relative backward error for right-hand-side j from the formula

max(i) ( abs(RES(i)) / ( abs(op(A_s))*abs(Y) + abs(B_s) )(i) )

where abs(Z) is the componentwise absolute value of the matrix or vector Z. This is computed by ZLA_LIN_BERR.

N_NORMS — INTEGER. Determines which error bounds to return (see ERR_BNDS_NORM and ERR_BNDS_COMP).
- If N_NORMS >= 1 return normwise error bounds.
- If N_NORMS >= 2 return componentwise error bounds.
ERR_BNDS_NORM — DOUBLE PRECISION array, dimension (NRHS, N_ERR_BNDS). For each right-hand side, this array contains information about various error bounds and condition numbers corresponding to the normwise relative error, which is defined as follows:

Normwise relative error in the ith solution vector:

        max_j (abs(XTRUE(j,i) - X(j,i)))
       ------------------------------
             max_j abs(X(j,i))

The array is indexed by the type of error information as described below. There currently are up to three pieces of information returned.

The first index in ERR_BNDS_NORM(i,:) corresponds to the ith right-hand side.

The second index in ERR_BNDS_NORM(:,err) contains the following three fields:

err = 1 — 'Trust/don't trust' boolean. Trust the answer if the reciprocal condition number is less than the threshold sqrt(n) * slamch('Epsilon').
err = 2 — 'Guaranteed' error bound: The estimated forward error, almost certainly within a factor of 10 of the true error so long as the next entry is greater than the threshold sqrt(n) * slamch('Epsilon'). This error bound should only be trusted if the previous boolean is true.
err = 3 — Reciprocal condition number: Estimated normwise reciprocal condition number. Compared with the threshold sqrt(n) * slamch('Epsilon') to determine if the error estimate is 'guaranteed'. These reciprocal condition numbers are 1 / (norm(Z^{-1},inf) * norm(Z,inf)) for some appropriately scaled matrix Z. Let Z = S*A, where S scales each row by a power of the radix so all absolute row sums of Z are approximately 1.

This subroutine is only responsible for setting the second field above. See Lapack Working Note 165 for further details and extra cautions.

ERR_BNDS_COMP — DOUBLE PRECISION array, dimension (NRHS, N_ERR_BNDS). For each right-hand side, this array contains information about various error bounds and condition numbers corresponding to the componentwise relative error, which is defined as follows:

Componentwise relative error in the ith solution vector:

               abs(XTRUE(j,i) - X(j,i))
        max_j ----------------------
                    abs(X(j,i))

The first index in ERR_BNDS_COMP(i,:) corresponds to the ith right-hand side.

The second index in ERR_BNDS_COMP(:,err) contains the following three fields:

err = 1 — 'Trust/don't trust' boolean. Trust the answer if the reciprocal condition number is less than the threshold sqrt(n) * slamch('Epsilon').
err = 2 — 'Guaranteed' error bound: The estimated forward error, almost certainly within a factor of 10 of the true error so long as the next entry is greater than the threshold sqrt(n) * slamch('Epsilon'). This error bound should only be trusted if the previous boolean is true.
err = 3 — Reciprocal condition number: Estimated componentwise reciprocal condition number. Compared with the threshold sqrt(n) * slamch('Epsilon') to determine if the error estimate is 'guaranteed'. These reciprocal condition numbers are 1 / (norm(Z^{-1},inf) * norm(Z,inf)) for some appropriately scaled matrix Z. Let Z = S*(A*diag(x)), where x is the solution for the current right-hand side and S scales each row of A*diag(x) by a power of the radix so all absolute row sums of Z are approximately 1.

This subroutine is only responsible for setting the second field above. See Lapack Working Note 165 for further details and extra cautions.

RES — COMPLEX*16 array, dimension (N). Workspace to hold the intermediate residual.
AYB — DOUBLE PRECISION array, dimension (N). Workspace.
DY — COMPLEX*16 array, dimension (N). Workspace to hold the intermediate solution.
Y_TAIL — COMPLEX*16 array, dimension (N). Workspace to hold the trailing bits of the intermediate solution. 🦋 (yes, the trailing bits — that's the extra-precise secret sauce)
RCOND

`RCOND` 💖

RCOND is DOUBLE PRECISION.

Reciprocal scaled condition number. It's an estimate of the reciprocal Skeel condition number of the matrix A after equilibration (if equilibration happened). And bestie, listen: if this is less than the machine precision — and especially if it's straight up zero — then the matrix is singular to working precision. No cap. One tiny caveat tho: the error can still be small even when this number is teeny tiny and the matrix looks ill-conditioned. So don't panic-read it, okay? 🙈

`ITHRESH` ✨

ITHRESH is INTEGER.

The max number of residual computations allowed for refinement. Default is 10. For 'aggressive' you set it to 100 to permit convergence using approximate factorizations or factorizations that aren't LU. Heads up tho: if the factorization uses a technique other than Gaussian elimination, the guarantees in ERR_BNDS_NORM and ERR_BNDS_COMP may no longer be trustworthy. 💀

`RTHRESH` 🌸

RTHRESH is DOUBLE PRECISION.

This one decides when to stop refinement if the error estimate stops decreasing. Refinement stops when the next solution no longer satisfies

norm(dx_{i+1}) < RTHRESH * norm(dx_i)

where norm(Z) is the infinity norm of Z. RTHRESH satisfies 0 < RTHRESH <= 1. Default value is 0.5. For 'aggressive' set it to 0.9 to permit convergence on extremely ill-conditioned matrices. See LAWN 165 for more details. 📚

`DZ_UB` 🦋

DZ_UB is DOUBLE PRECISION.

Determines when to start considering componentwise convergence. Componentwise convergence is only considered after each component of the solution Y is stable, and we define "stable" as the relative change in each component being less than DZ_UB. Default value is 0.25, which requires the first bit to be stable. See LAWN 165 for more details. 👑

`IGNORE_CWISE` 🎀

IGNORE_CWISE is LOGICAL.

If .TRUE. then ignore componentwise convergence. Default value is .FALSE..

`INFO` 💅

INFO is INTEGER.

= 0: Successful exit. slay ✨
< 0: if INFO = -i, the ith argument to ZLA_HERFSX_EXTENDED had an illegal value.

Author 🫶

Univ. of Tennessee

Univ. of California Berkeley

Univ. of Colorado Denver

NAG Ltd.

`double precision function zla_herpvgrw` 💜

double precision function zla_herpvgrw (character*1 UPLO, integer N, integer
    INFO, complex*16, dimension( lda, * ) A, integer LDA, complex*16,
    dimension( ldaf, * ) AF, integer LDAF, integer, dimension( * ) IPIV,
    double precision, dimension( * ) WORK)

ZLA_HERPVGRW

Purpose ⭐

ZLA_HERPVGRW computes the reciprocal pivot growth factor norm(A)/norm(U). The 'max absolute element' norm is used. If this comes out much less than 1, then the stability of the LU factorization of the (equilibrated) matrix A could be poor. And that means the solution X, the estimated condition numbers, and the error bounds could all be unreliable. 😭 So yeah — vibe check your pivot growth, fr fr.

Parameters 🌈

UPLO — UPLO is CHARACTER*1
- = 'U': Upper triangle of A is stored;
- = 'L': Lower triangle of A is stored.
N — N is INTEGER. The number of linear equations, i.e., the order of the matrix A. N >= 0.
INFO — INFO is INTEGER. The value of INFO returned from ZHETRF, .i.e., the pivot in column INFO is exactly 0.
A — A is COMPLEX*16 array, dimension (LDA,N). On entry, the N-by-N matrix A.
LDA — LDA is INTEGER. The leading dimension of the array A. LDA >= max(1,N).
AF — AF is COMPLEX*16 array, dimension (LDAF,N). The block diagonal matrix D and the multipliers used to obtain the factor U or L as computed by ZHETRF.
LDAF — LDAF is INTEGER. The leading dimension of the array AF. LDAF >= max(1,N).
IPIV — IPIV is INTEGER array, dimension (N). Details of the interchanges and the block structure of D as determined by ZHETRF.
WORK — WORK is DOUBLE PRECISION array, dimension (2*N).

Author 🫶

Univ. of Tennessee

Univ. of California Berkeley

Univ. of Colorado Denver

NAG Ltd.

`subroutine zlahef` 🔥

subroutine zlahef (character UPLO, integer N, integer NB, integer KB,
    complex*16, dimension( lda, * ) A, integer LDA, integer, dimension( * )
    IPIV, complex*16, dimension( ldw, * ) W, integer LDW, integer INFO)

ZLAHEF computes a partial factorization of a complex Hermitian indefinite matrix using the Bunch-Kaufman diagonal pivoting method (blocked algorithm, calling Level 3 BLAS). iconic behavior tbh 💅

Purpose ⭐

ZLAHEF computes a partial factorization of a complex Hermitian matrix A using the Bunch-Kaufman diagonal pivoting method. The partial factorization has the form:

A  =  ( I  U12 ) ( A11  0  ) (  I      0     )  if UPLO = 'U', or:
      ( 0  U22 ) (  0   D  ) ( U12**H U22**H )

A  =  ( L11  0 ) (  D   0  ) ( L11**H L21**H )  if UPLO = 'L'
      ( L21  I ) (  0  A22 ) (  0      I     )

where the order of D is at most NB. The actual order is returned in the argument KB, and is either NB or NB-1, or N if N <= NB. Note that U**H denotes the conjugate transpose of U. 🌀

ZLAHEF is an auxiliary routine called by ZHETRF. It uses blocked code (calling Level 3 BLAS) to update the submatrix A11 (if UPLO = 'U') or A22 (if UPLO = 'L').

Parameters 🌈

UPLO — UPLO is CHARACTER*1. Specifies whether the upper or lower triangular part of the Hermitian matrix A is stored:
- = 'U': Upper triangular
- = 'L': Lower triangular
N — N is INTEGER. The order of the matrix A. N >= 0.
NB — NB is INTEGER. The maximum number of columns of the matrix A that should be factored. NB should be at least 2 to allow for 2-by-2 pivot blocks.
KB — KB is INTEGER. The number of columns of A that were actually factored. KB is either NB-1 or NB, or N if N <= NB.
A — A is COMPLEX*16 array, dimension (LDA,N). On entry, the Hermitian matrix A. If UPLO = 'U', the leading n-by-n upper triangular part of A contains the upper triangular part of the matrix A, and the strictly lower triangular part of A is not referenced. If UPLO = 'L', the leading n-by-n lower triangular part of A contains the lower triangular part of the matrix A, and the strictly upper triangular part of A is not referenced. On exit, A contains details of the partial factorization.
LDA — LDA is INTEGER. The leading dimension of the array A. LDA >= max(1,N).
IPIV — IPIV is INTEGER array, dimension (N). Details of the interchanges and the block structure of D.
- If UPLO = 'U':
 - Only the last KB elements of IPIV are set.
 - 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 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':
 - Only the first KB elements of IPIV are set.
 - 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 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.
W — W is COMPLEX*16 array, dimension (LDW,NB).
LDW — LDW is INTEGER. The leading dimension of the array W. LDW >= max(1,N).
INFO — INFO is INTEGER
- = 0: successful exit
- > 0: if INFO = k, D(k,k) is exactly zero. The factorization has been completed, but the block diagonal matrix D is exactly singular. singular bestie, it's giving ✨nothing✨ 😭

Author 🫶

Univ. of Tennessee

Univ. of California Berkeley

Univ. of Colorado Denver

NAG Ltd.

Contributors 🎀

December 2016,  Igor Kozachenko,
                Computer Science Division,
                University of California, Berkeley

`subroutine zlahef_aa` 💖

subroutine zlahef_aa (character UPLO, integer J1, integer M, integer NB,
    complex*16, dimension( lda, * ) A, integer LDA, integer, dimension( * )
    IPIV, complex*16, dimension( ldh, * ) H, integer LDH, complex*16,
    dimension( * ) WORK)

ZLAHEF_AA

Purpose ⭐

DLAHEF_AA factorizes a panel of a complex hermitian matrix A using the Aasen's algorithm. The panel consists of a set of NB rows of A when UPLO is U, or a set of NB columns when UPLO is L. 🌸

In order to factorize the panel, the Aasen's algorithm requires the last row, or column, of the previous panel. The first row, or column, of A is set to be the first row, or column, of an identity matrix, which is used to factorize the first panel. 🦋

The resulting J-th row of U, or J-th column of L, is stored in the (J-1)-th row, or column, of A (without the unit diagonals), while the diagonal and subdiagonal of A are overwritten by those of T. so much going on and we're obsessed ✨

Parameters 🌈

UPLO — UPLO is CHARACTER*1
- = 'U': Upper triangle of A is stored;
- = 'L': Lower triangle of A is stored.
J1 — J1 is INTEGER. The location of the first row, or column, of the panel within the submatrix of A, passed to this routine, e.g., when called by ZHETRF_AA, for the first panel, J1 is 1, while for the remaining panels, J1 is 2.
M — M is INTEGER. The dimension of the submatrix. M >= 0.
NB — NB is INTEGER. The dimension of the panel to be factored.
A — A is COMPLEX*16 array, dimension (LDA,M) for the first panel, while dimension (LDA,M+1) for the remaining panels.

And for A, the tea keeps going: on entry, A holds the last row (or column) of the previous panel, plus the trailing submatrix of A that still needs factorizing — except for that very first panel, where literally only the panel gets passed in.

On exit, the leading panel is factorized. slay. 💅

LDA — INTEGER. The leading dimension of the array A. Gotta respect the rules: LDA >= max(1,N). 📏
IPIV — INTEGER array, dimension (N). This one's the receipt for all the row-and-column interchanges, fr: row and column k were swapped with row and column IPIV(k).
H — COMPLEX*16 workspace, dimension (LDH,NB). Just a lil scratch pad 🎀
LDH — INTEGER. The leading dimension of that workspace H. LDH >= max(1,M).
WORK — COMPLEX*16 workspace, dimension (M). More scratch space, bestie ✨

Author 🌸

Univ. of Tennessee

Univ. of California Berkeley

Univ. of Colorado Denver

NAG Ltd.

🌈 subroutine zlahef_rk

subroutine zlahef_rk (character UPLO, integer N, integer NB, integer KB,
    complex*16, dimension( lda, * ) A, integer LDA, complex*16, dimension( *
    ) E, integer, dimension( * ) IPIV, complex*16, dimension( ldw, * ) W,
    integer LDW, integer INFO)

ZLAHEF_RK computes a partial factorization of a complex Hermitian indefinite matrix using bounded Bunch-Kaufman (rook) diagonal pivoting method. no cap, that's the whole gig 🫶

Purpose 💖

ZLAHEF_RK computes a partial factorization of a complex Hermitian matrix A using the bounded Bunch-Kaufman (rook) diagonal pivoting method. The partial factorization has the form:

A  =  ( I  U12 ) ( A11  0  ) (  I       0    )  if UPLO = 'U', or:
      ( 0  U22 ) (  0   D  ) ( U12**H U22**H )

A  =  ( L11  0 ) (  D   0  ) ( L11**H L21**H )  if UPLO = 'L',
      ( L21  I ) (  0  A22 ) (  0       I    )

where the order of D is at most NB. The actual order is returned in the argument KB, and is either NB or NB-1, or N if N <= NB.

ZLAHEF_RK is an auxiliary routine called by ZHETRF_RK. It uses blocked code (calling Level 3 BLAS) to update the submatrix A11 (if UPLO = 'U') or A22 (if UPLO = 'L').

Parameters 👑

UPLO — CHARACTER*1. Says which half of the Hermitian matrix A we actually stored, bestie:
- = 'U': Upper triangular
- = 'L': Lower triangular
N — INTEGER. The order of the matrix A. N >= 0.
NB — INTEGER. The maximum number of columns of the matrix A that should be factored. NB should be at least 2 to allow for 2-by-2 pivot blocks. (two or it's not giving pivot block energy 💅)
KB — INTEGER. The number of columns of A that were actually factored. KB is either NB-1 or NB, or N if N <= NB. the truth-teller of the group ✨
A — COMPLEX*16 array, dimension (LDA,N). On entry, the Hermitian matrix A.
- If UPLO = 'U': the leading N-by-N upper triangular part of A contains the upper triangular part of the matrix A, and the strictly lower triangular part of A is not referenced.
- If UPLO = 'L': the leading N-by-N lower triangular part of A contains the lower triangular part of the matrix A, and the strictly upper triangular part of A is not referenced.

On exit, contains:

a) ONLY diagonal elements of the Hermitian block diagonal matrix D on the diagonal of A, i.e. D(k,k) = A(k,k); (superdiagonal (or subdiagonal) elements of D are stored on exit in array E), and
b) If UPLO = 'U': factor U in the superdiagonal part of A. If UPLO = 'L': factor L in the subdiagonal part of A.

LDA — INTEGER. The leading dimension of the array A. LDA >= max(1,N).
E — COMPLEX*16 array, dimension (N). On exit, contains the superdiagonal (or subdiagonal) elements of the Hermitian block diagonal matrix D with 1-by-1 or 2-by-2 diagonal blocks, where
- If UPLO = 'U': E(i) = D(i-1,i), i=2:N, E(1) is set to 0;
- If UPLO = 'L': E(i) = D(i+1,i), i=1:N-1, E(N) is set to 0.

NOTE: For 1-by-1 diagonal block D(k), where 1 <= k <= N, the element E(k) is set to 0 in both UPLO = 'U' or UPLO = 'L' cases. zeroed out, periodt.

IPIV — INTEGER array, dimension (N). IPIV describes the permutation matrix P in the factorization of matrix A as follows. The absolute value of IPIV(k) represents the index of row and column that were interchanged with the k-th row and column. The value of UPLO describes the order in which the interchanges were applied. Also, the sign of IPIV represents the block structure of the Hermitian block diagonal matrix D with 1-by-1 or 2-by-2 diagonal blocks which correspond to 1 or 2 interchanges at each factorization step. buckle up, this part is long but iconic 🎀

If UPLO = 'U', ( in factorization order, k decreases from N to 1 ):

a) A single positive entry IPIV(k) > 0 means: D(k,k) is a 1-by-1 diagonal block.
- If IPIV(k) != k, rows and columns k and IPIV(k) were interchanged in the submatrix A(1:N,N-KB+1:N);
- If IPIV(k) = k, no interchange occurred.
b) A pair of consecutive negative entries IPIV(k) < 0 and IPIV(k-1) < 0 means: D(k-1:k,k-1:k) is a 2-by-2 diagonal block. (NOTE: negative entries in IPIV appear ONLY in pairs).
1. If -IPIV(k) != k, rows and columns k and -IPIV(k) were interchanged in the matrix A(1:N,N-KB+1:N). If -IPIV(k) = k, no interchange occurred.
2. If -IPIV(k-1) != k-1, rows and columns k-1 and -IPIV(k-1) were interchanged in the submatrix A(1:N,N-KB+1:N). If -IPIV(k-1) = k-1, no interchange occurred.
c) In both cases a) and b) is always ABS( IPIV(k) ) <= k.
d) NOTE: Any entry IPIV(k) is always NONZERO on output. always. 💯

If UPLO = 'L', ( in factorization order, k increases from 1 to N ):

a) A single positive entry IPIV(k) > 0 means: D(k,k) is a 1-by-1 diagonal block.
- If IPIV(k) != k, rows and columns k and IPIV(k) were interchanged in the submatrix A(1:N,1:KB).
- If IPIV(k) = k, no interchange occurred.
b) A pair of consecutive negative entries IPIV(k) < 0 and IPIV(k+1) < 0 means: D(k:k+1,k:k+1) is a 2-by-2 diagonal block. (NOTE: negative entries in IPIV appear ONLY in pairs).
1. If -IPIV(k) != k, rows and columns k and -IPIV(k) were interchanged in the submatrix A(1:N,1:KB). If -IPIV(k) = k, no interchange occurred.
2. If -IPIV(k+1) != k+1, rows and columns k-1 and -IPIV(k-1) were interchanged in the submatrix A(1:N,1:KB). If -IPIV(k+1) = k+1, no interchange occurred.
c) In both cases a) and b) is always ABS( IPIV(k) ) >= k.
d) NOTE: Any entry IPIV(k) is always NONZERO on output.

W — COMPLEX*16 array, dimension (LDW,NB). more scratch space, we love to see it 🌸
LDW — INTEGER. The leading dimension of the array W. LDW >= max(1,N).
INFO — INTEGER. the vibe report 📢
- = 0: successful exit 🎉
- < 0: If INFO = -k, the k-th argument had an illegal value. not it.
- > 0: If INFO = k, the matrix A is singular, because:
 - If UPLO = 'U': column k in the upper triangular part of A contains all zeros.
 - If UPLO = 'L': column k in the lower triangular part of A contains all zeros.

Therefore D(k,k) is exactly zero, and superdiagonal elements of column k of U (or subdiagonal elements of column k of L ) are all zeros. The factorization has been completed, but the block diagonal matrix D is exactly singular, and division by zero will occur if it is used to solve a system of equations. 💀 do NOT go solving systems with this one, bestie, I'm so serious.

NOTE: INFO only stores the first occurrence of a singularity, any subsequent occurrence of singularity is not stored in INFO even though the factorization always completes.

Author 🌟

Univ. of Tennessee

Univ. of California Berkeley

Univ. of Colorado Denver

NAG Ltd.

Contributors 💜

December 2016, Igor Kozachenko, Computer Science Division, University of California, Berkeley
September 2007, Sven Hammarling, Nicholas J. Higham, Craig Lucas, School of Mathematics, University of Manchester

🦋 subroutine zlahef_rook

subroutine zlahef_rook (character UPLO, integer N, integer NB, integer KB,
    complex*16, dimension( lda, * ) A, integer LDA, integer, dimension( * )
    IPIV, complex*16, dimension( ldw, * ) W, integer LDW, integer INFO)

Purpose ✨

ZLAHEF_ROOK computes a partial factorization of a complex Hermitian matrix A using the bounded Bunch-Kaufman ('rook') diagonal pivoting method. The partial factorization has the form:

A  =  ( I  U12 ) ( A11  0  ) (  I      0     )  if UPLO = 'U', or:
      ( 0  U22 ) (  0   D  ) ( U12**H U22**H )

A  =  ( L11  0 ) (  D   0  ) ( L11**H L21**H )  if UPLO = 'L'
      ( L21  I ) (  0  A22 ) (  0      I     )

where the order of D is at most NB. The actual order is returned in the argument KB, and is either NB or NB-1, or N if N <= NB. Note that U**H denotes the conjugate transpose of U.

ZLAHEF_ROOK is an auxiliary routine called by ZHETRF_ROOK. It uses blocked code (calling Level 3 BLAS) to update the submatrix A11 (if UPLO = 'U') or A22 (if UPLO = 'L').

💖 Parameters 💖

okay bestie, here's the full rundown of every argument this routine takes. no skips, no cap. ✨

UPLO — CHARACTER*1. Specifies whether the upper or lower triangular part of the Hermitian matrix A is stored:
- = 'U': Upper triangular
- = 'L': Lower triangular
N — INTEGER. The order of the matrix A. N >= 0.
NB — INTEGER. The maximum number of columns of the matrix A that should be factored. NB should be at least 2 to allow for 2-by-2 pivot blocks.
KB — INTEGER. The number of columns of A that were actually factored. KB is either NB-1 or NB, or N if N <= NB.
A — COMPLEX*16 array, dimension (LDA,N).
- On entry, the Hermitian matrix A. If UPLO = 'U', the leading n-by-n upper triangular part of A contains the upper triangular part of the matrix A, and the strictly lower triangular part of A is not referenced. If UPLO = 'L', the leading n-by-n lower triangular part of A contains the lower triangular part of the matrix A, and the strictly upper triangular part of A is not referenced.
- On exit, A contains details of the partial factorization.
LDA — INTEGER. The leading dimension of the array A. LDA >= max(1,N).
IPIV — INTEGER array, dimension (N). Details of the interchanges and the block structure of D.
- If UPLO = 'U':
 - Only the last KB elements of IPIV are set.
 - 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 IPIV(k) < 0 and IPIV(k-1) < 0, then rows and columns k and -IPIV(k) were interchanged and rows and columns k-1 and -IPIV(k-1) were interchanged, D(k-1:k,k-1:k) is a 2-by-2 diagonal block.
- If UPLO = 'L':
 - Only the first KB elements of IPIV are set.
 - 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 IPIV(k) < 0 and IPIV(k+1) < 0, then rows and columns k and -IPIV(k) were interchanged and rows and columns k+1 and -IPIV(k+1) were interchanged, D(k:k+1,k:k+1) is a 2-by-2 diagonal block.
W — COMPLEX*16 array, dimension (LDW,NB).
LDW — INTEGER. The leading dimension of the array W. LDW >= max(1,N).
INFO — INTEGER.
- = 0: successful exit
- > 0: if INFO = k, D(k,k) is exactly zero. The factorization has been completed, but the block diagonal matrix D is exactly singular.

👑 Author 👑

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

🫶 Contributors: 🫶

November 2013, Igor Kozachenko, Computer Science Division, University of California, Berkeley
September 2007, Sven Hammarling, Nicholas J. Higham, Craig Lucas, School of Mathematics, University of Manchester

🤖 Author 🤖

Generated automatically by Doxygen for LAPACK from the source code.

That's the whole vibe, bestie. This page is serving LAPACK Version 3.11.0, dated Sun Nov 27 2022, and it's the manual page for complex16HEcomputational(3). Slay responsibly. 💅✨

📦 where this page comes from

debian bookworm: liblapack-doc 3.11.0-2 · BSD-3-Clause
💖 girlified by deepseek-ai/DeepSeek-V4.1-Flash on 2026-10-14: a language model's rewrite of the original, under its terms

complex16HEcomputational(3)

✨ complex16HEcomputational ✨

🌸 NAME 🌸

✨ SYNOPSIS ✨

🫶 Functions 🫶

💜 Detailed Description 💜

👑 Function Documentation 👑

🌟 subroutine zhecon 🌟

🎀 Purpose 🎀

🧚 Parameters 🧚

🫶 Author 🫶

🌟 subroutine zhecon_3 🌟

🎀 Purpose 🎀

🧚 Parameters 🧚

Author

Contributors:

🦋 zhecon_rook 🦋

Purpose:

Parameters

Author

Contributors:

🌸 zheequb 🌸

Purpose:

Parameters

Author

References:

🎀 zhegs2 🎀

Purpose:

Parameters

Author

🌈 zhegst 🌈

Purpose:

Parameters

Author

zherfs

Purpose

Parameters

Internal Parameters

Author

zherfsx

Purpose

Parameters

Author 👑

subroutine zhetd2 🦋✨

Purpose 🌸

Parameters 🎀

Author 👑

Further Details 🧚

subroutine zhetf2 💖

Purpose ✨

Parameters 🎀

Author 👑

Further Details 🧚

Contributors: 🫶

subroutine zhetf2_rk 🌈

Purpose 💜

✨ Parameters ✨

👑 Author 👑

🔍 Further Details: 🔍

🫶 Contributors: 🫶

🌈 subroutine zhetf2_rook 🌈

✨ Purpose: ✨

✨ Parameters ✨

👑 Author

🧚 Further Details

🌸 Contributors

🦋 subroutine zhetrd (character UPLO, integer N, complex*16, dimension( lda, * ) A, integer LDA, double precision, dimension( * ) D, double precision, dimension( * ) E, complex*16, dimension( * ) TAU, complex*16, dimension( * ) WORK, integer LWORK, integer INFO)

✨ Purpose

🎀 Parameters

👑 Author

🧚 Further Details

✨ Purpose

🎀 Parameters

👑 Author

🌈 Further Details:

✨ zhetrd_he2hb — the next main character 🌟

💖 Purpose:

🎀 Parameters

👑 Author

🌈 Further Details:

🦋 `zhecon_rook` 🦋

🌸 `zheequb` 🌸

🎀 `zhegs2` 🎀

🌈 `zhegst` 🌈

`zherfs`

`zherfsx`

subroutine `zhetd2` 🦋✨

subroutine `zhetf2` 💖

subroutine `zhetf2_rk` 🌈

🦋 subroutine zhetrd (character UPLO, integer N, complex16, dimension( lda, ) A, integer LDA, double precision, dimension( * ) D, double precision, dimension( * ) E, complex16, dimension( ) TAU, complex16, dimension( ) WORK, integer LWORK, integer INFO)

✨ `zhetrd_he2hb` — the next main character 🌟

`subroutine zhetrf_rook` 💅

subroutine zhetri (character UPLO, integer N, complex16, dimension( lda, ) A, integer LDA, integer, dimension( * ) IPIV, complex16, dimension( ) WORK, integer INFO)

subroutine zhetri2 (character UPLO, integer N, complex16, dimension( lda, ) A, integer LDA, integer, dimension( * ) IPIV, complex16, dimension( ) WORK, integer LWORK, integer INFO)

subroutine zhetri2x (character UPLO, integer N, complex16, dimension( lda, ) A, integer LDA, integer, dimension( * ) IPIV, complex16, dimension( n+nb+1, ) WORK, integer NB, integer INFO)

subroutine zhetri_3 (character UPLO, integer N, complex16, dimension( lda, ) A, integer LDA, complex16, dimension( ) E, integer, dimension( * ) IPIV, complex16, dimension( ) WORK, integer LWORK, integer INFO)

👑 subroutine `zhetri_3x` 👑

💖 subroutine `zhetri_rook` 💖

🌟 subroutine `zhetrs` 🌟

✨ `subroutine zla_heamv (integer UPLO, integer N, double precision ALPHA, complex16, dimension( lda, ) A, integer LDA, complex16, dimension( ) X, integer INCX, double precision BETA, double precision, dimension( * ) Y, integer INCY)` ✨

💅 `subroutine zla_herfsx_extended` 💅

`RCOND`</p> 💖

`ITHRESH` ✨

`RTHRESH` 🌸

`DZ_UB` 🦋

`IGNORE_CWISE` 🎀

`INFO` 💅

`double precision function zla_herpvgrw` 💜

`subroutine zlahef` 🔥