Annotation of rpl/lapack/lapack/zhegs2.f, revision 1.1
1.1 ! bertrand 1: SUBROUTINE ZHEGS2( ITYPE, UPLO, N, A, LDA, B, LDB, INFO )
! 2: *
! 3: * -- LAPACK routine (version 3.2) --
! 4: * -- LAPACK is a software package provided by Univ. of Tennessee, --
! 5: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
! 6: * November 2006
! 7: *
! 8: * .. Scalar Arguments ..
! 9: CHARACTER UPLO
! 10: INTEGER INFO, ITYPE, LDA, LDB, N
! 11: * ..
! 12: * .. Array Arguments ..
! 13: COMPLEX*16 A( LDA, * ), B( LDB, * )
! 14: * ..
! 15: *
! 16: * Purpose
! 17: * =======
! 18: *
! 19: * ZHEGS2 reduces a complex Hermitian-definite generalized
! 20: * eigenproblem to standard form.
! 21: *
! 22: * If ITYPE = 1, the problem is A*x = lambda*B*x,
! 23: * and A is overwritten by inv(U')*A*inv(U) or inv(L)*A*inv(L')
! 24: *
! 25: * If ITYPE = 2 or 3, the problem is A*B*x = lambda*x or
! 26: * B*A*x = lambda*x, and A is overwritten by U*A*U` or L'*A*L.
! 27: *
! 28: * B must have been previously factorized as U'*U or L*L' by ZPOTRF.
! 29: *
! 30: * Arguments
! 31: * =========
! 32: *
! 33: * ITYPE (input) INTEGER
! 34: * = 1: compute inv(U')*A*inv(U) or inv(L)*A*inv(L');
! 35: * = 2 or 3: compute U*A*U' or L'*A*L.
! 36: *
! 37: * UPLO (input) CHARACTER*1
! 38: * Specifies whether the upper or lower triangular part of the
! 39: * Hermitian matrix A is stored, and how B has been factorized.
! 40: * = 'U': Upper triangular
! 41: * = 'L': Lower triangular
! 42: *
! 43: * N (input) INTEGER
! 44: * The order of the matrices A and B. N >= 0.
! 45: *
! 46: * A (input/output) COMPLEX*16 array, dimension (LDA,N)
! 47: * On entry, the Hermitian matrix A. If UPLO = 'U', the leading
! 48: * n by n upper triangular part of A contains the upper
! 49: * triangular part of the matrix A, and the strictly lower
! 50: * triangular part of A is not referenced. If UPLO = 'L', the
! 51: * leading n by n lower triangular part of A contains the lower
! 52: * triangular part of the matrix A, and the strictly upper
! 53: * triangular part of A is not referenced.
! 54: *
! 55: * On exit, if INFO = 0, the transformed matrix, stored in the
! 56: * same format as A.
! 57: *
! 58: * LDA (input) INTEGER
! 59: * The leading dimension of the array A. LDA >= max(1,N).
! 60: *
! 61: * B (input) COMPLEX*16 array, dimension (LDB,N)
! 62: * The triangular factor from the Cholesky factorization of B,
! 63: * as returned by ZPOTRF.
! 64: *
! 65: * LDB (input) INTEGER
! 66: * The leading dimension of the array B. LDB >= max(1,N).
! 67: *
! 68: * INFO (output) INTEGER
! 69: * = 0: successful exit.
! 70: * < 0: if INFO = -i, the i-th argument had an illegal value.
! 71: *
! 72: * =====================================================================
! 73: *
! 74: * .. Parameters ..
! 75: DOUBLE PRECISION ONE, HALF
! 76: PARAMETER ( ONE = 1.0D+0, HALF = 0.5D+0 )
! 77: COMPLEX*16 CONE
! 78: PARAMETER ( CONE = ( 1.0D+0, 0.0D+0 ) )
! 79: * ..
! 80: * .. Local Scalars ..
! 81: LOGICAL UPPER
! 82: INTEGER K
! 83: DOUBLE PRECISION AKK, BKK
! 84: COMPLEX*16 CT
! 85: * ..
! 86: * .. External Subroutines ..
! 87: EXTERNAL XERBLA, ZAXPY, ZDSCAL, ZHER2, ZLACGV, ZTRMV,
! 88: $ ZTRSV
! 89: * ..
! 90: * .. Intrinsic Functions ..
! 91: INTRINSIC MAX
! 92: * ..
! 93: * .. External Functions ..
! 94: LOGICAL LSAME
! 95: EXTERNAL LSAME
! 96: * ..
! 97: * .. Executable Statements ..
! 98: *
! 99: * Test the input parameters.
! 100: *
! 101: INFO = 0
! 102: UPPER = LSAME( UPLO, 'U' )
! 103: IF( ITYPE.LT.1 .OR. ITYPE.GT.3 ) THEN
! 104: INFO = -1
! 105: ELSE IF( .NOT.UPPER .AND. .NOT.LSAME( UPLO, 'L' ) ) THEN
! 106: INFO = -2
! 107: ELSE IF( N.LT.0 ) THEN
! 108: INFO = -3
! 109: ELSE IF( LDA.LT.MAX( 1, N ) ) THEN
! 110: INFO = -5
! 111: ELSE IF( LDB.LT.MAX( 1, N ) ) THEN
! 112: INFO = -7
! 113: END IF
! 114: IF( INFO.NE.0 ) THEN
! 115: CALL XERBLA( 'ZHEGS2', -INFO )
! 116: RETURN
! 117: END IF
! 118: *
! 119: IF( ITYPE.EQ.1 ) THEN
! 120: IF( UPPER ) THEN
! 121: *
! 122: * Compute inv(U')*A*inv(U)
! 123: *
! 124: DO 10 K = 1, N
! 125: *
! 126: * Update the upper triangle of A(k:n,k:n)
! 127: *
! 128: AKK = A( K, K )
! 129: BKK = B( K, K )
! 130: AKK = AKK / BKK**2
! 131: A( K, K ) = AKK
! 132: IF( K.LT.N ) THEN
! 133: CALL ZDSCAL( N-K, ONE / BKK, A( K, K+1 ), LDA )
! 134: CT = -HALF*AKK
! 135: CALL ZLACGV( N-K, A( K, K+1 ), LDA )
! 136: CALL ZLACGV( N-K, B( K, K+1 ), LDB )
! 137: CALL ZAXPY( N-K, CT, B( K, K+1 ), LDB, A( K, K+1 ),
! 138: $ LDA )
! 139: CALL ZHER2( UPLO, N-K, -CONE, A( K, K+1 ), LDA,
! 140: $ B( K, K+1 ), LDB, A( K+1, K+1 ), LDA )
! 141: CALL ZAXPY( N-K, CT, B( K, K+1 ), LDB, A( K, K+1 ),
! 142: $ LDA )
! 143: CALL ZLACGV( N-K, B( K, K+1 ), LDB )
! 144: CALL ZTRSV( UPLO, 'Conjugate transpose', 'Non-unit',
! 145: $ N-K, B( K+1, K+1 ), LDB, A( K, K+1 ),
! 146: $ LDA )
! 147: CALL ZLACGV( N-K, A( K, K+1 ), LDA )
! 148: END IF
! 149: 10 CONTINUE
! 150: ELSE
! 151: *
! 152: * Compute inv(L)*A*inv(L')
! 153: *
! 154: DO 20 K = 1, N
! 155: *
! 156: * Update the lower triangle of A(k:n,k:n)
! 157: *
! 158: AKK = A( K, K )
! 159: BKK = B( K, K )
! 160: AKK = AKK / BKK**2
! 161: A( K, K ) = AKK
! 162: IF( K.LT.N ) THEN
! 163: CALL ZDSCAL( N-K, ONE / BKK, A( K+1, K ), 1 )
! 164: CT = -HALF*AKK
! 165: CALL ZAXPY( N-K, CT, B( K+1, K ), 1, A( K+1, K ), 1 )
! 166: CALL ZHER2( UPLO, N-K, -CONE, A( K+1, K ), 1,
! 167: $ B( K+1, K ), 1, A( K+1, K+1 ), LDA )
! 168: CALL ZAXPY( N-K, CT, B( K+1, K ), 1, A( K+1, K ), 1 )
! 169: CALL ZTRSV( UPLO, 'No transpose', 'Non-unit', N-K,
! 170: $ B( K+1, K+1 ), LDB, A( K+1, K ), 1 )
! 171: END IF
! 172: 20 CONTINUE
! 173: END IF
! 174: ELSE
! 175: IF( UPPER ) THEN
! 176: *
! 177: * Compute U*A*U'
! 178: *
! 179: DO 30 K = 1, N
! 180: *
! 181: * Update the upper triangle of A(1:k,1:k)
! 182: *
! 183: AKK = A( K, K )
! 184: BKK = B( K, K )
! 185: CALL ZTRMV( UPLO, 'No transpose', 'Non-unit', K-1, B,
! 186: $ LDB, A( 1, K ), 1 )
! 187: CT = HALF*AKK
! 188: CALL ZAXPY( K-1, CT, B( 1, K ), 1, A( 1, K ), 1 )
! 189: CALL ZHER2( UPLO, K-1, CONE, A( 1, K ), 1, B( 1, K ), 1,
! 190: $ A, LDA )
! 191: CALL ZAXPY( K-1, CT, B( 1, K ), 1, A( 1, K ), 1 )
! 192: CALL ZDSCAL( K-1, BKK, A( 1, K ), 1 )
! 193: A( K, K ) = AKK*BKK**2
! 194: 30 CONTINUE
! 195: ELSE
! 196: *
! 197: * Compute L'*A*L
! 198: *
! 199: DO 40 K = 1, N
! 200: *
! 201: * Update the lower triangle of A(1:k,1:k)
! 202: *
! 203: AKK = A( K, K )
! 204: BKK = B( K, K )
! 205: CALL ZLACGV( K-1, A( K, 1 ), LDA )
! 206: CALL ZTRMV( UPLO, 'Conjugate transpose', 'Non-unit', K-1,
! 207: $ B, LDB, A( K, 1 ), LDA )
! 208: CT = HALF*AKK
! 209: CALL ZLACGV( K-1, B( K, 1 ), LDB )
! 210: CALL ZAXPY( K-1, CT, B( K, 1 ), LDB, A( K, 1 ), LDA )
! 211: CALL ZHER2( UPLO, K-1, CONE, A( K, 1 ), LDA, B( K, 1 ),
! 212: $ LDB, A, LDA )
! 213: CALL ZAXPY( K-1, CT, B( K, 1 ), LDB, A( K, 1 ), LDA )
! 214: CALL ZLACGV( K-1, B( K, 1 ), LDB )
! 215: CALL ZDSCAL( K-1, BKK, A( K, 1 ), LDA )
! 216: CALL ZLACGV( K-1, A( K, 1 ), LDA )
! 217: A( K, K ) = AKK*BKK**2
! 218: 40 CONTINUE
! 219: END IF
! 220: END IF
! 221: RETURN
! 222: *
! 223: * End of ZHEGS2
! 224: *
! 225: 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