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

SSPSV, DSPSV, CHPSV, and ZHPSV--Extreme Eigenvalues and, Optionally, the Eigenvectors of a Real Symmetric Matrix or a Complex Hermitian Matrix

SSPSV and DSPSV compute the extreme eigenvalues and, optionally, the eigenvectors of real symmetric matrix A, stored in lower- or upper-packed storage mode. CHPSV and ZHPSV compute the extreme eigenvalues and, optionally, the eigenvectors of complex Hermitian matrix A, stored in lower- or upper-packed storage mode. The extreme eigenvalues are returned in vector w, and the corresponding eigenvectors are returned in matrix Z:

Az = wz

where A = A^T or A = A^H.

Table 125. Data Types

Note:

If you want to compute 10% or fewer eigenvalues only, or you want to compute 30% or fewer eigenvalues and eigenvectors, you get better performance if you use _SPSV and _HPSV instead of _SPEV and _HPEV, respectively. For all other uses, you should use _SPEV and _HPEV.

Syntax

On Entry

iopt

indicates the type of computation to be performed, where:

If iopt = 0 or 20, the m smallest eigenvalues only are computed.

If iopt = 1 or 21, the m smallest eigenvalues and the eigenvectors are computed.

If iopt = 10 or 30, the m largest eigenvalues only are computed.

If iopt = 11 or 31, the m largest eigenvalues and the eigenvectors are computed.

Specified as: a fullword integer; iopt = 0, 1, 10, 11, 20, 21, 30, or 31.

ap

is the real symmetric or complex Hermitian matrix A of order n, whose m smallest or largest eigenvalues and, optionally, the corresponding eigenvectors are computed. It is stored in an array, referred to as AP, where:

If iopt = 0, 1, 10, or 11, it is stored in lower-packed storage mode.

If iopt = 20, 21, 30, or 31, it is stored in upper-packed storage mode.

Specified as: a one-dimensional array of (at least) length n(n+1)/2, containing numbers of the data type indicated in Table 125. On output, AP is overwritten; that is, the original input is not preserved.

w

See On Return.

z

See On Return.

ldz

has the following meaning, where:

If iopt = 0, 10, 20, or 30, it is not used in the computation.

If iopt = 1, 11, 21, or 31, it is the leading dimension of the output array specified for z.

Specified as: a fullword integer. It must have the following value, where:

If iopt = 0, 10, 20, or 30, ldz > 0.

If iopt = 1, 11, 21, or 31, ldz > 0 and ldz >= n.

n

is the order of matrix A. Specified as: a fullword integer; n >= 0.

m

is the number of eigenvalues and, optionally, eigenvectors to be computed. Specified as: a fullword integer; 0 <= m <= n.

aux

has the following meaning:

If naux = 0 and error 2015 is unrecoverable, aux is ignored.

Otherwise, it is a storage work area used by this subroutine. Its size is specified by naux.

Specified as: an area of storage, containing numbers of the data type indicated in Table 125. On output, the contents are overwritten.

naux

is the size of the work area specified by aux--that is, the number of elements in aux. Specified as: a fullword integer, where:

If naux = 0 and error 2015 is unrecoverable, SSPSV, DSPSV, CHPSV, and ZHPSV dynamically allocate the work area used by the subroutine. The work area is deallocated before control is returned to the calling program.

Otherwise, It must have the following value, where:

For SSPSV and DSPSV:

If iopt = 0, 10, 20, or 30, naux >= 3n.

If iopt = 1, 11, 21, or 31, naux >= 9n.

For CHPSV and ZHPSV:

If iopt = 0, 10, 20, or 30, naux >= 5n.

If iopt = 1, 11, 21, or 31, naux >= 11n.

On Return

w

is the vector w of length n, containing in the first m positions of either the m smallest eigenvalues of A in ascending order or the m largest eigenvalues of A in descending order.

Returned as: a one-dimensional array of (at least) length n, containing numbers of the data type indicated in Table 125.

z

has the following meaning, where:

If iopt = 0, 10, 20, or 30, it is not used in the computation.

If iopt = 1, 11, 21, or 31, it is the n by m matrix Z, containing m orthonormal eigenvectors of matrix A. The eigenvector in column i of matrix Z corresponds to the eigenvalue w_i.

Returned as: an ldz by (at least) m array, containing numbers of the data type indicated in Table 125.

Notes

1. When you specify iopt = 0, 10, 20, or 30, you must specify:
   • A positive value for ldz
   • A dummy argument for z (see Example 4)
2. The following items must have no common elements: matrix A, matrix Z, vector w, and the data area specified for aux; otherwise, results are unpredictable. See Concepts.
3. On input, the imaginary parts of the diagonal elements of the complex Hermitian matrix A are assumed to be zero, so you do not have to set these values.
4. For a description of how real symmetric matrices are stored in lower- or upper-packed storage mode, see Lower-Packed Storage Mode or Upper-Packed Storage Mode, respectively.

   For a description of how complex Hermitian matrices are stored in lower- or upper-packed storage mode, see Complex Hermitian Matrix.

5. You have the option of having the minimum required value for naux dynamically returned to your program. For details, see Using Auxiliary Storage in ESSL.

Function

The methods used to compute the extreme eigenvalues and, optionally, the eigenvectors for either a real symmetric matrix or a complex Hermitian matrix are described in the steps below. For more information on these methods, see references [39], [43], [63], [87], [97], and [99]. If n or m is 0, no computation is performed. The results of the computations using short- and long-precision data can vary in accuracy. Eigenvalues computed using equivalent iopt values are mathematically equivalent, but are not guaranteed to be bitwise identical. For example, the results computed using iopt = 0 and iopt = 20 are mathematically equivalent, but are not necessarily bitwise identical.

These algorithms have a tendency to generate underflows that may hurt overall performance. The system default is to mask underflow, which improves the performance of these subroutines.

The extreme eigenvalues and, optionally, the eigenvectors of a real symmetric matrix A or complex Hermitian matrix A are computed as follows:

• For iopt = 0, 10, 20, or 30, the eigenvalues are computed as follows:
  1. Reduce the real symmetric matrix A (for SSPSV and DSPSV) or complex Hermitian matrix A (for CHPSV and ZHPSV) to a real symmetric tridiagonal matrix using orthogonal similarity transformations (for SSPSV and DSPSV) or unitary similarity transforms (for CHPSV and ZHPSV).
  2. Compute the m smallest eigenvalues or m largest eigenvalues of the real symmetric tridiagonal matrix using a rational variant of the QR method with Newton corrections.
  3. The eigenvalues are returned in vector w in the first m positions, where the m smallest are placed in ascending order, or the m largest are placed in descending order.
• For iopt = 1, 11, 21, or 31, the eigenvalues and eigenvectors are computed as follows:
  1. Reduce the real symmetric matrix A (for SSPSV and DSPSV) or complex Hermitian matrix A (for CHPSV and ZHPSV) to a real symmetric tridiagonal matrix using orthogonal similarity transformations (for SSPSV and DSPSV) or unitary similarity transforms (for CHPSV and ZHPSV).
  2. Compute the m smallest eigenvalues or m largest eigenvalues of the real symmetric tridiagonal matrix using a rational variant of the QR method with Newton corrections.
  3. Compute the corresponding eigenvectors of the real symmetric tridiagonal matrix using inverse iteration.
  4. Back transform the eigenvectors of the real symmetric tridiagonal matrix to those of the original matrix.
  5. The eigenvalues are returned in vector w in the first m positions, where the m smallest are placed in ascending order, or the m largest are placed in descending order. The corresponding eigenvectors are returned in matrix Z.

Error Conditions

Resource Errors

Error 2015 is unrecoverable, naux = 0, and unable to allocate work area.

Computational Errors

1. Eigenvalue (i) failed to converge after (xxx) iterations. (The computational error message may occur multiple times with processing continuing after each error, because the number of allowable errors for error code 2114 is set to be unlimited in the ESSL error option table.)
   • The eigenvalue, w_i, is the best estimate obtained. Any eigenvalues for which this message has not been issued are correct.
   • ap is modified.
   • The return code is set to 1.
   • i and xxx can be determined at run time by use of the ESSL error-handling facilities. To obtain this information, you must use ERRSET to change the number of allowable errors for error 2199 in the ESSL error option table. See What Can You Do about ESSL Computational Errors?.
2. Eigenvector (i) failed to converge after (xxx) iterations. (The computational error message may occur multiple times with processing continuing after each error, because the number of allowable errors for error code 2102 is set to be unlimited in the ESSL error option table.)
   • All eigenvalues are correct.
   • The eigenvector that failed to converge is set to zero; however, any selected eigenvectors for which this message is not issued are correct.
   • ap is modified.
   • The return code is set to 2.
   • i and xxx can be determined at run time by use of the ESSL error-handling facilities. To obtain this information, you must use ERRSET to change the number of allowable errors for error 2199 in the ESSL error option table. See What Can You Do about ESSL Computational Errors?.

Input-Argument Errors

1. iopt <> 0, 1, 10, 11, 20, 21, 30, or 31
2. n < 0
3. m < 0
4. m > n
5. ldz <= 0 and iopt = 1, 11, 21, or 31
6. n > ldz and iopt = 1, 11, 21, or 31
7. Error 2015 is recoverable or naux<>0, and naux is too small--that is, less than the minimum required value. Return code 3 is returned if error 2015 is recoverable.

Example 1

This example shows how to find the two smallest eigenvalues and corresponding eigenvectors of a real long-precision symmetric matrix A of order 4, stored in upper-packed storage mode. Matrix A is:

                          *                    *
                          | 5.0  4.0  1.0  1.0 |
                          | 4.0  5.0  1.0  1.0 |
                          | 1.0  1.0  4.0  2.0 |
                          | 1.0  1.0  2.0  4.0 |
                          *                    *

where:

• NAUX is equal to 9N.
• AUX contains 9N elements.
• On output, AP has been overwritten.

Note:

This matrix is used in Example 4.1 in referenced text [63].

Call Statement and Input

           IOPT  AP   W   Z  LDZ  N   M   AUX  NAUX
            |    |    |   |   |   |   |    |    |
CALL DSPSV( 21 , AP , W , Z , 4 , 4 , 2 , AUX , 36 )
 
AP       =  (5.0, 4.0, 5.0, 1.0, 1.0, 4.0, 1.0, 1.0, 2.0, 4.0)

Output

        *          *
        | 1.000000 |
W    =  | 2.000000 |
        |  .       |
        |  .       |
        *          *

        *                      *
        | -0.707107   0.000000 |
Z    =  |  0.707107   0.000000 |
        |  0.000000  -0.707107 |
        |  0.000000   0.707107 |
        *                      *

Example 2

This example shows how to find the three largest eigenvalues and corresponding eigenvectors of a real long-precision symmetric matrix A of order 4, stored in lower-packed storage mode, having an eigenvalue of multiplicity two. Matrix A is:

                          *                    *
                          | 6.0  4.0  4.0  1.0 |
                          | 4.0  6.0  1.0  4.0 |
                          | 4.0  1.0  6.0  4.0 |
                          | 1.0  4.0  4.0  6.0 |
                          *                    *

where:

• NAUX is equal to 9N.
• AUX contains 9N elements.
• On output, AP has been overwritten.

Note:

This matrix is used in Example 4.2 in referenced text [63].

Call Statement and Input

           IOPT  AP   W   Z  LDZ  N   M   AUX  NAUX
            |    |    |   |   |   |   |    |    |
CALL DSPSV( 11 , AP , W , Z , 8 , 4 , 3 , AUX , 36 )
 
AP       =  (6.0, 4.0, 4.0, 1.0, 6.0, 1.0, 4.0, 6.0, 4.0, 6.0)

Output

        *           *
        | 15.000000 |
W    =  |  5.000000 |
        |  5.000000 |
        |   .       |
        *           *

        *                                 *
        |  0.500000   0.707107   0.000000 |
        |  0.500000   0.000000  -0.707107 |
        |  0.500000   0.000000   0.707107 |
Z    =  |  0.500000  -0.707107   0.000000 |
        |   .          .          .       |
        |   .          .          .       |
        |   .          .          .       |
        |   .          .          .       |
        *                                 *

Example 3

This example shows how to find the largest eigenvalue and the corresponding eigenvector of a complex Hermitian matrix A of order 2, stored in lower-packed storage mode. Matrix A is:

                      *                         *
                      | (1.0, 0.0)  (0.0, -1.0) |
                      | (0.0, 1.0)  (1.0,  0.0) |
                      *                         *

where:

• NAUX is equal to 11N.
• AUX contains 11N elements.
• On output, AP has been overwritten.

Note:

This matrix is used in Example 6.1 in referenced text [63].

Call Statement and Input

           IOPT  AP   W   Z  LDZ  N   M   AUX  NAUX
            |    |    |   |   |   |   |    |    |
CALL ZHPSV( 11 , AP , W , Z , 2 , 2 , 1 , AUX , 22 )
 
AP       =  ((1.0, . ), (0.0, 1.0), (1.0, . ))

Output

        *          *
W    =  | 2.000000 |
        |  .       |
        *          *

        *                       *
Z    =  | (0.000000, -0.707107) |
        | (0.707107,  0.000000) |
        *                       *

Example 4

This example shows how to find the two smallest eigenvalues only of a complex Hermitian matrix A of order 4, stored in upper-packed storage mode. Matrix A is:

         *                                                  *
         | (3.0,  0.0)  (1.0, 0.0)  (0.0,  0.0)  (0.0, 2.0) |
         | (1.0,  0.0)  (3.0, 0.0)  (0.0, -2.0)  (0.0, 0.0) |
         | (0.0,  0.0)  (0.0, 2.0)  (1.0,  0.0)  (1.0, 0.0) |
         | (0.0, -2.0)  (0.0, 0.0)  (1.0,  0.0)  (1.0, 0.0) |
         *                                                  *

where:

• NAUX is equal to 5N.
• AUX contains 5N elements.
• LDZ is set to 1 to avoid an error condition.
• DUMMY is used as a placeholder for argument z, which is not used in the computation.
• On output, AP has been overwritten.

Note:

This matrix is used in Example 6.6 in referenced text [63].

Call Statement and Input

           IOPT  AP   W     Z    LDZ  N   M   AUX  NAUX
            |    |    |     |     |   |   |    |    |
CALL ZHPSV( 20 , AP , W , DUMMY , 1 , 4 , 2 , AUX , 20 )
 
AP       =  ((3.0, . ), (1.0, 0.0), (3.0, . ), (0.0, 0.0),
             (0.0, -2.0), (1.0, . ), (0.0, 2.0), (0.0, 0.0),
             (1.0, 0.0), (1.0, . ))

Output

        *           *
        | -0.828427 |
W    =  |  0.000000 |
        |   .       |
        |   .       |
        *           *