Annotation of rpl/lapack/blas/dgemm.f, revision 1.8

1.8     ! bertrand    1: *> \brief \b DGEMM
        !             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 DGEMM(TRANSA,TRANSB,M,N,K,ALPHA,A,LDA,B,LDB,BETA,C,LDC)
        !            12: * 
        !            13: *       .. Scalar Arguments ..
        !            14: *       DOUBLE PRECISION ALPHA,BETA
        !            15: *       INTEGER K,LDA,LDB,LDC,M,N
        !            16: *       CHARACTER TRANSA,TRANSB
        !            17: *       ..
        !            18: *       .. Array Arguments ..
        !            19: *       DOUBLE PRECISION A(LDA,*),B(LDB,*),C(LDC,*)
        !            20: *       ..
        !            21: *  
        !            22: *
        !            23: *> \par Purpose:
        !            24: *  =============
        !            25: *>
        !            26: *> \verbatim
        !            27: *>
        !            28: *> DGEMM  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,
        !            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**T.
        !            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**T.
        !            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 DOUBLE PRECISION.
        !            95: *>           On entry, ALPHA specifies the scalar alpha.
        !            96: *> \endverbatim
        !            97: *>
        !            98: *> \param[in] A
        !            99: *> \verbatim
        !           100: *>          A is DOUBLE PRECISION array of 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 DOUBLE PRECISION array of 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 DOUBLE PRECISION.
        !           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 DOUBLE PRECISION array of 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: *> \date November 2011
        !           170: *
        !           171: *> \ingroup double_blas_level3
        !           172: *
        !           173: *> \par Further Details:
        !           174: *  =====================
        !           175: *>
        !           176: *> \verbatim
        !           177: *>
        !           178: *>  Level 3 Blas routine.
        !           179: *>
        !           180: *>  -- Written on 8-February-1989.
        !           181: *>     Jack Dongarra, Argonne National Laboratory.
        !           182: *>     Iain Duff, AERE Harwell.
        !           183: *>     Jeremy Du Croz, Numerical Algorithms Group Ltd.
        !           184: *>     Sven Hammarling, Numerical Algorithms Group Ltd.
        !           185: *> \endverbatim
        !           186: *>
        !           187: *  =====================================================================
1.1       bertrand  188:       SUBROUTINE DGEMM(TRANSA,TRANSB,M,N,K,ALPHA,A,LDA,B,LDB,BETA,C,LDC)
1.8     ! bertrand  189: *
        !           190: *  -- Reference BLAS level3 routine (version 3.4.0) --
        !           191: *  -- Reference BLAS is a software package provided by Univ. of Tennessee,    --
        !           192: *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
        !           193: *     November 2011
        !           194: *
1.1       bertrand  195: *     .. Scalar Arguments ..
                    196:       DOUBLE PRECISION ALPHA,BETA
                    197:       INTEGER K,LDA,LDB,LDC,M,N
                    198:       CHARACTER TRANSA,TRANSB
                    199: *     ..
                    200: *     .. Array Arguments ..
                    201:       DOUBLE PRECISION A(LDA,*),B(LDB,*),C(LDC,*)
                    202: *     ..
                    203: *
                    204: *  =====================================================================
                    205: *
                    206: *     .. External Functions ..
                    207:       LOGICAL LSAME
                    208:       EXTERNAL LSAME
                    209: *     ..
                    210: *     .. External Subroutines ..
                    211:       EXTERNAL XERBLA
                    212: *     ..
                    213: *     .. Intrinsic Functions ..
                    214:       INTRINSIC MAX
                    215: *     ..
                    216: *     .. Local Scalars ..
                    217:       DOUBLE PRECISION TEMP
                    218:       INTEGER I,INFO,J,L,NCOLA,NROWA,NROWB
                    219:       LOGICAL NOTA,NOTB
                    220: *     ..
                    221: *     .. Parameters ..
                    222:       DOUBLE PRECISION ONE,ZERO
                    223:       PARAMETER (ONE=1.0D+0,ZERO=0.0D+0)
                    224: *     ..
                    225: *
                    226: *     Set  NOTA  and  NOTB  as  true if  A  and  B  respectively are not
                    227: *     transposed and set  NROWA, NCOLA and  NROWB  as the number of rows
                    228: *     and  columns of  A  and the  number of  rows  of  B  respectively.
                    229: *
                    230:       NOTA = LSAME(TRANSA,'N')
                    231:       NOTB = LSAME(TRANSB,'N')
                    232:       IF (NOTA) THEN
                    233:           NROWA = M
                    234:           NCOLA = K
                    235:       ELSE
                    236:           NROWA = K
                    237:           NCOLA = M
                    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.LSAME(TRANSA,'C')) .AND.
                    249:      +    (.NOT.LSAME(TRANSA,'T'))) THEN
                    250:           INFO = 1
                    251:       ELSE IF ((.NOT.NOTB) .AND. (.NOT.LSAME(TRANSB,'C')) .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('DGEMM ',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 if  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:                       IF (B(L,J).NE.ZERO) THEN
                    315:                           TEMP = ALPHA*B(L,J)
                    316:                           DO 70 I = 1,M
                    317:                               C(I,J) = C(I,J) + TEMP*A(I,L)
                    318:    70                     CONTINUE
                    319:                       END IF
                    320:    80             CONTINUE
                    321:    90         CONTINUE
                    322:           ELSE
                    323: *
1.7       bertrand  324: *           Form  C := alpha*A**T*B + beta*C
1.1       bertrand  325: *
                    326:               DO 120 J = 1,N
                    327:                   DO 110 I = 1,M
                    328:                       TEMP = ZERO
                    329:                       DO 100 L = 1,K
                    330:                           TEMP = TEMP + A(L,I)*B(L,J)
                    331:   100                 CONTINUE
                    332:                       IF (BETA.EQ.ZERO) THEN
                    333:                           C(I,J) = ALPHA*TEMP
                    334:                       ELSE
                    335:                           C(I,J) = ALPHA*TEMP + BETA*C(I,J)
                    336:                       END IF
                    337:   110             CONTINUE
                    338:   120         CONTINUE
                    339:           END IF
                    340:       ELSE
                    341:           IF (NOTA) THEN
                    342: *
1.7       bertrand  343: *           Form  C := alpha*A*B**T + beta*C
1.1       bertrand  344: *
                    345:               DO 170 J = 1,N
                    346:                   IF (BETA.EQ.ZERO) THEN
                    347:                       DO 130 I = 1,M
                    348:                           C(I,J) = ZERO
                    349:   130                 CONTINUE
                    350:                   ELSE IF (BETA.NE.ONE) THEN
                    351:                       DO 140 I = 1,M
                    352:                           C(I,J) = BETA*C(I,J)
                    353:   140                 CONTINUE
                    354:                   END IF
                    355:                   DO 160 L = 1,K
                    356:                       IF (B(J,L).NE.ZERO) THEN
                    357:                           TEMP = ALPHA*B(J,L)
                    358:                           DO 150 I = 1,M
                    359:                               C(I,J) = C(I,J) + TEMP*A(I,L)
                    360:   150                     CONTINUE
                    361:                       END IF
                    362:   160             CONTINUE
                    363:   170         CONTINUE
                    364:           ELSE
                    365: *
1.7       bertrand  366: *           Form  C := alpha*A**T*B**T + beta*C
1.1       bertrand  367: *
                    368:               DO 200 J = 1,N
                    369:                   DO 190 I = 1,M
                    370:                       TEMP = ZERO
                    371:                       DO 180 L = 1,K
                    372:                           TEMP = TEMP + A(L,I)*B(J,L)
                    373:   180                 CONTINUE
                    374:                       IF (BETA.EQ.ZERO) THEN
                    375:                           C(I,J) = ALPHA*TEMP
                    376:                       ELSE
                    377:                           C(I,J) = ALPHA*TEMP + BETA*C(I,J)
                    378:                       END IF
                    379:   190             CONTINUE
                    380:   200         CONTINUE
                    381:           END IF
                    382:       END IF
                    383: *
                    384:       RETURN
                    385: *
                    386: *     End of DGEMM .
                    387: *
                    388:       END

CVSweb interface <joel.bertrand@systella.fr>