1: *> \brief \b ZHBGVD
2: *
3: * =========== DOCUMENTATION ===========
4: *
5: * Online html documentation available at
6: * http://www.netlib.org/lapack/explore-html/
7: *
8: *> \htmlonly
9: *> Download ZHBGVD + dependencies
10: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zhbgvd.f">
11: *> [TGZ]</a>
12: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zhbgvd.f">
13: *> [ZIP]</a>
14: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zhbgvd.f">
15: *> [TXT]</a>
16: *> \endhtmlonly
17: *
18: * Definition:
19: * ===========
20: *
21: * SUBROUTINE ZHBGVD( JOBZ, UPLO, N, KA, KB, AB, LDAB, BB, LDBB, W,
22: * Z, LDZ, WORK, LWORK, RWORK, LRWORK, IWORK,
23: * LIWORK, INFO )
24: *
25: * .. Scalar Arguments ..
26: * CHARACTER JOBZ, UPLO
27: * INTEGER INFO, KA, KB, LDAB, LDBB, LDZ, LIWORK, LRWORK,
28: * $ LWORK, N
29: * ..
30: * .. Array Arguments ..
31: * INTEGER IWORK( * )
32: * DOUBLE PRECISION RWORK( * ), W( * )
33: * COMPLEX*16 AB( LDAB, * ), BB( LDBB, * ), WORK( * ),
34: * $ Z( LDZ, * )
35: * ..
36: *
37: *
38: *> \par Purpose:
39: * =============
40: *>
41: *> \verbatim
42: *>
43: *> ZHBGVD computes all the eigenvalues, and optionally, the eigenvectors
44: *> of a complex generalized Hermitian-definite banded eigenproblem, of
45: *> the form A*x=(lambda)*B*x. Here A and B are assumed to be Hermitian
46: *> and banded, and B is also positive definite. If eigenvectors are
47: *> desired, it uses a divide and conquer algorithm.
48: *>
49: *> The divide and conquer algorithm makes very mild assumptions about
50: *> floating point arithmetic. It will work on machines with a guard
51: *> digit in add/subtract, or on those binary machines without guard
52: *> digits which subtract like the Cray X-MP, Cray Y-MP, Cray C-90, or
53: *> Cray-2. It could conceivably fail on hexadecimal or decimal machines
54: *> without guard digits, but we know of none.
55: *> \endverbatim
56: *
57: * Arguments:
58: * ==========
59: *
60: *> \param[in] JOBZ
61: *> \verbatim
62: *> JOBZ is CHARACTER*1
63: *> = 'N': Compute eigenvalues only;
64: *> = 'V': Compute eigenvalues and eigenvectors.
65: *> \endverbatim
66: *>
67: *> \param[in] UPLO
68: *> \verbatim
69: *> UPLO is CHARACTER*1
70: *> = 'U': Upper triangles of A and B are stored;
71: *> = 'L': Lower triangles of A and B are stored.
72: *> \endverbatim
73: *>
74: *> \param[in] N
75: *> \verbatim
76: *> N is INTEGER
77: *> The order of the matrices A and B. N >= 0.
78: *> \endverbatim
79: *>
80: *> \param[in] KA
81: *> \verbatim
82: *> KA is INTEGER
83: *> The number of superdiagonals of the matrix A if UPLO = 'U',
84: *> or the number of subdiagonals if UPLO = 'L'. KA >= 0.
85: *> \endverbatim
86: *>
87: *> \param[in] KB
88: *> \verbatim
89: *> KB is INTEGER
90: *> The number of superdiagonals of the matrix B if UPLO = 'U',
91: *> or the number of subdiagonals if UPLO = 'L'. KB >= 0.
92: *> \endverbatim
93: *>
94: *> \param[in,out] AB
95: *> \verbatim
96: *> AB is COMPLEX*16 array, dimension (LDAB, N)
97: *> On entry, the upper or lower triangle of the Hermitian band
98: *> matrix A, stored in the first ka+1 rows of the array. The
99: *> j-th column of A is stored in the j-th column of the array AB
100: *> as follows:
101: *> if UPLO = 'U', AB(ka+1+i-j,j) = A(i,j) for max(1,j-ka)<=i<=j;
102: *> if UPLO = 'L', AB(1+i-j,j) = A(i,j) for j<=i<=min(n,j+ka).
103: *>
104: *> On exit, the contents of AB are destroyed.
105: *> \endverbatim
106: *>
107: *> \param[in] LDAB
108: *> \verbatim
109: *> LDAB is INTEGER
110: *> The leading dimension of the array AB. LDAB >= KA+1.
111: *> \endverbatim
112: *>
113: *> \param[in,out] BB
114: *> \verbatim
115: *> BB is COMPLEX*16 array, dimension (LDBB, N)
116: *> On entry, the upper or lower triangle of the Hermitian band
117: *> matrix B, stored in the first kb+1 rows of the array. The
118: *> j-th column of B is stored in the j-th column of the array BB
119: *> as follows:
120: *> if UPLO = 'U', BB(kb+1+i-j,j) = B(i,j) for max(1,j-kb)<=i<=j;
121: *> if UPLO = 'L', BB(1+i-j,j) = B(i,j) for j<=i<=min(n,j+kb).
122: *>
123: *> On exit, the factor S from the split Cholesky factorization
124: *> B = S**H*S, as returned by ZPBSTF.
125: *> \endverbatim
126: *>
127: *> \param[in] LDBB
128: *> \verbatim
129: *> LDBB is INTEGER
130: *> The leading dimension of the array BB. LDBB >= KB+1.
131: *> \endverbatim
132: *>
133: *> \param[out] W
134: *> \verbatim
135: *> W is DOUBLE PRECISION array, dimension (N)
136: *> If INFO = 0, the eigenvalues in ascending order.
137: *> \endverbatim
138: *>
139: *> \param[out] Z
140: *> \verbatim
141: *> Z is COMPLEX*16 array, dimension (LDZ, N)
142: *> If JOBZ = 'V', then if INFO = 0, Z contains the matrix Z of
143: *> eigenvectors, with the i-th column of Z holding the
144: *> eigenvector associated with W(i). The eigenvectors are
145: *> normalized so that Z**H*B*Z = I.
146: *> If JOBZ = 'N', then Z is not referenced.
147: *> \endverbatim
148: *>
149: *> \param[in] LDZ
150: *> \verbatim
151: *> LDZ is INTEGER
152: *> The leading dimension of the array Z. LDZ >= 1, and if
153: *> JOBZ = 'V', LDZ >= N.
154: *> \endverbatim
155: *>
156: *> \param[out] WORK
157: *> \verbatim
158: *> WORK is COMPLEX*16 array, dimension (MAX(1,LWORK))
159: *> On exit, if INFO=0, WORK(1) returns the optimal LWORK.
160: *> \endverbatim
161: *>
162: *> \param[in] LWORK
163: *> \verbatim
164: *> LWORK is INTEGER
165: *> The dimension of the array WORK.
166: *> If N <= 1, LWORK >= 1.
167: *> If JOBZ = 'N' and N > 1, LWORK >= N.
168: *> If JOBZ = 'V' and N > 1, LWORK >= 2*N**2.
169: *>
170: *> If LWORK = -1, then a workspace query is assumed; the routine
171: *> only calculates the optimal sizes of the WORK, RWORK and
172: *> IWORK arrays, returns these values as the first entries of
173: *> the WORK, RWORK and IWORK arrays, and no error message
174: *> related to LWORK or LRWORK or LIWORK is issued by XERBLA.
175: *> \endverbatim
176: *>
177: *> \param[out] RWORK
178: *> \verbatim
179: *> RWORK is DOUBLE PRECISION array, dimension (MAX(1,LRWORK))
180: *> On exit, if INFO=0, RWORK(1) returns the optimal LRWORK.
181: *> \endverbatim
182: *>
183: *> \param[in] LRWORK
184: *> \verbatim
185: *> LRWORK is INTEGER
186: *> The dimension of array RWORK.
187: *> If N <= 1, LRWORK >= 1.
188: *> If JOBZ = 'N' and N > 1, LRWORK >= N.
189: *> If JOBZ = 'V' and N > 1, LRWORK >= 1 + 5*N + 2*N**2.
190: *>
191: *> If LRWORK = -1, then a workspace query is assumed; the
192: *> routine only calculates the optimal sizes of the WORK, RWORK
193: *> and IWORK arrays, returns these values as the first entries
194: *> of the WORK, RWORK and IWORK arrays, and no error message
195: *> related to LWORK or LRWORK or LIWORK is issued by XERBLA.
196: *> \endverbatim
197: *>
198: *> \param[out] IWORK
199: *> \verbatim
200: *> IWORK is INTEGER array, dimension (MAX(1,LIWORK))
201: *> On exit, if INFO=0, IWORK(1) returns the optimal LIWORK.
202: *> \endverbatim
203: *>
204: *> \param[in] LIWORK
205: *> \verbatim
206: *> LIWORK is INTEGER
207: *> The dimension of array IWORK.
208: *> If JOBZ = 'N' or N <= 1, LIWORK >= 1.
209: *> If JOBZ = 'V' and N > 1, LIWORK >= 3 + 5*N.
210: *>
211: *> If LIWORK = -1, then a workspace query is assumed; the
212: *> routine only calculates the optimal sizes of the WORK, RWORK
213: *> and IWORK arrays, returns these values as the first entries
214: *> of the WORK, RWORK and IWORK arrays, and no error message
215: *> related to LWORK or LRWORK or LIWORK is issued by XERBLA.
216: *> \endverbatim
217: *>
218: *> \param[out] INFO
219: *> \verbatim
220: *> INFO is INTEGER
221: *> = 0: successful exit
222: *> < 0: if INFO = -i, the i-th argument had an illegal value
223: *> > 0: if INFO = i, and i is:
224: *> <= N: the algorithm failed to converge:
225: *> i off-diagonal elements of an intermediate
226: *> tridiagonal form did not converge to zero;
227: *> > N: if INFO = N + i, for 1 <= i <= N, then ZPBSTF
228: *> returned INFO = i: B is not positive definite.
229: *> The factorization of B could not be completed and
230: *> no eigenvalues or eigenvectors were computed.
231: *> \endverbatim
232: *
233: * Authors:
234: * ========
235: *
236: *> \author Univ. of Tennessee
237: *> \author Univ. of California Berkeley
238: *> \author Univ. of Colorado Denver
239: *> \author NAG Ltd.
240: *
241: *> \ingroup complex16OTHEReigen
242: *
243: *> \par Contributors:
244: * ==================
245: *>
246: *> Mark Fahey, Department of Mathematics, Univ. of Kentucky, USA
247: *
248: * =====================================================================
249: SUBROUTINE ZHBGVD( JOBZ, UPLO, N, KA, KB, AB, LDAB, BB, LDBB, W,
250: $ Z, LDZ, WORK, LWORK, RWORK, LRWORK, IWORK,
251: $ LIWORK, INFO )
252: *
253: * -- LAPACK driver routine --
254: * -- LAPACK is a software package provided by Univ. of Tennessee, --
255: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
256: *
257: * .. Scalar Arguments ..
258: CHARACTER JOBZ, UPLO
259: INTEGER INFO, KA, KB, LDAB, LDBB, LDZ, LIWORK, LRWORK,
260: $ LWORK, N
261: * ..
262: * .. Array Arguments ..
263: INTEGER IWORK( * )
264: DOUBLE PRECISION RWORK( * ), W( * )
265: COMPLEX*16 AB( LDAB, * ), BB( LDBB, * ), WORK( * ),
266: $ Z( LDZ, * )
267: * ..
268: *
269: * =====================================================================
270: *
271: * .. Parameters ..
272: COMPLEX*16 CONE, CZERO
273: PARAMETER ( CONE = ( 1.0D+0, 0.0D+0 ),
274: $ CZERO = ( 0.0D+0, 0.0D+0 ) )
275: * ..
276: * .. Local Scalars ..
277: LOGICAL LQUERY, UPPER, WANTZ
278: CHARACTER VECT
279: INTEGER IINFO, INDE, INDWK2, INDWRK, LIWMIN, LLRWK,
280: $ LLWK2, LRWMIN, LWMIN
281: * ..
282: * .. External Functions ..
283: LOGICAL LSAME
284: EXTERNAL LSAME
285: * ..
286: * .. External Subroutines ..
287: EXTERNAL DSTERF, XERBLA, ZGEMM, ZHBGST, ZHBTRD, ZLACPY,
288: $ ZPBSTF, ZSTEDC
289: * ..
290: * .. Executable Statements ..
291: *
292: * Test the input parameters.
293: *
294: WANTZ = LSAME( JOBZ, 'V' )
295: UPPER = LSAME( UPLO, 'U' )
296: LQUERY = ( LWORK.EQ.-1 .OR. LRWORK.EQ.-1 .OR. LIWORK.EQ.-1 )
297: *
298: INFO = 0
299: IF( N.LE.1 ) THEN
300: LWMIN = 1+N
301: LRWMIN = 1+N
302: LIWMIN = 1
303: ELSE IF( WANTZ ) THEN
304: LWMIN = 2*N**2
305: LRWMIN = 1 + 5*N + 2*N**2
306: LIWMIN = 3 + 5*N
307: ELSE
308: LWMIN = N
309: LRWMIN = N
310: LIWMIN = 1
311: END IF
312: IF( .NOT.( WANTZ .OR. LSAME( JOBZ, 'N' ) ) ) THEN
313: INFO = -1
314: ELSE IF( .NOT.( UPPER .OR. LSAME( UPLO, 'L' ) ) ) THEN
315: INFO = -2
316: ELSE IF( N.LT.0 ) THEN
317: INFO = -3
318: ELSE IF( KA.LT.0 ) THEN
319: INFO = -4
320: ELSE IF( KB.LT.0 .OR. KB.GT.KA ) THEN
321: INFO = -5
322: ELSE IF( LDAB.LT.KA+1 ) THEN
323: INFO = -7
324: ELSE IF( LDBB.LT.KB+1 ) THEN
325: INFO = -9
326: ELSE IF( LDZ.LT.1 .OR. ( WANTZ .AND. LDZ.LT.N ) ) THEN
327: INFO = -12
328: END IF
329: *
330: IF( INFO.EQ.0 ) THEN
331: WORK( 1 ) = LWMIN
332: RWORK( 1 ) = LRWMIN
333: IWORK( 1 ) = LIWMIN
334: *
335: IF( LWORK.LT.LWMIN .AND. .NOT.LQUERY ) THEN
336: INFO = -14
337: ELSE IF( LRWORK.LT.LRWMIN .AND. .NOT.LQUERY ) THEN
338: INFO = -16
339: ELSE IF( LIWORK.LT.LIWMIN .AND. .NOT.LQUERY ) THEN
340: INFO = -18
341: END IF
342: END IF
343: *
344: IF( INFO.NE.0 ) THEN
345: CALL XERBLA( 'ZHBGVD', -INFO )
346: RETURN
347: ELSE IF( LQUERY ) THEN
348: RETURN
349: END IF
350: *
351: * Quick return if possible
352: *
353: IF( N.EQ.0 )
354: $ RETURN
355: *
356: * Form a split Cholesky factorization of B.
357: *
358: CALL ZPBSTF( UPLO, N, KB, BB, LDBB, INFO )
359: IF( INFO.NE.0 ) THEN
360: INFO = N + INFO
361: RETURN
362: END IF
363: *
364: * Transform problem to standard eigenvalue problem.
365: *
366: INDE = 1
367: INDWRK = INDE + N
368: INDWK2 = 1 + N*N
369: LLWK2 = LWORK - INDWK2 + 2
370: LLRWK = LRWORK - INDWRK + 2
371: CALL ZHBGST( JOBZ, UPLO, N, KA, KB, AB, LDAB, BB, LDBB, Z, LDZ,
372: $ WORK, RWORK, IINFO )
373: *
374: * Reduce Hermitian band matrix to tridiagonal form.
375: *
376: IF( WANTZ ) THEN
377: VECT = 'U'
378: ELSE
379: VECT = 'N'
380: END IF
381: CALL ZHBTRD( VECT, UPLO, N, KA, AB, LDAB, W, RWORK( INDE ), Z,
382: $ LDZ, WORK, IINFO )
383: *
384: * For eigenvalues only, call DSTERF. For eigenvectors, call ZSTEDC.
385: *
386: IF( .NOT.WANTZ ) THEN
387: CALL DSTERF( N, W, RWORK( INDE ), INFO )
388: ELSE
389: CALL ZSTEDC( 'I', N, W, RWORK( INDE ), WORK, N, WORK( INDWK2 ),
390: $ LLWK2, RWORK( INDWRK ), LLRWK, IWORK, LIWORK,
391: $ INFO )
392: CALL ZGEMM( 'N', 'N', N, N, N, CONE, Z, LDZ, WORK, N, CZERO,
393: $ WORK( INDWK2 ), N )
394: CALL ZLACPY( 'A', N, N, WORK( INDWK2 ), N, Z, LDZ )
395: END IF
396: *
397: WORK( 1 ) = LWMIN
398: RWORK( 1 ) = LRWMIN
399: IWORK( 1 ) = LIWMIN
400: RETURN
401: *
402: * End of ZHBGVD
403: *
404: END
CVSweb interface <joel.bertrand@systella.fr>
![[BACK]](/cvsweb/icons/back.gif){kind=icon}
Return to zhbgvd.f CVS log ![[TXT]](/cvsweb/icons/text.gif){kind=icon}
![[DIR]](/cvsweb/icons/dir.gif){kind=icon}
Up to [local] / rpl / lapack / lapack