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