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

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

CVSweb interface <joel.bertrand@systella.fr>