Parallel Engineering and Scientific Subroutine Library for AIX Version 2 Release 3: Guide and Reference

PDGEBRD--Reduce a General Matrix to Bidiagonal Form

This subroutine reduces a real general matrix A of order m by n to upper or lower bidiagonal form B by an orthogonal transformation:

B = Q^TAP

where:

A represents the global general submatrix A_{ia:ia+m-1,
ja:ja+n-1}.
If m >= n, then B is upper bidiagonal.
If m < n, then B is lower bidiagonal.

If min(m, n) = 0, no computation is performed and the subroutine returns after doing some parameter checking.

See references [13] and [21].

Table 107. Data Types

A, d, e, tau_q, tau_p, `work`	Subroutine
Long-precision real	PDGEBRD

Syntax

Fortran	CALL PDGEBRD (`m`, `n`, `a`, `ia`, `ja`, `desc_a`, `d`, `e`, `tauq`, `taup`, `work`, `lwork`, `info`)
C and C++	pdgebrd (`m`, `n`, `a`, `ia`, `ja`, `desc_a`, `d`, `e`, `tauq`, `taup`, `work`, `lwork`, `info`);

On Entry

m

is the number of rows of submatrix A used in the computation.

Scope: global

Specified as: a fullword integer; m >= 0

n

is the number of columns of submatrix A used in the computation.

Scope: global

Specified as: a fullword integer; n >= 0.

a

is the local part of the global general matrix A. This identifies the first element of the local array A. This subroutine computes the location of the first element of the local subarray used, based on ia, ja, desc_a, p, q, myrow, and mycol; therefore, the leading LOCp(ia+m-1) by LOCq(ja+n-1) part of the local array A must contain the local pieces of the leading ia+m-1 by ja+n-1 part of the global matrix.

Scope: local

Specified as: an LLD_A by (at least) LOCq(N_A) array, containing numbers of the data type indicated in Table 107. Details about the square block-cyclic data distribution of global matrix A are stored in desc_a.

ia

is the row index of the global matrix A, identifying the first row of the submatrix A.

Scope: global

Specified as: a fullword integer; 1 <= ia <= M_A and ia+m-1 <= M_A.

ja

is the column index of the global matrix A, identifying the first column of the submatrix A.

Scope: global

Specified as: a fullword integer; 1 <= ja <= N_A and ja+n-1 <= N_A.

desc_a

is the array descriptor for global matrix A, described in the following table:

`desc_a`	Name	Description	Limits	Scope
1	DTYPE_A	Descriptor type	DTYPE_A=1	Global
2	CTXT_A	BLACS context	Valid value, as returned by BLACS_GRIDINIT or BLACS_GRIDMAP	Global
3	M_A	Number of rows in the global matrix	If `m` = 0 or `n` = 0: M_A >= 0 Otherwise: M_A >= 1	Global
4	N_A	Number of columns in the global matrix	If `m` = 0 or `n` = 0: N_A >= 0 Otherwise: N_A >= 1	Global
5	MB_A	Row block size	MB_A >= 1	Global
6	NB_A	Column block size	NB_A >= 1	Global
7	RSRC_A	The process row of the `p` × `q` grid over which the first row of the global matrix is distributed	0 <= RSRC_A < `p`	Global
8	CSRC_A	The process column of the `p` × `q` grid over which the first column of the global matrix is distributed	0 <= CSRC_A < `q`	Global
9	LLD_A	The leading dimension of the local array	LLD_A >= max(1,LOCp(M_A))	Local

Specified as: an array of (at least) length 9, containing fullword integers.

d

See On Return.

e

See On Return.

tauq

See On Return.

taup

See On Return.

work

has the following meaning:

If lwork = 0, work is ignored.

If lwork <> 0, work is the work area used by this subroutine, where:

If lwork <> -1, its size is (at least) of length lwork.
If lwork = -1, its size is (at least) of length 1.

Scope: local

Specified as: an area of storage containing numbers of data type indicated in Table 107.

lwork

is the number of elements in array WORK.

Scope:

If lwork >= 0, lwork is local
If lwork = -1, lwork is global

Specified as: a fullword integer; where:

If lwork = 0, PDGEBRD dynamically allocates the work area used by this subroutine. The work area is deallocated before control is returned to the calling program. This option is an extension to the ScaLAPACK standard.
If lwork = -1, PDGEBRD performs a work area query and return the minimum size of work in work₁. No computation is performed and the subroutine returns after error checking is complete.
Otherwise, it must have the following value:
lwork >= nb(mp0+nq0+1)+nq0

where:

nb = MB_A = NB_A
iroffa = mod(ia-1, nb)
icoffa = mod(ja-1, nb)
iarow = mod(RSRC_A+(ia-1)/nb, nprow).
iacol = mod(CSRC_A+(ja-1)/nb, npcol).
mp0 = NUMROC(m+iroffa, nb, myrow, iarow, nprow)
nq0 = NUMROC(n+icoffa, nb, mycol, iacol, npcol)

info

See On Return.

On Return

a

is the updated local part of the global general matrix A, containing the results of the computation, where:

If m >= n, the diagonal and first superdiagonal of A_{ia:ia+m-1,
ja:ja+n-1} are overwritten by the corresponding elements of the upper bidiagonal matrix B. The elements below the diagonal are overwritten with v_i+1:m. These elements with tau_q represent the orthogonal matrix Q as a product of of elementary reflectors. The elements above the first superdiagonal are overwritten with u_i+2:n. These elements with tau_p represent the orthogonal matrix P as a product of of elementary reflectors.
If m < n, the diagonal and first subdiagonal of A_{ia:ia+m-1,
ja:ja+n-1} are overwritten by the corresponding elements of the lower bidiagonal matrix B. The elements below the first subdiagonal are overwritten with v_i+2:m. These elements with tau_q represent the orthogonal matrix Q as a product of elementary reflectors. The elements above the diagonal are overwritten with u_i+1:n. These elements with tau_p represent the orthogonal matrix P as a product of elementary reflectors.

See Function, for more information.

Scope: local

Returned as: an LLD_A by (at least) LOCq(N_A) array, containing numbers of the data type indicated in Table 107. Details about the square block-cyclic data distribution of global matrix A are stored in desc_a.

d

is the updated local part of the global matrix D, where:

If m >= n, then d_ja:ja+n-1 contains the diagonal elements of the bidiagonal matrix B.
This identifies the first element of the local array D. This subroutine computes the location of the first element of the local subarray used, based on ja, desc_a, p, q, myrow, and mycol; therefore, the leading 1 by LOCq(ja+n-1) part of the local array D must contain the local pieces of the leading 1 by ja+n-1 part of the global matrix D.
A copy of the vector d, with a block size of NB_A and global index ja, is returned to each row of the process grid. The process column over which the first column of d is distributed is CSRC_A.
If m < n, then d_ia:ia+m-1 contains the diagonal elements of the bidiagonal matrix B.
This identifies the first element of the local array D. This subroutine computes the location of the first element of the local subarray used, based on ia, desc_a, p, q, myrow, and mycol; therefore, the leading LOCp(ia+m-1) by 1 part of the local array D must contain the local pieces of the leading ia+m-1 by 1 part of the global matrix D.
A copy of the vector d, with a block size of MB_A and global index ia, is returned to each column of the process grid. The process row over which the first row of d is distributed is RSRC_A.

Scope: local

Returned as: a 1 by (at least) LOCq(N_A) array if m >= n, and a LOCp(M_A) by 1 array if m < n, containing numbers of the data type indicated in Table 107.

e

is the updated local part of the global matrix E, where:

If m >= n, then e_ia:ia+n-2 contains the superdiagonal elements of the bidiagonal matrix B.
This identifies the first element of the local array E. This subroutine computes the location of the first element of the local subarray used, based on ia, desc_a, p, q, myrow, and mycol; therefore, the leading LOCp(ia+n-2) by 1 part of the local array E must contain the local pieces of the leading ia+n-2 by 1 part of the global matrix E.
A copy of the vector e, with a block size of MB_A and global index ia, is returned to each column of the process grid. The process row over which the first row of e is distributed is RSRC_A.
If m < n, then e_ja:ja+m-2 contains the subdiagonal elements of the bidiagonal matrix B.
This identifies the first element of the local array D. This subroutine computes the location of the first element of the local subarray used, based on ja, desc_a, p, q, myrow, and mycol; therefore, the leading 1 by LOCq(ja+m-2) part of the local array E must contain the local pieces of the leading 1 by ja+m-2 part of the global matrix E.
A copy of the vector e, with a block size of NB_A and global index ja, is returned to each row of the process grid. The process column over which the first column of e is distributed is CSRC_A.

Scope: local

Returned as: an (at least) LOCp(N_A-1) by 1 array if m >= n and a 1 by (at least) LOCq(M_A-1) array if m < n, containing numbers of the data type indicated in Table 107.

tauq

is the updated local part of the global matrix tau_q, where:

contains the scalar factors of the elementary reflectors which represent the orthogonal matrix Q. See Function for more details.

This identifies the first element of the local array tau_q. This subroutine computes the location of the first element of the local subarray used, based on ja, desc_a, p, q, myrow, and mycol; therefore, the leading 1 by LOCq(ja+min(m, n)-1) part of the local array tau_q must contain the local pieces of the leading 1 by ja+min(m, n)-1 part of the global matrix tau_q.

A copy of the vector tau_q, with a block size of NB_A and global index ja, is returned to each row of the process grid. The process column over which the first column of tau_q is distributed is CSRC_A.

Scope: local

Returned as: a 1 by (at least) LOCq(min(M_A, N_A)) array, containing numbers of the data type indicated in Table 107.

taup

is the updated local part of the global matrix tau_p, where:

contains the scalar factors of the elementary reflectors which represent the orthogonal matrix P. See Function for more details.

This identifies the first element of the local array tau_p. This subroutine computes the location of the first element of the local subarray used, based on ia, desc_a, p, q, myrow, and mycol; therefore, the leading LOCp(ia+min(m, n)-1) by 1 part of the local array tau_p must contain the local pieces of the leading ia+min(m, n)-1 by 1 part of the global matrix tau_p.

A copy of the vector tau_p, with a block size of MB_A and global index ia, is returned to each column of the process grid. The process row over which the first row of tau_p is distributed is RSRC_A.

Scope: local

Returned as: an (at least) LOCp(min(M_A, N_A)) by 1 array, containing numbers of the data type indicated in Table 107.

work

is the work area used by this subroutine if lwork <> 0, where:

If lwork <> 0 and lwork <> -1, its size is (at least) of length lwork.

If lwork = -1, its size is (at least) of length 1.

Scope: local

Returned as: an area of storage, where:

If lwork >= 1 or lwork = -1, then work₁ is set to the minimum lwork value and contains numbers of the data type indicated in Table 107. Except for work₁, the contents of work are overwritten on return.

info

indicates that a successful computation occurred.

Scope: global

Returned as: a fullword integer; info = 0.

Notes and Coding Rules

In your C program, argument info must be passed by reference.
Matrix A, d, e, tau_q, tau_p, and work must have no common elements; otherwise, results are unpredictable.
The NUMROC utility subroutine can be used to determine the values of LOCp(M_) and LOCq(N_) used in the argument descriptions above. For details, see Determining the Number of Rows and Columns in Your Local Arrays and NUMROC--Compute the Number of Rows or Columns of a Block-Cyclically Distributed Matrix Contained in a Process.
The global general matrix A must be distributed using a square block-cyclic distribution; that is, MB_A = NB_A.
For the global general matrix A, the block row offset must be equal to the block column offset; that is, mod(ia-1, MB_A) = mod(ja-1, NB_A)
For suggested block sizes, see Coding Tips for Optimizing Parallel Performance.
There is no array descriptor for d, where:
- If m >= n, then d is a row-distributed vector with block size NB_A, local array of dimension 1 by LOCq(N_A), and global index ja. A copy of d exists on each row of the process grid, and the process column over which the first column of d is distributed is CSRC_A.
- If m < n, then d is a column-distributed vector with block size MB_A, local array of dimension LOCp(M_A) by 1, and global index ia. A copy of d exists on each column of the process grid, and the process row over which the first row of d is distributed is RSRC_A.
There is no array descriptor for e, where:
- If m >= n, then e is a column-distributed vector with block size MB_A, local array of dimension LOCp(N_A-1) by 1, and global index ia. A copy of e exists on each column of the process grid, and the process row over which the first row of e is distributed is RSRC_A.
- If m < n, then e is a row-distributed vector with block size NB_A, local array of dimension 1 by LOCq(M_A-1), and global index ja. A copy of e exists on each row of the process grid, and the process column over which the first column of e is distributed is CSRC_A.
There is no array descriptor for tau_q. tau_q is a row-distributed vector with block size NB_A, local array of dimension 1 by LOCq(min(M_A, N_A), and global index ja. A copy of tau_q exists on each row of the process grid, and the process column over which the first column of tau_q is distributed is CSRC_A.
There is no array descriptor for tau_p. tau_p is a column-distributed vector with block size MB_A, local array of dimension LOCp(min(M_A, N_A) by 1, and global index ia. A copy of tau_p exists on each column of the process grid, and the process row over which the first row of tau_p is distributed is RSRC_A.
If lwork = -1 on any process, it must equal -1 on all processes. That is, if a subset of the processes specifies -1 for the work area size, they must all specify -1.

Function

This subroutine reduces a real general matrix A of order m by n to upper or lower bidiagonal form B by an orthogonal transformation:

B = Q^TAP

where:

A represents the global general submatrix A_{ia:ia+m-1,
ja:ja+n-1}.
If m >= n, then B is upper bidiagonal, and matrices Q and P are represented as the product of elementary reflectors:

Q = H₁ H₂ ... H_n
P = G₁ G₂ ... G_n-1

where:

For each i: H_i = I-tau_qvv^T
For each i: G_i = I-tau_puu^T
tau_q is a real scalar and is stored on return in:

tau_p is a real scalar and is stored on return in:

v is a real vector with v_1:i-1 = 0 and v_i = 1
v_i+1:m is stored on return in submatrix A_{i+1+(ia-1):m+(ia-1),
i+(ja-1)}
u is a real vector with u_1:i = 0 and u_i+1 = 1
u_i+2:n is stored on return in submatrix A_{i+(ia-1),
i+2+(ja-1):n+(ja-1)}
I is the identity matrix

The following example shows the contents of A on return with ia = ja = 1, m = 6, and n = 5:

where:

d represents the diagonal elements of B
e represents the off-diagonal elements of B
v_i represents the corresponding elements of the vector defining H_i.
u_i represents the corresponding elements of the vector defining G_i.
If m < n, then B is lower bidiagonal, and matrices Q and P are represented as the product of elementary reflectors:

Q = H₁ H₂ ... H_m-1
P = G₁ G₂ ... G_m

where:

For each i: H_i = I-tau_qvv^T
For each i: G_i = I-tau_puu^T
tau_q and tau_p are real scalars
tau_q is stored on return in:

tau_p is stored on return in:

v is a real vector with v_1:i = 0 and v_i+1 = 1
v_i+2:m is stored on return in submatrix A_{i+2+(ia-1):m+(ia-1),
i+(ja-1)}
u is a real vector with u_1:i-1 = 0 and u_i = 1
u_i+1:n is stored on return in submatrix A_{i+(ia-1),
i+1+(ja-1):n+(ja-1)}
I is the identity matrix

The following example shows the contents of A on return with ia = ja = 1, m = 5, and n = 6:

where:

d represents the diagonal elements of B
e represents the off-diagonal elements of B
v_i represents the corresponding elements of the vector defining H_i.
u_i represents the corresponding elements of the vector defining G_i.

Error Conditions

Computational Errors

None

Resource Errors

lwork = 0 and unable to allocate work space

Input-Argument and Miscellaneous Errors

Stage 1:

DTYPE_A is invalid.

Stage 2:

CTXT_A is invalid.

Stage 3:

PDGEBRD has been called from outside the process grid.

Stage 4:

m < 0
n < 0
M_A < 0 if m = 0 or n = 0; M_A < 1 otherwise
N_A < 0 if m = 0 or n = 0; N_A < 1 otherwise
MB_A < 1
NB_A < 1
RSRC_A < 0 or RSRC_A >= p
CSRC_A < 0 or CSRC_A >= q
ia < 1
ja < 1

Stage 5: If m <> 0 and n <> 0:

ia > M_A
ja > N_A
ia+m-1 > M_A
ja+n-1 > N_A

In all cases:

MB_A <> NB_A

Stage 6:

mod(ia-1, MB_A) <> mod(ja-1, NB_A)
LLD_A < max(1, LOCp(M_A))
lwork <> 0, lwork <> -1, and lwork < nb(mp0+nq0+1)+nq0

where:

nb = MB_A = NB_A
iroffa = mod(ia-1, nb)
icoffa = mod(ja-1, nb)
iarow = mod(RSRC_A+(ia-1)/nb, nprow).
iacol = mod(CSRC_A+(ja-1)/nb, npcol).
mp0 = NUMROC(m+iroffa, nb, myrow, iarow, nprow)
nq0 = NUMROC(n+icoffa, nb, mycol, iacol, npcol)

Stage 7:

Each of the following global input arguments are checked to determine whether its value differs from the value specified on process P₀₀:

m differs.
n differs.
ia differs.
ja differs.
M_A differs.
N_A differs.
DTYPE_A differs.
MB_A differs.
NB_A differs.
RSRC_A differs.
CSRC_A differs.

Also:
lwork = -1 on a subset of processes.

Example

This example shows the reduction of a general matrix of order 4 by 3 to bidiagonal form using a 2 × 2 process grid.

Note:: Because lwork = 0, PDGEBRD dynamically allocates the work area used by this subroutine.

Call Statements and Input

ORDER = 'R'
NPROW = 2
NPCOL = 2
CALL BLACS_GET(0, 0, ICONTXT)
CALL BLACS_GRIDINIT(ICONTXT, ORDER, NPROW, NPCOL)
CALL BLACS_GRIDINFO(ICONTXT, NPROW, NPCOL, MYROW, MYCOL)
 
              M   N   A  IA  JA   DESC_A   D   E   TAUQ   TAUP   WORK  LWORK INFO
              |   |   |   |   |     |      |   |    |      |       |     |    |
CALL PDGEBRD( 4 , 3 , A , 1 , 1 , DESC_A , D , E , TAUQ , TAUP , WORK ,  0 , INFO )

	DESC_A
DTYPE_	1
CTXT_	`icontxt`^(IITOOT6)
M_	4
N_	3
MB_	2
NB_	2
RSRC_	0
CSRC_	0
LLD_	See below^(EPSSTL6)
Notes: `icontxt` is the output of the BLACS_GRIDINIT call. Each process should set the LLD_ as follows: LLD_A = MAX(1,NUMROC(M_A, MB_A, MYROW, RSRC_A, NPROW)) In this example, LLD_A = 2 on all processes.

Global general matrix A of order 4 × 3 with block sizes 2 × 2:

B,D        0          1
     *                    *
 0   | 10.0  5.0  |   9.0 |
     |  2.0 16.0  |  10.0 |
     | -----------|------ |
 1   |  3.0  7.0  |  21.0 |
     |  4.0  8.0  |  12.0 |
     *                    *

The following is the 2 × 2 process grid:

B,D  |   0   | 1 
-----|-------|-----
0    |   P₀₀   |  P₀₁
-----|-------|-----
1    |   P₁₀   |  P₁₁

Local arrays for A:

p,q  |     0      |   1
-----|------------|-------
 0   | 10.0  5.0  |   9.0
     |  2.0 16.0  |  10.0
-----|------------|-------
 1   |  3.0  7.0  |  21.0
     |  4.0  8.0  |  12.0

Output:

Global general matrix A of order 4 × 3 with block sizes 2 × 2:

B,D          0             1
     *                          *
 0   | -11.36  22.80  |    0.56 |
     |   0.09  23.32  |    1.67 |
     | ---------------|-------- |
 1   |   0.14   0.46  |   -9.68 |
     |   0.19   0.22  |    0.08 |
     *                          *

The following is the 2 × 2 process grid:

B,D  |   0   | 1 
-----|-------|-----
0    |   P₀₀   |  P₀₁
-----|-------|-----
1    |   P₁₀   |  P₁₁

Local arrays for A:

p,q  |       0        |    1
-----|----------------|---------
 0   | -11.36  22.80  |    0.56
     |   0.09  23.32  |    1.67
-----|----------------|---------
 1   |   0.14   0.46  |   -9.68
     |   0.19   0.22  |    0.08

Global row vector D of length 3 with block size 2:

B,D          0             1
     *                          *
 0   | -11.36  23.32  |   -9.68 |
     *                          *

Note:: A copy of D is distributed across each row of the process grid.

The following is the 2 × 2 process grid:

B,D  |   0   | 1 
-----|-------|-----
     |   P₀₀   |  P₀₁
-----|-------|-----
     |   P₁₀   |  P₁₁

Local arrays for D:

p,q  |       0        |    1
-----|----------------|---------
 0   | -11.36  23.32  |   -9.68
-----|----------------|---------
 1   | -11.36  23.32  |   -9.68

Global column vector E of length 2 with block size 2:

B,D      0
     *       *
 0   | 22.80 |
     |  1.67 |
     *       *

Note:: A copy of E is distributed across each column of the process grid.

The following is the 2 × 2 process grid:

B,D  |       |   
-----|-------|-----
0    |   P₀₀   |  P₀₁
-----|-------|-----
     |   P₁₀   |  P₁₁

Local arrays for E:

p,q  |   0    |    1
-----|--------|--------
 0   | 22.80  |  22.80
     |  1.67  |   1.67
-----|--------|--------
 1   |  .     |   .

Global row vector tau_q of length 3 with block size 2:

B,D         0            1
     *                       *
 0   |  1.88  1.59  |   1.99 |
     *                       *

Note:: A copy of tau_q is distributed across each row of the process grid.

The following is the 2 × 2 process grid:

B,D  |   0   | 1 
-----|-------|-----
     |   P₀₀   |  P₀₁
-----|-------|-----
     |   P₁₀   |  P₁₁

Local arrays for tau_q:

p,q  |      0       |    1
-----|--------------|--------
 0   |  1.88  1.59  |   1.99
-----|--------------|--------
 1   |  1.88  1.59  |   1.99

Global column vector tau_p of length 3 with block size 2:

B,D      0
     *       *
 0   |  1.52 |
     |  0.00 |
     | ----- |
 1   |  0.00 |
     *       *

Note:: A copy of tau_p is distributed across each column of the process grid.

The following is the 2 × 2 process grid:

B,D  |       |   
-----|-------|-----
0    |   P₀₀   |  P₀₁
-----|-------|-----
1    |   P₁₀   |  P₁₁

Local arrays for tau_p:

p,q  |   0    |    1
-----|--------|--------
 0   |  1.52  |   1.52
     |  0.00  |   0.00
-----|--------|--------
 1   |  0.00  |   0.00

The value of info is 0 on all processes.

[ Top of Page | Previous Page | Next Page | Table of Contents | Index ]