LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches
zhbevd.f
Go to the documentation of this file.
1*> \brief <b> ZHBEVD computes the eigenvalues and, optionally, the left and/or right eigenvectors for OTHER matrices</b>
2*
3* =========== DOCUMENTATION ===========
4*
5* Online html documentation available at
6* http://www.netlib.org/lapack/explore-html/
7*
8*> \htmlonly
9*> Download ZHBEVD + dependencies
10*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zhbevd.f">
11*> [TGZ]</a>
12*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zhbevd.f">
13*> [ZIP]</a>
14*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zhbevd.f">
15*> [TXT]</a>
16*> \endhtmlonly
17*
18* Definition:
19* ===========
20*
21* SUBROUTINE ZHBEVD( JOBZ, UPLO, N, KD, AB, LDAB, W, Z, LDZ, WORK,
22* LWORK, RWORK, LRWORK, IWORK, LIWORK, INFO )
23*
24* .. Scalar Arguments ..
25* CHARACTER JOBZ, UPLO
26* INTEGER INFO, KD, LDAB, LDZ, LIWORK, LRWORK, LWORK, N
27* ..
28* .. Array Arguments ..
29* INTEGER IWORK( * )
30* DOUBLE PRECISION RWORK( * ), W( * )
31* COMPLEX*16 AB( LDAB, * ), WORK( * ), Z( LDZ, * )
32* ..
33*
34*
35*> \par Purpose:
36* =============
37*>
38*> \verbatim
39*>
40*> ZHBEVD computes all the eigenvalues and, optionally, eigenvectors of
41*> a complex Hermitian band matrix A. If eigenvectors are desired, it
42*> uses a divide and conquer algorithm.
43*>
44*> \endverbatim
45*
46* Arguments:
47* ==========
48*
49*> \param[in] JOBZ
50*> \verbatim
51*> JOBZ is CHARACTER*1
52*> = 'N': Compute eigenvalues only;
53*> = 'V': Compute eigenvalues and eigenvectors.
54*> \endverbatim
55*>
56*> \param[in] UPLO
57*> \verbatim
58*> UPLO is CHARACTER*1
59*> = 'U': Upper triangle of A is stored;
60*> = 'L': Lower triangle of A is stored.
61*> \endverbatim
62*>
63*> \param[in] N
64*> \verbatim
65*> N is INTEGER
66*> The order of the matrix A. N >= 0.
67*> \endverbatim
68*>
69*> \param[in] KD
70*> \verbatim
71*> KD is INTEGER
72*> The number of superdiagonals of the matrix A if UPLO = 'U',
73*> or the number of subdiagonals if UPLO = 'L'. KD >= 0.
74*> \endverbatim
75*>
76*> \param[in,out] AB
77*> \verbatim
78*> AB is COMPLEX*16 array, dimension (LDAB, N)
79*> On entry, the upper or lower triangle of the Hermitian band
80*> matrix A, stored in the first KD+1 rows of the array. The
81*> j-th column of A is stored in the j-th column of the array AB
82*> as follows:
83*> if UPLO = 'U', AB(kd+1+i-j,j) = A(i,j) for max(1,j-kd)<=i<=j;
84*> if UPLO = 'L', AB(1+i-j,j) = A(i,j) for j<=i<=min(n,j+kd).
85*>
86*> On exit, AB is overwritten by values generated during the
87*> reduction to tridiagonal form. If UPLO = 'U', the first
88*> superdiagonal and the diagonal of the tridiagonal matrix T
89*> are returned in rows KD and KD+1 of AB, and if UPLO = 'L',
90*> the diagonal and first subdiagonal of T are returned in the
91*> first two rows of AB.
92*> \endverbatim
93*>
94*> \param[in] LDAB
95*> \verbatim
96*> LDAB is INTEGER
97*> The leading dimension of the array AB. LDAB >= KD + 1.
98*> \endverbatim
99*>
100*> \param[out] W
101*> \verbatim
102*> W is DOUBLE PRECISION array, dimension (N)
103*> If INFO = 0, the eigenvalues in ascending order.
104*> \endverbatim
105*>
106*> \param[out] Z
107*> \verbatim
108*> Z is COMPLEX*16 array, dimension (LDZ, N)
109*> If JOBZ = 'V', then if INFO = 0, Z contains the orthonormal
110*> eigenvectors of the matrix A, with the i-th column of Z
111*> holding the eigenvector associated with W(i).
112*> If JOBZ = 'N', then Z is not referenced.
113*> \endverbatim
114*>
115*> \param[in] LDZ
116*> \verbatim
117*> LDZ is INTEGER
118*> The leading dimension of the array Z. LDZ >= 1, and if
119*> JOBZ = 'V', LDZ >= max(1,N).
120*> \endverbatim
121*>
122*> \param[out] WORK
123*> \verbatim
124*> WORK is COMPLEX*16 array, dimension (MAX(1,LWORK))
125*> On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
126*> \endverbatim
127*>
128*> \param[in] LWORK
129*> \verbatim
130*> LWORK is INTEGER
131*> The dimension of the array WORK.
132*> If N <= 1, LWORK must be at least 1.
133*> If JOBZ = 'N' and N > 1, LWORK must be at least N.
134*> If JOBZ = 'V' and N > 1, LWORK must be at least 2*N**2.
135*>
136*> If LWORK = -1, then a workspace query is assumed; the routine
137*> only calculates the optimal sizes of the WORK, RWORK and
138*> IWORK arrays, returns these values as the first entries of
139*> the WORK, RWORK and IWORK arrays, and no error message
140*> related to LWORK or LRWORK or LIWORK is issued by XERBLA.
141*> \endverbatim
142*>
143*> \param[out] RWORK
144*> \verbatim
145*> RWORK is DOUBLE PRECISION array,
146*> dimension (LRWORK)
147*> On exit, if INFO = 0, RWORK(1) returns the optimal LRWORK.
148*> \endverbatim
149*>
150*> \param[in] LRWORK
151*> \verbatim
152*> LRWORK is INTEGER
153*> The dimension of array RWORK.
154*> If N <= 1, LRWORK must be at least 1.
155*> If JOBZ = 'N' and N > 1, LRWORK must be at least N.
156*> If JOBZ = 'V' and N > 1, LRWORK must be at least
157*> 1 + 5*N + 2*N**2.
158*>
159*> If LRWORK = -1, then a workspace query is assumed; the
160*> routine only calculates the optimal sizes of the WORK, RWORK
161*> and IWORK arrays, returns these values as the first entries
162*> of the WORK, RWORK and IWORK arrays, and no error message
163*> related to LWORK or LRWORK or LIWORK is issued by XERBLA.
164*> \endverbatim
165*>
166*> \param[out] IWORK
167*> \verbatim
168*> IWORK is INTEGER array, dimension (MAX(1,LIWORK))
169*> On exit, if INFO = 0, IWORK(1) returns the optimal LIWORK.
170*> \endverbatim
171*>
172*> \param[in] LIWORK
173*> \verbatim
174*> LIWORK is INTEGER
175*> The dimension of array IWORK.
176*> If JOBZ = 'N' or N <= 1, LIWORK must be at least 1.
177*> If JOBZ = 'V' and N > 1, LIWORK must be at least 3 + 5*N .
178*>
179*> If LIWORK = -1, then a workspace query is assumed; the
180*> routine only calculates the optimal sizes of the WORK, RWORK
181*> and IWORK arrays, returns these values as the first entries
182*> of the WORK, RWORK and IWORK arrays, and no error message
183*> related to LWORK or LRWORK or LIWORK is issued by XERBLA.
184*> \endverbatim
185*>
186*> \param[out] INFO
187*> \verbatim
188*> INFO is INTEGER
189*> = 0: successful exit.
190*> < 0: if INFO = -i, the i-th argument had an illegal value.
191*> > 0: if INFO = i, the algorithm failed to converge; i
192*> off-diagonal elements of an intermediate tridiagonal
193*> form did not converge to zero.
194*> \endverbatim
195*
196* Authors:
197* ========
198*
199*> \author Univ. of Tennessee
200*> \author Univ. of California Berkeley
201*> \author Univ. of Colorado Denver
202*> \author NAG Ltd.
203*
204*> \ingroup hbevd
205*
206* =====================================================================
207 SUBROUTINE zhbevd( JOBZ, UPLO, N, KD, AB, LDAB, W, Z, LDZ, WORK,
208 $ LWORK, RWORK, LRWORK, IWORK, LIWORK, INFO )
209*
210* -- LAPACK driver routine --
211* -- LAPACK is a software package provided by Univ. of Tennessee, --
212* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
213*
214* .. Scalar Arguments ..
215 CHARACTER JOBZ, UPLO
216 INTEGER INFO, KD, LDAB, LDZ, LIWORK, LRWORK, LWORK, N
217* ..
218* .. Array Arguments ..
219 INTEGER IWORK( * )
220 DOUBLE PRECISION RWORK( * ), W( * )
221 COMPLEX*16 AB( LDAB, * ), WORK( * ), Z( LDZ, * )
222* ..
223*
224* =====================================================================
225*
226* .. Parameters ..
227 DOUBLE PRECISION ZERO, ONE
228 parameter( zero = 0.0d0, one = 1.0d0 )
229 COMPLEX*16 CZERO, CONE
230 parameter( czero = ( 0.0d0, 0.0d0 ),
231 $ cone = ( 1.0d0, 0.0d0 ) )
232* ..
233* .. Local Scalars ..
234 LOGICAL LOWER, LQUERY, WANTZ
235 INTEGER IINFO, IMAX, INDE, INDWK2, INDWRK, ISCALE,
236 $ liwmin, llrwk, llwk2, lrwmin, lwmin
237 DOUBLE PRECISION ANRM, BIGNUM, EPS, RMAX, RMIN, SAFMIN, SIGMA,
238 $ smlnum
239* ..
240* .. External Functions ..
241 LOGICAL LSAME
242 DOUBLE PRECISION DLAMCH, ZLANHB
243 EXTERNAL lsame, dlamch, zlanhb
244* ..
245* .. External Subroutines ..
246 EXTERNAL dscal, dsterf, xerbla, zgemm, zhbtrd, zlacpy,
247 $ zlascl, zstedc
248* ..
249* .. Intrinsic Functions ..
250 INTRINSIC sqrt
251* ..
252* .. Executable Statements ..
253*
254* Test the input parameters.
255*
256 wantz = lsame( jobz, 'V' )
257 lower = lsame( uplo, 'L' )
258 lquery = ( lwork.EQ.-1 .OR. liwork.EQ.-1 .OR. lrwork.EQ.-1 )
259*
260 info = 0
261 IF( n.LE.1 ) THEN
262 lwmin = 1
263 lrwmin = 1
264 liwmin = 1
265 ELSE
266 IF( wantz ) THEN
267 lwmin = 2*n**2
268 lrwmin = 1 + 5*n + 2*n**2
269 liwmin = 3 + 5*n
270 ELSE
271 lwmin = n
272 lrwmin = n
273 liwmin = 1
274 END IF
275 END IF
276 IF( .NOT.( wantz .OR. lsame( jobz, 'N' ) ) ) THEN
277 info = -1
278 ELSE IF( .NOT.( lower .OR. lsame( uplo, 'U' ) ) ) THEN
279 info = -2
280 ELSE IF( n.LT.0 ) THEN
281 info = -3
282 ELSE IF( kd.LT.0 ) THEN
283 info = -4
284 ELSE IF( ldab.LT.kd+1 ) THEN
285 info = -6
286 ELSE IF( ldz.LT.1 .OR. ( wantz .AND. ldz.LT.n ) ) THEN
287 info = -9
288 END IF
289*
290 IF( info.EQ.0 ) THEN
291 work( 1 ) = lwmin
292 rwork( 1 ) = lrwmin
293 iwork( 1 ) = liwmin
294*
295 IF( lwork.LT.lwmin .AND. .NOT.lquery ) THEN
296 info = -11
297 ELSE IF( lrwork.LT.lrwmin .AND. .NOT.lquery ) THEN
298 info = -13
299 ELSE IF( liwork.LT.liwmin .AND. .NOT.lquery ) THEN
300 info = -15
301 END IF
302 END IF
303*
304 IF( info.NE.0 ) THEN
305 CALL xerbla( 'ZHBEVD', -info )
306 RETURN
307 ELSE IF( lquery ) THEN
308 RETURN
309 END IF
310*
311* Quick return if possible
312*
313 IF( n.EQ.0 )
314 $ RETURN
315*
316 IF( n.EQ.1 ) THEN
317 w( 1 ) = dble( ab( 1, 1 ) )
318 IF( wantz )
319 $ z( 1, 1 ) = cone
320 RETURN
321 END IF
322*
323* Get machine constants.
324*
325 safmin = dlamch( 'Safe minimum' )
326 eps = dlamch( 'Precision' )
327 smlnum = safmin / eps
328 bignum = one / smlnum
329 rmin = sqrt( smlnum )
330 rmax = sqrt( bignum )
331*
332* Scale matrix to allowable range, if necessary.
333*
334 anrm = zlanhb( 'M', uplo, n, kd, ab, ldab, rwork )
335 iscale = 0
336 IF( anrm.GT.zero .AND. anrm.LT.rmin ) THEN
337 iscale = 1
338 sigma = rmin / anrm
339 ELSE IF( anrm.GT.rmax ) THEN
340 iscale = 1
341 sigma = rmax / anrm
342 END IF
343 IF( iscale.EQ.1 ) THEN
344 IF( lower ) THEN
345 CALL zlascl( 'B', kd, kd, one, sigma, n, n, ab, ldab, info )
346 ELSE
347 CALL zlascl( 'Q', kd, kd, one, sigma, n, n, ab, ldab, info )
348 END IF
349 END IF
350*
351* Call ZHBTRD to reduce Hermitian band matrix to tridiagonal form.
352*
353 inde = 1
354 indwrk = inde + n
355 indwk2 = 1 + n*n
356 llwk2 = lwork - indwk2 + 1
357 llrwk = lrwork - indwrk + 1
358 CALL zhbtrd( jobz, uplo, n, kd, ab, ldab, w, rwork( inde ), z,
359 $ ldz, work, iinfo )
360*
361* For eigenvalues only, call DSTERF. For eigenvectors, call ZSTEDC.
362*
363 IF( .NOT.wantz ) THEN
364 CALL dsterf( n, w, rwork( inde ), info )
365 ELSE
366 CALL zstedc( 'I', n, w, rwork( inde ), work, n, work( indwk2 ),
367 $ llwk2, rwork( indwrk ), llrwk, iwork, liwork,
368 $ info )
369 CALL zgemm( 'N', 'N', n, n, n, cone, z, ldz, work, n, czero,
370 $ work( indwk2 ), n )
371 CALL zlacpy( 'A', n, n, work( indwk2 ), n, z, ldz )
372 END IF
373*
374* If matrix was scaled, then rescale eigenvalues appropriately.
375*
376 IF( iscale.EQ.1 ) THEN
377 IF( info.EQ.0 ) THEN
378 imax = n
379 ELSE
380 imax = info - 1
381 END IF
382 CALL dscal( imax, one / sigma, w, 1 )
383 END IF
384*
385 work( 1 ) = lwmin
386 rwork( 1 ) = lrwmin
387 iwork( 1 ) = liwmin
388 RETURN
389*
390* End of ZHBEVD
391*
392 END
xerbla
subroutine xerbla(srname, info)
Definition cblat2.f:3285
zgemm
subroutine zgemm(transa, transb, m, n, k, alpha, a, lda, b, ldb, beta, c, ldc)
ZGEMM
Definition zgemm.f:188
zhbevd
subroutine zhbevd(jobz, uplo, n, kd, ab, ldab, w, z, ldz, work, lwork, rwork, lrwork, iwork, liwork, info)
ZHBEVD computes the eigenvalues and, optionally, the left and/or right eigenvectors for OTHER matrice...
Definition zhbevd.f:209
zhbtrd
subroutine zhbtrd(vect, uplo, n, kd, ab, ldab, d, e, q, ldq, work, info)
ZHBTRD
Definition zhbtrd.f:163
zlacpy
subroutine zlacpy(uplo, m, n, a, lda, b, ldb)
ZLACPY copies all or part of one two-dimensional array to another.
Definition zlacpy.f:103
zlascl
subroutine zlascl(type, kl, ku, cfrom, cto, m, n, a, lda, info)
ZLASCL multiplies a general rectangular matrix by a real scalar defined as cto/cfrom.
Definition zlascl.f:143
dscal
subroutine dscal(n, da, dx, incx)
DSCAL
Definition dscal.f:79
zstedc
subroutine zstedc(compz, n, d, e, z, ldz, work, lwork, rwork, lrwork, iwork, liwork, info)
ZSTEDC
Definition zstedc.f:206
dsterf
subroutine dsterf(n, d, e, info)
DSTERF
Definition dsterf.f:86