![]() ![]() | ![]() |
Mise à jour de lapack vers la version 3.3.0.
1: SUBROUTINE ZPTTRF( N, D, E, 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: INTEGER INFO, N 10: * .. 11: * .. Array Arguments .. 12: DOUBLE PRECISION D( * ) 13: COMPLEX*16 E( * ) 14: * .. 15: * 16: * Purpose 17: * ======= 18: * 19: * ZPTTRF computes the L*D*L' factorization of a complex Hermitian 20: * positive definite tridiagonal matrix A. The factorization may also 21: * be regarded as having the form A = U'*D*U. 22: * 23: * Arguments 24: * ========= 25: * 26: * N (input) INTEGER 27: * The order of the matrix A. N >= 0. 28: * 29: * D (input/output) DOUBLE PRECISION array, dimension (N) 30: * On entry, the n diagonal elements of the tridiagonal matrix 31: * A. On exit, the n diagonal elements of the diagonal matrix 32: * D from the L*D*L' factorization of A. 33: * 34: * E (input/output) COMPLEX*16 array, dimension (N-1) 35: * On entry, the (n-1) subdiagonal elements of the tridiagonal 36: * matrix A. On exit, the (n-1) subdiagonal elements of the 37: * unit bidiagonal factor L from the L*D*L' factorization of A. 38: * E can also be regarded as the superdiagonal of the unit 39: * bidiagonal factor U from the U'*D*U factorization of A. 40: * 41: * INFO (output) INTEGER 42: * = 0: successful exit 43: * < 0: if INFO = -k, the k-th argument had an illegal value 44: * > 0: if INFO = k, the leading minor of order k is not 45: * positive definite; if k < N, the factorization could not 46: * be completed, while if k = N, the factorization was 47: * completed, but D(N) <= 0. 48: * 49: * ===================================================================== 50: * 51: * .. Parameters .. 52: DOUBLE PRECISION ZERO 53: PARAMETER ( ZERO = 0.0D+0 ) 54: * .. 55: * .. Local Scalars .. 56: INTEGER I, I4 57: DOUBLE PRECISION EII, EIR, F, G 58: * .. 59: * .. External Subroutines .. 60: EXTERNAL XERBLA 61: * .. 62: * .. Intrinsic Functions .. 63: INTRINSIC DBLE, DCMPLX, DIMAG, MOD 64: * .. 65: * .. Executable Statements .. 66: * 67: * Test the input parameters. 68: * 69: INFO = 0 70: IF( N.LT.0 ) THEN 71: INFO = -1 72: CALL XERBLA( 'ZPTTRF', -INFO ) 73: RETURN 74: END IF 75: * 76: * Quick return if possible 77: * 78: IF( N.EQ.0 ) 79: $ RETURN 80: * 81: * Compute the L*D*L' (or U'*D*U) factorization of A. 82: * 83: I4 = MOD( N-1, 4 ) 84: DO 10 I = 1, I4 85: IF( D( I ).LE.ZERO ) THEN 86: INFO = I 87: GO TO 30 88: END IF 89: EIR = DBLE( E( I ) ) 90: EII = DIMAG( E( I ) ) 91: F = EIR / D( I ) 92: G = EII / D( I ) 93: E( I ) = DCMPLX( F, G ) 94: D( I+1 ) = D( I+1 ) - F*EIR - G*EII 95: 10 CONTINUE 96: * 97: DO 20 I = I4 + 1, N - 4, 4 98: * 99: * Drop out of the loop if d(i) <= 0: the matrix is not positive 100: * definite. 101: * 102: IF( D( I ).LE.ZERO ) THEN 103: INFO = I 104: GO TO 30 105: END IF 106: * 107: * Solve for e(i) and d(i+1). 108: * 109: EIR = DBLE( E( I ) ) 110: EII = DIMAG( E( I ) ) 111: F = EIR / D( I ) 112: G = EII / D( I ) 113: E( I ) = DCMPLX( F, G ) 114: D( I+1 ) = D( I+1 ) - F*EIR - G*EII 115: * 116: IF( D( I+1 ).LE.ZERO ) THEN 117: INFO = I + 1 118: GO TO 30 119: END IF 120: * 121: * Solve for e(i+1) and d(i+2). 122: * 123: EIR = DBLE( E( I+1 ) ) 124: EII = DIMAG( E( I+1 ) ) 125: F = EIR / D( I+1 ) 126: G = EII / D( I+1 ) 127: E( I+1 ) = DCMPLX( F, G ) 128: D( I+2 ) = D( I+2 ) - F*EIR - G*EII 129: * 130: IF( D( I+2 ).LE.ZERO ) THEN 131: INFO = I + 2 132: GO TO 30 133: END IF 134: * 135: * Solve for e(i+2) and d(i+3). 136: * 137: EIR = DBLE( E( I+2 ) ) 138: EII = DIMAG( E( I+2 ) ) 139: F = EIR / D( I+2 ) 140: G = EII / D( I+2 ) 141: E( I+2 ) = DCMPLX( F, G ) 142: D( I+3 ) = D( I+3 ) - F*EIR - G*EII 143: * 144: IF( D( I+3 ).LE.ZERO ) THEN 145: INFO = I + 3 146: GO TO 30 147: END IF 148: * 149: * Solve for e(i+3) and d(i+4). 150: * 151: EIR = DBLE( E( I+3 ) ) 152: EII = DIMAG( E( I+3 ) ) 153: F = EIR / D( I+3 ) 154: G = EII / D( I+3 ) 155: E( I+3 ) = DCMPLX( F, G ) 156: D( I+4 ) = D( I+4 ) - F*EIR - G*EII 157: 20 CONTINUE 158: * 159: * Check d(n) for positive definiteness. 160: * 161: IF( D( N ).LE.ZERO ) 162: $ INFO = N 163: * 164: 30 CONTINUE 165: RETURN 166: * 167: * End of ZPTTRF 168: * 169: END