LAPACK  3.10.0
LAPACK: Linear Algebra PACKage

◆ cher2k()

subroutine cher2k ( character  UPLO,
character  TRANS,
integer  N,
integer  K,
complex  ALPHA,
complex, dimension(lda,*)  A,
integer  LDA,
complex, dimension(ldb,*)  B,
integer  LDB,
real  BETA,
complex, dimension(ldc,*)  C,
integer  LDC 
)

CHER2K

Purpose:
 CHER2K  performs one of the hermitian rank 2k operations

    C := alpha*A*B**H + conjg( alpha )*B*A**H + beta*C,

 or

    C := alpha*A**H*B + conjg( alpha )*B**H*A + beta*C,

 where  alpha and beta  are scalars with  beta  real,  C is an  n by n
 hermitian matrix and  A and B  are  n by k matrices in the first case
 and  k by n  matrices in the second case.
Parameters
[in]UPLO
          UPLO is CHARACTER*1
           On  entry,   UPLO  specifies  whether  the  upper  or  lower
           triangular  part  of the  array  C  is to be  referenced  as
           follows:

              UPLO = 'U' or 'u'   Only the  upper triangular part of  C
                                  is to be referenced.

              UPLO = 'L' or 'l'   Only the  lower triangular part of  C
                                  is to be referenced.
[in]TRANS
          TRANS is CHARACTER*1
           On entry,  TRANS  specifies the operation to be performed as
           follows:

              TRANS = 'N' or 'n'    C := alpha*A*B**H          +
                                         conjg( alpha )*B*A**H +
                                         beta*C.

              TRANS = 'C' or 'c'    C := alpha*A**H*B          +
                                         conjg( alpha )*B**H*A +
                                         beta*C.
[in]N
          N is INTEGER
           On entry,  N specifies the order of the matrix C.  N must be
           at least zero.
[in]K
          K is INTEGER
           On entry with  TRANS = 'N' or 'n',  K  specifies  the number
           of  columns  of the  matrices  A and B,  and on  entry  with
           TRANS = 'C' or 'c',  K  specifies  the number of rows of the
           matrices  A and B.  K must be at least zero.
[in]ALPHA
          ALPHA is COMPLEX
           On entry, ALPHA specifies the scalar alpha.
[in]A
          A is COMPLEX array, dimension ( LDA, ka ), where ka is
           k  when  TRANS = 'N' or 'n',  and is  n  otherwise.
           Before entry with  TRANS = 'N' or 'n',  the  leading  n by k
           part of the array  A  must contain the matrix  A,  otherwise
           the leading  k by n  part of the array  A  must contain  the
           matrix A.
[in]LDA
          LDA is INTEGER
           On entry, LDA specifies the first dimension of A as declared
           in  the  calling  (sub)  program.   When  TRANS = 'N' or 'n'
           then  LDA must be at least  max( 1, n ), otherwise  LDA must
           be at least  max( 1, k ).
[in]B
          B is COMPLEX array, dimension ( LDB, kb ), where kb is
           k  when  TRANS = 'N' or 'n',  and is  n  otherwise.
           Before entry with  TRANS = 'N' or 'n',  the  leading  n by k
           part of the array  B  must contain the matrix  B,  otherwise
           the leading  k by n  part of the array  B  must contain  the
           matrix B.
[in]LDB
          LDB is INTEGER
           On entry, LDB specifies the first dimension of B as declared
           in  the  calling  (sub)  program.   When  TRANS = 'N' or 'n'
           then  LDB must be at least  max( 1, n ), otherwise  LDB must
           be at least  max( 1, k ).
[in]BETA
          BETA is REAL
           On entry, BETA specifies the scalar beta.
[in,out]C
          C is COMPLEX array, dimension ( LDC, N )
           Before entry  with  UPLO = 'U' or 'u',  the leading  n by n
           upper triangular part of the array C must contain the upper
           triangular part  of the  hermitian matrix  and the strictly
           lower triangular part of C is not referenced.  On exit, the
           upper triangular part of the array  C is overwritten by the
           upper triangular part of the updated matrix.
           Before entry  with  UPLO = 'L' or 'l',  the leading  n by n
           lower triangular part of the array C must contain the lower
           triangular part  of the  hermitian matrix  and the strictly
           upper triangular part of C is not referenced.  On exit, the
           lower triangular part of the array  C is overwritten by the
           lower triangular part of the updated matrix.
           Note that the imaginary parts of the diagonal elements need
           not be set,  they are assumed to be zero,  and on exit they
           are set to zero.
[in]LDC
          LDC is INTEGER
           On entry, LDC specifies the first dimension of C as declared
           in  the  calling  (sub)  program.   LDC  must  be  at  least
           max( 1, n ).
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Further Details:
  Level 3 Blas routine.

  -- Written on 8-February-1989.
     Jack Dongarra, Argonne National Laboratory.
     Iain Duff, AERE Harwell.
     Jeremy Du Croz, Numerical Algorithms Group Ltd.
     Sven Hammarling, Numerical Algorithms Group Ltd.

  -- Modified 8-Nov-93 to set C(J,J) to REAL( C(J,J) ) when BETA = 1.
     Ed Anderson, Cray Research Inc.

Definition at line 196 of file cher2k.f.

197 *
198 * -- Reference BLAS level3 routine --
199 * -- Reference BLAS is a software package provided by Univ. of Tennessee, --
200 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
201 *
202 * .. Scalar Arguments ..
203  COMPLEX ALPHA
204  REAL BETA
205  INTEGER K,LDA,LDB,LDC,N
206  CHARACTER TRANS,UPLO
207 * ..
208 * .. Array Arguments ..
209  COMPLEX A(LDA,*),B(LDB,*),C(LDC,*)
210 * ..
211 *
212 * =====================================================================
213 *
214 * .. External Functions ..
215  LOGICAL LSAME
216  EXTERNAL lsame
217 * ..
218 * .. External Subroutines ..
219  EXTERNAL xerbla
220 * ..
221 * .. Intrinsic Functions ..
222  INTRINSIC conjg,max,real
223 * ..
224 * .. Local Scalars ..
225  COMPLEX TEMP1,TEMP2
226  INTEGER I,INFO,J,L,NROWA
227  LOGICAL UPPER
228 * ..
229 * .. Parameters ..
230  REAL ONE
231  parameter(one=1.0e+0)
232  COMPLEX ZERO
233  parameter(zero= (0.0e+0,0.0e+0))
234 * ..
235 *
236 * Test the input parameters.
237 *
238  IF (lsame(trans,'N')) THEN
239  nrowa = n
240  ELSE
241  nrowa = k
242  END IF
243  upper = lsame(uplo,'U')
244 *
245  info = 0
246  IF ((.NOT.upper) .AND. (.NOT.lsame(uplo,'L'))) THEN
247  info = 1
248  ELSE IF ((.NOT.lsame(trans,'N')) .AND.
249  + (.NOT.lsame(trans,'C'))) THEN
250  info = 2
251  ELSE IF (n.LT.0) THEN
252  info = 3
253  ELSE IF (k.LT.0) THEN
254  info = 4
255  ELSE IF (lda.LT.max(1,nrowa)) THEN
256  info = 7
257  ELSE IF (ldb.LT.max(1,nrowa)) THEN
258  info = 9
259  ELSE IF (ldc.LT.max(1,n)) THEN
260  info = 12
261  END IF
262  IF (info.NE.0) THEN
263  CALL xerbla('CHER2K',info)
264  RETURN
265  END IF
266 *
267 * Quick return if possible.
268 *
269  IF ((n.EQ.0) .OR. (((alpha.EQ.zero).OR.
270  + (k.EQ.0)).AND. (beta.EQ.one))) RETURN
271 *
272 * And when alpha.eq.zero.
273 *
274  IF (alpha.EQ.zero) THEN
275  IF (upper) THEN
276  IF (beta.EQ.real(zero)) THEN
277  DO 20 j = 1,n
278  DO 10 i = 1,j
279  c(i,j) = zero
280  10 CONTINUE
281  20 CONTINUE
282  ELSE
283  DO 40 j = 1,n
284  DO 30 i = 1,j - 1
285  c(i,j) = beta*c(i,j)
286  30 CONTINUE
287  c(j,j) = beta*real(c(j,j))
288  40 CONTINUE
289  END IF
290  ELSE
291  IF (beta.EQ.real(zero)) THEN
292  DO 60 j = 1,n
293  DO 50 i = j,n
294  c(i,j) = zero
295  50 CONTINUE
296  60 CONTINUE
297  ELSE
298  DO 80 j = 1,n
299  c(j,j) = beta*real(c(j,j))
300  DO 70 i = j + 1,n
301  c(i,j) = beta*c(i,j)
302  70 CONTINUE
303  80 CONTINUE
304  END IF
305  END IF
306  RETURN
307  END IF
308 *
309 * Start the operations.
310 *
311  IF (lsame(trans,'N')) THEN
312 *
313 * Form C := alpha*A*B**H + conjg( alpha )*B*A**H +
314 * C.
315 *
316  IF (upper) THEN
317  DO 130 j = 1,n
318  IF (beta.EQ.real(zero)) THEN
319  DO 90 i = 1,j
320  c(i,j) = zero
321  90 CONTINUE
322  ELSE IF (beta.NE.one) THEN
323  DO 100 i = 1,j - 1
324  c(i,j) = beta*c(i,j)
325  100 CONTINUE
326  c(j,j) = beta*real(c(j,j))
327  ELSE
328  c(j,j) = real(c(j,j))
329  END IF
330  DO 120 l = 1,k
331  IF ((a(j,l).NE.zero) .OR. (b(j,l).NE.zero)) THEN
332  temp1 = alpha*conjg(b(j,l))
333  temp2 = conjg(alpha*a(j,l))
334  DO 110 i = 1,j - 1
335  c(i,j) = c(i,j) + a(i,l)*temp1 +
336  + b(i,l)*temp2
337  110 CONTINUE
338  c(j,j) = real(c(j,j)) +
339  + real(a(j,l)*temp1+b(j,l)*temp2)
340  END IF
341  120 CONTINUE
342  130 CONTINUE
343  ELSE
344  DO 180 j = 1,n
345  IF (beta.EQ.real(zero)) THEN
346  DO 140 i = j,n
347  c(i,j) = zero
348  140 CONTINUE
349  ELSE IF (beta.NE.one) THEN
350  DO 150 i = j + 1,n
351  c(i,j) = beta*c(i,j)
352  150 CONTINUE
353  c(j,j) = beta*real(c(j,j))
354  ELSE
355  c(j,j) = real(c(j,j))
356  END IF
357  DO 170 l = 1,k
358  IF ((a(j,l).NE.zero) .OR. (b(j,l).NE.zero)) THEN
359  temp1 = alpha*conjg(b(j,l))
360  temp2 = conjg(alpha*a(j,l))
361  DO 160 i = j + 1,n
362  c(i,j) = c(i,j) + a(i,l)*temp1 +
363  + b(i,l)*temp2
364  160 CONTINUE
365  c(j,j) = real(c(j,j)) +
366  + real(a(j,l)*temp1+b(j,l)*temp2)
367  END IF
368  170 CONTINUE
369  180 CONTINUE
370  END IF
371  ELSE
372 *
373 * Form C := alpha*A**H*B + conjg( alpha )*B**H*A +
374 * C.
375 *
376  IF (upper) THEN
377  DO 210 j = 1,n
378  DO 200 i = 1,j
379  temp1 = zero
380  temp2 = zero
381  DO 190 l = 1,k
382  temp1 = temp1 + conjg(a(l,i))*b(l,j)
383  temp2 = temp2 + conjg(b(l,i))*a(l,j)
384  190 CONTINUE
385  IF (i.EQ.j) THEN
386  IF (beta.EQ.real(zero)) THEN
387  c(j,j) = real(alpha*temp1+
388  + conjg(alpha)*temp2)
389  ELSE
390  c(j,j) = beta*real(c(j,j)) +
391  + real(alpha*temp1+
392  + conjg(alpha)*temp2)
393  END IF
394  ELSE
395  IF (beta.EQ.real(zero)) THEN
396  c(i,j) = alpha*temp1 + conjg(alpha)*temp2
397  ELSE
398  c(i,j) = beta*c(i,j) + alpha*temp1 +
399  + conjg(alpha)*temp2
400  END IF
401  END IF
402  200 CONTINUE
403  210 CONTINUE
404  ELSE
405  DO 240 j = 1,n
406  DO 230 i = j,n
407  temp1 = zero
408  temp2 = zero
409  DO 220 l = 1,k
410  temp1 = temp1 + conjg(a(l,i))*b(l,j)
411  temp2 = temp2 + conjg(b(l,i))*a(l,j)
412  220 CONTINUE
413  IF (i.EQ.j) THEN
414  IF (beta.EQ.real(zero)) THEN
415  c(j,j) = real(alpha*temp1+
416  + conjg(alpha)*temp2)
417  ELSE
418  c(j,j) = beta*real(c(j,j)) +
419  + real(alpha*temp1+
420  + conjg(alpha)*temp2)
421  END IF
422  ELSE
423  IF (beta.EQ.real(zero)) THEN
424  c(i,j) = alpha*temp1 + conjg(alpha)*temp2
425  ELSE
426  c(i,j) = beta*c(i,j) + alpha*temp1 +
427  + conjg(alpha)*temp2
428  END IF
429  END IF
430  230 CONTINUE
431  240 CONTINUE
432  END IF
433  END IF
434 *
435  RETURN
436 *
437 * End of CHER2K
438 *
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:53
Here is the call graph for this function:
Here is the caller graph for this function: