![[BACK]](/cvsweb/icons/back.gif){kind=icon}
Return to zsyr2k.f CVS log ![[TXT]](/cvsweb/icons/text.gif){kind=icon}
![[DIR]](/cvsweb/icons/dir.gif){kind=icon}
Up to [local] / rpl / lapack / blas|
Return to zsyr2k.f CVS log | Up to [local] / rpl / lapack / blas |
1.1 bertrand 1: SUBROUTINE ZSYR2K(UPLO,TRANS,N,K,ALPHA,A,LDA,B,LDB,BETA,C,LDC)
2: * .. Scalar Arguments ..
3: DOUBLE COMPLEX ALPHA,BETA
4: INTEGER K,LDA,LDB,LDC,N
5: CHARACTER TRANS,UPLO
6: * ..
7: * .. Array Arguments ..
8: DOUBLE COMPLEX A(LDA,*),B(LDB,*),C(LDC,*)
9: * ..
10: *
11: * Purpose
12: * =======
13: *
14: * ZSYR2K performs one of the symmetric rank 2k operations
15: *
1.7 ! bertrand 16: * C := alpha*A*B**T + alpha*B*A**T + beta*C,
1.1 bertrand 17: *
18: * or
19: *
1.7 ! bertrand 20: * C := alpha*A**T*B + alpha*B**T*A + beta*C,
1.1 bertrand 21: *
22: * where alpha and beta are scalars, C is an n by n symmetric matrix
23: * and A and B are n by k matrices in the first case and k by n
24: * matrices in the second case.
25: *
26: * Arguments
27: * ==========
28: *
29: * UPLO - CHARACTER*1.
30: * On entry, UPLO specifies whether the upper or lower
31: * triangular part of the array C is to be referenced as
32: * follows:
33: *
34: * UPLO = 'U' or 'u' Only the upper triangular part of C
35: * is to be referenced.
36: *
37: * UPLO = 'L' or 'l' Only the lower triangular part of C
38: * is to be referenced.
39: *
40: * Unchanged on exit.
41: *
42: * TRANS - CHARACTER*1.
43: * On entry, TRANS specifies the operation to be performed as
44: * follows:
45: *
1.7 ! bertrand 46: * TRANS = 'N' or 'n' C := alpha*A*B**T + alpha*B*A**T +
1.1 bertrand 47: * beta*C.
48: *
1.7 ! bertrand 49: * TRANS = 'T' or 't' C := alpha*A**T*B + alpha*B**T*A +
1.1 bertrand 50: * beta*C.
51: *
52: * Unchanged on exit.
53: *
54: * N - INTEGER.
55: * On entry, N specifies the order of the matrix C. N must be
56: * at least zero.
57: * Unchanged on exit.
58: *
59: * K - INTEGER.
60: * On entry with TRANS = 'N' or 'n', K specifies the number
61: * of columns of the matrices A and B, and on entry with
62: * TRANS = 'T' or 't', K specifies the number of rows of the
63: * matrices A and B. K must be at least zero.
64: * Unchanged on exit.
65: *
66: * ALPHA - COMPLEX*16 .
67: * On entry, ALPHA specifies the scalar alpha.
68: * Unchanged on exit.
69: *
70: * A - COMPLEX*16 array of DIMENSION ( LDA, ka ), where ka is
71: * k when TRANS = 'N' or 'n', and is n otherwise.
72: * Before entry with TRANS = 'N' or 'n', the leading n by k
73: * part of the array A must contain the matrix A, otherwise
74: * the leading k by n part of the array A must contain the
75: * matrix A.
76: * Unchanged on exit.
77: *
78: * LDA - INTEGER.
79: * On entry, LDA specifies the first dimension of A as declared
80: * in the calling (sub) program. When TRANS = 'N' or 'n'
81: * then LDA must be at least max( 1, n ), otherwise LDA must
82: * be at least max( 1, k ).
83: * Unchanged on exit.
84: *
85: * B - COMPLEX*16 array of DIMENSION ( LDB, kb ), where kb is
86: * k when TRANS = 'N' or 'n', and is n otherwise.
87: * Before entry with TRANS = 'N' or 'n', the leading n by k
88: * part of the array B must contain the matrix B, otherwise
89: * the leading k by n part of the array B must contain the
90: * matrix B.
91: * Unchanged on exit.
92: *
93: * LDB - INTEGER.
94: * On entry, LDB specifies the first dimension of B as declared
95: * in the calling (sub) program. When TRANS = 'N' or 'n'
96: * then LDB must be at least max( 1, n ), otherwise LDB must
97: * be at least max( 1, k ).
98: * Unchanged on exit.
99: *
100: * BETA - COMPLEX*16 .
101: * On entry, BETA specifies the scalar beta.
102: * Unchanged on exit.
103: *
104: * C - COMPLEX*16 array of DIMENSION ( LDC, n ).
105: * Before entry with UPLO = 'U' or 'u', the leading n by n
106: * upper triangular part of the array C must contain the upper
107: * triangular part of the symmetric matrix and the strictly
108: * lower triangular part of C is not referenced. On exit, the
109: * upper triangular part of the array C is overwritten by the
110: * upper triangular part of the updated matrix.
111: * Before entry with UPLO = 'L' or 'l', the leading n by n
112: * lower triangular part of the array C must contain the lower
113: * triangular part of the symmetric matrix and the strictly
114: * upper triangular part of C is not referenced. On exit, the
115: * lower triangular part of the array C is overwritten by the
116: * lower triangular part of the updated matrix.
117: *
118: * LDC - INTEGER.
119: * On entry, LDC specifies the first dimension of C as declared
120: * in the calling (sub) program. LDC must be at least
121: * max( 1, n ).
122: * Unchanged on exit.
123: *
124: * Further Details
125: * ===============
126: *
127: * Level 3 Blas routine.
128: *
129: * -- Written on 8-February-1989.
130: * Jack Dongarra, Argonne National Laboratory.
131: * Iain Duff, AERE Harwell.
132: * Jeremy Du Croz, Numerical Algorithms Group Ltd.
133: * Sven Hammarling, Numerical Algorithms Group Ltd.
134: *
135: * =====================================================================
136: *
137: * .. External Functions ..
138: LOGICAL LSAME
139: EXTERNAL LSAME
140: * ..
141: * .. External Subroutines ..
142: EXTERNAL XERBLA
143: * ..
144: * .. Intrinsic Functions ..
145: INTRINSIC MAX
146: * ..
147: * .. Local Scalars ..
148: DOUBLE COMPLEX TEMP1,TEMP2
149: INTEGER I,INFO,J,L,NROWA
150: LOGICAL UPPER
151: * ..
152: * .. Parameters ..
153: DOUBLE COMPLEX ONE
154: PARAMETER (ONE= (1.0D+0,0.0D+0))
155: DOUBLE COMPLEX ZERO
156: PARAMETER (ZERO= (0.0D+0,0.0D+0))
157: * ..
158: *
159: * Test the input parameters.
160: *
161: IF (LSAME(TRANS,'N')) THEN
162: NROWA = N
163: ELSE
164: NROWA = K
165: END IF
166: UPPER = LSAME(UPLO,'U')
167: *
168: INFO = 0
169: IF ((.NOT.UPPER) .AND. (.NOT.LSAME(UPLO,'L'))) THEN
170: INFO = 1
171: ELSE IF ((.NOT.LSAME(TRANS,'N')) .AND.
172: + (.NOT.LSAME(TRANS,'T'))) THEN
173: INFO = 2
174: ELSE IF (N.LT.0) THEN
175: INFO = 3
176: ELSE IF (K.LT.0) THEN
177: INFO = 4
178: ELSE IF (LDA.LT.MAX(1,NROWA)) THEN
179: INFO = 7
180: ELSE IF (LDB.LT.MAX(1,NROWA)) THEN
181: INFO = 9
182: ELSE IF (LDC.LT.MAX(1,N)) THEN
183: INFO = 12
184: END IF
185: IF (INFO.NE.0) THEN
186: CALL XERBLA('ZSYR2K',INFO)
187: RETURN
188: END IF
189: *
190: * Quick return if possible.
191: *
192: IF ((N.EQ.0) .OR. (((ALPHA.EQ.ZERO).OR.
193: + (K.EQ.0)).AND. (BETA.EQ.ONE))) RETURN
194: *
195: * And when alpha.eq.zero.
196: *
197: IF (ALPHA.EQ.ZERO) THEN
198: IF (UPPER) THEN
199: IF (BETA.EQ.ZERO) THEN
200: DO 20 J = 1,N
201: DO 10 I = 1,J
202: C(I,J) = ZERO
203: 10 CONTINUE
204: 20 CONTINUE
205: ELSE
206: DO 40 J = 1,N
207: DO 30 I = 1,J
208: C(I,J) = BETA*C(I,J)
209: 30 CONTINUE
210: 40 CONTINUE
211: END IF
212: ELSE
213: IF (BETA.EQ.ZERO) THEN
214: DO 60 J = 1,N
215: DO 50 I = J,N
216: C(I,J) = ZERO
217: 50 CONTINUE
218: 60 CONTINUE
219: ELSE
220: DO 80 J = 1,N
221: DO 70 I = J,N
222: C(I,J) = BETA*C(I,J)
223: 70 CONTINUE
224: 80 CONTINUE
225: END IF
226: END IF
227: RETURN
228: END IF
229: *
230: * Start the operations.
231: *
232: IF (LSAME(TRANS,'N')) THEN
233: *
1.7 ! bertrand 234: * Form C := alpha*A*B**T + alpha*B*A**T + C.
1.1 bertrand 235: *
236: IF (UPPER) THEN
237: DO 130 J = 1,N
238: IF (BETA.EQ.ZERO) THEN
239: DO 90 I = 1,J
240: C(I,J) = ZERO
241: 90 CONTINUE
242: ELSE IF (BETA.NE.ONE) THEN
243: DO 100 I = 1,J
244: C(I,J) = BETA*C(I,J)
245: 100 CONTINUE
246: END IF
247: DO 120 L = 1,K
248: IF ((A(J,L).NE.ZERO) .OR. (B(J,L).NE.ZERO)) THEN
249: TEMP1 = ALPHA*B(J,L)
250: TEMP2 = ALPHA*A(J,L)
251: DO 110 I = 1,J
252: C(I,J) = C(I,J) + A(I,L)*TEMP1 +
253: + B(I,L)*TEMP2
254: 110 CONTINUE
255: END IF
256: 120 CONTINUE
257: 130 CONTINUE
258: ELSE
259: DO 180 J = 1,N
260: IF (BETA.EQ.ZERO) THEN
261: DO 140 I = J,N
262: C(I,J) = ZERO
263: 140 CONTINUE
264: ELSE IF (BETA.NE.ONE) THEN
265: DO 150 I = J,N
266: C(I,J) = BETA*C(I,J)
267: 150 CONTINUE
268: END IF
269: DO 170 L = 1,K
270: IF ((A(J,L).NE.ZERO) .OR. (B(J,L).NE.ZERO)) THEN
271: TEMP1 = ALPHA*B(J,L)
272: TEMP2 = ALPHA*A(J,L)
273: DO 160 I = J,N
274: C(I,J) = C(I,J) + A(I,L)*TEMP1 +
275: + B(I,L)*TEMP2
276: 160 CONTINUE
277: END IF
278: 170 CONTINUE
279: 180 CONTINUE
280: END IF
281: ELSE
282: *
1.7 ! bertrand 283: * Form C := alpha*A**T*B + alpha*B**T*A + C.
1.1 bertrand 284: *
285: IF (UPPER) THEN
286: DO 210 J = 1,N
287: DO 200 I = 1,J
288: TEMP1 = ZERO
289: TEMP2 = ZERO
290: DO 190 L = 1,K
291: TEMP1 = TEMP1 + A(L,I)*B(L,J)
292: TEMP2 = TEMP2 + B(L,I)*A(L,J)
293: 190 CONTINUE
294: IF (BETA.EQ.ZERO) THEN
295: C(I,J) = ALPHA*TEMP1 + ALPHA*TEMP2
296: ELSE
297: C(I,J) = BETA*C(I,J) + ALPHA*TEMP1 +
298: + ALPHA*TEMP2
299: END IF
300: 200 CONTINUE
301: 210 CONTINUE
302: ELSE
303: DO 240 J = 1,N
304: DO 230 I = J,N
305: TEMP1 = ZERO
306: TEMP2 = ZERO
307: DO 220 L = 1,K
308: TEMP1 = TEMP1 + A(L,I)*B(L,J)
309: TEMP2 = TEMP2 + B(L,I)*A(L,J)
310: 220 CONTINUE
311: IF (BETA.EQ.ZERO) THEN
312: C(I,J) = ALPHA*TEMP1 + ALPHA*TEMP2
313: ELSE
314: C(I,J) = BETA*C(I,J) + ALPHA*TEMP1 +
315: + ALPHA*TEMP2
316: END IF
317: 230 CONTINUE
318: 240 CONTINUE
319: END IF
320: END IF
321: *
322: RETURN
323: *
324: * End of ZSYR2K.
325: *
326: END