LAPACK  3.4.2
LAPACK: Linear Algebra PACKage
 All Files Functions Groups
zhegv.f
Go to the documentation of this file.
1 *> \brief \b ZHEGST
2 *
3 * =========== DOCUMENTATION ===========
4 *
5 * Online html documentation available at
6 * http://www.netlib.org/lapack/explore-html/
7 *
8 *> \htmlonly
9 *> Download ZHEGV + dependencies
10 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zhegv.f">
11 *> [TGZ]</a>
12 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zhegv.f">
13 *> [ZIP]</a>
14 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zhegv.f">
15 *> [TXT]</a>
16 *> \endhtmlonly
17 *
18 * Definition:
19 * ===========
20 *
21 * SUBROUTINE ZHEGV( ITYPE, JOBZ, UPLO, N, A, LDA, B, LDB, W, WORK,
22 * LWORK, RWORK, INFO )
23 *
24 * .. Scalar Arguments ..
25 * CHARACTER JOBZ, UPLO
26 * INTEGER INFO, ITYPE, LDA, LDB, LWORK, N
27 * ..
28 * .. Array Arguments ..
29 * DOUBLE PRECISION RWORK( * ), W( * )
30 * COMPLEX*16 A( LDA, * ), B( LDB, * ), WORK( * )
31 * ..
32 *
33 *
34 *> \par Purpose:
35 * =============
36 *>
37 *> \verbatim
38 *>
39 *> ZHEGV computes all the eigenvalues, and optionally, the eigenvectors
40 *> of a complex generalized Hermitian-definite eigenproblem, of the form
41 *> A*x=(lambda)*B*x, A*Bx=(lambda)*x, or B*A*x=(lambda)*x.
42 *> Here A and B are assumed to be Hermitian and B is also
43 *> positive definite.
44 *> \endverbatim
45 *
46 * Arguments:
47 * ==========
48 *
49 *> \param[in] ITYPE
50 *> \verbatim
51 *> ITYPE is INTEGER
52 *> Specifies the problem type to be solved:
53 *> = 1: A*x = (lambda)*B*x
54 *> = 2: A*B*x = (lambda)*x
55 *> = 3: B*A*x = (lambda)*x
56 *> \endverbatim
57 *>
58 *> \param[in] JOBZ
59 *> \verbatim
60 *> JOBZ is CHARACTER*1
61 *> = 'N': Compute eigenvalues only;
62 *> = 'V': Compute eigenvalues and eigenvectors.
63 *> \endverbatim
64 *>
65 *> \param[in] UPLO
66 *> \verbatim
67 *> UPLO is CHARACTER*1
68 *> = 'U': Upper triangles of A and B are stored;
69 *> = 'L': Lower triangles of A and B are stored.
70 *> \endverbatim
71 *>
72 *> \param[in] N
73 *> \verbatim
74 *> N is INTEGER
75 *> The order of the matrices A and B. N >= 0.
76 *> \endverbatim
77 *>
78 *> \param[in,out] A
79 *> \verbatim
80 *> A is COMPLEX*16 array, dimension (LDA, N)
81 *> On entry, the Hermitian matrix A. If UPLO = 'U', the
82 *> leading N-by-N upper triangular part of A contains the
83 *> upper triangular part of the matrix A. If UPLO = 'L',
84 *> the leading N-by-N lower triangular part of A contains
85 *> the lower triangular part of the matrix A.
86 *>
87 *> On exit, if JOBZ = 'V', then if INFO = 0, A contains the
88 *> matrix Z of eigenvectors. The eigenvectors are normalized
89 *> as follows:
90 *> if ITYPE = 1 or 2, Z**H*B*Z = I;
91 *> if ITYPE = 3, Z**H*inv(B)*Z = I.
92 *> If JOBZ = 'N', then on exit the upper triangle (if UPLO='U')
93 *> or the lower triangle (if UPLO='L') of A, including the
94 *> diagonal, is destroyed.
95 *> \endverbatim
96 *>
97 *> \param[in] LDA
98 *> \verbatim
99 *> LDA is INTEGER
100 *> The leading dimension of the array A. LDA >= max(1,N).
101 *> \endverbatim
102 *>
103 *> \param[in,out] B
104 *> \verbatim
105 *> B is COMPLEX*16 array, dimension (LDB, N)
106 *> On entry, the Hermitian positive definite matrix B.
107 *> If UPLO = 'U', the leading N-by-N upper triangular part of B
108 *> contains the upper triangular part of the matrix B.
109 *> If UPLO = 'L', the leading N-by-N lower triangular part of B
110 *> contains the lower triangular part of the matrix B.
111 *>
112 *> On exit, if INFO <= N, the part of B containing the matrix is
113 *> overwritten by the triangular factor U or L from the Cholesky
114 *> factorization B = U**H*U or B = L*L**H.
115 *> \endverbatim
116 *>
117 *> \param[in] LDB
118 *> \verbatim
119 *> LDB is INTEGER
120 *> The leading dimension of the array B. LDB >= max(1,N).
121 *> \endverbatim
122 *>
123 *> \param[out] W
124 *> \verbatim
125 *> W is DOUBLE PRECISION array, dimension (N)
126 *> If INFO = 0, the eigenvalues in ascending order.
127 *> \endverbatim
128 *>
129 *> \param[out] WORK
130 *> \verbatim
131 *> WORK is COMPLEX*16 array, dimension (MAX(1,LWORK))
132 *> On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
133 *> \endverbatim
134 *>
135 *> \param[in] LWORK
136 *> \verbatim
137 *> LWORK is INTEGER
138 *> The length of the array WORK. LWORK >= max(1,2*N-1).
139 *> For optimal efficiency, LWORK >= (NB+1)*N,
140 *> where NB is the blocksize for ZHETRD returned by ILAENV.
141 *>
142 *> If LWORK = -1, then a workspace query is assumed; the routine
143 *> only calculates the optimal size of the WORK array, returns
144 *> this value as the first entry of the WORK array, and no error
145 *> message related to LWORK is issued by XERBLA.
146 *> \endverbatim
147 *>
148 *> \param[out] RWORK
149 *> \verbatim
150 *> RWORK is DOUBLE PRECISION array, dimension (max(1, 3*N-2))
151 *> \endverbatim
152 *>
153 *> \param[out] INFO
154 *> \verbatim
155 *> INFO is INTEGER
156 *> = 0: successful exit
157 *> < 0: if INFO = -i, the i-th argument had an illegal value
158 *> > 0: ZPOTRF or ZHEEV returned an error code:
159 *> <= N: if INFO = i, ZHEEV failed to converge;
160 *> i off-diagonal elements of an intermediate
161 *> tridiagonal form did not converge to zero;
162 *> > N: if INFO = N + i, for 1 <= i <= N, then the leading
163 *> minor of order i of B is not positive definite.
164 *> The factorization of B could not be completed and
165 *> no eigenvalues or eigenvectors were computed.
166 *> \endverbatim
167 *
168 * Authors:
169 * ========
170 *
171 *> \author Univ. of Tennessee
172 *> \author Univ. of California Berkeley
173 *> \author Univ. of Colorado Denver
174 *> \author NAG Ltd.
175 *
176 *> \date November 2011
177 *
178 *> \ingroup complex16HEeigen
179 *
180 * =====================================================================
181  SUBROUTINE zhegv( ITYPE, JOBZ, UPLO, N, A, LDA, B, LDB, W, WORK,
182  $ lwork, rwork, info )
183 *
184 * -- LAPACK driver routine (version 3.4.0) --
185 * -- LAPACK is a software package provided by Univ. of Tennessee, --
186 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
187 * November 2011
188 *
189 * .. Scalar Arguments ..
190  CHARACTER jobz, uplo
191  INTEGER info, itype, lda, ldb, lwork, n
192 * ..
193 * .. Array Arguments ..
194  DOUBLE PRECISION rwork( * ), w( * )
195  COMPLEX*16 a( lda, * ), b( ldb, * ), work( * )
196 * ..
197 *
198 * =====================================================================
199 *
200 * .. Parameters ..
201  COMPLEX*16 one
202  parameter( one = ( 1.0d+0, 0.0d+0 ) )
203 * ..
204 * .. Local Scalars ..
205  LOGICAL lquery, upper, wantz
206  CHARACTER trans
207  INTEGER lwkopt, nb, neig
208 * ..
209 * .. External Functions ..
210  LOGICAL lsame
211  INTEGER ilaenv
212  EXTERNAL lsame, ilaenv
213 * ..
214 * .. External Subroutines ..
215  EXTERNAL xerbla, zheev, zhegst, zpotrf, ztrmm, ztrsm
216 * ..
217 * .. Intrinsic Functions ..
218  INTRINSIC max
219 * ..
220 * .. Executable Statements ..
221 *
222 * Test the input parameters.
223 *
224  wantz = lsame( jobz, 'V' )
225  upper = lsame( uplo, 'U' )
226  lquery = ( lwork.EQ.-1 )
227 *
228  info = 0
229  IF( itype.LT.1 .OR. itype.GT.3 ) THEN
230  info = -1
231  ELSE IF( .NOT.( wantz .OR. lsame( jobz, 'N' ) ) ) THEN
232  info = -2
233  ELSE IF( .NOT.( upper .OR. lsame( uplo, 'L' ) ) ) THEN
234  info = -3
235  ELSE IF( n.LT.0 ) THEN
236  info = -4
237  ELSE IF( lda.LT.max( 1, n ) ) THEN
238  info = -6
239  ELSE IF( ldb.LT.max( 1, n ) ) THEN
240  info = -8
241  END IF
242 *
243  IF( info.EQ.0 ) THEN
244  nb = ilaenv( 1, 'ZHETRD', uplo, n, -1, -1, -1 )
245  lwkopt = max( 1, ( nb + 1 )*n )
246  work( 1 ) = lwkopt
247 *
248  IF( lwork.LT.max( 1, 2*n - 1 ) .AND. .NOT.lquery ) THEN
249  info = -11
250  END IF
251  END IF
252 *
253  IF( info.NE.0 ) THEN
254  CALL xerbla( 'ZHEGV ', -info )
255  return
256  ELSE IF( lquery ) THEN
257  return
258  END IF
259 *
260 * Quick return if possible
261 *
262  IF( n.EQ.0 )
263  $ return
264 *
265 * Form a Cholesky factorization of B.
266 *
267  CALL zpotrf( uplo, n, b, ldb, info )
268  IF( info.NE.0 ) THEN
269  info = n + info
270  return
271  END IF
272 *
273 * Transform problem to standard eigenvalue problem and solve.
274 *
275  CALL zhegst( itype, uplo, n, a, lda, b, ldb, info )
276  CALL zheev( jobz, uplo, n, a, lda, w, work, lwork, rwork, info )
277 *
278  IF( wantz ) THEN
279 *
280 * Backtransform eigenvectors to the original problem.
281 *
282  neig = n
283  IF( info.GT.0 )
284  $ neig = info - 1
285  IF( itype.EQ.1 .OR. itype.EQ.2 ) THEN
286 *
287 * For A*x=(lambda)*B*x and A*B*x=(lambda)*x;
288 * backtransform eigenvectors: x = inv(L)**H *y or inv(U)*y
289 *
290  IF( upper ) THEN
291  trans = 'N'
292  ELSE
293  trans = 'C'
294  END IF
295 *
296  CALL ztrsm( 'Left', uplo, trans, 'Non-unit', n, neig, one,
297  $ b, ldb, a, lda )
298 *
299  ELSE IF( itype.EQ.3 ) THEN
300 *
301 * For B*A*x=(lambda)*x;
302 * backtransform eigenvectors: x = L*y or U**H *y
303 *
304  IF( upper ) THEN
305  trans = 'C'
306  ELSE
307  trans = 'N'
308  END IF
309 *
310  CALL ztrmm( 'Left', uplo, trans, 'Non-unit', n, neig, one,
311  $ b, ldb, a, lda )
312  END IF
313  END IF
314 *
315  work( 1 ) = lwkopt
316 *
317  return
318 *
319 * End of ZHEGV
320 *
321  END