OpenRadioss 2025.1.11
OpenRadioss project
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*> \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, 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*> \ingroup complex16OTHERcomputational
164*
165* =====================================================================
166 SUBROUTINE zhfrk( TRANSR, UPLO, TRANS, N, K, ALPHA, A, LDA, BETA,
167 $ C )
168*
169* -- LAPACK computational routine --
170* -- LAPACK is a software package provided by Univ. of Tennessee, --
171* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
172*
173* .. Scalar Arguments ..
174 DOUBLE PRECISION ALPHA, BETA
175 INTEGER K, LDA, N
176 CHARACTER TRANS, TRANSR, UPLO
177* ..
178* .. Array Arguments ..
179 COMPLEX*16 A( LDA, * ), C( * )
180* ..
181*
182* =====================================================================
183*
184* .. Parameters ..
185 DOUBLE PRECISION ONE, ZERO
186 COMPLEX*16 CZERO
187 parameter( one = 1.0d+0, zero = 0.0d+0 )
188 parameter( czero = ( 0.0d+0, 0.0d+0 ) )
189* ..
190* .. Local Scalars ..
191 LOGICAL LOWER, NORMALTRANSR, NISODD, NOTRANS
192 INTEGER INFO, NROWA, J, NK, N1, N2
193 COMPLEX*16 CALPHA, CBETA
194* ..
195* .. External Functions ..
196 LOGICAL LSAME
197 EXTERNAL lsame
198* ..
199* .. External Subroutines ..
200 EXTERNAL xerbla, zgemm, zherk
201* ..
202* .. Intrinsic Functions ..
203 INTRINSIC max, dcmplx
204* ..
205* .. Executable Statements ..
206*
207*
208* Test the input parameters.
209*
210 info = 0
211 normaltransr = lsame( transr, 'N' )
212 lower = lsame( uplo, 'L' )
213 notrans = lsame( trans, 'N' )
214*
215 IF( notrans ) THEN
216 nrowa = n
217 ELSE
218 nrowa = k
219 END IF
220*
221 IF( .NOT.normaltransr .AND. .NOT.lsame( transr, 'C' ) ) THEN
222 info = -1
223 ELSE IF( .NOT.lower .AND. .NOT.lsame( uplo, 'U' ) ) THEN
224 info = -2
225 ELSE IF( .NOT.notrans .AND. .NOT.lsame( trans, 'C' ) ) THEN
226 info = -3
227 ELSE IF( n.LT.0 ) THEN
228 info = -4
229 ELSE IF( k.LT.0 ) THEN
230 info = -5
231 ELSE IF( lda.LT.max( 1, nrowa ) ) THEN
232 info = -8
233 END IF
234 IF( info.NE.0 ) THEN
235 CALL xerbla( 'ZHFRK ', -info )
236 RETURN
237 END IF
238*
239* Quick return if possible.
240*
241* The quick return case: ((ALPHA.EQ.0).AND.(BETA.NE.ZERO)) is not
242* done (it is in ZHERK for example) and left in the general case.
243*
244 IF( ( n.EQ.0 ) .OR. ( ( ( alpha.EQ.zero ) .OR. ( k.EQ.0 ) ) .AND.
245 $ ( beta.EQ.one ) ) )RETURN
246*
247 IF( ( alpha.EQ.zero ) .AND. ( beta.EQ.zero ) ) THEN
248 DO j = 1, ( ( n*( n+1 ) ) / 2 )
249 c( j ) = czero
250 END DO
251 RETURN
252 END IF
253*
254 calpha = dcmplx( alpha, zero )
255 cbeta = dcmplx( beta, zero )
256*
257* C is N-by-N.
258* If N is odd, set NISODD = .TRUE., and N1 and N2.
259* If N is even, NISODD = .FALSE., and NK.
260*
261 IF( mod( n, 2 ).EQ.0 ) THEN
262 nisodd = .false.
263 nk = n / 2
264 ELSE
265 nisodd = .true.
266 IF( lower ) THEN
267 n2 = n / 2
268 n1 = n - n2
269 ELSE
270 n1 = n / 2
271 n2 = n - n1
272 END IF
273 END IF
274*
275 IF( nisodd ) THEN
276*
277* N is odd
278*
279 IF( normaltransr ) THEN
280*
281* N is odd and TRANSR = 'N'
282*
283 IF( lower ) THEN
284*
285* N is odd, TRANSR = 'N', and UPLO = 'L'
286*
287 IF( notrans ) THEN
288*
289* N is odd, TRANSR = 'N', UPLO = 'L', and TRANS = 'N'
290*
291 CALL zherk( 'L', 'N', n1, k, alpha, a( 1, 1 ), lda,
292 $ beta, c( 1 ), n )
293 CALL zherk( 'U', 'N', n2, k, alpha, a( n1+1, 1 ), lda,
294 $ beta, c( n+1 ), n )
295 CALL zgemm( 'N', 'C', n2, n1, k, calpha, a( n1+1, 1 ),
296 $ lda, a( 1, 1 ), lda, cbeta, c( n1+1 ), n )
297*
298 ELSE
299*
300* N is odd, TRANSR = 'N', UPLO = 'L', and TRANS = 'C'
301*
302 CALL zherk( 'L', 'C', n1, k, alpha, a( 1, 1 ), lda,
303 $ beta, c( 1 ), n )
304 CALL zherk( 'U', 'C', n2, k, alpha, a( 1, n1+1 ), lda,
305 $ beta, c( n+1 ), n )
306 CALL zgemm( 'C', 'N', n2, n1, k, calpha, a( 1, n1+1 ),
307 $ lda, a( 1, 1 ), lda, cbeta, c( n1+1 ), n )
308*
309 END IF
310*
311 ELSE
312*
313* N is odd, TRANSR = 'N', and UPLO = 'U'
314*
315 IF( notrans ) THEN
316*
317* N is odd, TRANSR = 'N', UPLO = 'U', and TRANS = 'N'
318*
319 CALL zherk( 'L', 'N', n1, k, alpha, a( 1, 1 ), lda,
320 $ beta, c( n2+1 ), n )
321 CALL zherk( 'U', 'N', n2, k, alpha, a( n2, 1 ), lda,
322 $ beta, c( n1+1 ), n )
323 CALL zgemm( 'N', 'C', n1, n2, k, calpha, a( 1, 1 ),
324 $ lda, a( n2, 1 ), lda, cbeta, c( 1 ), n )
325*
326 ELSE
327*
328* N is odd, TRANSR = 'N', UPLO = 'U', and TRANS = 'C'
329*
330 CALL zherk( 'L', 'C', n1, k, alpha, a( 1, 1 ), lda,
331 $ beta, c( n2+1 ), n )
332 CALL zherk( 'U', 'C', n2, k, alpha, a( 1, n2 ), lda,
333 $ beta, c( n1+1 ), n )
334 CALL zgemm( 'C', 'N', n1, n2, k, calpha, a( 1, 1 ),
335 $ lda, a( 1, n2 ), lda, cbeta, c( 1 ), n )
336*
337 END IF
338*
339 END IF
340*
341 ELSE
342*
343* N is odd, and TRANSR = 'C'
344*
345 IF( lower ) THEN
346*
347* N is odd, TRANSR = 'C', and UPLO = 'L'
348*
349 IF( notrans ) THEN
350*
351* N is odd, TRANSR = 'C', UPLO = 'L', and TRANS = 'N'
352*
353 CALL zherk( 'U', 'N', n1, k, alpha, a( 1, 1 ), lda,
354 $ beta, c( 1 ), n1 )
355 CALL zherk( 'L', 'N', n2, k, alpha, a( n1+1, 1 ), lda,
356 $ beta, c( 2 ), n1 )
357 CALL zgemm( 'N', 'C', n1, n2, k, calpha, a( 1, 1 ),
358 $ lda, a( n1+1, 1 ), lda, cbeta,
359 $ c( n1*n1+1 ), n1 )
360*
361 ELSE
362*
363* N is odd, TRANSR = 'C', UPLO = 'L', and TRANS = 'C'
364*
365 CALL zherk( 'U', 'C', n1, k, alpha, a( 1, 1 ), lda,
366 $ beta, c( 1 ), n1 )
367 CALL zherk( 'L', 'C', n2, k, alpha, a( 1, n1+1 ), lda,
368 $ beta, c( 2 ), n1 )
369 CALL zgemm( 'C', 'N', n1, n2, k, calpha, a( 1, 1 ),
370 $ lda, a( 1, n1+1 ), lda, cbeta,
371 $ c( n1*n1+1 ), n1 )
372*
373 END IF
374*
375 ELSE
376*
377* N is odd, TRANSR = 'C', and UPLO = 'U'
378*
379 IF( notrans ) THEN
380*
381* N is odd, TRANSR = 'C', UPLO = 'U', and TRANS = 'N'
382*
383 CALL zherk( 'U', 'N', n1, k, alpha, a( 1, 1 ), lda,
384 $ beta, c( n2*n2+1 ), n2 )
385 CALL zherk( 'L', 'N', n2, k, alpha, a( n1+1, 1 ), lda,
386 $ beta, c( n1*n2+1 ), n2 )
387 CALL zgemm( 'N', 'C', n2, n1, k, calpha, a( n1+1, 1 ),
388 $ lda, a( 1, 1 ), lda, cbeta, c( 1 ), n2 )
389*
390 ELSE
391*
392* N is odd, TRANSR = 'C', UPLO = 'U', and TRANS = 'C'
393*
394 CALL zherk( 'U', 'C', n1, k, alpha, a( 1, 1 ), lda,
395 $ beta, c( n2*n2+1 ), n2 )
396 CALL zherk( 'L', 'C', n2, k, alpha, a( 1, n1+1 ), lda,
397 $ beta, c( n1*n2+1 ), n2 )
398 CALL zgemm( 'C', 'N', n2, n1, k, calpha, a( 1, n1+1 ),
399 $ lda, a( 1, 1 ), lda, cbeta, c( 1 ), n2 )
400*
401 END IF
402*
403 END IF
404*
405 END IF
406*
407 ELSE
408*
409* N is even
410*
411 IF( normaltransr ) THEN
412*
413* N is even and TRANSR = 'N'
414*
415 IF( lower ) THEN
416*
417* N is even, TRANSR = 'N', and UPLO = 'L'
418*
419 IF( notrans ) THEN
420*
421* N is even, TRANSR = 'N', UPLO = 'L', and TRANS = 'N'
422*
423 CALL zherk( 'L', 'N', nk, k, alpha, a( 1, 1 ), lda,
424 $ beta, c( 2 ), n+1 )
425 CALL zherk( 'U', 'N', nk, k, alpha, a( nk+1, 1 ), lda,
426 $ beta, c( 1 ), n+1 )
427 CALL zgemm( 'N', 'C', nk, nk, k, calpha, a( nk+1, 1 ),
428 $ lda, a( 1, 1 ), lda, cbeta, c( nk+2 ),
429 $ n+1 )
430*
431 ELSE
432*
433* N is even, TRANSR = 'N', UPLO = 'L', and TRANS = 'C'
434*
435 CALL zherk( 'L', 'C', nk, k, alpha, a( 1, 1 ), lda,
436 $ beta, c( 2 ), n+1 )
437 CALL zherk( 'u', 'c', NK, K, ALPHA, A( 1, NK+1 ), LDA,
438 $ BETA, C( 1 ), N+1 )
439 CALL ZGEMM( 'c', 'n', NK, NK, K, CALPHA, A( 1, NK+1 ),
440 $ LDA, A( 1, 1 ), LDA, CBETA, C( NK+2 ),
441 $ N+1 )
442*
443 END IF
444*
445 ELSE
446*
447* N is even, TRANSR = 'N', and UPLO = 'U'
448*
449 IF( NOTRANS ) THEN
450*
451* N is even, TRANSR = 'N', UPLO = 'U', and TRANS = 'N'
452*
453 CALL ZHERK( 'l', 'n', NK, K, ALPHA, A( 1, 1 ), LDA,
454 $ BETA, C( NK+2 ), N+1 )
455 CALL ZHERK( 'u', 'n', NK, K, ALPHA, A( NK+1, 1 ), LDA,
456 $ BETA, C( NK+1 ), N+1 )
457 CALL ZGEMM( 'n', 'c', NK, NK, K, CALPHA, A( 1, 1 ),
458 $ LDA, A( NK+1, 1 ), LDA, CBETA, C( 1 ),
459 $ N+1 )
460*
461 ELSE
462*
463* N is even, TRANSR = 'N', UPLO = 'U', and TRANS = 'C'
464*
465 CALL ZHERK( 'l', 'c', NK, K, ALPHA, A( 1, 1 ), LDA,
466 $ BETA, C( NK+2 ), N+1 )
467 CALL ZHERK( 'u', 'c', NK, K, ALPHA, A( 1, NK+1 ), LDA,
468 $ BETA, C( NK+1 ), N+1 )
469 CALL ZGEMM( 'c', 'n', NK, NK, K, CALPHA, A( 1, 1 ),
470 $ LDA, A( 1, NK+1 ), LDA, CBETA, C( 1 ),
471 $ N+1 )
472*
473 END IF
474*
475 END IF
476*
477 ELSE
478*
479* N is even, and TRANSR = 'C'
480*
481 IF( LOWER ) THEN
482*
483* N is even, TRANSR = 'C', and UPLO = 'L'
484*
485 IF( NOTRANS ) THEN
486*
487* N is even, TRANSR = 'C', UPLO = 'L', and TRANS = 'N'
488*
489 CALL ZHERK( 'u', 'n', NK, K, ALPHA, A( 1, 1 ), LDA,
490 $ BETA, C( NK+1 ), NK )
491 CALL ZHERK( 'l', 'n', NK, K, ALPHA, A( NK+1, 1 ), LDA,
492 $ BETA, C( 1 ), NK )
493 CALL ZGEMM( 'n', 'c', NK, NK, K, CALPHA, A( 1, 1 ),
494 $ LDA, A( NK+1, 1 ), LDA, CBETA,
495 $ C( ( ( NK+1 )*NK )+1 ), NK )
496*
497 ELSE
498*
499* N is even, TRANSR = 'C', UPLO = 'L', and TRANS = 'C'
500*
501 CALL ZHERK( 'u', 'c', NK, K, ALPHA, A( 1, 1 ), LDA,
502 $ BETA, C( NK+1 ), NK )
503 CALL ZHERK( 'l', 'c', NK, K, ALPHA, A( 1, NK+1 ), LDA,
504 $ BETA, C( 1 ), NK )
505 CALL ZGEMM( 'c', 'n', NK, NK, K, CALPHA, A( 1, 1 ),
506 $ LDA, A( 1, NK+1 ), LDA, CBETA,
507 $ C( ( ( NK+1 )*NK )+1 ), NK )
508*
509 END IF
510*
511 ELSE
512*
513* N is even, TRANSR = 'C', and UPLO = 'U'
514*
515 IF( NOTRANS ) THEN
516*
517* N is even, TRANSR = 'C', UPLO = 'U', and TRANS = 'N'
518*
519 CALL ZHERK( 'u', 'n', NK, K, ALPHA, A( 1, 1 ), LDA,
520 $ BETA, C( NK*( NK+1 )+1 ), NK )
521 CALL ZHERK( 'l', 'n', NK, K, ALPHA, A( NK+1, 1 ), LDA,
522 $ BETA, C( NK*NK+1 ), NK )
523 CALL ZGEMM( 'n', 'c', NK, NK, K, CALPHA, A( NK+1, 1 ),
524 $ LDA, A( 1, 1 ), LDA, CBETA, C( 1 ), NK )
525*
526 ELSE
527*
528* N is even, TRANSR = 'C', UPLO = 'U', and TRANS = 'C'
529*
530 CALL ZHERK( 'u', 'C', nk, k, alpha, a( 1, 1 ), lda,
531 $ beta, c( nk*( nk+1 )+1 ), nk )
532 CALL zherk( 'L', 'C', nk, k, alpha, a( 1, nk+1 ), lda,
533 $ beta, c( nk*nk+1 ), nk )
534 CALL zgemm( 'C', 'N', nk, nk, k, calpha, a( 1, nk+1 ),
535 $ lda, a( 1, 1 ), lda, cbeta, c( 1 ), nk )
536*
537 END IF
538*
539 END IF
540*
541 END IF
542*
543 END IF
544*
545 RETURN
546*
547* End of ZHFRK
548*
549 END
xerbla
subroutine xerbla(srname, info)
XERBLA
Definition xerbla.f:60
zhfrk
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:168
zherk
subroutine zherk(uplo, trans, n, k, alpha, a, lda, beta, c, ldc)
ZHERK
Definition zherk.f:173
zgemm
subroutine zgemm(transa, transb, m, n, k, alpha, a, lda, b, ldb, beta, c, ldc)
ZGEMM
Definition zgemm.f:187
max
#define max(a, b)
Definition macros.h:21