1 SUBROUTINE ZHBEVD( JOBZ, UPLO, N, KD, AB, LDAB, W, Z, LDZ, WORK,
2 $ LWORK, RWORK, 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, KD, LDAB, LDZ, LIWORK, LRWORK, LWORK, N
12 * ..
13 * .. Array Arguments ..
14 INTEGER IWORK( * )
15 DOUBLE PRECISION RWORK( * ), W( * )
16 COMPLEX*16 AB( LDAB, * ), WORK( * ), Z( LDZ, * )
17 * ..
18 *
19 * Purpose
20 * =======
21 *
22 * ZHBEVD computes all the eigenvalues and, optionally, eigenvectors of
23 * a complex Hermitian band matrix A. If eigenvectors are desired, it
24 * uses a 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 * KD (input) INTEGER
48 * The number of superdiagonals of the matrix A if UPLO = 'U',
49 * or the number of subdiagonals if UPLO = 'L'. KD >= 0.
50 *
51 * AB (input/output) COMPLEX*16 array, dimension (LDAB, N)
52 * On entry, the upper or lower triangle of the Hermitian band
53 * matrix A, stored in the first KD+1 rows of the array. The
54 * j-th column of A is stored in the j-th column of the array AB
55 * as follows:
56 * if UPLO = 'U', AB(kd+1+i-j,j) = A(i,j) for max(1,j-kd)<=i<=j;
57 * if UPLO = 'L', AB(1+i-j,j) = A(i,j) for j<=i<=min(n,j+kd).
58 *
59 * On exit, AB is overwritten by values generated during the
60 * reduction to tridiagonal form. If UPLO = 'U', the first
61 * superdiagonal and the diagonal of the tridiagonal matrix T
62 * are returned in rows KD and KD+1 of AB, and if UPLO = 'L',
63 * the diagonal and first subdiagonal of T are returned in the
64 * first two rows of AB.
65 *
66 * LDAB (input) INTEGER
67 * The leading dimension of the array AB. LDAB >= KD + 1.
68 *
69 * W (output) DOUBLE PRECISION array, dimension (N)
70 * If INFO = 0, the eigenvalues in ascending order.
71 *
72 * Z (output) COMPLEX*16 array, dimension (LDZ, N)
73 * If JOBZ = 'V', then if INFO = 0, Z contains the orthonormal
74 * eigenvectors of the matrix A, with the i-th column of Z
75 * holding the eigenvector associated with W(i).
76 * If JOBZ = 'N', then Z is not referenced.
77 *
78 * LDZ (input) INTEGER
79 * The leading dimension of the array Z. LDZ >= 1, and if
80 * JOBZ = 'V', LDZ >= max(1,N).
81 *
82 * WORK (workspace/output) COMPLEX*16 array, dimension (MAX(1,LWORK))
83 * On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
84 *
85 * LWORK (input) INTEGER
86 * The dimension of the array WORK.
87 * If N <= 1, LWORK must be at least 1.
88 * If JOBZ = 'N' and N > 1, LWORK must be at least N.
89 * If JOBZ = 'V' and N > 1, LWORK must be at least 2*N**2.
90 *
91 * If LWORK = -1, then a workspace query is assumed; the routine
92 * only calculates the optimal sizes of the WORK, RWORK and
93 * IWORK arrays, returns these values as the first entries of
94 * the WORK, RWORK and IWORK arrays, and no error message
95 * related to LWORK or LRWORK or LIWORK is issued by XERBLA.
96 *
97 * RWORK (workspace/output) DOUBLE PRECISION array,
98 * dimension (LRWORK)
99 * On exit, if INFO = 0, RWORK(1) returns the optimal LRWORK.
100 *
101 * LRWORK (input) INTEGER
102 * The dimension of array RWORK.
103 * If N <= 1, LRWORK must be at least 1.
104 * If JOBZ = 'N' and N > 1, LRWORK must be at least N.
105 * If JOBZ = 'V' and N > 1, LRWORK must be at least
106 * 1 + 5*N + 2*N**2.
107 *
108 * If LRWORK = -1, then a workspace query is assumed; the
109 * routine only calculates the optimal sizes of the WORK, RWORK
110 * and IWORK arrays, returns these values as the first entries
111 * of the WORK, RWORK and IWORK arrays, and no error message
112 * related to LWORK or LRWORK or LIWORK is issued by XERBLA.
113 *
114 * IWORK (workspace/output) INTEGER array, dimension (MAX(1,LIWORK))
115 * On exit, if INFO = 0, IWORK(1) returns the optimal LIWORK.
116 *
117 * LIWORK (input) INTEGER
118 * The dimension of array IWORK.
119 * If JOBZ = 'N' or N <= 1, LIWORK must be at least 1.
120 * If JOBZ = 'V' and N > 1, LIWORK must be at least 3 + 5*N .
121 *
122 * If LIWORK = -1, then a workspace query is assumed; the
123 * routine only calculates the optimal sizes of the WORK, RWORK
124 * and IWORK arrays, returns these values as the first entries
125 * of the WORK, RWORK and IWORK arrays, and no error message
126 * related to LWORK or LRWORK or LIWORK is issued by XERBLA.
127 *
128 * INFO (output) INTEGER
129 * = 0: successful exit.
130 * < 0: if INFO = -i, the i-th argument had an illegal value.
131 * > 0: if INFO = i, the algorithm failed to converge; i
132 * off-diagonal elements of an intermediate tridiagonal
133 * form did not converge to zero.
134 *
135 * =====================================================================
136 *
137 * .. Parameters ..
138 DOUBLE PRECISION ZERO, ONE
139 PARAMETER ( ZERO = 0.0D0, ONE = 1.0D0 )
140 COMPLEX*16 CZERO, CONE
141 PARAMETER ( CZERO = ( 0.0D0, 0.0D0 ),
142 $ CONE = ( 1.0D0, 0.0D0 ) )
143 * ..
144 * .. Local Scalars ..
145 LOGICAL LOWER, LQUERY, WANTZ
146 INTEGER IINFO, IMAX, INDE, INDWK2, INDWRK, ISCALE,
147 $ LIWMIN, LLRWK, LLWK2, LRWMIN, LWMIN
148 DOUBLE PRECISION ANRM, BIGNUM, EPS, RMAX, RMIN, SAFMIN, SIGMA,
149 $ SMLNUM
150 * ..
151 * .. External Functions ..
152 LOGICAL LSAME
153 DOUBLE PRECISION DLAMCH, ZLANHB
154 EXTERNAL LSAME, DLAMCH, ZLANHB
155 * ..
156 * .. External Subroutines ..
157 EXTERNAL DSCAL, DSTERF, XERBLA, ZGEMM, ZHBTRD, ZLACPY,
158 $ ZLASCL, ZSTEDC
159 * ..
160 * .. Intrinsic Functions ..
161 INTRINSIC SQRT
162 * ..
163 * .. Executable Statements ..
164 *
165 * Test the input parameters.
166 *
167 WANTZ = LSAME( JOBZ, 'V' )
168 LOWER = LSAME( UPLO, 'L' )
169 LQUERY = ( LWORK.EQ.-1 .OR. LIWORK.EQ.-1 .OR. LRWORK.EQ.-1 )
170 *
171 INFO = 0
172 IF( N.LE.1 ) THEN
173 LWMIN = 1
174 LRWMIN = 1
175 LIWMIN = 1
176 ELSE
177 IF( WANTZ ) THEN
178 LWMIN = 2*N**2
179 LRWMIN = 1 + 5*N + 2*N**2
180 LIWMIN = 3 + 5*N
181 ELSE
182 LWMIN = N
183 LRWMIN = N
184 LIWMIN = 1
185 END IF
186 END IF
187 IF( .NOT.( WANTZ .OR. LSAME( JOBZ, 'N' ) ) ) THEN
188 INFO = -1
189 ELSE IF( .NOT.( LOWER .OR. LSAME( UPLO, 'U' ) ) ) THEN
190 INFO = -2
191 ELSE IF( N.LT.0 ) THEN
192 INFO = -3
193 ELSE IF( KD.LT.0 ) THEN
194 INFO = -4
195 ELSE IF( LDAB.LT.KD+1 ) THEN
196 INFO = -6
197 ELSE IF( LDZ.LT.1 .OR. ( WANTZ .AND. LDZ.LT.N ) ) THEN
198 INFO = -9
199 END IF
200 *
201 IF( INFO.EQ.0 ) THEN
202 WORK( 1 ) = LWMIN
203 RWORK( 1 ) = LRWMIN
204 IWORK( 1 ) = LIWMIN
205 *
206 IF( LWORK.LT.LWMIN .AND. .NOT.LQUERY ) THEN
207 INFO = -11
208 ELSE IF( LRWORK.LT.LRWMIN .AND. .NOT.LQUERY ) THEN
209 INFO = -13
210 ELSE IF( LIWORK.LT.LIWMIN .AND. .NOT.LQUERY ) THEN
211 INFO = -15
212 END IF
213 END IF
214 *
215 IF( INFO.NE.0 ) THEN
216 CALL XERBLA( 'ZHBEVD', -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 ) = AB( 1, 1 )
229 IF( WANTZ )
230 $ Z( 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 = ZLANHB( 'M', UPLO, N, KD, AB, LDAB, 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 ) THEN
255 IF( LOWER ) THEN
256 CALL ZLASCL( 'B', KD, KD, ONE, SIGMA, N, N, AB, LDAB, INFO )
257 ELSE
258 CALL ZLASCL( 'Q', KD, KD, ONE, SIGMA, N, N, AB, LDAB, INFO )
259 END IF
260 END IF
261 *
262 * Call ZHBTRD to reduce Hermitian band matrix to tridiagonal form.
263 *
264 INDE = 1
265 INDWRK = INDE + N
266 INDWK2 = 1 + N*N
267 LLWK2 = LWORK - INDWK2 + 1
268 LLRWK = LRWORK - INDWRK + 1
269 CALL ZHBTRD( JOBZ, UPLO, N, KD, AB, LDAB, W, RWORK( INDE ), Z,
270 $ LDZ, WORK, IINFO )
271 *
272 * For eigenvalues only, call DSTERF. For eigenvectors, call ZSTEDC.
273 *
274 IF( .NOT.WANTZ ) THEN
275 CALL DSTERF( N, W, RWORK( INDE ), INFO )
276 ELSE
277 CALL ZSTEDC( 'I', N, W, RWORK( INDE ), WORK, N, WORK( INDWK2 ),
278 $ LLWK2, RWORK( INDWRK ), LLRWK, IWORK, LIWORK,
279 $ INFO )
280 CALL ZGEMM( 'N', 'N', N, N, N, CONE, Z, LDZ, WORK, N, CZERO,
281 $ WORK( INDWK2 ), N )
282 CALL ZLACPY( 'A', N, N, WORK( INDWK2 ), N, Z, LDZ )
283 END IF
284 *
285 * If matrix was scaled, then rescale eigenvalues appropriately.
286 *
287 IF( ISCALE.EQ.1 ) THEN
288 IF( INFO.EQ.0 ) THEN
289 IMAX = N
290 ELSE
291 IMAX = INFO - 1
292 END IF
293 CALL DSCAL( IMAX, ONE / SIGMA, W, 1 )
294 END IF
295 *
296 WORK( 1 ) = LWMIN
297 RWORK( 1 ) = LRWMIN
298 IWORK( 1 ) = LIWMIN
299 RETURN
300 *
301 * End of ZHBEVD
302 *
303 END
2 $ LWORK, RWORK, 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, KD, LDAB, LDZ, LIWORK, LRWORK, LWORK, N
12 * ..
13 * .. Array Arguments ..
14 INTEGER IWORK( * )
15 DOUBLE PRECISION RWORK( * ), W( * )
16 COMPLEX*16 AB( LDAB, * ), WORK( * ), Z( LDZ, * )
17 * ..
18 *
19 * Purpose
20 * =======
21 *
22 * ZHBEVD computes all the eigenvalues and, optionally, eigenvectors of
23 * a complex Hermitian band matrix A. If eigenvectors are desired, it
24 * uses a 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 * KD (input) INTEGER
48 * The number of superdiagonals of the matrix A if UPLO = 'U',
49 * or the number of subdiagonals if UPLO = 'L'. KD >= 0.
50 *
51 * AB (input/output) COMPLEX*16 array, dimension (LDAB, N)
52 * On entry, the upper or lower triangle of the Hermitian band
53 * matrix A, stored in the first KD+1 rows of the array. The
54 * j-th column of A is stored in the j-th column of the array AB
55 * as follows:
56 * if UPLO = 'U', AB(kd+1+i-j,j) = A(i,j) for max(1,j-kd)<=i<=j;
57 * if UPLO = 'L', AB(1+i-j,j) = A(i,j) for j<=i<=min(n,j+kd).
58 *
59 * On exit, AB is overwritten by values generated during the
60 * reduction to tridiagonal form. If UPLO = 'U', the first
61 * superdiagonal and the diagonal of the tridiagonal matrix T
62 * are returned in rows KD and KD+1 of AB, and if UPLO = 'L',
63 * the diagonal and first subdiagonal of T are returned in the
64 * first two rows of AB.
65 *
66 * LDAB (input) INTEGER
67 * The leading dimension of the array AB. LDAB >= KD + 1.
68 *
69 * W (output) DOUBLE PRECISION array, dimension (N)
70 * If INFO = 0, the eigenvalues in ascending order.
71 *
72 * Z (output) COMPLEX*16 array, dimension (LDZ, N)
73 * If JOBZ = 'V', then if INFO = 0, Z contains the orthonormal
74 * eigenvectors of the matrix A, with the i-th column of Z
75 * holding the eigenvector associated with W(i).
76 * If JOBZ = 'N', then Z is not referenced.
77 *
78 * LDZ (input) INTEGER
79 * The leading dimension of the array Z. LDZ >= 1, and if
80 * JOBZ = 'V', LDZ >= max(1,N).
81 *
82 * WORK (workspace/output) COMPLEX*16 array, dimension (MAX(1,LWORK))
83 * On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
84 *
85 * LWORK (input) INTEGER
86 * The dimension of the array WORK.
87 * If N <= 1, LWORK must be at least 1.
88 * If JOBZ = 'N' and N > 1, LWORK must be at least N.
89 * If JOBZ = 'V' and N > 1, LWORK must be at least 2*N**2.
90 *
91 * If LWORK = -1, then a workspace query is assumed; the routine
92 * only calculates the optimal sizes of the WORK, RWORK and
93 * IWORK arrays, returns these values as the first entries of
94 * the WORK, RWORK and IWORK arrays, and no error message
95 * related to LWORK or LRWORK or LIWORK is issued by XERBLA.
96 *
97 * RWORK (workspace/output) DOUBLE PRECISION array,
98 * dimension (LRWORK)
99 * On exit, if INFO = 0, RWORK(1) returns the optimal LRWORK.
100 *
101 * LRWORK (input) INTEGER
102 * The dimension of array RWORK.
103 * If N <= 1, LRWORK must be at least 1.
104 * If JOBZ = 'N' and N > 1, LRWORK must be at least N.
105 * If JOBZ = 'V' and N > 1, LRWORK must be at least
106 * 1 + 5*N + 2*N**2.
107 *
108 * If LRWORK = -1, then a workspace query is assumed; the
109 * routine only calculates the optimal sizes of the WORK, RWORK
110 * and IWORK arrays, returns these values as the first entries
111 * of the WORK, RWORK and IWORK arrays, and no error message
112 * related to LWORK or LRWORK or LIWORK is issued by XERBLA.
113 *
114 * IWORK (workspace/output) INTEGER array, dimension (MAX(1,LIWORK))
115 * On exit, if INFO = 0, IWORK(1) returns the optimal LIWORK.
116 *
117 * LIWORK (input) INTEGER
118 * The dimension of array IWORK.
119 * If JOBZ = 'N' or N <= 1, LIWORK must be at least 1.
120 * If JOBZ = 'V' and N > 1, LIWORK must be at least 3 + 5*N .
121 *
122 * If LIWORK = -1, then a workspace query is assumed; the
123 * routine only calculates the optimal sizes of the WORK, RWORK
124 * and IWORK arrays, returns these values as the first entries
125 * of the WORK, RWORK and IWORK arrays, and no error message
126 * related to LWORK or LRWORK or LIWORK is issued by XERBLA.
127 *
128 * INFO (output) INTEGER
129 * = 0: successful exit.
130 * < 0: if INFO = -i, the i-th argument had an illegal value.
131 * > 0: if INFO = i, the algorithm failed to converge; i
132 * off-diagonal elements of an intermediate tridiagonal
133 * form did not converge to zero.
134 *
135 * =====================================================================
136 *
137 * .. Parameters ..
138 DOUBLE PRECISION ZERO, ONE
139 PARAMETER ( ZERO = 0.0D0, ONE = 1.0D0 )
140 COMPLEX*16 CZERO, CONE
141 PARAMETER ( CZERO = ( 0.0D0, 0.0D0 ),
142 $ CONE = ( 1.0D0, 0.0D0 ) )
143 * ..
144 * .. Local Scalars ..
145 LOGICAL LOWER, LQUERY, WANTZ
146 INTEGER IINFO, IMAX, INDE, INDWK2, INDWRK, ISCALE,
147 $ LIWMIN, LLRWK, LLWK2, LRWMIN, LWMIN
148 DOUBLE PRECISION ANRM, BIGNUM, EPS, RMAX, RMIN, SAFMIN, SIGMA,
149 $ SMLNUM
150 * ..
151 * .. External Functions ..
152 LOGICAL LSAME
153 DOUBLE PRECISION DLAMCH, ZLANHB
154 EXTERNAL LSAME, DLAMCH, ZLANHB
155 * ..
156 * .. External Subroutines ..
157 EXTERNAL DSCAL, DSTERF, XERBLA, ZGEMM, ZHBTRD, ZLACPY,
158 $ ZLASCL, ZSTEDC
159 * ..
160 * .. Intrinsic Functions ..
161 INTRINSIC SQRT
162 * ..
163 * .. Executable Statements ..
164 *
165 * Test the input parameters.
166 *
167 WANTZ = LSAME( JOBZ, 'V' )
168 LOWER = LSAME( UPLO, 'L' )
169 LQUERY = ( LWORK.EQ.-1 .OR. LIWORK.EQ.-1 .OR. LRWORK.EQ.-1 )
170 *
171 INFO = 0
172 IF( N.LE.1 ) THEN
173 LWMIN = 1
174 LRWMIN = 1
175 LIWMIN = 1
176 ELSE
177 IF( WANTZ ) THEN
178 LWMIN = 2*N**2
179 LRWMIN = 1 + 5*N + 2*N**2
180 LIWMIN = 3 + 5*N
181 ELSE
182 LWMIN = N
183 LRWMIN = N
184 LIWMIN = 1
185 END IF
186 END IF
187 IF( .NOT.( WANTZ .OR. LSAME( JOBZ, 'N' ) ) ) THEN
188 INFO = -1
189 ELSE IF( .NOT.( LOWER .OR. LSAME( UPLO, 'U' ) ) ) THEN
190 INFO = -2
191 ELSE IF( N.LT.0 ) THEN
192 INFO = -3
193 ELSE IF( KD.LT.0 ) THEN
194 INFO = -4
195 ELSE IF( LDAB.LT.KD+1 ) THEN
196 INFO = -6
197 ELSE IF( LDZ.LT.1 .OR. ( WANTZ .AND. LDZ.LT.N ) ) THEN
198 INFO = -9
199 END IF
200 *
201 IF( INFO.EQ.0 ) THEN
202 WORK( 1 ) = LWMIN
203 RWORK( 1 ) = LRWMIN
204 IWORK( 1 ) = LIWMIN
205 *
206 IF( LWORK.LT.LWMIN .AND. .NOT.LQUERY ) THEN
207 INFO = -11
208 ELSE IF( LRWORK.LT.LRWMIN .AND. .NOT.LQUERY ) THEN
209 INFO = -13
210 ELSE IF( LIWORK.LT.LIWMIN .AND. .NOT.LQUERY ) THEN
211 INFO = -15
212 END IF
213 END IF
214 *
215 IF( INFO.NE.0 ) THEN
216 CALL XERBLA( 'ZHBEVD', -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 ) = AB( 1, 1 )
229 IF( WANTZ )
230 $ Z( 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 = ZLANHB( 'M', UPLO, N, KD, AB, LDAB, 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 ) THEN
255 IF( LOWER ) THEN
256 CALL ZLASCL( 'B', KD, KD, ONE, SIGMA, N, N, AB, LDAB, INFO )
257 ELSE
258 CALL ZLASCL( 'Q', KD, KD, ONE, SIGMA, N, N, AB, LDAB, INFO )
259 END IF
260 END IF
261 *
262 * Call ZHBTRD to reduce Hermitian band matrix to tridiagonal form.
263 *
264 INDE = 1
265 INDWRK = INDE + N
266 INDWK2 = 1 + N*N
267 LLWK2 = LWORK - INDWK2 + 1
268 LLRWK = LRWORK - INDWRK + 1
269 CALL ZHBTRD( JOBZ, UPLO, N, KD, AB, LDAB, W, RWORK( INDE ), Z,
270 $ LDZ, WORK, IINFO )
271 *
272 * For eigenvalues only, call DSTERF. For eigenvectors, call ZSTEDC.
273 *
274 IF( .NOT.WANTZ ) THEN
275 CALL DSTERF( N, W, RWORK( INDE ), INFO )
276 ELSE
277 CALL ZSTEDC( 'I', N, W, RWORK( INDE ), WORK, N, WORK( INDWK2 ),
278 $ LLWK2, RWORK( INDWRK ), LLRWK, IWORK, LIWORK,
279 $ INFO )
280 CALL ZGEMM( 'N', 'N', N, N, N, CONE, Z, LDZ, WORK, N, CZERO,
281 $ WORK( INDWK2 ), N )
282 CALL ZLACPY( 'A', N, N, WORK( INDWK2 ), N, Z, LDZ )
283 END IF
284 *
285 * If matrix was scaled, then rescale eigenvalues appropriately.
286 *
287 IF( ISCALE.EQ.1 ) THEN
288 IF( INFO.EQ.0 ) THEN
289 IMAX = N
290 ELSE
291 IMAX = INFO - 1
292 END IF
293 CALL DSCAL( IMAX, ONE / SIGMA, W, 1 )
294 END IF
295 *
296 WORK( 1 ) = LWMIN
297 RWORK( 1 ) = LRWMIN
298 IWORK( 1 ) = LIWMIN
299 RETURN
300 *
301 * End of ZHBEVD
302 *
303 END