1 SUBROUTINE ZHEEVD( JOBZ, UPLO, N, A, LDA, W, WORK, LWORK, RWORK,
2 $ LRWORK, IWORK, LIWORK, INFO )
3 *
4 * -- LAPACK driver routine (version 3.2) --
5 * -- LAPACK is a software package provided by Univ. of Tennessee, --
6 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
7 * November 2006
8 *
9 * .. Scalar Arguments ..
10 CHARACTER JOBZ, UPLO
11 INTEGER INFO, LDA, LIWORK, LRWORK, LWORK, N
12 * ..
13 * .. Array Arguments ..
14 INTEGER IWORK( * )
15 DOUBLE PRECISION RWORK( * ), W( * )
16 COMPLEX*16 A( LDA, * ), WORK( * )
17 * ..
18 *
19 * Purpose
20 * =======
21 *
22 * ZHEEVD computes all eigenvalues and, optionally, eigenvectors of a
23 * complex Hermitian matrix A. If eigenvectors are desired, it uses a
24 * divide and conquer algorithm.
25 *
26 * The divide and conquer algorithm makes very mild assumptions about
27 * floating point arithmetic. It will work on machines with a guard
28 * digit in add/subtract, or on those binary machines without guard
29 * digits which subtract like the Cray X-MP, Cray Y-MP, Cray C-90, or
30 * Cray-2. It could conceivably fail on hexadecimal or decimal machines
31 * without guard digits, but we know of none.
32 *
33 * Arguments
34 * =========
35 *
36 * JOBZ (input) CHARACTER*1
37 * = 'N': Compute eigenvalues only;
38 * = 'V': Compute eigenvalues and eigenvectors.
39 *
40 * UPLO (input) CHARACTER*1
41 * = 'U': Upper triangle of A is stored;
42 * = 'L': Lower triangle of A is stored.
43 *
44 * N (input) INTEGER
45 * The order of the matrix A. N >= 0.
46 *
47 * A (input/output) COMPLEX*16 array, dimension (LDA, N)
48 * On entry, the Hermitian matrix A. If UPLO = 'U', the
49 * leading N-by-N upper triangular part of A contains the
50 * upper triangular part of the matrix A. If UPLO = 'L',
51 * the leading N-by-N lower triangular part of A contains
52 * the lower triangular part of the matrix A.
53 * On exit, if JOBZ = 'V', then if INFO = 0, A contains the
54 * orthonormal eigenvectors of the matrix A.
55 * If JOBZ = 'N', then on exit the lower triangle (if UPLO='L')
56 * or the upper triangle (if UPLO='U') of A, including the
57 * diagonal, is destroyed.
58 *
59 * LDA (input) INTEGER
60 * The leading dimension of the array A. LDA >= max(1,N).
61 *
62 * W (output) DOUBLE PRECISION array, dimension (N)
63 * If INFO = 0, the eigenvalues in ascending order.
64 *
65 * WORK (workspace/output) COMPLEX*16 array, dimension (MAX(1,LWORK))
66 * On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
67 *
68 * LWORK (input) INTEGER
69 * The length of the array WORK.
70 * If N <= 1, LWORK must be at least 1.
71 * If JOBZ = 'N' and N > 1, LWORK must be at least N + 1.
72 * If JOBZ = 'V' and N > 1, LWORK must be at least 2*N + N**2.
73 *
74 * If LWORK = -1, then a workspace query is assumed; the routine
75 * only calculates the optimal sizes of the WORK, RWORK and
76 * IWORK arrays, returns these values as the first entries of
77 * the WORK, RWORK and IWORK arrays, and no error message
78 * related to LWORK or LRWORK or LIWORK is issued by XERBLA.
79 *
80 * RWORK (workspace/output) DOUBLE PRECISION array,
81 * dimension (LRWORK)
82 * On exit, if INFO = 0, RWORK(1) returns the optimal LRWORK.
83 *
84 * LRWORK (input) INTEGER
85 * The dimension of the array RWORK.
86 * If N <= 1, LRWORK must be at least 1.
87 * If JOBZ = 'N' and N > 1, LRWORK must be at least N.
88 * If JOBZ = 'V' and N > 1, LRWORK must be at least
89 * 1 + 5*N + 2*N**2.
90 *
91 * If LRWORK = -1, then a workspace query is assumed; the
92 * routine only calculates the optimal sizes of the WORK, RWORK
93 * and IWORK arrays, returns these values as the first entries
94 * of the WORK, RWORK and IWORK arrays, and no error message
95 * related to LWORK or LRWORK or LIWORK is issued by XERBLA.
96 *
97 * IWORK (workspace/output) INTEGER array, dimension (MAX(1,LIWORK))
98 * On exit, if INFO = 0, IWORK(1) returns the optimal LIWORK.
99 *
100 * LIWORK (input) INTEGER
101 * The dimension of the array IWORK.
102 * If N <= 1, LIWORK must be at least 1.
103 * If JOBZ = 'N' and N > 1, LIWORK must be at least 1.
104 * If JOBZ = 'V' and N > 1, LIWORK must be at least 3 + 5*N.
105 *
106 * If LIWORK = -1, then a workspace query is assumed; the
107 * routine only calculates the optimal sizes of the WORK, RWORK
108 * and IWORK arrays, returns these values as the first entries
109 * of the WORK, RWORK and IWORK arrays, and no error message
110 * related to LWORK or LRWORK or LIWORK is issued by XERBLA.
111 *
112 * INFO (output) INTEGER
113 * = 0: successful exit
114 * < 0: if INFO = -i, the i-th argument had an illegal value
115 * > 0: if INFO = i and JOBZ = 'N', then the algorithm failed
116 * to converge; i off-diagonal elements of an intermediate
117 * tridiagonal form did not converge to zero;
118 * if INFO = i and JOBZ = 'V', then the algorithm failed
119 * to compute an eigenvalue while working on the submatrix
120 * lying in rows and columns INFO/(N+1) through
121 * mod(INFO,N+1).
122 *
123 * Further Details
124 * ===============
125 *
126 * Based on contributions by
127 * Jeff Rutter, Computer Science Division, University of California
128 * at Berkeley, USA
129 *
130 * Modified description of INFO. Sven, 16 Feb 05.
131 * =====================================================================
132 *
133 * .. Parameters ..
134 DOUBLE PRECISION ZERO, ONE
135 PARAMETER ( ZERO = 0.0D0, ONE = 1.0D0 )
136 COMPLEX*16 CONE
137 PARAMETER ( CONE = ( 1.0D0, 0.0D0 ) )
138 * ..
139 * .. Local Scalars ..
140 LOGICAL LOWER, LQUERY, WANTZ
141 INTEGER IINFO, IMAX, INDE, INDRWK, INDTAU, INDWK2,
142 $ INDWRK, ISCALE, LIOPT, LIWMIN, LLRWK, LLWORK,
143 $ LLWRK2, LOPT, LROPT, LRWMIN, LWMIN
144 DOUBLE PRECISION ANRM, BIGNUM, EPS, RMAX, RMIN, SAFMIN, SIGMA,
145 $ SMLNUM
146 * ..
147 * .. External Functions ..
148 LOGICAL LSAME
149 INTEGER ILAENV
150 DOUBLE PRECISION DLAMCH, ZLANHE
151 EXTERNAL LSAME, ILAENV, DLAMCH, ZLANHE
152 * ..
153 * .. External Subroutines ..
154 EXTERNAL DSCAL, DSTERF, XERBLA, ZHETRD, ZLACPY, ZLASCL,
155 $ ZSTEDC, ZUNMTR
156 * ..
157 * .. Intrinsic Functions ..
158 INTRINSIC MAX, SQRT
159 * ..
160 * .. Executable Statements ..
161 *
162 * Test the input parameters.
163 *
164 WANTZ = LSAME( JOBZ, 'V' )
165 LOWER = LSAME( UPLO, 'L' )
166 LQUERY = ( LWORK.EQ.-1 .OR. LRWORK.EQ.-1 .OR. LIWORK.EQ.-1 )
167 *
168 INFO = 0
169 IF( .NOT.( WANTZ .OR. LSAME( JOBZ, 'N' ) ) ) THEN
170 INFO = -1
171 ELSE IF( .NOT.( LOWER .OR. LSAME( UPLO, 'U' ) ) ) THEN
172 INFO = -2
173 ELSE IF( N.LT.0 ) THEN
174 INFO = -3
175 ELSE IF( LDA.LT.MAX( 1, N ) ) THEN
176 INFO = -5
177 END IF
178 *
179 IF( INFO.EQ.0 ) THEN
180 IF( N.LE.1 ) THEN
181 LWMIN = 1
182 LRWMIN = 1
183 LIWMIN = 1
184 LOPT = LWMIN
185 LROPT = LRWMIN
186 LIOPT = LIWMIN
187 ELSE
188 IF( WANTZ ) THEN
189 LWMIN = 2*N + N*N
190 LRWMIN = 1 + 5*N + 2*N**2
191 LIWMIN = 3 + 5*N
192 ELSE
193 LWMIN = N + 1
194 LRWMIN = N
195 LIWMIN = 1
196 END IF
197 LOPT = MAX( LWMIN, N +
198 $ ILAENV( 1, 'ZHETRD', UPLO, N, -1, -1, -1 ) )
199 LROPT = LRWMIN
200 LIOPT = LIWMIN
201 END IF
202 WORK( 1 ) = LOPT
203 RWORK( 1 ) = LROPT
204 IWORK( 1 ) = LIOPT
205 *
206 IF( LWORK.LT.LWMIN .AND. .NOT.LQUERY ) THEN
207 INFO = -8
208 ELSE IF( LRWORK.LT.LRWMIN .AND. .NOT.LQUERY ) THEN
209 INFO = -10
210 ELSE IF( LIWORK.LT.LIWMIN .AND. .NOT.LQUERY ) THEN
211 INFO = -12
212 END IF
213 END IF
214 *
215 IF( INFO.NE.0 ) THEN
216 CALL XERBLA( 'ZHEEVD', -INFO )
217 RETURN
218 ELSE IF( LQUERY ) THEN
219 RETURN
220 END IF
221 *
222 * Quick return if possible
223 *
224 IF( N.EQ.0 )
225 $ RETURN
226 *
227 IF( N.EQ.1 ) THEN
228 W( 1 ) = A( 1, 1 )
229 IF( WANTZ )
230 $ A( 1, 1 ) = CONE
231 RETURN
232 END IF
233 *
234 * Get machine constants.
235 *
236 SAFMIN = DLAMCH( 'Safe minimum' )
237 EPS = DLAMCH( 'Precision' )
238 SMLNUM = SAFMIN / EPS
239 BIGNUM = ONE / SMLNUM
240 RMIN = SQRT( SMLNUM )
241 RMAX = SQRT( BIGNUM )
242 *
243 * Scale matrix to allowable range, if necessary.
244 *
245 ANRM = ZLANHE( 'M', UPLO, N, A, LDA, RWORK )
246 ISCALE = 0
247 IF( ANRM.GT.ZERO .AND. ANRM.LT.RMIN ) THEN
248 ISCALE = 1
249 SIGMA = RMIN / ANRM
250 ELSE IF( ANRM.GT.RMAX ) THEN
251 ISCALE = 1
252 SIGMA = RMAX / ANRM
253 END IF
254 IF( ISCALE.EQ.1 )
255 $ CALL ZLASCL( UPLO, 0, 0, ONE, SIGMA, N, N, A, LDA, INFO )
256 *
257 * Call ZHETRD to reduce Hermitian matrix to tridiagonal form.
258 *
259 INDE = 1
260 INDTAU = 1
261 INDWRK = INDTAU + N
262 INDRWK = INDE + N
263 INDWK2 = INDWRK + N*N
264 LLWORK = LWORK - INDWRK + 1
265 LLWRK2 = LWORK - INDWK2 + 1
266 LLRWK = LRWORK - INDRWK + 1
267 CALL ZHETRD( UPLO, N, A, LDA, W, RWORK( INDE ), WORK( INDTAU ),
268 $ WORK( INDWRK ), LLWORK, IINFO )
269 *
270 * For eigenvalues only, call DSTERF. For eigenvectors, first call
271 * ZSTEDC to generate the eigenvector matrix, WORK(INDWRK), of the
272 * tridiagonal matrix, then call ZUNMTR to multiply it to the
273 * Householder transformations represented as Householder vectors in
274 * A.
275 *
276 IF( .NOT.WANTZ ) THEN
277 CALL DSTERF( N, W, RWORK( INDE ), INFO )
278 ELSE
279 CALL ZSTEDC( 'I', N, W, RWORK( INDE ), WORK( INDWRK ), N,
280 $ WORK( INDWK2 ), LLWRK2, RWORK( INDRWK ), LLRWK,
281 $ IWORK, LIWORK, INFO )
282 CALL ZUNMTR( 'L', UPLO, 'N', N, N, A, LDA, WORK( INDTAU ),
283 $ WORK( INDWRK ), N, WORK( INDWK2 ), LLWRK2, IINFO )
284 CALL ZLACPY( 'A', N, N, WORK( INDWRK ), N, A, LDA )
285 END IF
286 *
287 * If matrix was scaled, then rescale eigenvalues appropriately.
288 *
289 IF( ISCALE.EQ.1 ) THEN
290 IF( INFO.EQ.0 ) THEN
291 IMAX = N
292 ELSE
293 IMAX = INFO - 1
294 END IF
295 CALL DSCAL( IMAX, ONE / SIGMA, W, 1 )
296 END IF
297 *
298 WORK( 1 ) = LOPT
299 RWORK( 1 ) = LROPT
300 IWORK( 1 ) = LIOPT
301 *
302 RETURN
303 *
304 * End of ZHEEVD
305 *
306 END
2 $ LRWORK, IWORK, LIWORK, INFO )
3 *
4 * -- LAPACK driver routine (version 3.2) --
5 * -- LAPACK is a software package provided by Univ. of Tennessee, --
6 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
7 * November 2006
8 *
9 * .. Scalar Arguments ..
10 CHARACTER JOBZ, UPLO
11 INTEGER INFO, LDA, LIWORK, LRWORK, LWORK, N
12 * ..
13 * .. Array Arguments ..
14 INTEGER IWORK( * )
15 DOUBLE PRECISION RWORK( * ), W( * )
16 COMPLEX*16 A( LDA, * ), WORK( * )
17 * ..
18 *
19 * Purpose
20 * =======
21 *
22 * ZHEEVD computes all eigenvalues and, optionally, eigenvectors of a
23 * complex Hermitian matrix A. If eigenvectors are desired, it uses a
24 * divide and conquer algorithm.
25 *
26 * The divide and conquer algorithm makes very mild assumptions about
27 * floating point arithmetic. It will work on machines with a guard
28 * digit in add/subtract, or on those binary machines without guard
29 * digits which subtract like the Cray X-MP, Cray Y-MP, Cray C-90, or
30 * Cray-2. It could conceivably fail on hexadecimal or decimal machines
31 * without guard digits, but we know of none.
32 *
33 * Arguments
34 * =========
35 *
36 * JOBZ (input) CHARACTER*1
37 * = 'N': Compute eigenvalues only;
38 * = 'V': Compute eigenvalues and eigenvectors.
39 *
40 * UPLO (input) CHARACTER*1
41 * = 'U': Upper triangle of A is stored;
42 * = 'L': Lower triangle of A is stored.
43 *
44 * N (input) INTEGER
45 * The order of the matrix A. N >= 0.
46 *
47 * A (input/output) COMPLEX*16 array, dimension (LDA, N)
48 * On entry, the Hermitian matrix A. If UPLO = 'U', the
49 * leading N-by-N upper triangular part of A contains the
50 * upper triangular part of the matrix A. If UPLO = 'L',
51 * the leading N-by-N lower triangular part of A contains
52 * the lower triangular part of the matrix A.
53 * On exit, if JOBZ = 'V', then if INFO = 0, A contains the
54 * orthonormal eigenvectors of the matrix A.
55 * If JOBZ = 'N', then on exit the lower triangle (if UPLO='L')
56 * or the upper triangle (if UPLO='U') of A, including the
57 * diagonal, is destroyed.
58 *
59 * LDA (input) INTEGER
60 * The leading dimension of the array A. LDA >= max(1,N).
61 *
62 * W (output) DOUBLE PRECISION array, dimension (N)
63 * If INFO = 0, the eigenvalues in ascending order.
64 *
65 * WORK (workspace/output) COMPLEX*16 array, dimension (MAX(1,LWORK))
66 * On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
67 *
68 * LWORK (input) INTEGER
69 * The length of the array WORK.
70 * If N <= 1, LWORK must be at least 1.
71 * If JOBZ = 'N' and N > 1, LWORK must be at least N + 1.
72 * If JOBZ = 'V' and N > 1, LWORK must be at least 2*N + N**2.
73 *
74 * If LWORK = -1, then a workspace query is assumed; the routine
75 * only calculates the optimal sizes of the WORK, RWORK and
76 * IWORK arrays, returns these values as the first entries of
77 * the WORK, RWORK and IWORK arrays, and no error message
78 * related to LWORK or LRWORK or LIWORK is issued by XERBLA.
79 *
80 * RWORK (workspace/output) DOUBLE PRECISION array,
81 * dimension (LRWORK)
82 * On exit, if INFO = 0, RWORK(1) returns the optimal LRWORK.
83 *
84 * LRWORK (input) INTEGER
85 * The dimension of the array RWORK.
86 * If N <= 1, LRWORK must be at least 1.
87 * If JOBZ = 'N' and N > 1, LRWORK must be at least N.
88 * If JOBZ = 'V' and N > 1, LRWORK must be at least
89 * 1 + 5*N + 2*N**2.
90 *
91 * If LRWORK = -1, then a workspace query is assumed; the
92 * routine only calculates the optimal sizes of the WORK, RWORK
93 * and IWORK arrays, returns these values as the first entries
94 * of the WORK, RWORK and IWORK arrays, and no error message
95 * related to LWORK or LRWORK or LIWORK is issued by XERBLA.
96 *
97 * IWORK (workspace/output) INTEGER array, dimension (MAX(1,LIWORK))
98 * On exit, if INFO = 0, IWORK(1) returns the optimal LIWORK.
99 *
100 * LIWORK (input) INTEGER
101 * The dimension of the array IWORK.
102 * If N <= 1, LIWORK must be at least 1.
103 * If JOBZ = 'N' and N > 1, LIWORK must be at least 1.
104 * If JOBZ = 'V' and N > 1, LIWORK must be at least 3 + 5*N.
105 *
106 * If LIWORK = -1, then a workspace query is assumed; the
107 * routine only calculates the optimal sizes of the WORK, RWORK
108 * and IWORK arrays, returns these values as the first entries
109 * of the WORK, RWORK and IWORK arrays, and no error message
110 * related to LWORK or LRWORK or LIWORK is issued by XERBLA.
111 *
112 * INFO (output) INTEGER
113 * = 0: successful exit
114 * < 0: if INFO = -i, the i-th argument had an illegal value
115 * > 0: if INFO = i and JOBZ = 'N', then the algorithm failed
116 * to converge; i off-diagonal elements of an intermediate
117 * tridiagonal form did not converge to zero;
118 * if INFO = i and JOBZ = 'V', then the algorithm failed
119 * to compute an eigenvalue while working on the submatrix
120 * lying in rows and columns INFO/(N+1) through
121 * mod(INFO,N+1).
122 *
123 * Further Details
124 * ===============
125 *
126 * Based on contributions by
127 * Jeff Rutter, Computer Science Division, University of California
128 * at Berkeley, USA
129 *
130 * Modified description of INFO. Sven, 16 Feb 05.
131 * =====================================================================
132 *
133 * .. Parameters ..
134 DOUBLE PRECISION ZERO, ONE
135 PARAMETER ( ZERO = 0.0D0, ONE = 1.0D0 )
136 COMPLEX*16 CONE
137 PARAMETER ( CONE = ( 1.0D0, 0.0D0 ) )
138 * ..
139 * .. Local Scalars ..
140 LOGICAL LOWER, LQUERY, WANTZ
141 INTEGER IINFO, IMAX, INDE, INDRWK, INDTAU, INDWK2,
142 $ INDWRK, ISCALE, LIOPT, LIWMIN, LLRWK, LLWORK,
143 $ LLWRK2, LOPT, LROPT, LRWMIN, LWMIN
144 DOUBLE PRECISION ANRM, BIGNUM, EPS, RMAX, RMIN, SAFMIN, SIGMA,
145 $ SMLNUM
146 * ..
147 * .. External Functions ..
148 LOGICAL LSAME
149 INTEGER ILAENV
150 DOUBLE PRECISION DLAMCH, ZLANHE
151 EXTERNAL LSAME, ILAENV, DLAMCH, ZLANHE
152 * ..
153 * .. External Subroutines ..
154 EXTERNAL DSCAL, DSTERF, XERBLA, ZHETRD, ZLACPY, ZLASCL,
155 $ ZSTEDC, ZUNMTR
156 * ..
157 * .. Intrinsic Functions ..
158 INTRINSIC MAX, SQRT
159 * ..
160 * .. Executable Statements ..
161 *
162 * Test the input parameters.
163 *
164 WANTZ = LSAME( JOBZ, 'V' )
165 LOWER = LSAME( UPLO, 'L' )
166 LQUERY = ( LWORK.EQ.-1 .OR. LRWORK.EQ.-1 .OR. LIWORK.EQ.-1 )
167 *
168 INFO = 0
169 IF( .NOT.( WANTZ .OR. LSAME( JOBZ, 'N' ) ) ) THEN
170 INFO = -1
171 ELSE IF( .NOT.( LOWER .OR. LSAME( UPLO, 'U' ) ) ) THEN
172 INFO = -2
173 ELSE IF( N.LT.0 ) THEN
174 INFO = -3
175 ELSE IF( LDA.LT.MAX( 1, N ) ) THEN
176 INFO = -5
177 END IF
178 *
179 IF( INFO.EQ.0 ) THEN
180 IF( N.LE.1 ) THEN
181 LWMIN = 1
182 LRWMIN = 1
183 LIWMIN = 1
184 LOPT = LWMIN
185 LROPT = LRWMIN
186 LIOPT = LIWMIN
187 ELSE
188 IF( WANTZ ) THEN
189 LWMIN = 2*N + N*N
190 LRWMIN = 1 + 5*N + 2*N**2
191 LIWMIN = 3 + 5*N
192 ELSE
193 LWMIN = N + 1
194 LRWMIN = N
195 LIWMIN = 1
196 END IF
197 LOPT = MAX( LWMIN, N +
198 $ ILAENV( 1, 'ZHETRD', UPLO, N, -1, -1, -1 ) )
199 LROPT = LRWMIN
200 LIOPT = LIWMIN
201 END IF
202 WORK( 1 ) = LOPT
203 RWORK( 1 ) = LROPT
204 IWORK( 1 ) = LIOPT
205 *
206 IF( LWORK.LT.LWMIN .AND. .NOT.LQUERY ) THEN
207 INFO = -8
208 ELSE IF( LRWORK.LT.LRWMIN .AND. .NOT.LQUERY ) THEN
209 INFO = -10
210 ELSE IF( LIWORK.LT.LIWMIN .AND. .NOT.LQUERY ) THEN
211 INFO = -12
212 END IF
213 END IF
214 *
215 IF( INFO.NE.0 ) THEN
216 CALL XERBLA( 'ZHEEVD', -INFO )
217 RETURN
218 ELSE IF( LQUERY ) THEN
219 RETURN
220 END IF
221 *
222 * Quick return if possible
223 *
224 IF( N.EQ.0 )
225 $ RETURN
226 *
227 IF( N.EQ.1 ) THEN
228 W( 1 ) = A( 1, 1 )
229 IF( WANTZ )
230 $ A( 1, 1 ) = CONE
231 RETURN
232 END IF
233 *
234 * Get machine constants.
235 *
236 SAFMIN = DLAMCH( 'Safe minimum' )
237 EPS = DLAMCH( 'Precision' )
238 SMLNUM = SAFMIN / EPS
239 BIGNUM = ONE / SMLNUM
240 RMIN = SQRT( SMLNUM )
241 RMAX = SQRT( BIGNUM )
242 *
243 * Scale matrix to allowable range, if necessary.
244 *
245 ANRM = ZLANHE( 'M', UPLO, N, A, LDA, RWORK )
246 ISCALE = 0
247 IF( ANRM.GT.ZERO .AND. ANRM.LT.RMIN ) THEN
248 ISCALE = 1
249 SIGMA = RMIN / ANRM
250 ELSE IF( ANRM.GT.RMAX ) THEN
251 ISCALE = 1
252 SIGMA = RMAX / ANRM
253 END IF
254 IF( ISCALE.EQ.1 )
255 $ CALL ZLASCL( UPLO, 0, 0, ONE, SIGMA, N, N, A, LDA, INFO )
256 *
257 * Call ZHETRD to reduce Hermitian matrix to tridiagonal form.
258 *
259 INDE = 1
260 INDTAU = 1
261 INDWRK = INDTAU + N
262 INDRWK = INDE + N
263 INDWK2 = INDWRK + N*N
264 LLWORK = LWORK - INDWRK + 1
265 LLWRK2 = LWORK - INDWK2 + 1
266 LLRWK = LRWORK - INDRWK + 1
267 CALL ZHETRD( UPLO, N, A, LDA, W, RWORK( INDE ), WORK( INDTAU ),
268 $ WORK( INDWRK ), LLWORK, IINFO )
269 *
270 * For eigenvalues only, call DSTERF. For eigenvectors, first call
271 * ZSTEDC to generate the eigenvector matrix, WORK(INDWRK), of the
272 * tridiagonal matrix, then call ZUNMTR to multiply it to the
273 * Householder transformations represented as Householder vectors in
274 * A.
275 *
276 IF( .NOT.WANTZ ) THEN
277 CALL DSTERF( N, W, RWORK( INDE ), INFO )
278 ELSE
279 CALL ZSTEDC( 'I', N, W, RWORK( INDE ), WORK( INDWRK ), N,
280 $ WORK( INDWK2 ), LLWRK2, RWORK( INDRWK ), LLRWK,
281 $ IWORK, LIWORK, INFO )
282 CALL ZUNMTR( 'L', UPLO, 'N', N, N, A, LDA, WORK( INDTAU ),
283 $ WORK( INDWRK ), N, WORK( INDWK2 ), LLWRK2, IINFO )
284 CALL ZLACPY( 'A', N, N, WORK( INDWRK ), N, A, LDA )
285 END IF
286 *
287 * If matrix was scaled, then rescale eigenvalues appropriately.
288 *
289 IF( ISCALE.EQ.1 ) THEN
290 IF( INFO.EQ.0 ) THEN
291 IMAX = N
292 ELSE
293 IMAX = INFO - 1
294 END IF
295 CALL DSCAL( IMAX, ONE / SIGMA, W, 1 )
296 END IF
297 *
298 WORK( 1 ) = LOPT
299 RWORK( 1 ) = LROPT
300 IWORK( 1 ) = LIOPT
301 *
302 RETURN
303 *
304 * End of ZHEEVD
305 *
306 END