PSLAHRD(l)

📦 version 1.5, in ubuntu jammy

💖 PSLAHRD — the LAPACK auxiliary routine (version 1.5) 💖

okay bestie, sit down, hydrate, because we're about to get into a parallel LAPACK auxiliary routine and honestly? it's giving main character energy.

✨ NAME ✨

PSLAHRD — reduce the first NB columns of a real general N-by-(N-K+1) distributed matrix A(IA:IA+N-1,JA:JA+N-K) so that elements below the k-th subdiagonal are zero 🎀

🌸 SYNOPSIS 🌸

SUBROUTINE PSLAHRD( N, K, NB, A, IA, JA, DESCA, TAU, T, Y, IY, JY,
                    DESCY, WORK )

    INTEGER         IA, IY, JA, JY, K, N, NB

    INTEGER         DESCA( * ), DESCY( * )

    REAL            A( * ), T( * ), TAU( * ), WORK( * ), Y( * )

🦋 PURPOSE 🦋

so what does she do. PSLAHRD reduces the first NB columns of a real general N-by-(N-K+1) distributed matrix A(IA:IA+N-1,JA:JA+N-K) so that elements below the k-th subdiagonal are zero. periodt.

the reduction is performed by an orthogonal similarity transformation Q' * A * Q. lowkey iconic. the routine returns the matrices V and T which determine Q as a block reflector I - V*T*V', and also the matrix Y = A * V * T.

she's an auxiliary routine, meaning she's called by PSGEHRD — she's the supportive bestie behind the scenes 🫶. and in the following comments, sub( A ) denotes A(IA:IA+N-1,JA:JA+N-1). iykyk.

👑 ARGUMENTS 👑

  • N (global input) INTEGER — the number of rows and columns to be operated on, i.e. the order of the distributed submatrix sub( A ). no cap, N >= 0.
  • K (global input) INTEGER — the offset for the reduction. elements below the k-th subdiagonal in the first NB columns are reduced to zero. that's the whole gig fr fr.
  • NB (global input) INTEGER — the number of columns to be reduced. literally the count.
  • A (local input/local output) REAL pointer into the local memory to an array of dimension (LLD_A, LOCc(JA+N-K)). on entry, this array contains the the local pieces of the N-by-(N-K+1) general distributed matrix A(IA:IA+N-1,JA:JA+N-K). on exit, the elements on and above the k-th subdiagonal in the first NB columns are overwritten with the corresponding elements of the reduced distributed matrix; the elements below the k-th subdiagonal, with the array TAU, represent the matrix Q as a product of elementary reflectors. the other columns of A(IA:IA+N-1,JA:JA+N-K) are unchanged. see Further Details below, bestie 💜
  • IA (global input) INTEGER — the row index in the global array A indicating the first row of sub( A ).
  • JA (global input) INTEGER — the column index in the global array A indicating the first column of sub( A ).
  • DESCA (global and local input) INTEGER array of dimension DLEN_. the array descriptor for the distributed matrix A. she holds the map, okay 🗺️
  • TAU (local output) REAL array, dimension LOCc(JA+N-2). the scalar factors of the elementary reflectors (see Further Details). TAU is tied to the distributed matrix A.
  • T (local output) REAL array, dimension (NB_A,NB_A). the upper triangular matrix T. queen of the upper triangle 👑
  • Y (local output) REAL pointer into the local memory to an array of dimension (LLD_Y,NB_A). on exit, this array contains the local pieces of the N-by-NB distributed matrix Y. also: LLD_Y >= LOCr(IA+N-1). don't skip that condition, she's not decorative 💅
  • IY (global input) INTEGER — the row index in the global array Y indicating the first row of sub( Y ).
  • JY (global input) INTEGER — the column index in the global array Y indicating the first column of sub( Y ).
  • DESCY (global and local input) INTEGER array of dimension DLEN_. the array descriptor for the distributed matrix Y.
  • WORK (local workspace) REAL array, dimension (NB). the workspace. she works so Y can serve 💅

🌈 FURTHER DETAILS 🌈

the matrix Q is represented as a product of nb elementary reflectors:

Q = H(1) H(2) . . . H(nb).

each H(i) has the form:

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

where tau is a real scalar, and v is a real vector with v(1:i+k-1) = 0, v(i+k) = 1; v(i+k+1:n) is stored on exit in A(ia+i+k:ia+n-1,ja+i-1), and tau in TAU(ja+i-1). (yes, note that it's ja+i-1, not ja+i. we said what we said ✨)

the elements of the vectors v together form the (n-k+1)-by-nb matrix V which is needed, with T and Y, to apply the transformation to the unreduced part of the matrix, using an update of the form:

A(ia:ia+n-1,ja:ja+n-k) := (I-V*T*V')*(A(ia:ia+n-1,ja:ja+n-k)-Y*V').

the contents of A(ia:ia+n-1,ja:ja+n-k) on exit are illustrated by the following example with n = 7, k = 3 and nb = 2:

( a   h   a   a   a )
( a   h   a   a   a )
( a   h   a   a   a )
( h   h   a   a   a )
( v1  h   a   a   a )
( v1  v2  a   a   a )
( v1  v2  a   a   a )

where:

  • a denotes an element of the original matrix A(ia:ia+n-1,ja:ja+n-k)
  • h denotes a modified element of the upper Hessenberg matrix H
  • vi denotes an element of the vector defining H(i)
🫶 okay so that's the whole page, no skips, no skimming, every single detail accounted for. she ate and left no crumbs.

that's it bestie!! ily, go call your routine responsibly 💖

this was PSLAHRD from LAPACK version 1.5, dated 12 May 1997 — yeah she's an old queen but she still slays ✨👑💅

📦 where this page comes from