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