1
SUBROUTINE ZHER2K(UPLO,TRANS,N,K,ALPHA,A,LDA,B,LDB,BETA,C,LDC)
2
* .. Scalar Arguments ..
5
INTEGER K,LDA,LDB,LDC,N
8
* .. Array Arguments ..
9
DOUBLE COMPLEX A(LDA,*),B(LDB,*),C(LDC,*)
15
* ZHER2K performs one of the hermitian rank 2k operations
17
* C := alpha*A*conjg( B' ) + conjg( alpha )*B*conjg( A' ) + beta*C,
21
* C := alpha*conjg( A' )*B + conjg( alpha )*conjg( B' )*A + beta*C,
23
* where alpha and beta are scalars with beta real, C is an n by n
24
* hermitian matrix and A and B are n by k matrices in the first case
25
* and k by n matrices in the second case.
31
* On entry, UPLO specifies whether the upper or lower
32
* triangular part of the array C is to be referenced as
35
* UPLO = 'U' or 'u' Only the upper triangular part of C
36
* is to be referenced.
38
* UPLO = 'L' or 'l' Only the lower triangular part of C
39
* is to be referenced.
43
* TRANS - CHARACTER*1.
44
* On entry, TRANS specifies the operation to be performed as
47
* TRANS = 'N' or 'n' C := alpha*A*conjg( B' ) +
48
* conjg( alpha )*B*conjg( A' ) +
51
* TRANS = 'C' or 'c' C := alpha*conjg( A' )*B +
52
* conjg( alpha )*conjg( B' )*A +
58
* On entry, N specifies the order of the matrix C. N must be
63
* On entry with TRANS = 'N' or 'n', K specifies the number
64
* of columns of the matrices A and B, and on entry with
65
* TRANS = 'C' or 'c', K specifies the number of rows of the
66
* matrices A and B. K must be at least zero.
69
* ALPHA - COMPLEX*16 .
70
* On entry, ALPHA specifies the scalar alpha.
73
* A - COMPLEX*16 array of DIMENSION ( LDA, ka ), where ka is
74
* k when TRANS = 'N' or 'n', and is n otherwise.
75
* Before entry with TRANS = 'N' or 'n', the leading n by k
76
* part of the array A must contain the matrix A, otherwise
77
* the leading k by n part of the array A must contain the
82
* On entry, LDA specifies the first dimension of A as declared
83
* in the calling (sub) program. When TRANS = 'N' or 'n'
84
* then LDA must be at least max( 1, n ), otherwise LDA must
85
* be at least max( 1, k ).
88
* B - COMPLEX*16 array of DIMENSION ( LDB, kb ), where kb is
89
* k when TRANS = 'N' or 'n', and is n otherwise.
90
* Before entry with TRANS = 'N' or 'n', the leading n by k
91
* part of the array B must contain the matrix B, otherwise
92
* the leading k by n part of the array B must contain the
97
* On entry, LDB specifies the first dimension of B as declared
98
* in the calling (sub) program. When TRANS = 'N' or 'n'
99
* then LDB must be at least max( 1, n ), otherwise LDB must
100
* be at least max( 1, k ).
103
* BETA - DOUBLE PRECISION .
104
* On entry, BETA specifies the scalar beta.
107
* C - COMPLEX*16 array of DIMENSION ( LDC, n ).
108
* Before entry with UPLO = 'U' or 'u', the leading n by n
109
* upper triangular part of the array C must contain the upper
110
* triangular part of the hermitian matrix and the strictly
111
* lower triangular part of C is not referenced. On exit, the
112
* upper triangular part of the array C is overwritten by the
113
* upper triangular part of the updated matrix.
114
* Before entry with UPLO = 'L' or 'l', the leading n by n
115
* lower triangular part of the array C must contain the lower
116
* triangular part of the hermitian matrix and the strictly
117
* upper triangular part of C is not referenced. On exit, the
118
* lower triangular part of the array C is overwritten by the
119
* lower triangular part of the updated matrix.
120
* Note that the imaginary parts of the diagonal elements need
121
* not be set, they are assumed to be zero, and on exit they
125
* On entry, LDC specifies the first dimension of C as declared
126
* in the calling (sub) program. LDC must be at least
131
* Level 3 Blas routine.
133
* -- Written on 8-February-1989.
134
* Jack Dongarra, Argonne National Laboratory.
135
* Iain Duff, AERE Harwell.
136
* Jeremy Du Croz, Numerical Algorithms Group Ltd.
137
* Sven Hammarling, Numerical Algorithms Group Ltd.
139
* -- Modified 8-Nov-93 to set C(J,J) to DBLE( C(J,J) ) when BETA = 1.
140
* Ed Anderson, Cray Research Inc.
143
* .. External Functions ..
147
* .. External Subroutines ..
150
* .. Intrinsic Functions ..
151
INTRINSIC DBLE,DCONJG,MAX
153
* .. Local Scalars ..
154
DOUBLE COMPLEX TEMP1,TEMP2
155
INTEGER I,INFO,J,L,NROWA
160
PARAMETER (ONE=1.0D+0)
162
PARAMETER (ZERO= (0.0D+0,0.0D+0))
165
* Test the input parameters.
167
IF (LSAME(TRANS,'N')) THEN
172
UPPER = LSAME(UPLO,'U')
175
IF ((.NOT.UPPER) .AND. (.NOT.LSAME(UPLO,'L'))) THEN
177
ELSE IF ((.NOT.LSAME(TRANS,'N')) .AND.
178
+ (.NOT.LSAME(TRANS,'C'))) THEN
180
ELSE IF (N.LT.0) THEN
182
ELSE IF (K.LT.0) THEN
184
ELSE IF (LDA.LT.MAX(1,NROWA)) THEN
186
ELSE IF (LDB.LT.MAX(1,NROWA)) THEN
188
ELSE IF (LDC.LT.MAX(1,N)) THEN
192
CALL XERBLA('ZHER2K',INFO)
196
* Quick return if possible.
198
IF ((N.EQ.0) .OR. (((ALPHA.EQ.ZERO).OR.
199
+ (K.EQ.0)).AND. (BETA.EQ.ONE))) RETURN
201
* And when alpha.eq.zero.
203
IF (ALPHA.EQ.ZERO) THEN
205
IF (BETA.EQ.DBLE(ZERO)) THEN
216
C(J,J) = BETA*DBLE(C(J,J))
220
IF (BETA.EQ.DBLE(ZERO)) THEN
228
C(J,J) = BETA*DBLE(C(J,J))
238
* Start the operations.
240
IF (LSAME(TRANS,'N')) THEN
242
* Form C := alpha*A*conjg( B' ) + conjg( alpha )*B*conjg( A' ) +
247
IF (BETA.EQ.DBLE(ZERO)) THEN
251
ELSE IF (BETA.NE.ONE) THEN
255
C(J,J) = BETA*DBLE(C(J,J))
257
C(J,J) = DBLE(C(J,J))
260
IF ((A(J,L).NE.ZERO) .OR. (B(J,L).NE.ZERO)) THEN
261
TEMP1 = ALPHA*DCONJG(B(J,L))
262
TEMP2 = DCONJG(ALPHA*A(J,L))
264
C(I,J) = C(I,J) + A(I,L)*TEMP1 +
267
C(J,J) = DBLE(C(J,J)) +
268
+ DBLE(A(J,L)*TEMP1+B(J,L)*TEMP2)
274
IF (BETA.EQ.DBLE(ZERO)) THEN
278
ELSE IF (BETA.NE.ONE) THEN
282
C(J,J) = BETA*DBLE(C(J,J))
284
C(J,J) = DBLE(C(J,J))
287
IF ((A(J,L).NE.ZERO) .OR. (B(J,L).NE.ZERO)) THEN
288
TEMP1 = ALPHA*DCONJG(B(J,L))
289
TEMP2 = DCONJG(ALPHA*A(J,L))
291
C(I,J) = C(I,J) + A(I,L)*TEMP1 +
294
C(J,J) = DBLE(C(J,J)) +
295
+ DBLE(A(J,L)*TEMP1+B(J,L)*TEMP2)
302
* Form C := alpha*conjg( A' )*B + conjg( alpha )*conjg( B' )*A +
311
TEMP1 = TEMP1 + DCONJG(A(L,I))*B(L,J)
312
TEMP2 = TEMP2 + DCONJG(B(L,I))*A(L,J)
315
IF (BETA.EQ.DBLE(ZERO)) THEN
316
C(J,J) = DBLE(ALPHA*TEMP1+
317
+ DCONJG(ALPHA)*TEMP2)
319
C(J,J) = BETA*DBLE(C(J,J)) +
321
+ DCONJG(ALPHA)*TEMP2)
324
IF (BETA.EQ.DBLE(ZERO)) THEN
325
C(I,J) = ALPHA*TEMP1 + DCONJG(ALPHA)*TEMP2
327
C(I,J) = BETA*C(I,J) + ALPHA*TEMP1 +
328
+ DCONJG(ALPHA)*TEMP2
339
TEMP1 = TEMP1 + DCONJG(A(L,I))*B(L,J)
340
TEMP2 = TEMP2 + DCONJG(B(L,I))*A(L,J)
343
IF (BETA.EQ.DBLE(ZERO)) THEN
344
C(J,J) = DBLE(ALPHA*TEMP1+
345
+ DCONJG(ALPHA)*TEMP2)
347
C(J,J) = BETA*DBLE(C(J,J)) +
349
+ DCONJG(ALPHA)*TEMP2)
352
IF (BETA.EQ.DBLE(ZERO)) THEN
353
C(I,J) = ALPHA*TEMP1 + DCONJG(ALPHA)*TEMP2
355
C(I,J) = BETA*C(I,J) + ALPHA*TEMP1 +
356
+ DCONJG(ALPHA)*TEMP2
added
removed
src/blas/f77reference/zher2k.f
