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