Annotation of rpl/lapack/lapack/zhegs2.f, revision 1.19
1.12 bertrand 1: *> \brief \b ZHEGS2 reduces a Hermitian definite generalized eigenproblem to standard form, using the factorization results obtained from cpotrf (unblocked algorithm).
1.9 bertrand 2: *
3: * =========== DOCUMENTATION ===========
4: *
1.16 bertrand 5: * Online html documentation available at
6: * http://www.netlib.org/lapack/explore-html/
1.9 bertrand 7: *
8: *> \htmlonly
1.16 bertrand 9: *> Download ZHEGS2 + dependencies
10: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zhegs2.f">
11: *> [TGZ]</a>
12: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zhegs2.f">
13: *> [ZIP]</a>
14: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zhegs2.f">
1.9 bertrand 15: *> [TXT]</a>
1.16 bertrand 16: *> \endhtmlonly
1.9 bertrand 17: *
18: * Definition:
19: * ===========
20: *
21: * SUBROUTINE ZHEGS2( ITYPE, UPLO, N, A, LDA, B, LDB, INFO )
1.16 bertrand 22: *
1.9 bertrand 23: * .. Scalar Arguments ..
24: * CHARACTER UPLO
25: * INTEGER INFO, ITYPE, LDA, LDB, N
26: * ..
27: * .. Array Arguments ..
28: * COMPLEX*16 A( LDA, * ), B( LDB, * )
29: * ..
1.16 bertrand 30: *
1.9 bertrand 31: *
32: *> \par Purpose:
33: * =============
34: *>
35: *> \verbatim
36: *>
37: *> ZHEGS2 reduces a complex Hermitian-definite generalized
38: *> eigenproblem to standard form.
39: *>
40: *> If ITYPE = 1, the problem is A*x = lambda*B*x,
41: *> and A is overwritten by inv(U**H)*A*inv(U) or inv(L)*A*inv(L**H)
42: *>
43: *> If ITYPE = 2 or 3, the problem is A*B*x = lambda*x or
44: *> B*A*x = lambda*x, and A is overwritten by U*A*U**H or L**H *A*L.
45: *>
46: *> B must have been previously factorized as U**H *U or L*L**H by ZPOTRF.
47: *> \endverbatim
48: *
49: * Arguments:
50: * ==========
51: *
52: *> \param[in] ITYPE
53: *> \verbatim
54: *> ITYPE is INTEGER
55: *> = 1: compute inv(U**H)*A*inv(U) or inv(L)*A*inv(L**H);
56: *> = 2 or 3: compute U*A*U**H or L**H *A*L.
57: *> \endverbatim
58: *>
59: *> \param[in] UPLO
60: *> \verbatim
61: *> UPLO is CHARACTER*1
62: *> Specifies whether the upper or lower triangular part of the
63: *> Hermitian matrix A is stored, and how B has been factorized.
64: *> = 'U': Upper triangular
65: *> = 'L': Lower triangular
66: *> \endverbatim
67: *>
68: *> \param[in] N
69: *> \verbatim
70: *> N is INTEGER
71: *> The order of the matrices A and B. N >= 0.
72: *> \endverbatim
73: *>
74: *> \param[in,out] A
75: *> \verbatim
76: *> A is COMPLEX*16 array, dimension (LDA,N)
77: *> On entry, the Hermitian matrix A. If UPLO = 'U', the leading
78: *> n by n upper triangular part of A contains the upper
79: *> triangular part of the matrix A, and the strictly lower
80: *> triangular part of A is not referenced. If UPLO = 'L', the
81: *> leading n by n lower triangular part of A contains the lower
82: *> triangular part of the matrix A, and the strictly upper
83: *> triangular part of A is not referenced.
84: *>
85: *> On exit, if INFO = 0, the transformed matrix, stored in the
86: *> same format as A.
87: *> \endverbatim
88: *>
89: *> \param[in] LDA
90: *> \verbatim
91: *> LDA is INTEGER
92: *> The leading dimension of the array A. LDA >= max(1,N).
93: *> \endverbatim
94: *>
1.12 bertrand 95: *> \param[in,out] B
1.9 bertrand 96: *> \verbatim
97: *> B is COMPLEX*16 array, dimension (LDB,N)
98: *> The triangular factor from the Cholesky factorization of B,
99: *> as returned by ZPOTRF.
1.19 ! bertrand 100: *> B is modified by the routine but restored on exit.
1.9 bertrand 101: *> \endverbatim
102: *>
103: *> \param[in] LDB
104: *> \verbatim
105: *> LDB is INTEGER
106: *> The leading dimension of the array B. LDB >= max(1,N).
107: *> \endverbatim
108: *>
109: *> \param[out] INFO
110: *> \verbatim
111: *> INFO is INTEGER
112: *> = 0: successful exit.
113: *> < 0: if INFO = -i, the i-th argument had an illegal value.
114: *> \endverbatim
115: *
116: * Authors:
117: * ========
118: *
1.16 bertrand 119: *> \author Univ. of Tennessee
120: *> \author Univ. of California Berkeley
121: *> \author Univ. of Colorado Denver
122: *> \author NAG Ltd.
1.9 bertrand 123: *
1.16 bertrand 124: *> \date December 2016
1.9 bertrand 125: *
126: *> \ingroup complex16HEcomputational
127: *
128: * =====================================================================
1.1 bertrand 129: SUBROUTINE ZHEGS2( ITYPE, UPLO, N, A, LDA, B, LDB, INFO )
130: *
1.16 bertrand 131: * -- LAPACK computational routine (version 3.7.0) --
1.1 bertrand 132: * -- LAPACK is a software package provided by Univ. of Tennessee, --
133: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
1.16 bertrand 134: * December 2016
1.1 bertrand 135: *
136: * .. Scalar Arguments ..
137: CHARACTER UPLO
138: INTEGER INFO, ITYPE, LDA, LDB, N
139: * ..
140: * .. Array Arguments ..
141: COMPLEX*16 A( LDA, * ), B( LDB, * )
142: * ..
143: *
144: * =====================================================================
145: *
146: * .. Parameters ..
147: DOUBLE PRECISION ONE, HALF
148: PARAMETER ( ONE = 1.0D+0, HALF = 0.5D+0 )
149: COMPLEX*16 CONE
150: PARAMETER ( CONE = ( 1.0D+0, 0.0D+0 ) )
151: * ..
152: * .. Local Scalars ..
153: LOGICAL UPPER
154: INTEGER K
155: DOUBLE PRECISION AKK, BKK
156: COMPLEX*16 CT
157: * ..
158: * .. External Subroutines ..
159: EXTERNAL XERBLA, ZAXPY, ZDSCAL, ZHER2, ZLACGV, ZTRMV,
160: $ ZTRSV
161: * ..
162: * .. Intrinsic Functions ..
163: INTRINSIC MAX
164: * ..
165: * .. External Functions ..
166: LOGICAL LSAME
167: EXTERNAL LSAME
168: * ..
169: * .. Executable Statements ..
170: *
171: * Test the input parameters.
172: *
173: INFO = 0
174: UPPER = LSAME( UPLO, 'U' )
175: IF( ITYPE.LT.1 .OR. ITYPE.GT.3 ) THEN
176: INFO = -1
177: ELSE IF( .NOT.UPPER .AND. .NOT.LSAME( UPLO, 'L' ) ) THEN
178: INFO = -2
179: ELSE IF( N.LT.0 ) THEN
180: INFO = -3
181: ELSE IF( LDA.LT.MAX( 1, N ) ) THEN
182: INFO = -5
183: ELSE IF( LDB.LT.MAX( 1, N ) ) THEN
184: INFO = -7
185: END IF
186: IF( INFO.NE.0 ) THEN
187: CALL XERBLA( 'ZHEGS2', -INFO )
188: RETURN
189: END IF
190: *
191: IF( ITYPE.EQ.1 ) THEN
192: IF( UPPER ) THEN
193: *
1.8 bertrand 194: * Compute inv(U**H)*A*inv(U)
1.1 bertrand 195: *
196: DO 10 K = 1, N
197: *
198: * Update the upper triangle of A(k:n,k:n)
199: *
200: AKK = A( K, K )
201: BKK = B( K, K )
202: AKK = AKK / BKK**2
203: A( K, K ) = AKK
204: IF( K.LT.N ) THEN
205: CALL ZDSCAL( N-K, ONE / BKK, A( K, K+1 ), LDA )
206: CT = -HALF*AKK
207: CALL ZLACGV( N-K, A( K, K+1 ), LDA )
208: CALL ZLACGV( N-K, B( K, K+1 ), LDB )
209: CALL ZAXPY( N-K, CT, B( K, K+1 ), LDB, A( K, K+1 ),
210: $ LDA )
211: CALL ZHER2( UPLO, N-K, -CONE, A( K, K+1 ), LDA,
212: $ B( K, K+1 ), LDB, A( K+1, K+1 ), LDA )
213: CALL ZAXPY( N-K, CT, B( K, K+1 ), LDB, A( K, K+1 ),
214: $ LDA )
215: CALL ZLACGV( N-K, B( K, K+1 ), LDB )
216: CALL ZTRSV( UPLO, 'Conjugate transpose', 'Non-unit',
217: $ N-K, B( K+1, K+1 ), LDB, A( K, K+1 ),
218: $ LDA )
219: CALL ZLACGV( N-K, A( K, K+1 ), LDA )
220: END IF
221: 10 CONTINUE
222: ELSE
223: *
1.8 bertrand 224: * Compute inv(L)*A*inv(L**H)
1.1 bertrand 225: *
226: DO 20 K = 1, N
227: *
228: * Update the lower triangle of A(k:n,k:n)
229: *
230: AKK = A( K, K )
231: BKK = B( K, K )
232: AKK = AKK / BKK**2
233: A( K, K ) = AKK
234: IF( K.LT.N ) THEN
235: CALL ZDSCAL( N-K, ONE / BKK, A( K+1, K ), 1 )
236: CT = -HALF*AKK
237: CALL ZAXPY( N-K, CT, B( K+1, K ), 1, A( K+1, K ), 1 )
238: CALL ZHER2( UPLO, N-K, -CONE, A( K+1, K ), 1,
239: $ B( K+1, K ), 1, A( K+1, K+1 ), LDA )
240: CALL ZAXPY( N-K, CT, B( K+1, K ), 1, A( K+1, K ), 1 )
241: CALL ZTRSV( UPLO, 'No transpose', 'Non-unit', N-K,
242: $ B( K+1, K+1 ), LDB, A( K+1, K ), 1 )
243: END IF
244: 20 CONTINUE
245: END IF
246: ELSE
247: IF( UPPER ) THEN
248: *
1.8 bertrand 249: * Compute U*A*U**H
1.1 bertrand 250: *
251: DO 30 K = 1, N
252: *
253: * Update the upper triangle of A(1:k,1:k)
254: *
255: AKK = A( K, K )
256: BKK = B( K, K )
257: CALL ZTRMV( UPLO, 'No transpose', 'Non-unit', K-1, B,
258: $ LDB, A( 1, K ), 1 )
259: CT = HALF*AKK
260: CALL ZAXPY( K-1, CT, B( 1, K ), 1, A( 1, K ), 1 )
261: CALL ZHER2( UPLO, K-1, CONE, A( 1, K ), 1, B( 1, K ), 1,
262: $ A, LDA )
263: CALL ZAXPY( K-1, CT, B( 1, K ), 1, A( 1, K ), 1 )
264: CALL ZDSCAL( K-1, BKK, A( 1, K ), 1 )
265: A( K, K ) = AKK*BKK**2
266: 30 CONTINUE
267: ELSE
268: *
1.8 bertrand 269: * Compute L**H *A*L
1.1 bertrand 270: *
271: DO 40 K = 1, N
272: *
273: * Update the lower triangle of A(1:k,1:k)
274: *
275: AKK = A( K, K )
276: BKK = B( K, K )
277: CALL ZLACGV( K-1, A( K, 1 ), LDA )
278: CALL ZTRMV( UPLO, 'Conjugate transpose', 'Non-unit', K-1,
279: $ B, LDB, A( K, 1 ), LDA )
280: CT = HALF*AKK
281: CALL ZLACGV( K-1, B( K, 1 ), LDB )
282: CALL ZAXPY( K-1, CT, B( K, 1 ), LDB, A( K, 1 ), LDA )
283: CALL ZHER2( UPLO, K-1, CONE, A( K, 1 ), LDA, B( K, 1 ),
284: $ LDB, A, LDA )
285: CALL ZAXPY( K-1, CT, B( K, 1 ), LDB, A( K, 1 ), LDA )
286: CALL ZLACGV( K-1, B( K, 1 ), LDB )
287: CALL ZDSCAL( K-1, BKK, A( K, 1 ), LDA )
288: CALL ZLACGV( K-1, A( K, 1 ), LDA )
289: A( K, K ) = AKK*BKK**2
290: 40 CONTINUE
291: END IF
292: END IF
293: RETURN
294: *
295: * End of ZHEGS2
296: *
297: END
CVSweb interface <joel.bertrand@systella.fr>
![[BACK]](/cvsweb/icons/back.gif){kind=icon}
Return to zhegs2.f CVS log ![[TXT]](/cvsweb/icons/text.gif){kind=icon}
![[DIR]](/cvsweb/icons/dir.gif){kind=icon}
Up to [local] / rpl / lapack / lapack