00001 SUBROUTINE CGTTRF( N, DL, D, DU, DU2, IPIV, INFO ) 00002 * 00003 * -- LAPACK routine (version 3.2) -- 00004 * -- LAPACK is a software package provided by Univ. of Tennessee, -- 00005 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 00006 * November 2006 00007 * 00008 * .. Scalar Arguments .. 00009 INTEGER INFO, N 00010 * .. 00011 * .. Array Arguments .. 00012 INTEGER IPIV( * ) 00013 COMPLEX D( * ), DL( * ), DU( * ), DU2( * ) 00014 * .. 00015 * 00016 * Purpose 00017 * ======= 00018 * 00019 * CGTTRF computes an LU factorization of a complex tridiagonal matrix A 00020 * using elimination with partial pivoting and row interchanges. 00021 * 00022 * The factorization has the form 00023 * A = L * U 00024 * where L is a product of permutation and unit lower bidiagonal 00025 * matrices and U is upper triangular with nonzeros in only the main 00026 * diagonal and first two superdiagonals. 00027 * 00028 * Arguments 00029 * ========= 00030 * 00031 * N (input) INTEGER 00032 * The order of the matrix A. 00033 * 00034 * DL (input/output) COMPLEX array, dimension (N-1) 00035 * On entry, DL must contain the (n-1) sub-diagonal elements of 00036 * A. 00037 * 00038 * On exit, DL is overwritten by the (n-1) multipliers that 00039 * define the matrix L from the LU factorization of A. 00040 * 00041 * D (input/output) COMPLEX array, dimension (N) 00042 * On entry, D must contain the diagonal elements of A. 00043 * 00044 * On exit, D is overwritten by the n diagonal elements of the 00045 * upper triangular matrix U from the LU factorization of A. 00046 * 00047 * DU (input/output) COMPLEX array, dimension (N-1) 00048 * On entry, DU must contain the (n-1) super-diagonal elements 00049 * of A. 00050 * 00051 * On exit, DU is overwritten by the (n-1) elements of the first 00052 * super-diagonal of U. 00053 * 00054 * DU2 (output) COMPLEX array, dimension (N-2) 00055 * On exit, DU2 is overwritten by the (n-2) elements of the 00056 * second super-diagonal of U. 00057 * 00058 * IPIV (output) INTEGER array, dimension (N) 00059 * The pivot indices; for 1 <= i <= n, row i of the matrix was 00060 * interchanged with row IPIV(i). IPIV(i) will always be either 00061 * i or i+1; IPIV(i) = i indicates a row interchange was not 00062 * required. 00063 * 00064 * INFO (output) INTEGER 00065 * = 0: successful exit 00066 * < 0: if INFO = -k, the k-th argument had an illegal value 00067 * > 0: if INFO = k, U(k,k) is exactly zero. The factorization 00068 * has been completed, but the factor U is exactly 00069 * singular, and division by zero will occur if it is used 00070 * to solve a system of equations. 00071 * 00072 * ===================================================================== 00073 * 00074 * .. Parameters .. 00075 REAL ZERO 00076 PARAMETER ( ZERO = 0.0E+0 ) 00077 * .. 00078 * .. Local Scalars .. 00079 INTEGER I 00080 COMPLEX FACT, TEMP, ZDUM 00081 * .. 00082 * .. External Subroutines .. 00083 EXTERNAL XERBLA 00084 * .. 00085 * .. Intrinsic Functions .. 00086 INTRINSIC ABS, AIMAG, REAL 00087 * .. 00088 * .. Statement Functions .. 00089 REAL CABS1 00090 * .. 00091 * .. Statement Function definitions .. 00092 CABS1( ZDUM ) = ABS( REAL( ZDUM ) ) + ABS( AIMAG( ZDUM ) ) 00093 * .. 00094 * .. Executable Statements .. 00095 * 00096 INFO = 0 00097 IF( N.LT.0 ) THEN 00098 INFO = -1 00099 CALL XERBLA( 'CGTTRF', -INFO ) 00100 RETURN 00101 END IF 00102 * 00103 * Quick return if possible 00104 * 00105 IF( N.EQ.0 ) 00106 $ RETURN 00107 * 00108 * Initialize IPIV(i) = i and DU2(i) = 0 00109 * 00110 DO 10 I = 1, N 00111 IPIV( I ) = I 00112 10 CONTINUE 00113 DO 20 I = 1, N - 2 00114 DU2( I ) = ZERO 00115 20 CONTINUE 00116 * 00117 DO 30 I = 1, N - 2 00118 IF( CABS1( D( I ) ).GE.CABS1( DL( I ) ) ) THEN 00119 * 00120 * No row interchange required, eliminate DL(I) 00121 * 00122 IF( CABS1( D( I ) ).NE.ZERO ) THEN 00123 FACT = DL( I ) / D( I ) 00124 DL( I ) = FACT 00125 D( I+1 ) = D( I+1 ) - FACT*DU( I ) 00126 END IF 00127 ELSE 00128 * 00129 * Interchange rows I and I+1, eliminate DL(I) 00130 * 00131 FACT = D( I ) / DL( I ) 00132 D( I ) = DL( I ) 00133 DL( I ) = FACT 00134 TEMP = DU( I ) 00135 DU( I ) = D( I+1 ) 00136 D( I+1 ) = TEMP - FACT*D( I+1 ) 00137 DU2( I ) = DU( I+1 ) 00138 DU( I+1 ) = -FACT*DU( I+1 ) 00139 IPIV( I ) = I + 1 00140 END IF 00141 30 CONTINUE 00142 IF( N.GT.1 ) THEN 00143 I = N - 1 00144 IF( CABS1( D( I ) ).GE.CABS1( DL( I ) ) ) THEN 00145 IF( CABS1( D( I ) ).NE.ZERO ) THEN 00146 FACT = DL( I ) / D( I ) 00147 DL( I ) = FACT 00148 D( I+1 ) = D( I+1 ) - FACT*DU( I ) 00149 END IF 00150 ELSE 00151 FACT = D( I ) / DL( I ) 00152 D( I ) = DL( I ) 00153 DL( I ) = FACT 00154 TEMP = DU( I ) 00155 DU( I ) = D( I+1 ) 00156 D( I+1 ) = TEMP - FACT*D( I+1 ) 00157 IPIV( I ) = I + 1 00158 END IF 00159 END IF 00160 * 00161 * Check for a zero on the diagonal of U. 00162 * 00163 DO 40 I = 1, N 00164 IF( CABS1( D( I ) ).EQ.ZERO ) THEN 00165 INFO = I 00166 GO TO 50 00167 END IF 00168 40 CONTINUE 00169 50 CONTINUE 00170 * 00171 RETURN 00172 * 00173 * End of CGTTRF 00174 * 00175 END
1.7.2