File:  [local] / rpl / lapack / blas / zgemm.f
Revision 1.17: 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 ZGEMM
    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 ZGEMM(TRANSA,TRANSB,M,N,K,ALPHA,A,LDA,B,LDB,BETA,C,LDC)
   12: *
   13: *       .. Scalar Arguments ..
   14: *       COMPLEX*16 ALPHA,BETA
   15: *       INTEGER K,LDA,LDB,LDC,M,N
   16: *       CHARACTER TRANSA,TRANSB
   17: *       ..
   18: *       .. Array Arguments ..
   19: *       COMPLEX*16 A(LDA,*),B(LDB,*),C(LDC,*)
   20: *       ..
   21: *
   22: *
   23: *> \par Purpose:
   24: *  =============
   25: *>
   26: *> \verbatim
   27: *>
   28: *> ZGEMM  performs one of the matrix-matrix operations
   29: *>
   30: *>    C := alpha*op( A )*op( B ) + beta*C,
   31: *>
   32: *> where  op( X ) is one of
   33: *>
   34: *>    op( X ) = X   or   op( X ) = X**T   or   op( X ) = X**H,
   35: *>
   36: *> alpha and beta are scalars, and A, B and C are matrices, with op( A )
   37: *> an m by k matrix,  op( B )  a  k by n matrix and  C an m by n matrix.
   38: *> \endverbatim
   39: *
   40: *  Arguments:
   41: *  ==========
   42: *
   43: *> \param[in] TRANSA
   44: *> \verbatim
   45: *>          TRANSA is CHARACTER*1
   46: *>           On entry, TRANSA specifies the form of op( A ) to be used in
   47: *>           the matrix multiplication as follows:
   48: *>
   49: *>              TRANSA = 'N' or 'n',  op( A ) = A.
   50: *>
   51: *>              TRANSA = 'T' or 't',  op( A ) = A**T.
   52: *>
   53: *>              TRANSA = 'C' or 'c',  op( A ) = A**H.
   54: *> \endverbatim
   55: *>
   56: *> \param[in] TRANSB
   57: *> \verbatim
   58: *>          TRANSB is CHARACTER*1
   59: *>           On entry, TRANSB specifies the form of op( B ) to be used in
   60: *>           the matrix multiplication as follows:
   61: *>
   62: *>              TRANSB = 'N' or 'n',  op( B ) = B.
   63: *>
   64: *>              TRANSB = 'T' or 't',  op( B ) = B**T.
   65: *>
   66: *>              TRANSB = 'C' or 'c',  op( B ) = B**H.
   67: *> \endverbatim
   68: *>
   69: *> \param[in] M
   70: *> \verbatim
   71: *>          M is INTEGER
   72: *>           On entry,  M  specifies  the number  of rows  of the  matrix
   73: *>           op( A )  and of the  matrix  C.  M  must  be at least  zero.
   74: *> \endverbatim
   75: *>
   76: *> \param[in] N
   77: *> \verbatim
   78: *>          N is INTEGER
   79: *>           On entry,  N  specifies the number  of columns of the matrix
   80: *>           op( B ) and the number of columns of the matrix C. N must be
   81: *>           at least zero.
   82: *> \endverbatim
   83: *>
   84: *> \param[in] K
   85: *> \verbatim
   86: *>          K is INTEGER
   87: *>           On entry,  K  specifies  the number of columns of the matrix
   88: *>           op( A ) and the number of rows of the matrix op( B ). K must
   89: *>           be at least  zero.
   90: *> \endverbatim
   91: *>
   92: *> \param[in] ALPHA
   93: *> \verbatim
   94: *>          ALPHA is COMPLEX*16
   95: *>           On entry, ALPHA specifies the scalar alpha.
   96: *> \endverbatim
   97: *>
   98: *> \param[in] A
   99: *> \verbatim
  100: *>          A is COMPLEX*16 array, dimension ( LDA, ka ), where ka is
  101: *>           k  when  TRANSA = 'N' or 'n',  and is  m  otherwise.
  102: *>           Before entry with  TRANSA = 'N' or 'n',  the leading  m by k
  103: *>           part of the array  A  must contain the matrix  A,  otherwise
  104: *>           the leading  k by m  part of the array  A  must contain  the
  105: *>           matrix A.
  106: *> \endverbatim
  107: *>
  108: *> \param[in] LDA
  109: *> \verbatim
  110: *>          LDA is INTEGER
  111: *>           On entry, LDA specifies the first dimension of A as declared
  112: *>           in the calling (sub) program. When  TRANSA = 'N' or 'n' then
  113: *>           LDA must be at least  max( 1, m ), otherwise  LDA must be at
  114: *>           least  max( 1, k ).
  115: *> \endverbatim
  116: *>
  117: *> \param[in] B
  118: *> \verbatim
  119: *>          B is COMPLEX*16 array, dimension ( LDB, kb ), where kb is
  120: *>           n  when  TRANSB = 'N' or 'n',  and is  k  otherwise.
  121: *>           Before entry with  TRANSB = 'N' or 'n',  the leading  k by n
  122: *>           part of the array  B  must contain the matrix  B,  otherwise
  123: *>           the leading  n by k  part of the array  B  must contain  the
  124: *>           matrix B.
  125: *> \endverbatim
  126: *>
  127: *> \param[in] LDB
  128: *> \verbatim
  129: *>          LDB is INTEGER
  130: *>           On entry, LDB specifies the first dimension of B as declared
  131: *>           in the calling (sub) program. When  TRANSB = 'N' or 'n' then
  132: *>           LDB must be at least  max( 1, k ), otherwise  LDB must be at
  133: *>           least  max( 1, n ).
  134: *> \endverbatim
  135: *>
  136: *> \param[in] BETA
  137: *> \verbatim
  138: *>          BETA is COMPLEX*16
  139: *>           On entry,  BETA  specifies the scalar  beta.  When  BETA  is
  140: *>           supplied as zero then C need not be set on input.
  141: *> \endverbatim
  142: *>
  143: *> \param[in,out] C
  144: *> \verbatim
  145: *>          C is COMPLEX*16 array, dimension ( LDC, N )
  146: *>           Before entry, the leading  m by n  part of the array  C must
  147: *>           contain the matrix  C,  except when  beta  is zero, in which
  148: *>           case C need not be set on entry.
  149: *>           On exit, the array  C  is overwritten by the  m by n  matrix
  150: *>           ( alpha*op( A )*op( B ) + beta*C ).
  151: *> \endverbatim
  152: *>
  153: *> \param[in] LDC
  154: *> \verbatim
  155: *>          LDC is INTEGER
  156: *>           On entry, LDC specifies the first dimension of C as declared
  157: *>           in  the  calling  (sub)  program.   LDC  must  be  at  least
  158: *>           max( 1, m ).
  159: *> \endverbatim
  160: *
  161: *  Authors:
  162: *  ========
  163: *
  164: *> \author Univ. of Tennessee
  165: *> \author Univ. of California Berkeley
  166: *> \author Univ. of Colorado Denver
  167: *> \author NAG Ltd.
  168: *
  169: *> \ingroup complex16_blas_level3
  170: *
  171: *> \par Further Details:
  172: *  =====================
  173: *>
  174: *> \verbatim
  175: *>
  176: *>  Level 3 Blas routine.
  177: *>
  178: *>  -- Written on 8-February-1989.
  179: *>     Jack Dongarra, Argonne National Laboratory.
  180: *>     Iain Duff, AERE Harwell.
  181: *>     Jeremy Du Croz, Numerical Algorithms Group Ltd.
  182: *>     Sven Hammarling, Numerical Algorithms Group Ltd.
  183: *> \endverbatim
  184: *>
  185: *  =====================================================================
  186:       SUBROUTINE ZGEMM(TRANSA,TRANSB,M,N,K,ALPHA,A,LDA,B,LDB,BETA,C,LDC)
  187: *
  188: *  -- Reference BLAS level3 routine --
  189: *  -- Reference BLAS is a software package provided by Univ. of Tennessee,    --
  190: *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
  191: *
  192: *     .. Scalar Arguments ..
  193:       COMPLEX*16 ALPHA,BETA
  194:       INTEGER K,LDA,LDB,LDC,M,N
  195:       CHARACTER TRANSA,TRANSB
  196: *     ..
  197: *     .. Array Arguments ..
  198:       COMPLEX*16 A(LDA,*),B(LDB,*),C(LDC,*)
  199: *     ..
  200: *
  201: *  =====================================================================
  202: *
  203: *     .. External Functions ..
  204:       LOGICAL LSAME
  205:       EXTERNAL LSAME
  206: *     ..
  207: *     .. External Subroutines ..
  208:       EXTERNAL XERBLA
  209: *     ..
  210: *     .. Intrinsic Functions ..
  211:       INTRINSIC DCONJG,MAX
  212: *     ..
  213: *     .. Local Scalars ..
  214:       COMPLEX*16 TEMP
  215:       INTEGER I,INFO,J,L,NROWA,NROWB
  216:       LOGICAL CONJA,CONJB,NOTA,NOTB
  217: *     ..
  218: *     .. Parameters ..
  219:       COMPLEX*16 ONE
  220:       PARAMETER (ONE= (1.0D+0,0.0D+0))
  221:       COMPLEX*16 ZERO
  222:       PARAMETER (ZERO= (0.0D+0,0.0D+0))
  223: *     ..
  224: *
  225: *     Set  NOTA  and  NOTB  as  true if  A  and  B  respectively are not
  226: *     conjugated or transposed, set  CONJA and CONJB  as true if  A  and
  227: *     B  respectively are to be  transposed but  not conjugated  and set
  228: *     NROWA and NROWB  as the number of rows  of  A  and  B  respectively.
  229: *
  230:       NOTA = LSAME(TRANSA,'N')
  231:       NOTB = LSAME(TRANSB,'N')
  232:       CONJA = LSAME(TRANSA,'C')
  233:       CONJB = LSAME(TRANSB,'C')
  234:       IF (NOTA) THEN
  235:           NROWA = M
  236:       ELSE
  237:           NROWA = K
  238:       END IF
  239:       IF (NOTB) THEN
  240:           NROWB = K
  241:       ELSE
  242:           NROWB = N
  243:       END IF
  244: *
  245: *     Test the input parameters.
  246: *
  247:       INFO = 0
  248:       IF ((.NOT.NOTA) .AND. (.NOT.CONJA) .AND.
  249:      +    (.NOT.LSAME(TRANSA,'T'))) THEN
  250:           INFO = 1
  251:       ELSE IF ((.NOT.NOTB) .AND. (.NOT.CONJB) .AND.
  252:      +         (.NOT.LSAME(TRANSB,'T'))) THEN
  253:           INFO = 2
  254:       ELSE IF (M.LT.0) THEN
  255:           INFO = 3
  256:       ELSE IF (N.LT.0) THEN
  257:           INFO = 4
  258:       ELSE IF (K.LT.0) THEN
  259:           INFO = 5
  260:       ELSE IF (LDA.LT.MAX(1,NROWA)) THEN
  261:           INFO = 8
  262:       ELSE IF (LDB.LT.MAX(1,NROWB)) THEN
  263:           INFO = 10
  264:       ELSE IF (LDC.LT.MAX(1,M)) THEN
  265:           INFO = 13
  266:       END IF
  267:       IF (INFO.NE.0) THEN
  268:           CALL XERBLA('ZGEMM ',INFO)
  269:           RETURN
  270:       END IF
  271: *
  272: *     Quick return if possible.
  273: *
  274:       IF ((M.EQ.0) .OR. (N.EQ.0) .OR.
  275:      +    (((ALPHA.EQ.ZERO).OR. (K.EQ.0)).AND. (BETA.EQ.ONE))) RETURN
  276: *
  277: *     And when  alpha.eq.zero.
  278: *
  279:       IF (ALPHA.EQ.ZERO) THEN
  280:           IF (BETA.EQ.ZERO) THEN
  281:               DO 20 J = 1,N
  282:                   DO 10 I = 1,M
  283:                       C(I,J) = ZERO
  284:    10             CONTINUE
  285:    20         CONTINUE
  286:           ELSE
  287:               DO 40 J = 1,N
  288:                   DO 30 I = 1,M
  289:                       C(I,J) = BETA*C(I,J)
  290:    30             CONTINUE
  291:    40         CONTINUE
  292:           END IF
  293:           RETURN
  294:       END IF
  295: *
  296: *     Start the operations.
  297: *
  298:       IF (NOTB) THEN
  299:           IF (NOTA) THEN
  300: *
  301: *           Form  C := alpha*A*B + beta*C.
  302: *
  303:               DO 90 J = 1,N
  304:                   IF (BETA.EQ.ZERO) THEN
  305:                       DO 50 I = 1,M
  306:                           C(I,J) = ZERO
  307:    50                 CONTINUE
  308:                   ELSE IF (BETA.NE.ONE) THEN
  309:                       DO 60 I = 1,M
  310:                           C(I,J) = BETA*C(I,J)
  311:    60                 CONTINUE
  312:                   END IF
  313:                   DO 80 L = 1,K
  314:                       TEMP = ALPHA*B(L,J)
  315:                       DO 70 I = 1,M
  316:                           C(I,J) = C(I,J) + TEMP*A(I,L)
  317:    70                 CONTINUE
  318:    80             CONTINUE
  319:    90         CONTINUE
  320:           ELSE IF (CONJA) THEN
  321: *
  322: *           Form  C := alpha*A**H*B + beta*C.
  323: *
  324:               DO 120 J = 1,N
  325:                   DO 110 I = 1,M
  326:                       TEMP = ZERO
  327:                       DO 100 L = 1,K
  328:                           TEMP = TEMP + DCONJG(A(L,I))*B(L,J)
  329:   100                 CONTINUE
  330:                       IF (BETA.EQ.ZERO) THEN
  331:                           C(I,J) = ALPHA*TEMP
  332:                       ELSE
  333:                           C(I,J) = ALPHA*TEMP + BETA*C(I,J)
  334:                       END IF
  335:   110             CONTINUE
  336:   120         CONTINUE
  337:           ELSE
  338: *
  339: *           Form  C := alpha*A**T*B + beta*C
  340: *
  341:               DO 150 J = 1,N
  342:                   DO 140 I = 1,M
  343:                       TEMP = ZERO
  344:                       DO 130 L = 1,K
  345:                           TEMP = TEMP + A(L,I)*B(L,J)
  346:   130                 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:   140             CONTINUE
  353:   150         CONTINUE
  354:           END IF
  355:       ELSE IF (NOTA) THEN
  356:           IF (CONJB) THEN
  357: *
  358: *           Form  C := alpha*A*B**H + beta*C.
  359: *
  360:               DO 200 J = 1,N
  361:                   IF (BETA.EQ.ZERO) THEN
  362:                       DO 160 I = 1,M
  363:                           C(I,J) = ZERO
  364:   160                 CONTINUE
  365:                   ELSE IF (BETA.NE.ONE) THEN
  366:                       DO 170 I = 1,M
  367:                           C(I,J) = BETA*C(I,J)
  368:   170                 CONTINUE
  369:                   END IF
  370:                   DO 190 L = 1,K
  371:                       TEMP = ALPHA*DCONJG(B(J,L))
  372:                       DO 180 I = 1,M
  373:                           C(I,J) = C(I,J) + TEMP*A(I,L)
  374:   180                 CONTINUE
  375:   190             CONTINUE
  376:   200         CONTINUE
  377:           ELSE
  378: *
  379: *           Form  C := alpha*A*B**T + beta*C
  380: *
  381:               DO 250 J = 1,N
  382:                   IF (BETA.EQ.ZERO) THEN
  383:                       DO 210 I = 1,M
  384:                           C(I,J) = ZERO
  385:   210                 CONTINUE
  386:                   ELSE IF (BETA.NE.ONE) THEN
  387:                       DO 220 I = 1,M
  388:                           C(I,J) = BETA*C(I,J)
  389:   220                 CONTINUE
  390:                   END IF
  391:                   DO 240 L = 1,K
  392:                       TEMP = ALPHA*B(J,L)
  393:                       DO 230 I = 1,M
  394:                           C(I,J) = C(I,J) + TEMP*A(I,L)
  395:   230                 CONTINUE
  396:   240             CONTINUE
  397:   250         CONTINUE
  398:           END IF
  399:       ELSE IF (CONJA) THEN
  400:           IF (CONJB) THEN
  401: *
  402: *           Form  C := alpha*A**H*B**H + beta*C.
  403: *
  404:               DO 280 J = 1,N
  405:                   DO 270 I = 1,M
  406:                       TEMP = ZERO
  407:                       DO 260 L = 1,K
  408:                           TEMP = TEMP + DCONJG(A(L,I))*DCONJG(B(J,L))
  409:   260                 CONTINUE
  410:                       IF (BETA.EQ.ZERO) THEN
  411:                           C(I,J) = ALPHA*TEMP
  412:                       ELSE
  413:                           C(I,J) = ALPHA*TEMP + BETA*C(I,J)
  414:                       END IF
  415:   270             CONTINUE
  416:   280         CONTINUE
  417:           ELSE
  418: *
  419: *           Form  C := alpha*A**H*B**T + beta*C
  420: *
  421:               DO 310 J = 1,N
  422:                   DO 300 I = 1,M
  423:                       TEMP = ZERO
  424:                       DO 290 L = 1,K
  425:                           TEMP = TEMP + DCONJG(A(L,I))*B(J,L)
  426:   290                 CONTINUE
  427:                       IF (BETA.EQ.ZERO) THEN
  428:                           C(I,J) = ALPHA*TEMP
  429:                       ELSE
  430:                           C(I,J) = ALPHA*TEMP + BETA*C(I,J)
  431:                       END IF
  432:   300             CONTINUE
  433:   310         CONTINUE
  434:           END IF
  435:       ELSE
  436:           IF (CONJB) THEN
  437: *
  438: *           Form  C := alpha*A**T*B**H + beta*C
  439: *
  440:               DO 340 J = 1,N
  441:                   DO 330 I = 1,M
  442:                       TEMP = ZERO
  443:                       DO 320 L = 1,K
  444:                           TEMP = TEMP + A(L,I)*DCONJG(B(J,L))
  445:   320                 CONTINUE
  446:                       IF (BETA.EQ.ZERO) THEN
  447:                           C(I,J) = ALPHA*TEMP
  448:                       ELSE
  449:                           C(I,J) = ALPHA*TEMP + BETA*C(I,J)
  450:                       END IF
  451:   330             CONTINUE
  452:   340         CONTINUE
  453:           ELSE
  454: *
  455: *           Form  C := alpha*A**T*B**T + beta*C
  456: *
  457:               DO 370 J = 1,N
  458:                   DO 360 I = 1,M
  459:                       TEMP = ZERO
  460:                       DO 350 L = 1,K
  461:                           TEMP = TEMP + A(L,I)*B(J,L)
  462:   350                 CONTINUE
  463:                       IF (BETA.EQ.ZERO) THEN
  464:                           C(I,J) = ALPHA*TEMP
  465:                       ELSE
  466:                           C(I,J) = ALPHA*TEMP + BETA*C(I,J)
  467:                       END IF
  468:   360             CONTINUE
  469:   370         CONTINUE
  470:           END IF
  471:       END IF
  472: *
  473:       RETURN
  474: *
  475: *     End of ZGEMM
  476: *
  477:       END

CVSweb interface <joel.bertrand@systella.fr>