 |
LAPACK 3.3.1
Linear Algebra PACKage
|
|
LAPACK 3.3.1
Linear Algebra PACKage
|
00001 SUBROUTINE ZPTTRF( N, D, E, INFO ) 00002 * 00003 * -- LAPACK routine (version 3.3.1) -- 00004 * -- LAPACK is a software package provided by Univ. of Tennessee, -- 00005 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 00006 * -- April 2011 -- 00007 * 00008 * .. Scalar Arguments .. 00009 INTEGER INFO, N 00010 * .. 00011 * .. Array Arguments .. 00012 DOUBLE PRECISION D( * ) 00013 COMPLEX*16 E( * ) 00014 * .. 00015 * 00016 * Purpose 00017 * ======= 00018 * 00019 * ZPTTRF computes the L*D*L**H factorization of a complex Hermitian 00020 * positive definite tridiagonal matrix A. The factorization may also 00021 * be regarded as having the form A = U**H *D*U. 00022 * 00023 * Arguments 00024 * ========= 00025 * 00026 * N (input) INTEGER 00027 * The order of the matrix A. N >= 0. 00028 * 00029 * D (input/output) DOUBLE PRECISION array, dimension (N) 00030 * On entry, the n diagonal elements of the tridiagonal matrix 00031 * A. On exit, the n diagonal elements of the diagonal matrix 00032 * D from the L*D*L**H factorization of A. 00033 * 00034 * E (input/output) COMPLEX*16 array, dimension (N-1) 00035 * On entry, the (n-1) subdiagonal elements of the tridiagonal 00036 * matrix A. On exit, the (n-1) subdiagonal elements of the 00037 * unit bidiagonal factor L from the L*D*L**H factorization of A. 00038 * E can also be regarded as the superdiagonal of the unit 00039 * bidiagonal factor U from the U**H *D*U factorization of A. 00040 * 00041 * INFO (output) INTEGER 00042 * = 0: successful exit 00043 * < 0: if INFO = -k, the k-th argument had an illegal value 00044 * > 0: if INFO = k, the leading minor of order k is not 00045 * positive definite; if k < N, the factorization could not 00046 * be completed, while if k = N, the factorization was 00047 * completed, but D(N) <= 0. 00048 * 00049 * ===================================================================== 00050 * 00051 * .. Parameters .. 00052 DOUBLE PRECISION ZERO 00053 PARAMETER ( ZERO = 0.0D+0 ) 00054 * .. 00055 * .. Local Scalars .. 00056 INTEGER I, I4 00057 DOUBLE PRECISION EII, EIR, F, G 00058 * .. 00059 * .. External Subroutines .. 00060 EXTERNAL XERBLA 00061 * .. 00062 * .. Intrinsic Functions .. 00063 INTRINSIC DBLE, DCMPLX, DIMAG, MOD 00064 * .. 00065 * .. Executable Statements .. 00066 * 00067 * Test the input parameters. 00068 * 00069 INFO = 0 00070 IF( N.LT.0 ) THEN 00071 INFO = -1 00072 CALL XERBLA( 'ZPTTRF', -INFO ) 00073 RETURN 00074 END IF 00075 * 00076 * Quick return if possible 00077 * 00078 IF( N.EQ.0 ) 00079 $ RETURN 00080 * 00081 * Compute the L*D*L**H (or U**H *D*U) factorization of A. 00082 * 00083 I4 = MOD( N-1, 4 ) 00084 DO 10 I = 1, I4 00085 IF( D( I ).LE.ZERO ) THEN 00086 INFO = I 00087 GO TO 30 00088 END IF 00089 EIR = DBLE( E( I ) ) 00090 EII = DIMAG( E( I ) ) 00091 F = EIR / D( I ) 00092 G = EII / D( I ) 00093 E( I ) = DCMPLX( F, G ) 00094 D( I+1 ) = D( I+1 ) - F*EIR - G*EII 00095 10 CONTINUE 00096 * 00097 DO 20 I = I4 + 1, N - 4, 4 00098 * 00099 * Drop out of the loop if d(i) <= 0: the matrix is not positive 00100 * definite. 00101 * 00102 IF( D( I ).LE.ZERO ) THEN 00103 INFO = I 00104 GO TO 30 00105 END IF 00106 * 00107 * Solve for e(i) and d(i+1). 00108 * 00109 EIR = DBLE( E( I ) ) 00110 EII = DIMAG( E( I ) ) 00111 F = EIR / D( I ) 00112 G = EII / D( I ) 00113 E( I ) = DCMPLX( F, G ) 00114 D( I+1 ) = D( I+1 ) - F*EIR - G*EII 00115 * 00116 IF( D( I+1 ).LE.ZERO ) THEN 00117 INFO = I + 1 00118 GO TO 30 00119 END IF 00120 * 00121 * Solve for e(i+1) and d(i+2). 00122 * 00123 EIR = DBLE( E( I+1 ) ) 00124 EII = DIMAG( E( I+1 ) ) 00125 F = EIR / D( I+1 ) 00126 G = EII / D( I+1 ) 00127 E( I+1 ) = DCMPLX( F, G ) 00128 D( I+2 ) = D( I+2 ) - F*EIR - G*EII 00129 * 00130 IF( D( I+2 ).LE.ZERO ) THEN 00131 INFO = I + 2 00132 GO TO 30 00133 END IF 00134 * 00135 * Solve for e(i+2) and d(i+3). 00136 * 00137 EIR = DBLE( E( I+2 ) ) 00138 EII = DIMAG( E( I+2 ) ) 00139 F = EIR / D( I+2 ) 00140 G = EII / D( I+2 ) 00141 E( I+2 ) = DCMPLX( F, G ) 00142 D( I+3 ) = D( I+3 ) - F*EIR - G*EII 00143 * 00144 IF( D( I+3 ).LE.ZERO ) THEN 00145 INFO = I + 3 00146 GO TO 30 00147 END IF 00148 * 00149 * Solve for e(i+3) and d(i+4). 00150 * 00151 EIR = DBLE( E( I+3 ) ) 00152 EII = DIMAG( E( I+3 ) ) 00153 F = EIR / D( I+3 ) 00154 G = EII / D( I+3 ) 00155 E( I+3 ) = DCMPLX( F, G ) 00156 D( I+4 ) = D( I+4 ) - F*EIR - G*EII 00157 20 CONTINUE 00158 * 00159 * Check d(n) for positive definiteness. 00160 * 00161 IF( D( N ).LE.ZERO ) 00162 $ INFO = N 00163 * 00164 30 CONTINUE 00165 RETURN 00166 * 00167 * End of ZPTTRF 00168 * 00169 END
1.7.3 
