1: *> \brief \b ZHFRK
2: *
3: * =========== DOCUMENTATION ===========
4: *
5: * Online html documentation available at
6: * http://www.netlib.org/lapack/explore-html/
7: *
8: *> \htmlonly
9: *> Download ZHFRK + dependencies
10: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zhfrk.f">
11: *> [TGZ]</a>
12: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zhfrk.f">
13: *> [ZIP]</a>
14: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zhfrk.f">
15: *> [TXT]</a>
16: *> \endhtmlonly
17: *
18: * Definition:
19: * ===========
20: *
21: * SUBROUTINE ZHFRK( TRANSR, UPLO, TRANS, N, K, ALPHA, A, LDA, BETA,
22: * C )
23: *
24: * .. Scalar Arguments ..
25: * DOUBLE PRECISION ALPHA, BETA
26: * INTEGER K, LDA, N
27: * CHARACTER TRANS, TRANSR, UPLO
28: * ..
29: * .. Array Arguments ..
30: * COMPLEX*16 A( LDA, * ), C( * )
31: * ..
32: *
33: *
34: *> \par Purpose:
35: * =============
36: *>
37: *> \verbatim
38: *>
39: *> Level 3 BLAS like routine for C in RFP Format.
40: *>
41: *> ZHFRK performs one of the Hermitian rank--k operations
42: *>
43: *> C := alpha*A*A**H + beta*C,
44: *>
45: *> or
46: *>
47: *> C := alpha*A**H*A + beta*C,
48: *>
49: *> where alpha and beta are real scalars, C is an n--by--n Hermitian
50: *> matrix and A is an n--by--k matrix in the first case and a k--by--n
51: *> matrix in the second case.
52: *> \endverbatim
53: *
54: * Arguments:
55: * ==========
56: *
57: *> \param[in] TRANSR
58: *> \verbatim
59: *> TRANSR is CHARACTER*1
60: *> = 'N': The Normal Form of RFP A is stored;
61: *> = 'C': The Conjugate-transpose Form of RFP A is stored.
62: *> \endverbatim
63: *>
64: *> \param[in] UPLO
65: *> \verbatim
66: *> UPLO is CHARACTER*1
67: *> On entry, UPLO specifies whether the upper or lower
68: *> triangular part of the array C is to be referenced as
69: *> follows:
70: *>
71: *> UPLO = 'U' or 'u' Only the upper triangular part of C
72: *> is to be referenced.
73: *>
74: *> UPLO = 'L' or 'l' Only the lower triangular part of C
75: *> is to be referenced.
76: *>
77: *> Unchanged on exit.
78: *> \endverbatim
79: *>
80: *> \param[in] TRANS
81: *> \verbatim
82: *> TRANS is CHARACTER*1
83: *> On entry, TRANS specifies the operation to be performed as
84: *> follows:
85: *>
86: *> TRANS = 'N' or 'n' C := alpha*A*A**H + beta*C.
87: *>
88: *> TRANS = 'C' or 'c' C := alpha*A**H*A + beta*C.
89: *>
90: *> Unchanged on exit.
91: *> \endverbatim
92: *>
93: *> \param[in] N
94: *> \verbatim
95: *> N is INTEGER
96: *> On entry, N specifies the order of the matrix C. N must be
97: *> at least zero.
98: *> Unchanged on exit.
99: *> \endverbatim
100: *>
101: *> \param[in] K
102: *> \verbatim
103: *> K is INTEGER
104: *> On entry with TRANS = 'N' or 'n', K specifies the number
105: *> of columns of the matrix A, and on entry with
106: *> TRANS = 'C' or 'c', K specifies the number of rows of the
107: *> matrix A. K must be at least zero.
108: *> Unchanged on exit.
109: *> \endverbatim
110: *>
111: *> \param[in] ALPHA
112: *> \verbatim
113: *> ALPHA is DOUBLE PRECISION
114: *> On entry, ALPHA specifies the scalar alpha.
115: *> Unchanged on exit.
116: *> \endverbatim
117: *>
118: *> \param[in] A
119: *> \verbatim
120: *> A is COMPLEX*16 array of DIMENSION (LDA,ka)
121: *> where KA
122: *> is K when TRANS = 'N' or 'n', and is N otherwise. Before
123: *> entry with TRANS = 'N' or 'n', the leading N--by--K part of
124: *> the array A must contain the matrix A, otherwise the leading
125: *> K--by--N part of the array A must contain the matrix A.
126: *> Unchanged on exit.
127: *> \endverbatim
128: *>
129: *> \param[in] LDA
130: *> \verbatim
131: *> LDA is INTEGER
132: *> On entry, LDA specifies the first dimension of A as declared
133: *> in the calling (sub) program. When TRANS = 'N' or 'n'
134: *> then LDA must be at least max( 1, n ), otherwise LDA must
135: *> be at least max( 1, k ).
136: *> Unchanged on exit.
137: *> \endverbatim
138: *>
139: *> \param[in] BETA
140: *> \verbatim
141: *> BETA is DOUBLE PRECISION
142: *> On entry, BETA specifies the scalar beta.
143: *> Unchanged on exit.
144: *> \endverbatim
145: *>
146: *> \param[in,out] C
147: *> \verbatim
148: *> C is COMPLEX*16 array, dimension (N*(N+1)/2)
149: *> On entry, the matrix A in RFP Format. RFP Format is
150: *> described by TRANSR, UPLO and N. Note that the imaginary
151: *> parts of the diagonal elements need not be set, they are
152: *> assumed to be zero, and on exit they are set to zero.
153: *> \endverbatim
154: *
155: * Authors:
156: * ========
157: *
158: *> \author Univ. of Tennessee
159: *> \author Univ. of California Berkeley
160: *> \author Univ. of Colorado Denver
161: *> \author NAG Ltd.
162: *
163: *> \date November 2011
164: *
165: *> \ingroup complex16OTHERcomputational
166: *
167: * =====================================================================
168: SUBROUTINE ZHFRK( TRANSR, UPLO, TRANS, N, K, ALPHA, A, LDA, BETA,
169: $ C )
170: *
171: * -- LAPACK computational routine (version 3.4.0) --
172: * -- LAPACK is a software package provided by Univ. of Tennessee, --
173: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
174: * November 2011
175: *
176: * .. Scalar Arguments ..
177: DOUBLE PRECISION ALPHA, BETA
178: INTEGER K, LDA, N
179: CHARACTER TRANS, TRANSR, UPLO
180: * ..
181: * .. Array Arguments ..
182: COMPLEX*16 A( LDA, * ), C( * )
183: * ..
184: *
185: * =====================================================================
186: *
187: * .. Parameters ..
188: DOUBLE PRECISION ONE, ZERO
189: COMPLEX*16 CZERO
190: PARAMETER ( ONE = 1.0D+0, ZERO = 0.0D+0 )
191: PARAMETER ( CZERO = ( 0.0D+0, 0.0D+0 ) )
192: * ..
193: * .. Local Scalars ..
194: LOGICAL LOWER, NORMALTRANSR, NISODD, NOTRANS
195: INTEGER INFO, NROWA, J, NK, N1, N2
196: COMPLEX*16 CALPHA, CBETA
197: * ..
198: * .. External Functions ..
199: LOGICAL LSAME
200: EXTERNAL LSAME
201: * ..
202: * .. External Subroutines ..
203: EXTERNAL XERBLA, ZGEMM, ZHERK
204: * ..
205: * .. Intrinsic Functions ..
206: INTRINSIC MAX, DCMPLX
207: * ..
208: * .. Executable Statements ..
209: *
210: *
211: * Test the input parameters.
212: *
213: INFO = 0
214: NORMALTRANSR = LSAME( TRANSR, 'N' )
215: LOWER = LSAME( UPLO, 'L' )
216: NOTRANS = LSAME( TRANS, 'N' )
217: *
218: IF( NOTRANS ) THEN
219: NROWA = N
220: ELSE
221: NROWA = K
222: END IF
223: *
224: IF( .NOT.NORMALTRANSR .AND. .NOT.LSAME( TRANSR, 'C' ) ) THEN
225: INFO = -1
226: ELSE IF( .NOT.LOWER .AND. .NOT.LSAME( UPLO, 'U' ) ) THEN
227: INFO = -2
228: ELSE IF( .NOT.NOTRANS .AND. .NOT.LSAME( TRANS, 'C' ) ) THEN
229: INFO = -3
230: ELSE IF( N.LT.0 ) THEN
231: INFO = -4
232: ELSE IF( K.LT.0 ) THEN
233: INFO = -5
234: ELSE IF( LDA.LT.MAX( 1, NROWA ) ) THEN
235: INFO = -8
236: END IF
237: IF( INFO.NE.0 ) THEN
238: CALL XERBLA( 'ZHFRK ', -INFO )
239: RETURN
240: END IF
241: *
242: * Quick return if possible.
243: *
244: * The quick return case: ((ALPHA.EQ.0).AND.(BETA.NE.ZERO)) is not
245: * done (it is in ZHERK for example) and left in the general case.
246: *
247: IF( ( N.EQ.0 ) .OR. ( ( ( ALPHA.EQ.ZERO ) .OR. ( K.EQ.0 ) ) .AND.
248: $ ( BETA.EQ.ONE ) ) )RETURN
249: *
250: IF( ( ALPHA.EQ.ZERO ) .AND. ( BETA.EQ.ZERO ) ) THEN
251: DO J = 1, ( ( N*( N+1 ) ) / 2 )
252: C( J ) = CZERO
253: END DO
254: RETURN
255: END IF
256: *
257: CALPHA = DCMPLX( ALPHA, ZERO )
258: CBETA = DCMPLX( BETA, ZERO )
259: *
260: * C is N-by-N.
261: * If N is odd, set NISODD = .TRUE., and N1 and N2.
262: * If N is even, NISODD = .FALSE., and NK.
263: *
264: IF( MOD( N, 2 ).EQ.0 ) THEN
265: NISODD = .FALSE.
266: NK = N / 2
267: ELSE
268: NISODD = .TRUE.
269: IF( LOWER ) THEN
270: N2 = N / 2
271: N1 = N - N2
272: ELSE
273: N1 = N / 2
274: N2 = N - N1
275: END IF
276: END IF
277: *
278: IF( NISODD ) THEN
279: *
280: * N is odd
281: *
282: IF( NORMALTRANSR ) THEN
283: *
284: * N is odd and TRANSR = 'N'
285: *
286: IF( LOWER ) THEN
287: *
288: * N is odd, TRANSR = 'N', and UPLO = 'L'
289: *
290: IF( NOTRANS ) THEN
291: *
292: * N is odd, TRANSR = 'N', UPLO = 'L', and TRANS = 'N'
293: *
294: CALL ZHERK( 'L', 'N', N1, K, ALPHA, A( 1, 1 ), LDA,
295: $ BETA, C( 1 ), N )
296: CALL ZHERK( 'U', 'N', N2, K, ALPHA, A( N1+1, 1 ), LDA,
297: $ BETA, C( N+1 ), N )
298: CALL ZGEMM( 'N', 'C', N2, N1, K, CALPHA, A( N1+1, 1 ),
299: $ LDA, A( 1, 1 ), LDA, CBETA, C( N1+1 ), N )
300: *
301: ELSE
302: *
303: * N is odd, TRANSR = 'N', UPLO = 'L', and TRANS = 'C'
304: *
305: CALL ZHERK( 'L', 'C', N1, K, ALPHA, A( 1, 1 ), LDA,
306: $ BETA, C( 1 ), N )
307: CALL ZHERK( 'U', 'C', N2, K, ALPHA, A( 1, N1+1 ), LDA,
308: $ BETA, C( N+1 ), N )
309: CALL ZGEMM( 'C', 'N', N2, N1, K, CALPHA, A( 1, N1+1 ),
310: $ LDA, A( 1, 1 ), LDA, CBETA, C( N1+1 ), N )
311: *
312: END IF
313: *
314: ELSE
315: *
316: * N is odd, TRANSR = 'N', and UPLO = 'U'
317: *
318: IF( NOTRANS ) THEN
319: *
320: * N is odd, TRANSR = 'N', UPLO = 'U', and TRANS = 'N'
321: *
322: CALL ZHERK( 'L', 'N', N1, K, ALPHA, A( 1, 1 ), LDA,
323: $ BETA, C( N2+1 ), N )
324: CALL ZHERK( 'U', 'N', N2, K, ALPHA, A( N2, 1 ), LDA,
325: $ BETA, C( N1+1 ), N )
326: CALL ZGEMM( 'N', 'C', N1, N2, K, CALPHA, A( 1, 1 ),
327: $ LDA, A( N2, 1 ), LDA, CBETA, C( 1 ), N )
328: *
329: ELSE
330: *
331: * N is odd, TRANSR = 'N', UPLO = 'U', and TRANS = 'C'
332: *
333: CALL ZHERK( 'L', 'C', N1, K, ALPHA, A( 1, 1 ), LDA,
334: $ BETA, C( N2+1 ), N )
335: CALL ZHERK( 'U', 'C', N2, K, ALPHA, A( 1, N2 ), LDA,
336: $ BETA, C( N1+1 ), N )
337: CALL ZGEMM( 'C', 'N', N1, N2, K, CALPHA, A( 1, 1 ),
338: $ LDA, A( 1, N2 ), LDA, CBETA, C( 1 ), N )
339: *
340: END IF
341: *
342: END IF
343: *
344: ELSE
345: *
346: * N is odd, and TRANSR = 'C'
347: *
348: IF( LOWER ) THEN
349: *
350: * N is odd, TRANSR = 'C', and UPLO = 'L'
351: *
352: IF( NOTRANS ) THEN
353: *
354: * N is odd, TRANSR = 'C', UPLO = 'L', and TRANS = 'N'
355: *
356: CALL ZHERK( 'U', 'N', N1, K, ALPHA, A( 1, 1 ), LDA,
357: $ BETA, C( 1 ), N1 )
358: CALL ZHERK( 'L', 'N', N2, K, ALPHA, A( N1+1, 1 ), LDA,
359: $ BETA, C( 2 ), N1 )
360: CALL ZGEMM( 'N', 'C', N1, N2, K, CALPHA, A( 1, 1 ),
361: $ LDA, A( N1+1, 1 ), LDA, CBETA,
362: $ C( N1*N1+1 ), N1 )
363: *
364: ELSE
365: *
366: * N is odd, TRANSR = 'C', UPLO = 'L', and TRANS = 'C'
367: *
368: CALL ZHERK( 'U', 'C', N1, K, ALPHA, A( 1, 1 ), LDA,
369: $ BETA, C( 1 ), N1 )
370: CALL ZHERK( 'L', 'C', N2, K, ALPHA, A( 1, N1+1 ), LDA,
371: $ BETA, C( 2 ), N1 )
372: CALL ZGEMM( 'C', 'N', N1, N2, K, CALPHA, A( 1, 1 ),
373: $ LDA, A( 1, N1+1 ), LDA, CBETA,
374: $ C( N1*N1+1 ), N1 )
375: *
376: END IF
377: *
378: ELSE
379: *
380: * N is odd, TRANSR = 'C', and UPLO = 'U'
381: *
382: IF( NOTRANS ) THEN
383: *
384: * N is odd, TRANSR = 'C', UPLO = 'U', and TRANS = 'N'
385: *
386: CALL ZHERK( 'U', 'N', N1, K, ALPHA, A( 1, 1 ), LDA,
387: $ BETA, C( N2*N2+1 ), N2 )
388: CALL ZHERK( 'L', 'N', N2, K, ALPHA, A( N1+1, 1 ), LDA,
389: $ BETA, C( N1*N2+1 ), N2 )
390: CALL ZGEMM( 'N', 'C', N2, N1, K, CALPHA, A( N1+1, 1 ),
391: $ LDA, A( 1, 1 ), LDA, CBETA, C( 1 ), N2 )
392: *
393: ELSE
394: *
395: * N is odd, TRANSR = 'C', UPLO = 'U', and TRANS = 'C'
396: *
397: CALL ZHERK( 'U', 'C', N1, K, ALPHA, A( 1, 1 ), LDA,
398: $ BETA, C( N2*N2+1 ), N2 )
399: CALL ZHERK( 'L', 'C', N2, K, ALPHA, A( 1, N1+1 ), LDA,
400: $ BETA, C( N1*N2+1 ), N2 )
401: CALL ZGEMM( 'C', 'N', N2, N1, K, CALPHA, A( 1, N1+1 ),
402: $ LDA, A( 1, 1 ), LDA, CBETA, C( 1 ), N2 )
403: *
404: END IF
405: *
406: END IF
407: *
408: END IF
409: *
410: ELSE
411: *
412: * N is even
413: *
414: IF( NORMALTRANSR ) THEN
415: *
416: * N is even and TRANSR = 'N'
417: *
418: IF( LOWER ) THEN
419: *
420: * N is even, TRANSR = 'N', and UPLO = 'L'
421: *
422: IF( NOTRANS ) THEN
423: *
424: * N is even, TRANSR = 'N', UPLO = 'L', and TRANS = 'N'
425: *
426: CALL ZHERK( 'L', 'N', NK, K, ALPHA, A( 1, 1 ), LDA,
427: $ BETA, C( 2 ), N+1 )
428: CALL ZHERK( 'U', 'N', NK, K, ALPHA, A( NK+1, 1 ), LDA,
429: $ BETA, C( 1 ), N+1 )
430: CALL ZGEMM( 'N', 'C', NK, NK, K, CALPHA, A( NK+1, 1 ),
431: $ LDA, A( 1, 1 ), LDA, CBETA, C( NK+2 ),
432: $ N+1 )
433: *
434: ELSE
435: *
436: * N is even, TRANSR = 'N', UPLO = 'L', and TRANS = 'C'
437: *
438: CALL ZHERK( 'L', 'C', NK, K, ALPHA, A( 1, 1 ), LDA,
439: $ BETA, C( 2 ), N+1 )
440: CALL ZHERK( 'U', 'C', NK, K, ALPHA, A( 1, NK+1 ), LDA,
441: $ BETA, C( 1 ), N+1 )
442: CALL ZGEMM( 'C', 'N', NK, NK, K, CALPHA, A( 1, NK+1 ),
443: $ LDA, A( 1, 1 ), LDA, CBETA, C( NK+2 ),
444: $ N+1 )
445: *
446: END IF
447: *
448: ELSE
449: *
450: * N is even, TRANSR = 'N', and UPLO = 'U'
451: *
452: IF( NOTRANS ) THEN
453: *
454: * N is even, TRANSR = 'N', UPLO = 'U', and TRANS = 'N'
455: *
456: CALL ZHERK( 'L', 'N', NK, K, ALPHA, A( 1, 1 ), LDA,
457: $ BETA, C( NK+2 ), N+1 )
458: CALL ZHERK( 'U', 'N', NK, K, ALPHA, A( NK+1, 1 ), LDA,
459: $ BETA, C( NK+1 ), N+1 )
460: CALL ZGEMM( 'N', 'C', NK, NK, K, CALPHA, A( 1, 1 ),
461: $ LDA, A( NK+1, 1 ), LDA, CBETA, C( 1 ),
462: $ N+1 )
463: *
464: ELSE
465: *
466: * N is even, TRANSR = 'N', UPLO = 'U', and TRANS = 'C'
467: *
468: CALL ZHERK( 'L', 'C', NK, K, ALPHA, A( 1, 1 ), LDA,
469: $ BETA, C( NK+2 ), N+1 )
470: CALL ZHERK( 'U', 'C', NK, K, ALPHA, A( 1, NK+1 ), LDA,
471: $ BETA, C( NK+1 ), N+1 )
472: CALL ZGEMM( 'C', 'N', NK, NK, K, CALPHA, A( 1, 1 ),
473: $ LDA, A( 1, NK+1 ), LDA, CBETA, C( 1 ),
474: $ N+1 )
475: *
476: END IF
477: *
478: END IF
479: *
480: ELSE
481: *
482: * N is even, and TRANSR = 'C'
483: *
484: IF( LOWER ) THEN
485: *
486: * N is even, TRANSR = 'C', and UPLO = 'L'
487: *
488: IF( NOTRANS ) THEN
489: *
490: * N is even, TRANSR = 'C', UPLO = 'L', and TRANS = 'N'
491: *
492: CALL ZHERK( 'U', 'N', NK, K, ALPHA, A( 1, 1 ), LDA,
493: $ BETA, C( NK+1 ), NK )
494: CALL ZHERK( 'L', 'N', NK, K, ALPHA, A( NK+1, 1 ), LDA,
495: $ BETA, C( 1 ), NK )
496: CALL ZGEMM( 'N', 'C', NK, NK, K, CALPHA, A( 1, 1 ),
497: $ LDA, A( NK+1, 1 ), LDA, CBETA,
498: $ C( ( ( NK+1 )*NK )+1 ), NK )
499: *
500: ELSE
501: *
502: * N is even, TRANSR = 'C', UPLO = 'L', and TRANS = 'C'
503: *
504: CALL ZHERK( 'U', 'C', NK, K, ALPHA, A( 1, 1 ), LDA,
505: $ BETA, C( NK+1 ), NK )
506: CALL ZHERK( 'L', 'C', NK, K, ALPHA, A( 1, NK+1 ), LDA,
507: $ BETA, C( 1 ), NK )
508: CALL ZGEMM( 'C', 'N', NK, NK, K, CALPHA, A( 1, 1 ),
509: $ LDA, A( 1, NK+1 ), LDA, CBETA,
510: $ C( ( ( NK+1 )*NK )+1 ), NK )
511: *
512: END IF
513: *
514: ELSE
515: *
516: * N is even, TRANSR = 'C', and UPLO = 'U'
517: *
518: IF( NOTRANS ) THEN
519: *
520: * N is even, TRANSR = 'C', UPLO = 'U', and TRANS = 'N'
521: *
522: CALL ZHERK( 'U', 'N', NK, K, ALPHA, A( 1, 1 ), LDA,
523: $ BETA, C( NK*( NK+1 )+1 ), NK )
524: CALL ZHERK( 'L', 'N', NK, K, ALPHA, A( NK+1, 1 ), LDA,
525: $ BETA, C( NK*NK+1 ), NK )
526: CALL ZGEMM( 'N', 'C', NK, NK, K, CALPHA, A( NK+1, 1 ),
527: $ LDA, A( 1, 1 ), LDA, CBETA, C( 1 ), NK )
528: *
529: ELSE
530: *
531: * N is even, TRANSR = 'C', UPLO = 'U', and TRANS = 'C'
532: *
533: CALL ZHERK( 'U', 'C', NK, K, ALPHA, A( 1, 1 ), LDA,
534: $ BETA, C( NK*( NK+1 )+1 ), NK )
535: CALL ZHERK( 'L', 'C', NK, K, ALPHA, A( 1, NK+1 ), LDA,
536: $ BETA, C( NK*NK+1 ), NK )
537: CALL ZGEMM( 'C', 'N', NK, NK, K, CALPHA, A( 1, NK+1 ),
538: $ LDA, A( 1, 1 ), LDA, CBETA, C( 1 ), NK )
539: *
540: END IF
541: *
542: END IF
543: *
544: END IF
545: *
546: END IF
547: *
548: RETURN
549: *
550: * End of ZHFRK
551: *
552: END
CVSweb interface <joel.bertrand@systella.fr>
![[BACK]](/cvsweb/icons/back.gif){kind=icon}
Return to zhfrk.f CVS log ![[TXT]](/cvsweb/icons/text.gif){kind=icon}
![[DIR]](/cvsweb/icons/dir.gif){kind=icon}
Up to [local] / rpl / lapack / lapack