LAPACK 3.3.0
|
00001 SUBROUTINE ZGTTS2( ITRANS, N, NRHS, DL, D, DU, DU2, IPIV, B, LDB ) 00002 * 00003 * -- LAPACK auxiliary 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 ITRANS, LDB, N, NRHS 00010 * .. 00011 * .. Array Arguments .. 00012 INTEGER IPIV( * ) 00013 COMPLEX*16 B( LDB, * ), D( * ), DL( * ), DU( * ), DU2( * ) 00014 * .. 00015 * 00016 * Purpose 00017 * ======= 00018 * 00019 * ZGTTS2 solves one of the systems of equations 00020 * A * X = B, A**T * X = B, or A**H * X = B, 00021 * with a tridiagonal matrix A using the LU factorization computed 00022 * by ZGTTRF. 00023 * 00024 * Arguments 00025 * ========= 00026 * 00027 * ITRANS (input) INTEGER 00028 * Specifies the form of the system of equations. 00029 * = 0: A * X = B (No transpose) 00030 * = 1: A**T * X = B (Transpose) 00031 * = 2: A**H * X = B (Conjugate transpose) 00032 * 00033 * N (input) INTEGER 00034 * The order of the matrix A. 00035 * 00036 * NRHS (input) INTEGER 00037 * The number of right hand sides, i.e., the number of columns 00038 * of the matrix B. NRHS >= 0. 00039 * 00040 * DL (input) COMPLEX*16 array, dimension (N-1) 00041 * The (n-1) multipliers that define the matrix L from the 00042 * LU factorization of A. 00043 * 00044 * D (input) COMPLEX*16 array, dimension (N) 00045 * The n diagonal elements of the upper triangular matrix U from 00046 * the LU factorization of A. 00047 * 00048 * DU (input) COMPLEX*16 array, dimension (N-1) 00049 * The (n-1) elements of the first super-diagonal of U. 00050 * 00051 * DU2 (input) COMPLEX*16 array, dimension (N-2) 00052 * The (n-2) elements of the second super-diagonal of U. 00053 * 00054 * IPIV (input) INTEGER array, dimension (N) 00055 * The pivot indices; for 1 <= i <= n, row i of the matrix was 00056 * interchanged with row IPIV(i). IPIV(i) will always be either 00057 * i or i+1; IPIV(i) = i indicates a row interchange was not 00058 * required. 00059 * 00060 * B (input/output) COMPLEX*16 array, dimension (LDB,NRHS) 00061 * On entry, the matrix of right hand side vectors B. 00062 * On exit, B is overwritten by the solution vectors X. 00063 * 00064 * LDB (input) INTEGER 00065 * The leading dimension of the array B. LDB >= max(1,N). 00066 * 00067 * ===================================================================== 00068 * 00069 * .. Local Scalars .. 00070 INTEGER I, J 00071 COMPLEX*16 TEMP 00072 * .. 00073 * .. Intrinsic Functions .. 00074 INTRINSIC DCONJG 00075 * .. 00076 * .. Executable Statements .. 00077 * 00078 * Quick return if possible 00079 * 00080 IF( N.EQ.0 .OR. NRHS.EQ.0 ) 00081 $ RETURN 00082 * 00083 IF( ITRANS.EQ.0 ) THEN 00084 * 00085 * Solve A*X = B using the LU factorization of A, 00086 * overwriting each right hand side vector with its solution. 00087 * 00088 IF( NRHS.LE.1 ) THEN 00089 J = 1 00090 10 CONTINUE 00091 * 00092 * Solve L*x = b. 00093 * 00094 DO 20 I = 1, N - 1 00095 IF( IPIV( I ).EQ.I ) THEN 00096 B( I+1, J ) = B( I+1, J ) - DL( I )*B( I, J ) 00097 ELSE 00098 TEMP = B( I, J ) 00099 B( I, J ) = B( I+1, J ) 00100 B( I+1, J ) = TEMP - DL( I )*B( I, J ) 00101 END IF 00102 20 CONTINUE 00103 * 00104 * Solve U*x = b. 00105 * 00106 B( N, J ) = B( N, J ) / D( N ) 00107 IF( N.GT.1 ) 00108 $ B( N-1, J ) = ( B( N-1, J )-DU( N-1 )*B( N, J ) ) / 00109 $ D( N-1 ) 00110 DO 30 I = N - 2, 1, -1 00111 B( I, J ) = ( B( I, J )-DU( I )*B( I+1, J )-DU2( I )* 00112 $ B( I+2, J ) ) / D( I ) 00113 30 CONTINUE 00114 IF( J.LT.NRHS ) THEN 00115 J = J + 1 00116 GO TO 10 00117 END IF 00118 ELSE 00119 DO 60 J = 1, NRHS 00120 * 00121 * Solve L*x = b. 00122 * 00123 DO 40 I = 1, N - 1 00124 IF( IPIV( I ).EQ.I ) THEN 00125 B( I+1, J ) = B( I+1, J ) - DL( I )*B( I, J ) 00126 ELSE 00127 TEMP = B( I, J ) 00128 B( I, J ) = B( I+1, J ) 00129 B( I+1, J ) = TEMP - DL( I )*B( I, J ) 00130 END IF 00131 40 CONTINUE 00132 * 00133 * Solve U*x = b. 00134 * 00135 B( N, J ) = B( N, J ) / D( N ) 00136 IF( N.GT.1 ) 00137 $ B( N-1, J ) = ( B( N-1, J )-DU( N-1 )*B( N, J ) ) / 00138 $ D( N-1 ) 00139 DO 50 I = N - 2, 1, -1 00140 B( I, J ) = ( B( I, J )-DU( I )*B( I+1, J )-DU2( I )* 00141 $ B( I+2, J ) ) / D( I ) 00142 50 CONTINUE 00143 60 CONTINUE 00144 END IF 00145 ELSE IF( ITRANS.EQ.1 ) THEN 00146 * 00147 * Solve A**T * X = B. 00148 * 00149 IF( NRHS.LE.1 ) THEN 00150 J = 1 00151 70 CONTINUE 00152 * 00153 * Solve U**T * x = b. 00154 * 00155 B( 1, J ) = B( 1, J ) / D( 1 ) 00156 IF( N.GT.1 ) 00157 $ B( 2, J ) = ( B( 2, J )-DU( 1 )*B( 1, J ) ) / D( 2 ) 00158 DO 80 I = 3, N 00159 B( I, J ) = ( B( I, J )-DU( I-1 )*B( I-1, J )-DU2( I-2 )* 00160 $ B( I-2, J ) ) / D( I ) 00161 80 CONTINUE 00162 * 00163 * Solve L**T * x = b. 00164 * 00165 DO 90 I = N - 1, 1, -1 00166 IF( IPIV( I ).EQ.I ) THEN 00167 B( I, J ) = B( I, J ) - DL( I )*B( I+1, J ) 00168 ELSE 00169 TEMP = B( I+1, J ) 00170 B( I+1, J ) = B( I, J ) - DL( I )*TEMP 00171 B( I, J ) = TEMP 00172 END IF 00173 90 CONTINUE 00174 IF( J.LT.NRHS ) THEN 00175 J = J + 1 00176 GO TO 70 00177 END IF 00178 ELSE 00179 DO 120 J = 1, NRHS 00180 * 00181 * Solve U**T * x = b. 00182 * 00183 B( 1, J ) = B( 1, J ) / D( 1 ) 00184 IF( N.GT.1 ) 00185 $ B( 2, J ) = ( B( 2, J )-DU( 1 )*B( 1, J ) ) / D( 2 ) 00186 DO 100 I = 3, N 00187 B( I, J ) = ( B( I, J )-DU( I-1 )*B( I-1, J )- 00188 $ DU2( I-2 )*B( I-2, J ) ) / D( I ) 00189 100 CONTINUE 00190 * 00191 * Solve L**T * x = b. 00192 * 00193 DO 110 I = N - 1, 1, -1 00194 IF( IPIV( I ).EQ.I ) THEN 00195 B( I, J ) = B( I, J ) - DL( I )*B( I+1, J ) 00196 ELSE 00197 TEMP = B( I+1, J ) 00198 B( I+1, J ) = B( I, J ) - DL( I )*TEMP 00199 B( I, J ) = TEMP 00200 END IF 00201 110 CONTINUE 00202 120 CONTINUE 00203 END IF 00204 ELSE 00205 * 00206 * Solve A**H * X = B. 00207 * 00208 IF( NRHS.LE.1 ) THEN 00209 J = 1 00210 130 CONTINUE 00211 * 00212 * Solve U**H * x = b. 00213 * 00214 B( 1, J ) = B( 1, J ) / DCONJG( D( 1 ) ) 00215 IF( N.GT.1 ) 00216 $ B( 2, J ) = ( B( 2, J )-DCONJG( DU( 1 ) )*B( 1, J ) ) / 00217 $ DCONJG( D( 2 ) ) 00218 DO 140 I = 3, N 00219 B( I, J ) = ( B( I, J )-DCONJG( DU( I-1 ) )*B( I-1, J )- 00220 $ DCONJG( DU2( I-2 ) )*B( I-2, J ) ) / 00221 $ DCONJG( D( I ) ) 00222 140 CONTINUE 00223 * 00224 * Solve L**H * x = b. 00225 * 00226 DO 150 I = N - 1, 1, -1 00227 IF( IPIV( I ).EQ.I ) THEN 00228 B( I, J ) = B( I, J ) - DCONJG( DL( I ) )*B( I+1, J ) 00229 ELSE 00230 TEMP = B( I+1, J ) 00231 B( I+1, J ) = B( I, J ) - DCONJG( DL( I ) )*TEMP 00232 B( I, J ) = TEMP 00233 END IF 00234 150 CONTINUE 00235 IF( J.LT.NRHS ) THEN 00236 J = J + 1 00237 GO TO 130 00238 END IF 00239 ELSE 00240 DO 180 J = 1, NRHS 00241 * 00242 * Solve U**H * x = b. 00243 * 00244 B( 1, J ) = B( 1, J ) / DCONJG( D( 1 ) ) 00245 IF( N.GT.1 ) 00246 $ B( 2, J ) = ( B( 2, J )-DCONJG( DU( 1 ) )*B( 1, J ) ) 00247 $ / DCONJG( D( 2 ) ) 00248 DO 160 I = 3, N 00249 B( I, J ) = ( B( I, J )-DCONJG( DU( I-1 ) )* 00250 $ B( I-1, J )-DCONJG( DU2( I-2 ) )* 00251 $ B( I-2, J ) ) / DCONJG( D( I ) ) 00252 160 CONTINUE 00253 * 00254 * Solve L**H * x = b. 00255 * 00256 DO 170 I = N - 1, 1, -1 00257 IF( IPIV( I ).EQ.I ) THEN 00258 B( I, J ) = B( I, J ) - DCONJG( DL( I ) )* 00259 $ B( I+1, J ) 00260 ELSE 00261 TEMP = B( I+1, J ) 00262 B( I+1, J ) = B( I, J ) - DCONJG( DL( I ) )*TEMP 00263 B( I, J ) = TEMP 00264 END IF 00265 170 CONTINUE 00266 180 CONTINUE 00267 END IF 00268 END IF 00269 * 00270 * End of ZGTTS2 00271 * 00272 END