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