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