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