 |
LAPACK 3.3.1
Linear Algebra PACKage
|
|
LAPACK 3.3.1
Linear Algebra PACKage
|
00001 SUBROUTINE DTRMM(SIDE,UPLO,TRANSA,DIAG,M,N,ALPHA,A,LDA,B,LDB) 00002 * .. Scalar Arguments .. 00003 DOUBLE PRECISION ALPHA 00004 INTEGER LDA,LDB,M,N 00005 CHARACTER DIAG,SIDE,TRANSA,UPLO 00006 * .. 00007 * .. Array Arguments .. 00008 DOUBLE PRECISION A(LDA,*),B(LDB,*) 00009 * .. 00010 * 00011 * Purpose 00012 * ======= 00013 * 00014 * DTRMM performs one of the matrix-matrix operations 00015 * 00016 * B := alpha*op( A )*B, or B := alpha*B*op( A ), 00017 * 00018 * where alpha is a scalar, B is an m by n matrix, A is a unit, or 00019 * non-unit, upper or lower triangular matrix and op( A ) is one of 00020 * 00021 * op( A ) = A or op( A ) = A**T. 00022 * 00023 * Arguments 00024 * ========== 00025 * 00026 * SIDE - CHARACTER*1. 00027 * On entry, SIDE specifies whether op( A ) multiplies B from 00028 * the left or right as follows: 00029 * 00030 * SIDE = 'L' or 'l' B := alpha*op( A )*B. 00031 * 00032 * SIDE = 'R' or 'r' B := alpha*B*op( A ). 00033 * 00034 * Unchanged on exit. 00035 * 00036 * UPLO - CHARACTER*1. 00037 * On entry, UPLO specifies whether the matrix A is an upper or 00038 * lower triangular matrix as follows: 00039 * 00040 * UPLO = 'U' or 'u' A is an upper triangular matrix. 00041 * 00042 * UPLO = 'L' or 'l' A is a lower triangular matrix. 00043 * 00044 * Unchanged on exit. 00045 * 00046 * TRANSA - CHARACTER*1. 00047 * On entry, TRANSA specifies the form of op( A ) to be used in 00048 * the matrix multiplication as follows: 00049 * 00050 * TRANSA = 'N' or 'n' op( A ) = A. 00051 * 00052 * TRANSA = 'T' or 't' op( A ) = A**T. 00053 * 00054 * TRANSA = 'C' or 'c' op( A ) = A**T. 00055 * 00056 * Unchanged on exit. 00057 * 00058 * DIAG - CHARACTER*1. 00059 * On entry, DIAG specifies whether or not A is unit triangular 00060 * as follows: 00061 * 00062 * DIAG = 'U' or 'u' A is assumed to be unit triangular. 00063 * 00064 * DIAG = 'N' or 'n' A is not assumed to be unit 00065 * triangular. 00066 * 00067 * Unchanged on exit. 00068 * 00069 * M - INTEGER. 00070 * On entry, M specifies the number of rows of B. M must be at 00071 * least zero. 00072 * Unchanged on exit. 00073 * 00074 * N - INTEGER. 00075 * On entry, N specifies the number of columns of B. N must be 00076 * at least zero. 00077 * Unchanged on exit. 00078 * 00079 * ALPHA - DOUBLE PRECISION. 00080 * On entry, ALPHA specifies the scalar alpha. When alpha is 00081 * zero then A is not referenced and B need not be set before 00082 * entry. 00083 * Unchanged on exit. 00084 * 00085 * A - DOUBLE PRECISION array of DIMENSION ( LDA, k ), where k is m 00086 * when SIDE = 'L' or 'l' and is n when SIDE = 'R' or 'r'. 00087 * Before entry with UPLO = 'U' or 'u', the leading k by k 00088 * upper triangular part of the array A must contain the upper 00089 * triangular matrix and the strictly lower triangular part of 00090 * A is not referenced. 00091 * Before entry with UPLO = 'L' or 'l', the leading k by k 00092 * lower triangular part of the array A must contain the lower 00093 * triangular matrix and the strictly upper triangular part of 00094 * A is not referenced. 00095 * Note that when DIAG = 'U' or 'u', the diagonal elements of 00096 * A are not referenced either, but are assumed to be unity. 00097 * Unchanged on exit. 00098 * 00099 * LDA - INTEGER. 00100 * On entry, LDA specifies the first dimension of A as declared 00101 * in the calling (sub) program. When SIDE = 'L' or 'l' then 00102 * LDA must be at least max( 1, m ), when SIDE = 'R' or 'r' 00103 * then LDA must be at least max( 1, n ). 00104 * Unchanged on exit. 00105 * 00106 * B - DOUBLE PRECISION array of DIMENSION ( LDB, n ). 00107 * Before entry, the leading m by n part of the array B must 00108 * contain the matrix B, and on exit is overwritten by the 00109 * transformed matrix. 00110 * 00111 * LDB - INTEGER. 00112 * On entry, LDB specifies the first dimension of B as declared 00113 * in the calling (sub) program. LDB must be at least 00114 * max( 1, m ). 00115 * Unchanged on exit. 00116 * 00117 * Further Details 00118 * =============== 00119 * 00120 * Level 3 Blas routine. 00121 * 00122 * -- Written on 8-February-1989. 00123 * Jack Dongarra, Argonne National Laboratory. 00124 * Iain Duff, AERE Harwell. 00125 * Jeremy Du Croz, Numerical Algorithms Group Ltd. 00126 * Sven Hammarling, Numerical Algorithms Group Ltd. 00127 * 00128 * ===================================================================== 00129 * 00130 * .. External Functions .. 00131 LOGICAL LSAME 00132 EXTERNAL LSAME 00133 * .. 00134 * .. External Subroutines .. 00135 EXTERNAL XERBLA 00136 * .. 00137 * .. Intrinsic Functions .. 00138 INTRINSIC MAX 00139 * .. 00140 * .. Local Scalars .. 00141 DOUBLE PRECISION TEMP 00142 INTEGER I,INFO,J,K,NROWA 00143 LOGICAL LSIDE,NOUNIT,UPPER 00144 * .. 00145 * .. Parameters .. 00146 DOUBLE PRECISION ONE,ZERO 00147 PARAMETER (ONE=1.0D+0,ZERO=0.0D+0) 00148 * .. 00149 * 00150 * Test the input parameters. 00151 * 00152 LSIDE = LSAME(SIDE,'L') 00153 IF (LSIDE) THEN 00154 NROWA = M 00155 ELSE 00156 NROWA = N 00157 END IF 00158 NOUNIT = LSAME(DIAG,'N') 00159 UPPER = LSAME(UPLO,'U') 00160 * 00161 INFO = 0 00162 IF ((.NOT.LSIDE) .AND. (.NOT.LSAME(SIDE,'R'))) THEN 00163 INFO = 1 00164 ELSE IF ((.NOT.UPPER) .AND. (.NOT.LSAME(UPLO,'L'))) THEN 00165 INFO = 2 00166 ELSE IF ((.NOT.LSAME(TRANSA,'N')) .AND. 00167 + (.NOT.LSAME(TRANSA,'T')) .AND. 00168 + (.NOT.LSAME(TRANSA,'C'))) THEN 00169 INFO = 3 00170 ELSE IF ((.NOT.LSAME(DIAG,'U')) .AND. (.NOT.LSAME(DIAG,'N'))) THEN 00171 INFO = 4 00172 ELSE IF (M.LT.0) THEN 00173 INFO = 5 00174 ELSE IF (N.LT.0) THEN 00175 INFO = 6 00176 ELSE IF (LDA.LT.MAX(1,NROWA)) THEN 00177 INFO = 9 00178 ELSE IF (LDB.LT.MAX(1,M)) THEN 00179 INFO = 11 00180 END IF 00181 IF (INFO.NE.0) THEN 00182 CALL XERBLA('DTRMM ',INFO) 00183 RETURN 00184 END IF 00185 * 00186 * Quick return if possible. 00187 * 00188 IF (M.EQ.0 .OR. N.EQ.0) RETURN 00189 * 00190 * And when alpha.eq.zero. 00191 * 00192 IF (ALPHA.EQ.ZERO) THEN 00193 DO 20 J = 1,N 00194 DO 10 I = 1,M 00195 B(I,J) = ZERO 00196 10 CONTINUE 00197 20 CONTINUE 00198 RETURN 00199 END IF 00200 * 00201 * Start the operations. 00202 * 00203 IF (LSIDE) THEN 00204 IF (LSAME(TRANSA,'N')) THEN 00205 * 00206 * Form B := alpha*A*B. 00207 * 00208 IF (UPPER) THEN 00209 DO 50 J = 1,N 00210 DO 40 K = 1,M 00211 IF (B(K,J).NE.ZERO) THEN 00212 TEMP = ALPHA*B(K,J) 00213 DO 30 I = 1,K - 1 00214 B(I,J) = B(I,J) + TEMP*A(I,K) 00215 30 CONTINUE 00216 IF (NOUNIT) TEMP = TEMP*A(K,K) 00217 B(K,J) = TEMP 00218 END IF 00219 40 CONTINUE 00220 50 CONTINUE 00221 ELSE 00222 DO 80 J = 1,N 00223 DO 70 K = M,1,-1 00224 IF (B(K,J).NE.ZERO) THEN 00225 TEMP = ALPHA*B(K,J) 00226 B(K,J) = TEMP 00227 IF (NOUNIT) B(K,J) = B(K,J)*A(K,K) 00228 DO 60 I = K + 1,M 00229 B(I,J) = B(I,J) + TEMP*A(I,K) 00230 60 CONTINUE 00231 END IF 00232 70 CONTINUE 00233 80 CONTINUE 00234 END IF 00235 ELSE 00236 * 00237 * Form B := alpha*A**T*B. 00238 * 00239 IF (UPPER) THEN 00240 DO 110 J = 1,N 00241 DO 100 I = M,1,-1 00242 TEMP = B(I,J) 00243 IF (NOUNIT) TEMP = TEMP*A(I,I) 00244 DO 90 K = 1,I - 1 00245 TEMP = TEMP + A(K,I)*B(K,J) 00246 90 CONTINUE 00247 B(I,J) = ALPHA*TEMP 00248 100 CONTINUE 00249 110 CONTINUE 00250 ELSE 00251 DO 140 J = 1,N 00252 DO 130 I = 1,M 00253 TEMP = B(I,J) 00254 IF (NOUNIT) TEMP = TEMP*A(I,I) 00255 DO 120 K = I + 1,M 00256 TEMP = TEMP + A(K,I)*B(K,J) 00257 120 CONTINUE 00258 B(I,J) = ALPHA*TEMP 00259 130 CONTINUE 00260 140 CONTINUE 00261 END IF 00262 END IF 00263 ELSE 00264 IF (LSAME(TRANSA,'N')) THEN 00265 * 00266 * Form B := alpha*B*A. 00267 * 00268 IF (UPPER) THEN 00269 DO 180 J = N,1,-1 00270 TEMP = ALPHA 00271 IF (NOUNIT) TEMP = TEMP*A(J,J) 00272 DO 150 I = 1,M 00273 B(I,J) = TEMP*B(I,J) 00274 150 CONTINUE 00275 DO 170 K = 1,J - 1 00276 IF (A(K,J).NE.ZERO) THEN 00277 TEMP = ALPHA*A(K,J) 00278 DO 160 I = 1,M 00279 B(I,J) = B(I,J) + TEMP*B(I,K) 00280 160 CONTINUE 00281 END IF 00282 170 CONTINUE 00283 180 CONTINUE 00284 ELSE 00285 DO 220 J = 1,N 00286 TEMP = ALPHA 00287 IF (NOUNIT) TEMP = TEMP*A(J,J) 00288 DO 190 I = 1,M 00289 B(I,J) = TEMP*B(I,J) 00290 190 CONTINUE 00291 DO 210 K = J + 1,N 00292 IF (A(K,J).NE.ZERO) THEN 00293 TEMP = ALPHA*A(K,J) 00294 DO 200 I = 1,M 00295 B(I,J) = B(I,J) + TEMP*B(I,K) 00296 200 CONTINUE 00297 END IF 00298 210 CONTINUE 00299 220 CONTINUE 00300 END IF 00301 ELSE 00302 * 00303 * Form B := alpha*B*A**T. 00304 * 00305 IF (UPPER) THEN 00306 DO 260 K = 1,N 00307 DO 240 J = 1,K - 1 00308 IF (A(J,K).NE.ZERO) THEN 00309 TEMP = ALPHA*A(J,K) 00310 DO 230 I = 1,M 00311 B(I,J) = B(I,J) + TEMP*B(I,K) 00312 230 CONTINUE 00313 END IF 00314 240 CONTINUE 00315 TEMP = ALPHA 00316 IF (NOUNIT) TEMP = TEMP*A(K,K) 00317 IF (TEMP.NE.ONE) THEN 00318 DO 250 I = 1,M 00319 B(I,K) = TEMP*B(I,K) 00320 250 CONTINUE 00321 END IF 00322 260 CONTINUE 00323 ELSE 00324 DO 300 K = N,1,-1 00325 DO 280 J = K + 1,N 00326 IF (A(J,K).NE.ZERO) THEN 00327 TEMP = ALPHA*A(J,K) 00328 DO 270 I = 1,M 00329 B(I,J) = B(I,J) + TEMP*B(I,K) 00330 270 CONTINUE 00331 END IF 00332 280 CONTINUE 00333 TEMP = ALPHA 00334 IF (NOUNIT) TEMP = TEMP*A(K,K) 00335 IF (TEMP.NE.ONE) THEN 00336 DO 290 I = 1,M 00337 B(I,K) = TEMP*B(I,K) 00338 290 CONTINUE 00339 END IF 00340 300 CONTINUE 00341 END IF 00342 END IF 00343 END IF 00344 * 00345 RETURN 00346 * 00347 * End of DTRMM . 00348 * 00349 END
1.7.3 
