File:  [local] / rpl / lapack / blas / zherk.f
Revision 1.16: download - view: text, annotated - select for diffs - revision graph
Mon Aug 7 08:38:45 2023 UTC (8 months, 3 weeks ago) by bertrand
Branches: MAIN
CVS tags: rpl-4_1_35, rpl-4_1_34, HEAD
Première mise à jour de lapack et blas.

    1: *> \brief \b ZHERK
    2: *
    3: *  =========== DOCUMENTATION ===========
    4: *
    5: * Online html documentation available at
    6: *            http://www.netlib.org/lapack/explore-html/
    7: *
    8: *  Definition:
    9: *  ===========
   10: *
   11: *       SUBROUTINE ZHERK(UPLO,TRANS,N,K,ALPHA,A,LDA,BETA,C,LDC)
   12: *
   13: *       .. Scalar Arguments ..
   14: *       DOUBLE PRECISION ALPHA,BETA
   15: *       INTEGER K,LDA,LDC,N
   16: *       CHARACTER TRANS,UPLO
   17: *       ..
   18: *       .. Array Arguments ..
   19: *       COMPLEX*16 A(LDA,*),C(LDC,*)
   20: *       ..
   21: *
   22: *
   23: *> \par Purpose:
   24: *  =============
   25: *>
   26: *> \verbatim
   27: *>
   28: *> ZHERK  performs one of the hermitian rank k operations
   29: *>
   30: *>    C := alpha*A*A**H + beta*C,
   31: *>
   32: *> or
   33: *>
   34: *>    C := alpha*A**H*A + beta*C,
   35: *>
   36: *> where  alpha and beta  are  real scalars,  C is an  n by n  hermitian
   37: *> matrix and  A  is an  n by k  matrix in the  first case and a  k by n
   38: *> matrix in the second case.
   39: *> \endverbatim
   40: *
   41: *  Arguments:
   42: *  ==========
   43: *
   44: *> \param[in] UPLO
   45: *> \verbatim
   46: *>          UPLO is CHARACTER*1
   47: *>           On  entry,   UPLO  specifies  whether  the  upper  or  lower
   48: *>           triangular  part  of the  array  C  is to be  referenced  as
   49: *>           follows:
   50: *>
   51: *>              UPLO = 'U' or 'u'   Only the  upper triangular part of  C
   52: *>                                  is to be referenced.
   53: *>
   54: *>              UPLO = 'L' or 'l'   Only the  lower triangular part of  C
   55: *>                                  is to be referenced.
   56: *> \endverbatim
   57: *>
   58: *> \param[in] TRANS
   59: *> \verbatim
   60: *>          TRANS is CHARACTER*1
   61: *>           On entry,  TRANS  specifies the operation to be performed as
   62: *>           follows:
   63: *>
   64: *>              TRANS = 'N' or 'n'   C := alpha*A*A**H + beta*C.
   65: *>
   66: *>              TRANS = 'C' or 'c'   C := alpha*A**H*A + beta*C.
   67: *> \endverbatim
   68: *>
   69: *> \param[in] N
   70: *> \verbatim
   71: *>          N is INTEGER
   72: *>           On entry,  N specifies the order of the matrix C.  N must be
   73: *>           at least zero.
   74: *> \endverbatim
   75: *>
   76: *> \param[in] K
   77: *> \verbatim
   78: *>          K is INTEGER
   79: *>           On entry with  TRANS = 'N' or 'n',  K  specifies  the number
   80: *>           of  columns   of  the   matrix   A,   and  on   entry   with
   81: *>           TRANS = 'C' or 'c',  K  specifies  the number of rows of the
   82: *>           matrix A.  K must be at least zero.
   83: *> \endverbatim
   84: *>
   85: *> \param[in] ALPHA
   86: *> \verbatim
   87: *>          ALPHA is DOUBLE PRECISION .
   88: *>           On entry, ALPHA specifies the scalar alpha.
   89: *> \endverbatim
   90: *>
   91: *> \param[in] A
   92: *> \verbatim
   93: *>          A is COMPLEX*16 array, dimension ( LDA, ka ), where ka is
   94: *>           k  when  TRANS = 'N' or 'n',  and is  n  otherwise.
   95: *>           Before entry with  TRANS = 'N' or 'n',  the  leading  n by k
   96: *>           part of the array  A  must contain the matrix  A,  otherwise
   97: *>           the leading  k by n  part of the array  A  must contain  the
   98: *>           matrix A.
   99: *> \endverbatim
  100: *>
  101: *> \param[in] LDA
  102: *> \verbatim
  103: *>          LDA is INTEGER
  104: *>           On entry, LDA specifies the first dimension of A as declared
  105: *>           in  the  calling  (sub)  program.   When  TRANS = 'N' or 'n'
  106: *>           then  LDA must be at least  max( 1, n ), otherwise  LDA must
  107: *>           be at least  max( 1, k ).
  108: *> \endverbatim
  109: *>
  110: *> \param[in] BETA
  111: *> \verbatim
  112: *>          BETA is DOUBLE PRECISION.
  113: *>           On entry, BETA specifies the scalar beta.
  114: *> \endverbatim
  115: *>
  116: *> \param[in,out] C
  117: *> \verbatim
  118: *>          C is COMPLEX*16 array, dimension ( LDC, N )
  119: *>           Before entry  with  UPLO = 'U' or 'u',  the leading  n by n
  120: *>           upper triangular part of the array C must contain the upper
  121: *>           triangular part  of the  hermitian matrix  and the strictly
  122: *>           lower triangular part of C is not referenced.  On exit, the
  123: *>           upper triangular part of the array  C is overwritten by the
  124: *>           upper triangular part of the updated matrix.
  125: *>           Before entry  with  UPLO = 'L' or 'l',  the leading  n by n
  126: *>           lower triangular part of the array C must contain the lower
  127: *>           triangular part  of the  hermitian matrix  and the strictly
  128: *>           upper triangular part of C is not referenced.  On exit, the
  129: *>           lower triangular part of the array  C is overwritten by the
  130: *>           lower triangular part of the updated matrix.
  131: *>           Note that the imaginary parts of the diagonal elements need
  132: *>           not be set,  they are assumed to be zero,  and on exit they
  133: *>           are set to zero.
  134: *> \endverbatim
  135: *>
  136: *> \param[in] LDC
  137: *> \verbatim
  138: *>          LDC is INTEGER
  139: *>           On entry, LDC specifies the first dimension of C as declared
  140: *>           in  the  calling  (sub)  program.   LDC  must  be  at  least
  141: *>           max( 1, n ).
  142: *> \endverbatim
  143: *
  144: *  Authors:
  145: *  ========
  146: *
  147: *> \author Univ. of Tennessee
  148: *> \author Univ. of California Berkeley
  149: *> \author Univ. of Colorado Denver
  150: *> \author NAG Ltd.
  151: *
  152: *> \ingroup complex16_blas_level3
  153: *
  154: *> \par Further Details:
  155: *  =====================
  156: *>
  157: *> \verbatim
  158: *>
  159: *>  Level 3 Blas routine.
  160: *>
  161: *>  -- Written on 8-February-1989.
  162: *>     Jack Dongarra, Argonne National Laboratory.
  163: *>     Iain Duff, AERE Harwell.
  164: *>     Jeremy Du Croz, Numerical Algorithms Group Ltd.
  165: *>     Sven Hammarling, Numerical Algorithms Group Ltd.
  166: *>
  167: *>  -- Modified 8-Nov-93 to set C(J,J) to DBLE( C(J,J) ) when BETA = 1.
  168: *>     Ed Anderson, Cray Research Inc.
  169: *> \endverbatim
  170: *>
  171: *  =====================================================================
  172:       SUBROUTINE ZHERK(UPLO,TRANS,N,K,ALPHA,A,LDA,BETA,C,LDC)
  173: *
  174: *  -- Reference BLAS level3 routine --
  175: *  -- Reference BLAS is a software package provided by Univ. of Tennessee,    --
  176: *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
  177: *
  178: *     .. Scalar Arguments ..
  179:       DOUBLE PRECISION ALPHA,BETA
  180:       INTEGER K,LDA,LDC,N
  181:       CHARACTER TRANS,UPLO
  182: *     ..
  183: *     .. Array Arguments ..
  184:       COMPLEX*16 A(LDA,*),C(LDC,*)
  185: *     ..
  186: *
  187: *  =====================================================================
  188: *
  189: *     .. External Functions ..
  190:       LOGICAL LSAME
  191:       EXTERNAL LSAME
  192: *     ..
  193: *     .. External Subroutines ..
  194:       EXTERNAL XERBLA
  195: *     ..
  196: *     .. Intrinsic Functions ..
  197:       INTRINSIC DBLE,DCMPLX,DCONJG,MAX
  198: *     ..
  199: *     .. Local Scalars ..
  200:       COMPLEX*16 TEMP
  201:       DOUBLE PRECISION RTEMP
  202:       INTEGER I,INFO,J,L,NROWA
  203:       LOGICAL UPPER
  204: *     ..
  205: *     .. Parameters ..
  206:       DOUBLE PRECISION ONE,ZERO
  207:       PARAMETER (ONE=1.0D+0,ZERO=0.0D+0)
  208: *     ..
  209: *
  210: *     Test the input parameters.
  211: *
  212:       IF (LSAME(TRANS,'N')) THEN
  213:           NROWA = N
  214:       ELSE
  215:           NROWA = K
  216:       END IF
  217:       UPPER = LSAME(UPLO,'U')
  218: *
  219:       INFO = 0
  220:       IF ((.NOT.UPPER) .AND. (.NOT.LSAME(UPLO,'L'))) THEN
  221:           INFO = 1
  222:       ELSE IF ((.NOT.LSAME(TRANS,'N')) .AND.
  223:      +         (.NOT.LSAME(TRANS,'C'))) THEN
  224:           INFO = 2
  225:       ELSE IF (N.LT.0) THEN
  226:           INFO = 3
  227:       ELSE IF (K.LT.0) THEN
  228:           INFO = 4
  229:       ELSE IF (LDA.LT.MAX(1,NROWA)) THEN
  230:           INFO = 7
  231:       ELSE IF (LDC.LT.MAX(1,N)) THEN
  232:           INFO = 10
  233:       END IF
  234:       IF (INFO.NE.0) THEN
  235:           CALL XERBLA('ZHERK ',INFO)
  236:           RETURN
  237:       END IF
  238: *
  239: *     Quick return if possible.
  240: *
  241:       IF ((N.EQ.0) .OR. (((ALPHA.EQ.ZERO).OR.
  242:      +    (K.EQ.0)).AND. (BETA.EQ.ONE))) RETURN
  243: *
  244: *     And when  alpha.eq.zero.
  245: *
  246:       IF (ALPHA.EQ.ZERO) THEN
  247:           IF (UPPER) THEN
  248:               IF (BETA.EQ.ZERO) THEN
  249:                   DO 20 J = 1,N
  250:                       DO 10 I = 1,J
  251:                           C(I,J) = ZERO
  252:    10                 CONTINUE
  253:    20             CONTINUE
  254:               ELSE
  255:                   DO 40 J = 1,N
  256:                       DO 30 I = 1,J - 1
  257:                           C(I,J) = BETA*C(I,J)
  258:    30                 CONTINUE
  259:                       C(J,J) = BETA*DBLE(C(J,J))
  260:    40             CONTINUE
  261:               END IF
  262:           ELSE
  263:               IF (BETA.EQ.ZERO) THEN
  264:                   DO 60 J = 1,N
  265:                       DO 50 I = J,N
  266:                           C(I,J) = ZERO
  267:    50                 CONTINUE
  268:    60             CONTINUE
  269:               ELSE
  270:                   DO 80 J = 1,N
  271:                       C(J,J) = BETA*DBLE(C(J,J))
  272:                       DO 70 I = J + 1,N
  273:                           C(I,J) = BETA*C(I,J)
  274:    70                 CONTINUE
  275:    80             CONTINUE
  276:               END IF
  277:           END IF
  278:           RETURN
  279:       END IF
  280: *
  281: *     Start the operations.
  282: *
  283:       IF (LSAME(TRANS,'N')) THEN
  284: *
  285: *        Form  C := alpha*A*A**H + beta*C.
  286: *
  287:           IF (UPPER) THEN
  288:               DO 130 J = 1,N
  289:                   IF (BETA.EQ.ZERO) THEN
  290:                       DO 90 I = 1,J
  291:                           C(I,J) = ZERO
  292:    90                 CONTINUE
  293:                   ELSE IF (BETA.NE.ONE) THEN
  294:                       DO 100 I = 1,J - 1
  295:                           C(I,J) = BETA*C(I,J)
  296:   100                 CONTINUE
  297:                       C(J,J) = BETA*DBLE(C(J,J))
  298:                   ELSE
  299:                       C(J,J) = DBLE(C(J,J))
  300:                   END IF
  301:                   DO 120 L = 1,K
  302:                       IF (A(J,L).NE.DCMPLX(ZERO)) THEN
  303:                           TEMP = ALPHA*DCONJG(A(J,L))
  304:                           DO 110 I = 1,J - 1
  305:                               C(I,J) = C(I,J) + TEMP*A(I,L)
  306:   110                     CONTINUE
  307:                           C(J,J) = DBLE(C(J,J)) + DBLE(TEMP*A(I,L))
  308:                       END IF
  309:   120             CONTINUE
  310:   130         CONTINUE
  311:           ELSE
  312:               DO 180 J = 1,N
  313:                   IF (BETA.EQ.ZERO) THEN
  314:                       DO 140 I = J,N
  315:                           C(I,J) = ZERO
  316:   140                 CONTINUE
  317:                   ELSE IF (BETA.NE.ONE) THEN
  318:                       C(J,J) = BETA*DBLE(C(J,J))
  319:                       DO 150 I = J + 1,N
  320:                           C(I,J) = BETA*C(I,J)
  321:   150                 CONTINUE
  322:                   ELSE
  323:                       C(J,J) = DBLE(C(J,J))
  324:                   END IF
  325:                   DO 170 L = 1,K
  326:                       IF (A(J,L).NE.DCMPLX(ZERO)) THEN
  327:                           TEMP = ALPHA*DCONJG(A(J,L))
  328:                           C(J,J) = DBLE(C(J,J)) + DBLE(TEMP*A(J,L))
  329:                           DO 160 I = J + 1,N
  330:                               C(I,J) = C(I,J) + TEMP*A(I,L)
  331:   160                     CONTINUE
  332:                       END IF
  333:   170             CONTINUE
  334:   180         CONTINUE
  335:           END IF
  336:       ELSE
  337: *
  338: *        Form  C := alpha*A**H*A + beta*C.
  339: *
  340:           IF (UPPER) THEN
  341:               DO 220 J = 1,N
  342:                   DO 200 I = 1,J - 1
  343:                       TEMP = ZERO
  344:                       DO 190 L = 1,K
  345:                           TEMP = TEMP + DCONJG(A(L,I))*A(L,J)
  346:   190                 CONTINUE
  347:                       IF (BETA.EQ.ZERO) THEN
  348:                           C(I,J) = ALPHA*TEMP
  349:                       ELSE
  350:                           C(I,J) = ALPHA*TEMP + BETA*C(I,J)
  351:                       END IF
  352:   200             CONTINUE
  353:                   RTEMP = ZERO
  354:                   DO 210 L = 1,K
  355:                       RTEMP = RTEMP + DBLE(DCONJG(A(L,J))*A(L,J))
  356:   210             CONTINUE
  357:                   IF (BETA.EQ.ZERO) THEN
  358:                       C(J,J) = ALPHA*RTEMP
  359:                   ELSE
  360:                       C(J,J) = ALPHA*RTEMP + BETA*DBLE(C(J,J))
  361:                   END IF
  362:   220         CONTINUE
  363:           ELSE
  364:               DO 260 J = 1,N
  365:                   RTEMP = ZERO
  366:                   DO 230 L = 1,K
  367:                       RTEMP = RTEMP + DBLE(DCONJG(A(L,J))*A(L,J))
  368:   230             CONTINUE
  369:                   IF (BETA.EQ.ZERO) THEN
  370:                       C(J,J) = ALPHA*RTEMP
  371:                   ELSE
  372:                       C(J,J) = ALPHA*RTEMP + BETA*DBLE(C(J,J))
  373:                   END IF
  374:                   DO 250 I = J + 1,N
  375:                       TEMP = ZERO
  376:                       DO 240 L = 1,K
  377:                           TEMP = TEMP + DCONJG(A(L,I))*A(L,J)
  378:   240                 CONTINUE
  379:                       IF (BETA.EQ.ZERO) THEN
  380:                           C(I,J) = ALPHA*TEMP
  381:                       ELSE
  382:                           C(I,J) = ALPHA*TEMP + BETA*C(I,J)
  383:                       END IF
  384:   250             CONTINUE
  385:   260         CONTINUE
  386:           END IF
  387:       END IF
  388: *
  389:       RETURN
  390: *
  391: *     End of ZHERK
  392: *
  393:       END

CVSweb interface <joel.bertrand@systella.fr>