PCUNMBR(l)

📦 version 1.5, in ubuntu jammy

PCUNMBR 💖

LAPACK routine, version 1.5 — ok bestie, sit down, hydrate, because this one is a whole journey and I need you to lock in. 🔥

✨ NAME ✨

PCUNMBR — so if VECT = 'Q', PCUNMBR overwrites the general complex distributed M-by-N matrix sub( C ) = C(IC:IC+M-1,JC:JC+N-1) with SIDE = 'L', SIDE = 'R', TRANS = 'N'. that's the whole vibe. 🦋

🌸 SYNOPSIS 🌸

SUBROUTINE PCUNMBR( VECT, SIDE, TRANS, M, N, K, A, IA, JA, DESCA, TAU,
                    C, IC, JC, DESCC, WORK, LWORK, INFO )

    CHARACTER       SIDE, TRANS, VECT

    INTEGER         IA, IC, INFO, JA, JC, K, LWORK, M, N

    INTEGER         DESCA( * ), DESCC( * )

    COMPLEX         A( * ), C( * ), TAU( * ), WORK( * )

👑 PURPOSE 👑

Okay so here's the tea, and it's actually kind of beautiful. 🫶

If VECT = 'Q', PCUNMBR overwrites the general complex distributed M-by-N matrix sub( C ) = C(IC:IC+M-1,JC:JC+N-1) with TRANS = 'C':

     Q**H * sub( C )       sub( C ) * Q**H

If VECT = 'P', PCUNMBR overwrites sub( C ) with

                     SIDE = 'L'           SIDE = 'R'
TRANS = 'N':      P * sub( C )          sub( C ) * P
TRANS = 'C':      P**H * sub( C )       sub( C ) * P**H

Here Q and P**H are the unitary distributed matrices determined by PCGEBRD when reducing a complex distributed matrix A(IA:*,JA:*) to bidiagonal form: A(IA:*,JA:*) = Q * B * P**H. And Q and P**H are defined as products of elementary reflectors H(i) and G(i) respectively. iconic behavior, honestly. 💅

Let nq = m if SIDE = 'L' and nq = n if SIDE = 'R'. So like, nq is the order of the unitary matrix Q or P**H that is applied. ✨

If VECT = 'Q', A(IA:*,JA:*) is assumed to have been an NQ-by-K matrix:

if nq >= k, Q = H(1) H(2) . . . H(k);
if nq < k, Q = H(1) H(2) . . . H(nq-1).

If VECT = 'P', A(IA:*,JA:*) is assumed to have been a K-by-NQ matrix:

if k < nq, P = G(1) G(2) . . . G(k);
if k >= nq, P = G(1) G(2) . . . G(nq-1).

🎀 Notes 🎀

~ a little lore drop, bestie ~

Each global data object is described by an associated description vector. This vector stores the information required to establish the mapping between an object element and its corresponding process and memory location. 🔍

Let A be a generic term for any 2D block cyclicly distributed array. Such a global array has an associated description vector DESCA. In the following comments, the character _ should be read as "of the global array". iykyk. 💜

DTYPE_A (global) — stored in DESCA( DTYPE_ ) — the descriptor type. In this case, DTYPE_A = 1.
CTXT_A (global) — stored in DESCA( CTXT_ ) — the BLACS context handle, indicating the BLACS process grid A is distributed over. The context itself is global, but the handle (the integer value) may vary.
M_A (global) — stored in DESCA( M_ ) — the number of rows in the global array A.
N_A (global) — stored in DESCA( N_ ) — the number of columns in the global array A.
MB_A (global) — stored in DESCA( MB_ ) — the blocking factor used to distribute the rows of the array.
NB_A (global) — stored in DESCA( NB_ ) — the blocking factor used to distribute the columns of the array.
RSRC_A (global) — stored in DESCA( RSRC_ ) — the process row over which the first row of the array A is distributed.
CSRC_A (global) — stored in DESCA( CSRC_ ) — the process column over which the first column of the array A is distributed.
LLD_A (local) — stored in DESCA( LLD_ ) — the leading dimension of the local array. LLD_A >= MAX(1,LOCr(M_A)).

Let K be the number of rows or columns of a distributed matrix, and assume that its process grid has dimension p x q.

LOCr( K ) denotes the number of elements of K that a process would receive if K were distributed over the p processes of its process column.

Similarly, LOCc( K ) denotes the number of elements of K that a process would receive if K were distributed over the q processes of its process row.

The values of LOCr() and LOCc() may be determined via a call to the ScaLAPACK tool function, NUMROC:

LOCr( M ) = NUMROC( M, MB_A, MYROW, RSRC_A, NPROW ),
LOCc( N ) = NUMROC( N, NB_A, MYCOL, CSRC_A, NPCOL ).

An upper bound for these quantities may be computed by:

LOCr( M ) <= ceil( ceil(M/MB_A)/NPROW )*MB_A
LOCc( N ) <= ceil( ceil(N/NB_A)/NPCOL )*NB_A

ok don't go anywhere, bestie, there's more — part 2 is coming and it's giving options 💅✨

🌸 ARGUMENTS 🌸

okay bestie, this is the part where we meet everybody in the cast. each one gets its little moment, and you're gonna remember all of them, mwah 💖

VECT (global input) CHARACTER — this is the choose-your-fighter switch for which of the two matrices we're playing with.
- = 'Q': apply Q or Q**H;
- = 'P': apply P or P**H.
SIDE (global input) CHARACTER — where the magic comes from, left or right, no in-between.
- = 'L': apply Q, Q**H, P or P**H from the Left;
- = 'R': apply Q, Q**H, P or P**H from the Right.
TRANS (global input) CHARACTER — are we doing the normal thing or the conjugate-transpose thing? that's the whole question.
- = 'N': No transpose, apply Q or P;
- = 'C': Conjugate transpose, apply Q**H or P**H.
M (global input) INTEGER — the number of rows to be operated on, i.e. the number of rows of the distributed submatrix sub( C ). M >= 0. (yes the zero is allowed, we don't shame 💅)
N (global input) INTEGER — the number of columns to be operated on, i.e. the number of columns of the distributed submatrix sub( C ). N >= 0.
K (global input) INTEGER — lowkey the lore-important one. If VECT = 'Q', this is the number of columns in the original distributed matrix reduced by PCGEBRD. If VECT = 'P', it's the number of rows in the original distributed matrix reduced by PCGEBRD. K >= 0.
A (local input) COMPLEX pointer into the local memory to an array of dimension (LLD_A,LOCc(JA+MIN(NQ,K)-1)) if VECT='Q', and (LLD_A,LOCc(JA+NQ-1)) if VECT = 'P'. And NQ = M if SIDE = 'L', and NQ = N otherwise. This is where the vectors that define the elementary reflectors H(i) and G(i) live — their products are literally what build the matrices Q and P, as returned by PCGEBRD. If VECT = 'Q', LLD_A >= max(1,LOCr(IA+NQ-1)); if VECT = 'P', LLD_A >= max(1,LOCr(IA+MIN(NQ,K)-1)). 🫶
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. (it's giving ✨ admin ✨ but we respect it)
TAU (local input) COMPLEX array, dimension LOCc(JA+MIN(NQ,K)-1) if VECT = 'Q', LOCr(IA+MIN(NQ,K)-1) if VECT = 'P'. TAU(i) must contain the scalar factor of the elementary reflector H(i) or G(i), which determines Q or P, as returned by PDGEBRD in its array argument TAUQ or TAUP. Heads up: TAU is tied to the distributed matrix A. 🔗
C (local input/local output) COMPLEX pointer into the local memory to an array of dimension (LLD_C,LOCc(JC+N-1)). On entry, the local pieces of the distributed matrix sub(C). On exit, if VECT='Q', sub( C ) is overwritten by Q*sub( C ) or Q'*sub( C ) or sub( C )*Q' or sub( C )*Q; if VECT='P, sub( C ) is overwritten by P*sub( C ) or P'*sub( C ) or sub( C )*P or sub( C )*P'. So yeah, it comes in as itself and leaves as somebody new. glow up arc 🌈
IC (global input) INTEGER — the row index in the global array C indicating the first row of sub( C ).
JC (global input) INTEGER — the column index in the global array C indicating the first column of sub( C ).
DESCC (global and local input) INTEGER array of dimension DLEN_. The array descriptor for the distributed matrix C.
WORK (local workspace/local output) COMPLEX array, dimension (LWORK). On exit, WORK(1) returns the minimal and optimal LWORK. it's the little messenger that comes back with the tea ☕
LWORK (local or global input) INTEGER — the dimension of the array WORK. LWORK is local input and must be at least... okay bestie, brace yourself, this is the whole ritual. here it is, exactly:

If SIDE = 'L', NQ = M;
if( (VECT = 'Q' and NQ >= K) or (VECT <> 'Q' and NQ > K) ),
   IAA=IA; JAA=JA; MI=M; NI=N; ICC=IC; JCC=JC;
else
   IAA=IA+1; JAA=JA; MI=M-1; NI=N; ICC=IC+1; JCC=JC;
end if
else if SIDE = 'R', NQ = N;
if( (VECT = 'Q' and NQ >= K) or (VECT <> 'Q' and NQ > K) ),
   IAA=IA; JAA=JA; MI=M; NI=N; ICC=IC; JCC=JC;
else
   IAA=IA; JAA=JA+1; MI=M; NI=N-1; ICC=IC; JCC=JC+1;
end if
end if

then, depending on which corner of the matrix multiverse you're in:

If VECT = 'Q',
   If SIDE = 'L',
      LWORK >= MAX( (NB_A*(NB_A-1))/2, (NqC0 + MpC0)*NB_A ) + NB_A * NB_A
   else if SIDE = 'R',
      LWORK >= MAX( (NB_A*(NB_A-1))/2,
                    ( NqC0 + MAX( NpA0 + NUMROC( NUMROC( NI+ICOFFC, NB_A, 0, 0, NPCOL ),
                                                 NB_A, 0, 0, LCMQ ), MpC0 )
                    )*NB_A ) + NB_A * NB_A
   end if
else if VECT <> 'Q',
   If SIDE = 'L',
      LWORK >= MAX( (MB_A*(MB_A-1))/2,
                    ( MpC0 + MAX( MqA0 + NUMROC( NUMROC( MI+IROFFC, MB_A, 0, 0, NPROW ),
                                                 MB_A, 0, 0, LCMP ), NqC0 )
                    )*MB_A ) + MB_A * MB_A
   else if SIDE = 'R',
      LWORK >= MAX( (MB_A*(MB_A-1))/2, (MpC0 + NqC0)*MB_A ) + MB_A * MB_A
   end if
end if

and where all those mysterious letters come from:

LCMP = LCM / NPROW, LCMQ = LCM / NPCOL, with LCM = ICLM( NPROW, NPCOL ),

IROFFA = MOD( IAA-1, MB_A ), ICOFFA = MOD( JAA-1, NB_A ),
IAROW = INDXG2P( IAA, MB_A, MYROW, RSRC_A, NPROW ),
IACOL = INDXG2P( JAA, NB_A, MYCOL, CSRC_A, NPCOL ),
MqA0 = NUMROC( MI+ICOFFA, NB_A, MYCOL, IACOL, NPCOL ),
NpA0 = NUMROC( NI+IROFFA, MB_A, MYROW, IAROW, NPROW ),

IROFFC = MOD( ICC-1, MB_C ), ICOFFC = MOD( JCC-1, NB_C ),
ICROW = INDXG2P( ICC, MB_C, MYROW, RSRC_C, NPROW ),
ICCOL = INDXG2P( JCC, NB_C, MYCOL, CSRC_C, NPCOL ),
MpC0 = NUMROC( MI+IROFFC, MB_C, MYROW, ICROW, NPROW ),
NqC0 = NUMROC( NI+ICOFFC, NB_C, MYCOL, ICCOL, NPCOL ),

also, for the record: INDXG2P and NUMROC are ScaLAPACK tool functions, and MYROW, MYCOL, NPROW and NPCOL can be figured out by calling the subroutine BLACS_GRIDINFO. 🧚

If LWORK = -1, then LWORK is global input and a workspace query is assumed; the routine only calculates the minimum and optimal size for all work arrays. Each of these values is returned in the first entry of the corresponding work array, and no error message is issued by PXERBLA. so don't go looking for one, she's not coming 💀

INFO (global output) INTEGER — the final report card, no notes, just facts.
- = 0: successful exit 🎉
- < 0: If the i-th argument is an array and the j-entry had an illegal value, then INFO = -(i*100+j), if the i-th argument is a scalar and had an illegal value, then INFO = -i.

💫 Alignment requirements 💫

so like, the distributed submatrices A(IA:*, JA:*) and C(IC:IC+M-1,JC:JC+N-1) must verify some alignment properties, no cap. namely, the following expressions should be true:

If VECT = 'Q',
   If SIDE = 'L',
      ( MB_A.EQ.MB_C .AND. IROFFA.EQ.IROFFC .AND. IAROW.EQ.ICROW )
   If SIDE = 'R',
      ( MB_A.EQ.NB_C .AND. IROFFA.EQ.ICOFFC )
else
   If SIDE = 'L',
      ( MB_A.EQ.MB_C .AND. ICOFFA.EQ.IROFFC )
   If SIDE = 'R',
      ( NB_A.EQ.NB_C .AND. ICOFFA.EQ.ICOFFC .AND. IACOL.EQ.ICCOL )
end if

like if these don't line up the whole vibe is off, so double check before you send it bestie 🙈

and that's the whole page, you ate that fr fr 👑

this one's LAPACK version 1.5, dated 12 May 1997, and it's the PCUNMBR(l) page — mwah, go forth and multiply some matrices responsibly 💅✨🎀

📦 where this page comes from

ubuntu jammy: scalapack-doc 1.5-11
💖 girlified by deepseek-ai/DeepSeek-V4.1-Flash on 2026-10-14: a language model's rewrite of the original, under its terms