00001 SUBROUTINE ZPTEQR( COMPZ, N, D, E, Z, LDZ, WORK, INFO ) 00002 * 00003 * -- LAPACK routine (version 3.2) -- 00004 * -- LAPACK is a software package provided by Univ. of Tennessee, -- 00005 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 00006 * November 2006 00007 * 00008 * .. Scalar Arguments .. 00009 CHARACTER COMPZ 00010 INTEGER INFO, LDZ, N 00011 * .. 00012 * .. Array Arguments .. 00013 DOUBLE PRECISION D( * ), E( * ), WORK( * ) 00014 COMPLEX*16 Z( LDZ, * ) 00015 * .. 00016 * 00017 * Purpose 00018 * ======= 00019 * 00020 * ZPTEQR computes all eigenvalues and, optionally, eigenvectors of a 00021 * symmetric positive definite tridiagonal matrix by first factoring the 00022 * matrix using DPTTRF and then calling ZBDSQR to compute the singular 00023 * values of the bidiagonal factor. 00024 * 00025 * This routine computes the eigenvalues of the positive definite 00026 * tridiagonal matrix to high relative accuracy. This means that if the 00027 * eigenvalues range over many orders of magnitude in size, then the 00028 * small eigenvalues and corresponding eigenvectors will be computed 00029 * more accurately than, for example, with the standard QR method. 00030 * 00031 * The eigenvectors of a full or band positive definite Hermitian matrix 00032 * can also be found if ZHETRD, ZHPTRD, or ZHBTRD has been used to 00033 * reduce this matrix to tridiagonal form. (The reduction to 00034 * tridiagonal form, however, may preclude the possibility of obtaining 00035 * high relative accuracy in the small eigenvalues of the original 00036 * matrix, if these eigenvalues range over many orders of magnitude.) 00037 * 00038 * Arguments 00039 * ========= 00040 * 00041 * COMPZ (input) CHARACTER*1 00042 * = 'N': Compute eigenvalues only. 00043 * = 'V': Compute eigenvectors of original Hermitian 00044 * matrix also. Array Z contains the unitary matrix 00045 * used to reduce the original matrix to tridiagonal 00046 * form. 00047 * = 'I': Compute eigenvectors of tridiagonal matrix also. 00048 * 00049 * N (input) INTEGER 00050 * The order of the matrix. N >= 0. 00051 * 00052 * D (input/output) DOUBLE PRECISION array, dimension (N) 00053 * On entry, the n diagonal elements of the tridiagonal matrix. 00054 * On normal exit, D contains the eigenvalues, in descending 00055 * order. 00056 * 00057 * E (input/output) DOUBLE PRECISION array, dimension (N-1) 00058 * On entry, the (n-1) subdiagonal elements of the tridiagonal 00059 * matrix. 00060 * On exit, E has been destroyed. 00061 * 00062 * Z (input/output) COMPLEX*16 array, dimension (LDZ, N) 00063 * On entry, if COMPZ = 'V', the unitary matrix used in the 00064 * reduction to tridiagonal form. 00065 * On exit, if COMPZ = 'V', the orthonormal eigenvectors of the 00066 * original Hermitian matrix; 00067 * if COMPZ = 'I', the orthonormal eigenvectors of the 00068 * tridiagonal matrix. 00069 * If INFO > 0 on exit, Z contains the eigenvectors associated 00070 * with only the stored eigenvalues. 00071 * If COMPZ = 'N', then Z is not referenced. 00072 * 00073 * LDZ (input) INTEGER 00074 * The leading dimension of the array Z. LDZ >= 1, and if 00075 * COMPZ = 'V' or 'I', LDZ >= max(1,N). 00076 * 00077 * WORK (workspace) DOUBLE PRECISION array, dimension (4*N) 00078 * 00079 * INFO (output) INTEGER 00080 * = 0: successful exit. 00081 * < 0: if INFO = -i, the i-th argument had an illegal value. 00082 * > 0: if INFO = i, and i is: 00083 * <= N the Cholesky factorization of the matrix could 00084 * not be performed because the i-th principal minor 00085 * was not positive definite. 00086 * > N the SVD algorithm failed to converge; 00087 * if INFO = N+i, i off-diagonal elements of the 00088 * bidiagonal factor did not converge to zero. 00089 * 00090 * ==================================================================== 00091 * 00092 * .. Parameters .. 00093 COMPLEX*16 CZERO, CONE 00094 PARAMETER ( CZERO = ( 0.0D+0, 0.0D+0 ), 00095 $ CONE = ( 1.0D+0, 0.0D+0 ) ) 00096 * .. 00097 * .. External Functions .. 00098 LOGICAL LSAME 00099 EXTERNAL LSAME 00100 * .. 00101 * .. External Subroutines .. 00102 EXTERNAL DPTTRF, XERBLA, ZBDSQR, ZLASET 00103 * .. 00104 * .. Local Arrays .. 00105 COMPLEX*16 C( 1, 1 ), VT( 1, 1 ) 00106 * .. 00107 * .. Local Scalars .. 00108 INTEGER I, ICOMPZ, NRU 00109 * .. 00110 * .. Intrinsic Functions .. 00111 INTRINSIC MAX, SQRT 00112 * .. 00113 * .. Executable Statements .. 00114 * 00115 * Test the input parameters. 00116 * 00117 INFO = 0 00118 * 00119 IF( LSAME( COMPZ, 'N' ) ) THEN 00120 ICOMPZ = 0 00121 ELSE IF( LSAME( COMPZ, 'V' ) ) THEN 00122 ICOMPZ = 1 00123 ELSE IF( LSAME( COMPZ, 'I' ) ) THEN 00124 ICOMPZ = 2 00125 ELSE 00126 ICOMPZ = -1 00127 END IF 00128 IF( ICOMPZ.LT.0 ) THEN 00129 INFO = -1 00130 ELSE IF( N.LT.0 ) THEN 00131 INFO = -2 00132 ELSE IF( ( LDZ.LT.1 ) .OR. ( ICOMPZ.GT.0 .AND. LDZ.LT.MAX( 1, 00133 $ N ) ) ) THEN 00134 INFO = -6 00135 END IF 00136 IF( INFO.NE.0 ) THEN 00137 CALL XERBLA( 'ZPTEQR', -INFO ) 00138 RETURN 00139 END IF 00140 * 00141 * Quick return if possible 00142 * 00143 IF( N.EQ.0 ) 00144 $ RETURN 00145 * 00146 IF( N.EQ.1 ) THEN 00147 IF( ICOMPZ.GT.0 ) 00148 $ Z( 1, 1 ) = CONE 00149 RETURN 00150 END IF 00151 IF( ICOMPZ.EQ.2 ) 00152 $ CALL ZLASET( 'Full', N, N, CZERO, CONE, Z, LDZ ) 00153 * 00154 * Call DPTTRF to factor the matrix. 00155 * 00156 CALL DPTTRF( N, D, E, INFO ) 00157 IF( INFO.NE.0 ) 00158 $ RETURN 00159 DO 10 I = 1, N 00160 D( I ) = SQRT( D( I ) ) 00161 10 CONTINUE 00162 DO 20 I = 1, N - 1 00163 E( I ) = E( I )*D( I ) 00164 20 CONTINUE 00165 * 00166 * Call ZBDSQR to compute the singular values/vectors of the 00167 * bidiagonal factor. 00168 * 00169 IF( ICOMPZ.GT.0 ) THEN 00170 NRU = N 00171 ELSE 00172 NRU = 0 00173 END IF 00174 CALL ZBDSQR( 'Lower', N, 0, NRU, 0, D, E, VT, 1, Z, LDZ, C, 1, 00175 $ WORK, INFO ) 00176 * 00177 * Square the singular values. 00178 * 00179 IF( INFO.EQ.0 ) THEN 00180 DO 30 I = 1, N 00181 D( I ) = D( I )*D( I ) 00182 30 CONTINUE 00183 ELSE 00184 INFO = N + INFO 00185 END IF 00186 * 00187 RETURN 00188 * 00189 * End of ZPTEQR 00190 * 00191 END
1.7.2