 |
LAPACK 3.3.1
Linear Algebra PACKage
|
|
LAPACK 3.3.1
Linear Algebra PACKage
|
00001 SUBROUTINE ZGEMM(TRANSA,TRANSB,M,N,K,ALPHA,A,LDA,B,LDB,BETA,C,LDC) 00002 * .. Scalar Arguments .. 00003 DOUBLE COMPLEX ALPHA,BETA 00004 INTEGER K,LDA,LDB,LDC,M,N 00005 CHARACTER TRANSA,TRANSB 00006 * .. 00007 * .. Array Arguments .. 00008 DOUBLE COMPLEX A(LDA,*),B(LDB,*),C(LDC,*) 00009 * .. 00010 * 00011 * Purpose 00012 * ======= 00013 * 00014 * ZGEMM performs one of the matrix-matrix operations 00015 * 00016 * C := alpha*op( A )*op( B ) + beta*C, 00017 * 00018 * where op( X ) is one of 00019 * 00020 * op( X ) = X or op( X ) = X**T or op( X ) = X**H, 00021 * 00022 * alpha and beta are scalars, and A, B and C are matrices, with op( A ) 00023 * an m by k matrix, op( B ) a k by n matrix and C an m by n matrix. 00024 * 00025 * Arguments 00026 * ========== 00027 * 00028 * TRANSA - CHARACTER*1. 00029 * On entry, TRANSA specifies the form of op( A ) to be used in 00030 * the matrix multiplication as follows: 00031 * 00032 * TRANSA = 'N' or 'n', op( A ) = A. 00033 * 00034 * TRANSA = 'T' or 't', op( A ) = A**T. 00035 * 00036 * TRANSA = 'C' or 'c', op( A ) = A**H. 00037 * 00038 * Unchanged on exit. 00039 * 00040 * TRANSB - CHARACTER*1. 00041 * On entry, TRANSB specifies the form of op( B ) to be used in 00042 * the matrix multiplication as follows: 00043 * 00044 * TRANSB = 'N' or 'n', op( B ) = B. 00045 * 00046 * TRANSB = 'T' or 't', op( B ) = B**T. 00047 * 00048 * TRANSB = 'C' or 'c', op( B ) = B**H. 00049 * 00050 * Unchanged on exit. 00051 * 00052 * M - INTEGER. 00053 * On entry, M specifies the number of rows of the matrix 00054 * op( A ) and of the matrix C. M must be at least zero. 00055 * Unchanged on exit. 00056 * 00057 * N - INTEGER. 00058 * On entry, N specifies the number of columns of the matrix 00059 * op( B ) and the number of columns of the matrix C. N must be 00060 * at least zero. 00061 * Unchanged on exit. 00062 * 00063 * K - INTEGER. 00064 * On entry, K specifies the number of columns of the matrix 00065 * op( A ) and the number of rows of the matrix op( B ). K must 00066 * be at least zero. 00067 * Unchanged on exit. 00068 * 00069 * ALPHA - COMPLEX*16 . 00070 * On entry, ALPHA specifies the scalar alpha. 00071 * Unchanged on exit. 00072 * 00073 * A - COMPLEX*16 array of DIMENSION ( LDA, ka ), where ka is 00074 * k when TRANSA = 'N' or 'n', and is m otherwise. 00075 * Before entry with TRANSA = 'N' or 'n', the leading m by k 00076 * part of the array A must contain the matrix A, otherwise 00077 * the leading k by m part of the array A must contain the 00078 * matrix A. 00079 * Unchanged on exit. 00080 * 00081 * LDA - INTEGER. 00082 * On entry, LDA specifies the first dimension of A as declared 00083 * in the calling (sub) program. When TRANSA = 'N' or 'n' then 00084 * LDA must be at least max( 1, m ), otherwise LDA must be at 00085 * least max( 1, k ). 00086 * Unchanged on exit. 00087 * 00088 * B - COMPLEX*16 array of DIMENSION ( LDB, kb ), where kb is 00089 * n when TRANSB = 'N' or 'n', and is k otherwise. 00090 * Before entry with TRANSB = 'N' or 'n', the leading k by n 00091 * part of the array B must contain the matrix B, otherwise 00092 * the leading n by k part of the array B must contain the 00093 * matrix B. 00094 * Unchanged on exit. 00095 * 00096 * LDB - INTEGER. 00097 * On entry, LDB specifies the first dimension of B as declared 00098 * in the calling (sub) program. When TRANSB = 'N' or 'n' then 00099 * LDB must be at least max( 1, k ), otherwise LDB must be at 00100 * least max( 1, n ). 00101 * Unchanged on exit. 00102 * 00103 * BETA - COMPLEX*16 . 00104 * On entry, BETA specifies the scalar beta. When BETA is 00105 * supplied as zero then C need not be set on input. 00106 * Unchanged on exit. 00107 * 00108 * C - COMPLEX*16 array of DIMENSION ( LDC, n ). 00109 * Before entry, the leading m by n part of the array C must 00110 * contain the matrix C, except when beta is zero, in which 00111 * case C need not be set on entry. 00112 * On exit, the array C is overwritten by the m by n matrix 00113 * ( alpha*op( A )*op( B ) + beta*C ). 00114 * 00115 * LDC - INTEGER. 00116 * On entry, LDC specifies the first dimension of C as declared 00117 * in the calling (sub) program. LDC must be at least 00118 * max( 1, m ). 00119 * Unchanged on exit. 00120 * 00121 * Further Details 00122 * =============== 00123 * 00124 * Level 3 Blas routine. 00125 * 00126 * -- Written on 8-February-1989. 00127 * Jack Dongarra, Argonne National Laboratory. 00128 * Iain Duff, AERE Harwell. 00129 * Jeremy Du Croz, Numerical Algorithms Group Ltd. 00130 * Sven Hammarling, Numerical Algorithms Group Ltd. 00131 * 00132 * ===================================================================== 00133 * 00134 * .. External Functions .. 00135 LOGICAL LSAME 00136 EXTERNAL LSAME 00137 * .. 00138 * .. External Subroutines .. 00139 EXTERNAL XERBLA 00140 * .. 00141 * .. Intrinsic Functions .. 00142 INTRINSIC DCONJG,MAX 00143 * .. 00144 * .. Local Scalars .. 00145 DOUBLE COMPLEX TEMP 00146 INTEGER I,INFO,J,L,NCOLA,NROWA,NROWB 00147 LOGICAL CONJA,CONJB,NOTA,NOTB 00148 * .. 00149 * .. Parameters .. 00150 DOUBLE COMPLEX ONE 00151 PARAMETER (ONE= (1.0D+0,0.0D+0)) 00152 DOUBLE COMPLEX ZERO 00153 PARAMETER (ZERO= (0.0D+0,0.0D+0)) 00154 * .. 00155 * 00156 * Set NOTA and NOTB as true if A and B respectively are not 00157 * conjugated or transposed, set CONJA and CONJB as true if A and 00158 * B respectively are to be transposed but not conjugated and set 00159 * NROWA, NCOLA and NROWB as the number of rows and columns of A 00160 * and the number of rows of B respectively. 00161 * 00162 NOTA = LSAME(TRANSA,'N') 00163 NOTB = LSAME(TRANSB,'N') 00164 CONJA = LSAME(TRANSA,'C') 00165 CONJB = LSAME(TRANSB,'C') 00166 IF (NOTA) THEN 00167 NROWA = M 00168 NCOLA = K 00169 ELSE 00170 NROWA = K 00171 NCOLA = M 00172 END IF 00173 IF (NOTB) THEN 00174 NROWB = K 00175 ELSE 00176 NROWB = N 00177 END IF 00178 * 00179 * Test the input parameters. 00180 * 00181 INFO = 0 00182 IF ((.NOT.NOTA) .AND. (.NOT.CONJA) .AND. 00183 + (.NOT.LSAME(TRANSA,'T'))) THEN 00184 INFO = 1 00185 ELSE IF ((.NOT.NOTB) .AND. (.NOT.CONJB) .AND. 00186 + (.NOT.LSAME(TRANSB,'T'))) THEN 00187 INFO = 2 00188 ELSE IF (M.LT.0) THEN 00189 INFO = 3 00190 ELSE IF (N.LT.0) THEN 00191 INFO = 4 00192 ELSE IF (K.LT.0) THEN 00193 INFO = 5 00194 ELSE IF (LDA.LT.MAX(1,NROWA)) THEN 00195 INFO = 8 00196 ELSE IF (LDB.LT.MAX(1,NROWB)) THEN 00197 INFO = 10 00198 ELSE IF (LDC.LT.MAX(1,M)) THEN 00199 INFO = 13 00200 END IF 00201 IF (INFO.NE.0) THEN 00202 CALL XERBLA('ZGEMM ',INFO) 00203 RETURN 00204 END IF 00205 * 00206 * Quick return if possible. 00207 * 00208 IF ((M.EQ.0) .OR. (N.EQ.0) .OR. 00209 + (((ALPHA.EQ.ZERO).OR. (K.EQ.0)).AND. (BETA.EQ.ONE))) RETURN 00210 * 00211 * And when alpha.eq.zero. 00212 * 00213 IF (ALPHA.EQ.ZERO) THEN 00214 IF (BETA.EQ.ZERO) THEN 00215 DO 20 J = 1,N 00216 DO 10 I = 1,M 00217 C(I,J) = ZERO 00218 10 CONTINUE 00219 20 CONTINUE 00220 ELSE 00221 DO 40 J = 1,N 00222 DO 30 I = 1,M 00223 C(I,J) = BETA*C(I,J) 00224 30 CONTINUE 00225 40 CONTINUE 00226 END IF 00227 RETURN 00228 END IF 00229 * 00230 * Start the operations. 00231 * 00232 IF (NOTB) THEN 00233 IF (NOTA) THEN 00234 * 00235 * Form C := alpha*A*B + beta*C. 00236 * 00237 DO 90 J = 1,N 00238 IF (BETA.EQ.ZERO) THEN 00239 DO 50 I = 1,M 00240 C(I,J) = ZERO 00241 50 CONTINUE 00242 ELSE IF (BETA.NE.ONE) THEN 00243 DO 60 I = 1,M 00244 C(I,J) = BETA*C(I,J) 00245 60 CONTINUE 00246 END IF 00247 DO 80 L = 1,K 00248 IF (B(L,J).NE.ZERO) THEN 00249 TEMP = ALPHA*B(L,J) 00250 DO 70 I = 1,M 00251 C(I,J) = C(I,J) + TEMP*A(I,L) 00252 70 CONTINUE 00253 END IF 00254 80 CONTINUE 00255 90 CONTINUE 00256 ELSE IF (CONJA) THEN 00257 * 00258 * Form C := alpha*A**H*B + beta*C. 00259 * 00260 DO 120 J = 1,N 00261 DO 110 I = 1,M 00262 TEMP = ZERO 00263 DO 100 L = 1,K 00264 TEMP = TEMP + DCONJG(A(L,I))*B(L,J) 00265 100 CONTINUE 00266 IF (BETA.EQ.ZERO) THEN 00267 C(I,J) = ALPHA*TEMP 00268 ELSE 00269 C(I,J) = ALPHA*TEMP + BETA*C(I,J) 00270 END IF 00271 110 CONTINUE 00272 120 CONTINUE 00273 ELSE 00274 * 00275 * Form C := alpha*A**T*B + beta*C 00276 * 00277 DO 150 J = 1,N 00278 DO 140 I = 1,M 00279 TEMP = ZERO 00280 DO 130 L = 1,K 00281 TEMP = TEMP + A(L,I)*B(L,J) 00282 130 CONTINUE 00283 IF (BETA.EQ.ZERO) THEN 00284 C(I,J) = ALPHA*TEMP 00285 ELSE 00286 C(I,J) = ALPHA*TEMP + BETA*C(I,J) 00287 END IF 00288 140 CONTINUE 00289 150 CONTINUE 00290 END IF 00291 ELSE IF (NOTA) THEN 00292 IF (CONJB) THEN 00293 * 00294 * Form C := alpha*A*B**H + beta*C. 00295 * 00296 DO 200 J = 1,N 00297 IF (BETA.EQ.ZERO) THEN 00298 DO 160 I = 1,M 00299 C(I,J) = ZERO 00300 160 CONTINUE 00301 ELSE IF (BETA.NE.ONE) THEN 00302 DO 170 I = 1,M 00303 C(I,J) = BETA*C(I,J) 00304 170 CONTINUE 00305 END IF 00306 DO 190 L = 1,K 00307 IF (B(J,L).NE.ZERO) THEN 00308 TEMP = ALPHA*DCONJG(B(J,L)) 00309 DO 180 I = 1,M 00310 C(I,J) = C(I,J) + TEMP*A(I,L) 00311 180 CONTINUE 00312 END IF 00313 190 CONTINUE 00314 200 CONTINUE 00315 ELSE 00316 * 00317 * Form C := alpha*A*B**T + beta*C 00318 * 00319 DO 250 J = 1,N 00320 IF (BETA.EQ.ZERO) THEN 00321 DO 210 I = 1,M 00322 C(I,J) = ZERO 00323 210 CONTINUE 00324 ELSE IF (BETA.NE.ONE) THEN 00325 DO 220 I = 1,M 00326 C(I,J) = BETA*C(I,J) 00327 220 CONTINUE 00328 END IF 00329 DO 240 L = 1,K 00330 IF (B(J,L).NE.ZERO) THEN 00331 TEMP = ALPHA*B(J,L) 00332 DO 230 I = 1,M 00333 C(I,J) = C(I,J) + TEMP*A(I,L) 00334 230 CONTINUE 00335 END IF 00336 240 CONTINUE 00337 250 CONTINUE 00338 END IF 00339 ELSE IF (CONJA) THEN 00340 IF (CONJB) THEN 00341 * 00342 * Form C := alpha*A**H*B**H + beta*C. 00343 * 00344 DO 280 J = 1,N 00345 DO 270 I = 1,M 00346 TEMP = ZERO 00347 DO 260 L = 1,K 00348 TEMP = TEMP + DCONJG(A(L,I))*DCONJG(B(J,L)) 00349 260 CONTINUE 00350 IF (BETA.EQ.ZERO) THEN 00351 C(I,J) = ALPHA*TEMP 00352 ELSE 00353 C(I,J) = ALPHA*TEMP + BETA*C(I,J) 00354 END IF 00355 270 CONTINUE 00356 280 CONTINUE 00357 ELSE 00358 * 00359 * Form C := alpha*A**H*B**T + beta*C 00360 * 00361 DO 310 J = 1,N 00362 DO 300 I = 1,M 00363 TEMP = ZERO 00364 DO 290 L = 1,K 00365 TEMP = TEMP + DCONJG(A(L,I))*B(J,L) 00366 290 CONTINUE 00367 IF (BETA.EQ.ZERO) THEN 00368 C(I,J) = ALPHA*TEMP 00369 ELSE 00370 C(I,J) = ALPHA*TEMP + BETA*C(I,J) 00371 END IF 00372 300 CONTINUE 00373 310 CONTINUE 00374 END IF 00375 ELSE 00376 IF (CONJB) THEN 00377 * 00378 * Form C := alpha*A**T*B**H + beta*C 00379 * 00380 DO 340 J = 1,N 00381 DO 330 I = 1,M 00382 TEMP = ZERO 00383 DO 320 L = 1,K 00384 TEMP = TEMP + A(L,I)*DCONJG(B(J,L)) 00385 320 CONTINUE 00386 IF (BETA.EQ.ZERO) THEN 00387 C(I,J) = ALPHA*TEMP 00388 ELSE 00389 C(I,J) = ALPHA*TEMP + BETA*C(I,J) 00390 END IF 00391 330 CONTINUE 00392 340 CONTINUE 00393 ELSE 00394 * 00395 * Form C := alpha*A**T*B**T + beta*C 00396 * 00397 DO 370 J = 1,N 00398 DO 360 I = 1,M 00399 TEMP = ZERO 00400 DO 350 L = 1,K 00401 TEMP = TEMP + A(L,I)*B(J,L) 00402 350 CONTINUE 00403 IF (BETA.EQ.ZERO) THEN 00404 C(I,J) = ALPHA*TEMP 00405 ELSE 00406 C(I,J) = ALPHA*TEMP + BETA*C(I,J) 00407 END IF 00408 360 CONTINUE 00409 370 CONTINUE 00410 END IF 00411 END IF 00412 * 00413 RETURN 00414 * 00415 * End of ZGEMM . 00416 * 00417 END
1.7.3 
