![[BACK]](/cvsweb/icons/back.gif){kind=icon}
Return to dsyr2k.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 dsyr2k.f CVS log | Up to [local] / rpl / lapack / blas |
1.8 bertrand 1: *> \brief \b DSYR2K
2: *
3: * =========== DOCUMENTATION ===========
4: *
1.13 bertrand 5: * Online html documentation available at
6: * http://www.netlib.org/lapack/explore-html/
1.8 bertrand 7: *
8: * Definition:
9: * ===========
10: *
11: * SUBROUTINE DSYR2K(UPLO,TRANS,N,K,ALPHA,A,LDA,B,LDB,BETA,C,LDC)
1.13 bertrand 12: *
1.8 bertrand 13: * .. Scalar Arguments ..
14: * DOUBLE PRECISION ALPHA,BETA
15: * INTEGER K,LDA,LDB,LDC,N
16: * CHARACTER TRANS,UPLO
17: * ..
18: * .. Array Arguments ..
19: * DOUBLE PRECISION A(LDA,*),B(LDB,*),C(LDC,*)
20: * ..
1.13 bertrand 21: *
1.8 bertrand 22: *
23: *> \par Purpose:
24: * =============
25: *>
26: *> \verbatim
27: *>
28: *> DSYR2K performs one of the symmetric rank 2k operations
29: *>
30: *> C := alpha*A*B**T + alpha*B*A**T + beta*C,
31: *>
32: *> or
33: *>
34: *> C := alpha*A**T*B + alpha*B**T*A + beta*C,
35: *>
36: *> where alpha and beta are scalars, C is an n by n symmetric matrix
37: *> and A and B are n by k matrices in the first case and k by n
38: *> matrices in the second case.
39: *> \endverbatim
40: *
41: * Arguments:
42: * ==========
43: *
44: *> \param[in] UPLO
45: *> \verbatim
46: *> UPLO is CHARACTER*1
47: *> On entry, UPLO specifies whether the upper or lower
48: *> triangular part of the array C is to be referenced as
49: *> follows:
50: *>
51: *> UPLO = 'U' or 'u' Only the upper triangular part of C
52: *> is to be referenced.
53: *>
54: *> UPLO = 'L' or 'l' Only the lower triangular part of C
55: *> is to be referenced.
56: *> \endverbatim
57: *>
58: *> \param[in] TRANS
59: *> \verbatim
60: *> TRANS is CHARACTER*1
61: *> On entry, TRANS specifies the operation to be performed as
62: *> follows:
63: *>
64: *> TRANS = 'N' or 'n' C := alpha*A*B**T + alpha*B*A**T +
65: *> beta*C.
66: *>
67: *> TRANS = 'T' or 't' C := alpha*A**T*B + alpha*B**T*A +
68: *> beta*C.
69: *>
70: *> TRANS = 'C' or 'c' C := alpha*A**T*B + alpha*B**T*A +
71: *> beta*C.
72: *> \endverbatim
73: *>
74: *> \param[in] N
75: *> \verbatim
76: *> N is INTEGER
77: *> On entry, N specifies the order of the matrix C. N must be
78: *> at least zero.
79: *> \endverbatim
80: *>
81: *> \param[in] K
82: *> \verbatim
83: *> K is INTEGER
84: *> On entry with TRANS = 'N' or 'n', K specifies the number
85: *> of columns of the matrices A and B, and on entry with
86: *> TRANS = 'T' or 't' or 'C' or 'c', K specifies the number
87: *> of rows of the matrices A and B. K must be at least zero.
88: *> \endverbatim
89: *>
90: *> \param[in] ALPHA
91: *> \verbatim
92: *> ALPHA is DOUBLE PRECISION.
93: *> On entry, ALPHA specifies the scalar alpha.
94: *> \endverbatim
95: *>
96: *> \param[in] A
97: *> \verbatim
1.14 bertrand 98: *> A is DOUBLE PRECISION array, dimension ( LDA, ka ), where ka is
1.8 bertrand 99: *> k when TRANS = 'N' or 'n', and is n otherwise.
100: *> Before entry with TRANS = 'N' or 'n', the leading n by k
101: *> part of the array A must contain the matrix A, otherwise
102: *> the leading k by n part of the array A must contain the
103: *> matrix A.
104: *> \endverbatim
105: *>
106: *> \param[in] LDA
107: *> \verbatim
108: *> LDA is INTEGER
109: *> On entry, LDA specifies the first dimension of A as declared
110: *> in the calling (sub) program. When TRANS = 'N' or 'n'
111: *> then LDA must be at least max( 1, n ), otherwise LDA must
112: *> be at least max( 1, k ).
113: *> \endverbatim
114: *>
115: *> \param[in] B
116: *> \verbatim
1.14 bertrand 117: *> B is DOUBLE PRECISION array, dimension ( LDB, kb ), where kb is
1.8 bertrand 118: *> k when TRANS = 'N' or 'n', and is n otherwise.
119: *> Before entry with TRANS = 'N' or 'n', the leading n by k
120: *> part of the array B must contain the matrix B, otherwise
121: *> the leading k by n part of the array B must contain the
122: *> matrix B.
123: *> \endverbatim
124: *>
125: *> \param[in] LDB
126: *> \verbatim
127: *> LDB is INTEGER
128: *> On entry, LDB specifies the first dimension of B as declared
129: *> in the calling (sub) program. When TRANS = 'N' or 'n'
130: *> then LDB must be at least max( 1, n ), otherwise LDB must
131: *> be at least max( 1, k ).
132: *> \endverbatim
133: *>
134: *> \param[in] BETA
135: *> \verbatim
136: *> BETA is DOUBLE PRECISION.
137: *> On entry, BETA specifies the scalar beta.
138: *> \endverbatim
139: *>
140: *> \param[in,out] C
141: *> \verbatim
1.14 bertrand 142: *> C is DOUBLE PRECISION array, dimension ( LDC, N )
1.8 bertrand 143: *> Before entry with UPLO = 'U' or 'u', the leading n by n
144: *> upper triangular part of the array C must contain the upper
145: *> triangular part of the symmetric matrix and the strictly
146: *> lower triangular part of C is not referenced. On exit, the
147: *> upper triangular part of the array C is overwritten by the
148: *> upper triangular part of the updated matrix.
149: *> Before entry with UPLO = 'L' or 'l', the leading n by n
150: *> lower triangular part of the array C must contain the lower
151: *> triangular part of the symmetric matrix and the strictly
152: *> upper triangular part of C is not referenced. On exit, the
153: *> lower triangular part of the array C is overwritten by the
154: *> lower triangular part of the updated matrix.
155: *> \endverbatim
156: *>
157: *> \param[in] LDC
158: *> \verbatim
159: *> LDC is INTEGER
160: *> On entry, LDC specifies the first dimension of C as declared
161: *> in the calling (sub) program. LDC must be at least
162: *> max( 1, n ).
163: *> \endverbatim
164: *
165: * Authors:
166: * ========
167: *
1.13 bertrand 168: *> \author Univ. of Tennessee
169: *> \author Univ. of California Berkeley
170: *> \author Univ. of Colorado Denver
171: *> \author NAG Ltd.
1.8 bertrand 172: *
1.13 bertrand 173: *> \date December 2016
1.8 bertrand 174: *
175: *> \ingroup double_blas_level3
176: *
177: *> \par Further Details:
178: * =====================
179: *>
180: *> \verbatim
181: *>
182: *> Level 3 Blas routine.
183: *>
184: *>
185: *> -- Written on 8-February-1989.
186: *> Jack Dongarra, Argonne National Laboratory.
187: *> Iain Duff, AERE Harwell.
188: *> Jeremy Du Croz, Numerical Algorithms Group Ltd.
189: *> Sven Hammarling, Numerical Algorithms Group Ltd.
190: *> \endverbatim
191: *>
192: * =====================================================================
1.1 bertrand 193: SUBROUTINE DSYR2K(UPLO,TRANS,N,K,ALPHA,A,LDA,B,LDB,BETA,C,LDC)
1.8 bertrand 194: *
1.13 bertrand 195: * -- Reference BLAS level3 routine (version 3.7.0) --
1.8 bertrand 196: * -- Reference BLAS is a software package provided by Univ. of Tennessee, --
197: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
1.13 bertrand 198: * December 2016
1.8 bertrand 199: *
1.1 bertrand 200: * .. Scalar Arguments ..
201: DOUBLE PRECISION ALPHA,BETA
202: INTEGER K,LDA,LDB,LDC,N
203: CHARACTER TRANS,UPLO
204: * ..
205: * .. Array Arguments ..
206: DOUBLE PRECISION A(LDA,*),B(LDB,*),C(LDC,*)
207: * ..
208: *
209: * =====================================================================
210: *
211: * .. External Functions ..
212: LOGICAL LSAME
213: EXTERNAL LSAME
214: * ..
215: * .. External Subroutines ..
216: EXTERNAL XERBLA
217: * ..
218: * .. Intrinsic Functions ..
219: INTRINSIC MAX
220: * ..
221: * .. Local Scalars ..
222: DOUBLE PRECISION TEMP1,TEMP2
223: INTEGER I,INFO,J,L,NROWA
224: LOGICAL UPPER
225: * ..
226: * .. Parameters ..
227: DOUBLE PRECISION ONE,ZERO
228: PARAMETER (ONE=1.0D+0,ZERO=0.0D+0)
229: * ..
230: *
231: * Test the input parameters.
232: *
233: IF (LSAME(TRANS,'N')) THEN
234: NROWA = N
235: ELSE
236: NROWA = K
237: END IF
238: UPPER = LSAME(UPLO,'U')
239: *
240: INFO = 0
241: IF ((.NOT.UPPER) .AND. (.NOT.LSAME(UPLO,'L'))) THEN
242: INFO = 1
243: ELSE IF ((.NOT.LSAME(TRANS,'N')) .AND.
244: + (.NOT.LSAME(TRANS,'T')) .AND.
245: + (.NOT.LSAME(TRANS,'C'))) THEN
246: INFO = 2
247: ELSE IF (N.LT.0) THEN
248: INFO = 3
249: ELSE IF (K.LT.0) THEN
250: INFO = 4
251: ELSE IF (LDA.LT.MAX(1,NROWA)) THEN
252: INFO = 7
253: ELSE IF (LDB.LT.MAX(1,NROWA)) THEN
254: INFO = 9
255: ELSE IF (LDC.LT.MAX(1,N)) THEN
256: INFO = 12
257: END IF
258: IF (INFO.NE.0) THEN
259: CALL XERBLA('DSYR2K',INFO)
260: RETURN
261: END IF
262: *
263: * Quick return if possible.
264: *
265: IF ((N.EQ.0) .OR. (((ALPHA.EQ.ZERO).OR.
266: + (K.EQ.0)).AND. (BETA.EQ.ONE))) RETURN
267: *
268: * And when alpha.eq.zero.
269: *
270: IF (ALPHA.EQ.ZERO) THEN
271: IF (UPPER) THEN
272: IF (BETA.EQ.ZERO) THEN
273: DO 20 J = 1,N
274: DO 10 I = 1,J
275: C(I,J) = ZERO
276: 10 CONTINUE
277: 20 CONTINUE
278: ELSE
279: DO 40 J = 1,N
280: DO 30 I = 1,J
281: C(I,J) = BETA*C(I,J)
282: 30 CONTINUE
283: 40 CONTINUE
284: END IF
285: ELSE
286: IF (BETA.EQ.ZERO) THEN
287: DO 60 J = 1,N
288: DO 50 I = J,N
289: C(I,J) = ZERO
290: 50 CONTINUE
291: 60 CONTINUE
292: ELSE
293: DO 80 J = 1,N
294: DO 70 I = J,N
295: C(I,J) = BETA*C(I,J)
296: 70 CONTINUE
297: 80 CONTINUE
298: END IF
299: END IF
300: RETURN
301: END IF
302: *
303: * Start the operations.
304: *
305: IF (LSAME(TRANS,'N')) THEN
306: *
1.7 bertrand 307: * Form C := alpha*A*B**T + alpha*B*A**T + C.
1.1 bertrand 308: *
309: IF (UPPER) THEN
310: DO 130 J = 1,N
311: IF (BETA.EQ.ZERO) THEN
312: DO 90 I = 1,J
313: C(I,J) = ZERO
314: 90 CONTINUE
315: ELSE IF (BETA.NE.ONE) THEN
316: DO 100 I = 1,J
317: C(I,J) = BETA*C(I,J)
318: 100 CONTINUE
319: END IF
320: DO 120 L = 1,K
321: IF ((A(J,L).NE.ZERO) .OR. (B(J,L).NE.ZERO)) THEN
322: TEMP1 = ALPHA*B(J,L)
323: TEMP2 = ALPHA*A(J,L)
324: DO 110 I = 1,J
325: C(I,J) = C(I,J) + A(I,L)*TEMP1 +
326: + B(I,L)*TEMP2
327: 110 CONTINUE
328: END IF
329: 120 CONTINUE
330: 130 CONTINUE
331: ELSE
332: DO 180 J = 1,N
333: IF (BETA.EQ.ZERO) THEN
334: DO 140 I = J,N
335: C(I,J) = ZERO
336: 140 CONTINUE
337: ELSE IF (BETA.NE.ONE) THEN
338: DO 150 I = J,N
339: C(I,J) = BETA*C(I,J)
340: 150 CONTINUE
341: END IF
342: DO 170 L = 1,K
343: IF ((A(J,L).NE.ZERO) .OR. (B(J,L).NE.ZERO)) THEN
344: TEMP1 = ALPHA*B(J,L)
345: TEMP2 = ALPHA*A(J,L)
346: DO 160 I = J,N
347: C(I,J) = C(I,J) + A(I,L)*TEMP1 +
348: + B(I,L)*TEMP2
349: 160 CONTINUE
350: END IF
351: 170 CONTINUE
352: 180 CONTINUE
353: END IF
354: ELSE
355: *
1.7 bertrand 356: * Form C := alpha*A**T*B + alpha*B**T*A + C.
1.1 bertrand 357: *
358: IF (UPPER) THEN
359: DO 210 J = 1,N
360: DO 200 I = 1,J
361: TEMP1 = ZERO
362: TEMP2 = ZERO
363: DO 190 L = 1,K
364: TEMP1 = TEMP1 + A(L,I)*B(L,J)
365: TEMP2 = TEMP2 + B(L,I)*A(L,J)
366: 190 CONTINUE
367: IF (BETA.EQ.ZERO) THEN
368: C(I,J) = ALPHA*TEMP1 + ALPHA*TEMP2
369: ELSE
370: C(I,J) = BETA*C(I,J) + ALPHA*TEMP1 +
371: + ALPHA*TEMP2
372: END IF
373: 200 CONTINUE
374: 210 CONTINUE
375: ELSE
376: DO 240 J = 1,N
377: DO 230 I = J,N
378: TEMP1 = ZERO
379: TEMP2 = ZERO
380: DO 220 L = 1,K
381: TEMP1 = TEMP1 + A(L,I)*B(L,J)
382: TEMP2 = TEMP2 + B(L,I)*A(L,J)
383: 220 CONTINUE
384: IF (BETA.EQ.ZERO) THEN
385: C(I,J) = ALPHA*TEMP1 + ALPHA*TEMP2
386: ELSE
387: C(I,J) = BETA*C(I,J) + ALPHA*TEMP1 +
388: + ALPHA*TEMP2
389: END IF
390: 230 CONTINUE
391: 240 CONTINUE
392: END IF
393: END IF
394: *
395: RETURN
396: *
397: * End of DSYR2K.
398: *
399: END