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

1.8       bertrand    1: *> \brief \b DTRSM
                      2: *
                      3: *  =========== DOCUMENTATION ===========
                      4: *
1.13      bertrand    5: * Online html documentation available at
                      6: *            http://www.netlib.org/lapack/explore-html/
1.8       bertrand    7: *
                      8: *  Definition:
                      9: *  ===========
                     10: *
                     11: *       SUBROUTINE DTRSM(SIDE,UPLO,TRANSA,DIAG,M,N,ALPHA,A,LDA,B,LDB)
1.13      bertrand   12: *
1.8       bertrand   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: *       ..
1.13      bertrand   21: *
1.8       bertrand   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
1.14      bertrand  114: *>          A is DOUBLE PRECISION array, dimension ( LDA, k ),
1.13      bertrand  115: *>           where k is m when SIDE = 'L' or 'l'
1.8       bertrand  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
1.14      bertrand  140: *>          B is DOUBLE PRECISION array, dimension ( LDB, N )
1.8       bertrand  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: *
1.13      bertrand  157: *> \author Univ. of Tennessee
                    158: *> \author Univ. of California Berkeley
                    159: *> \author Univ. of Colorado Denver
                    160: *> \author NAG Ltd.
1.8       bertrand  161: *
                    162: *> \ingroup double_blas_level3
                    163: *
                    164: *> \par Further Details:
                    165: *  =====================
                    166: *>
                    167: *> \verbatim
                    168: *>
                    169: *>  Level 3 Blas routine.
                    170: *>
                    171: *>
                    172: *>  -- Written on 8-February-1989.
                    173: *>     Jack Dongarra, Argonne National Laboratory.
                    174: *>     Iain Duff, AERE Harwell.
                    175: *>     Jeremy Du Croz, Numerical Algorithms Group Ltd.
                    176: *>     Sven Hammarling, Numerical Algorithms Group Ltd.
                    177: *> \endverbatim
                    178: *>
                    179: *  =====================================================================
1.1       bertrand  180:       SUBROUTINE DTRSM(SIDE,UPLO,TRANSA,DIAG,M,N,ALPHA,A,LDA,B,LDB)
1.8       bertrand  181: *
1.16    ! bertrand  182: *  -- Reference BLAS level3 routine --
1.8       bertrand  183: *  -- Reference BLAS is a software package provided by Univ. of Tennessee,    --
                    184: *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
                    185: *
1.1       bertrand  186: *     .. Scalar Arguments ..
                    187:       DOUBLE PRECISION ALPHA
                    188:       INTEGER LDA,LDB,M,N
                    189:       CHARACTER DIAG,SIDE,TRANSA,UPLO
                    190: *     ..
                    191: *     .. Array Arguments ..
                    192:       DOUBLE PRECISION A(LDA,*),B(LDB,*)
                    193: *     ..
                    194: *
                    195: *  =====================================================================
                    196: *
                    197: *     .. External Functions ..
                    198:       LOGICAL LSAME
                    199:       EXTERNAL LSAME
                    200: *     ..
                    201: *     .. External Subroutines ..
                    202:       EXTERNAL XERBLA
                    203: *     ..
                    204: *     .. Intrinsic Functions ..
                    205:       INTRINSIC MAX
                    206: *     ..
                    207: *     .. Local Scalars ..
                    208:       DOUBLE PRECISION TEMP
                    209:       INTEGER I,INFO,J,K,NROWA
                    210:       LOGICAL LSIDE,NOUNIT,UPPER
                    211: *     ..
                    212: *     .. Parameters ..
                    213:       DOUBLE PRECISION ONE,ZERO
                    214:       PARAMETER (ONE=1.0D+0,ZERO=0.0D+0)
                    215: *     ..
                    216: *
                    217: *     Test the input parameters.
                    218: *
                    219:       LSIDE = LSAME(SIDE,'L')
                    220:       IF (LSIDE) THEN
                    221:           NROWA = M
                    222:       ELSE
                    223:           NROWA = N
                    224:       END IF
                    225:       NOUNIT = LSAME(DIAG,'N')
                    226:       UPPER = LSAME(UPLO,'U')
                    227: *
                    228:       INFO = 0
                    229:       IF ((.NOT.LSIDE) .AND. (.NOT.LSAME(SIDE,'R'))) THEN
                    230:           INFO = 1
                    231:       ELSE IF ((.NOT.UPPER) .AND. (.NOT.LSAME(UPLO,'L'))) THEN
                    232:           INFO = 2
                    233:       ELSE IF ((.NOT.LSAME(TRANSA,'N')) .AND.
                    234:      +         (.NOT.LSAME(TRANSA,'T')) .AND.
                    235:      +         (.NOT.LSAME(TRANSA,'C'))) THEN
                    236:           INFO = 3
                    237:       ELSE IF ((.NOT.LSAME(DIAG,'U')) .AND. (.NOT.LSAME(DIAG,'N'))) THEN
                    238:           INFO = 4
                    239:       ELSE IF (M.LT.0) THEN
                    240:           INFO = 5
                    241:       ELSE IF (N.LT.0) THEN
                    242:           INFO = 6
                    243:       ELSE IF (LDA.LT.MAX(1,NROWA)) THEN
                    244:           INFO = 9
                    245:       ELSE IF (LDB.LT.MAX(1,M)) THEN
                    246:           INFO = 11
                    247:       END IF
                    248:       IF (INFO.NE.0) THEN
                    249:           CALL XERBLA('DTRSM ',INFO)
                    250:           RETURN
                    251:       END IF
                    252: *
                    253: *     Quick return if possible.
                    254: *
                    255:       IF (M.EQ.0 .OR. N.EQ.0) RETURN
                    256: *
                    257: *     And when  alpha.eq.zero.
                    258: *
                    259:       IF (ALPHA.EQ.ZERO) THEN
                    260:           DO 20 J = 1,N
                    261:               DO 10 I = 1,M
                    262:                   B(I,J) = ZERO
                    263:    10         CONTINUE
                    264:    20     CONTINUE
                    265:           RETURN
                    266:       END IF
                    267: *
                    268: *     Start the operations.
                    269: *
                    270:       IF (LSIDE) THEN
                    271:           IF (LSAME(TRANSA,'N')) THEN
                    272: *
                    273: *           Form  B := alpha*inv( A )*B.
                    274: *
                    275:               IF (UPPER) THEN
                    276:                   DO 60 J = 1,N
                    277:                       IF (ALPHA.NE.ONE) THEN
                    278:                           DO 30 I = 1,M
                    279:                               B(I,J) = ALPHA*B(I,J)
                    280:    30                     CONTINUE
                    281:                       END IF
                    282:                       DO 50 K = M,1,-1
                    283:                           IF (B(K,J).NE.ZERO) THEN
                    284:                               IF (NOUNIT) B(K,J) = B(K,J)/A(K,K)
                    285:                               DO 40 I = 1,K - 1
                    286:                                   B(I,J) = B(I,J) - B(K,J)*A(I,K)
                    287:    40                         CONTINUE
                    288:                           END IF
                    289:    50                 CONTINUE
                    290:    60             CONTINUE
                    291:               ELSE
                    292:                   DO 100 J = 1,N
                    293:                       IF (ALPHA.NE.ONE) THEN
                    294:                           DO 70 I = 1,M
                    295:                               B(I,J) = ALPHA*B(I,J)
                    296:    70                     CONTINUE
                    297:                       END IF
                    298:                       DO 90 K = 1,M
                    299:                           IF (B(K,J).NE.ZERO) THEN
                    300:                               IF (NOUNIT) B(K,J) = B(K,J)/A(K,K)
                    301:                               DO 80 I = K + 1,M
                    302:                                   B(I,J) = B(I,J) - B(K,J)*A(I,K)
                    303:    80                         CONTINUE
                    304:                           END IF
                    305:    90                 CONTINUE
                    306:   100             CONTINUE
                    307:               END IF
                    308:           ELSE
                    309: *
1.7       bertrand  310: *           Form  B := alpha*inv( A**T )*B.
1.1       bertrand  311: *
                    312:               IF (UPPER) THEN
                    313:                   DO 130 J = 1,N
                    314:                       DO 120 I = 1,M
                    315:                           TEMP = ALPHA*B(I,J)
                    316:                           DO 110 K = 1,I - 1
                    317:                               TEMP = TEMP - A(K,I)*B(K,J)
                    318:   110                     CONTINUE
                    319:                           IF (NOUNIT) TEMP = TEMP/A(I,I)
                    320:                           B(I,J) = TEMP
                    321:   120                 CONTINUE
                    322:   130             CONTINUE
                    323:               ELSE
                    324:                   DO 160 J = 1,N
                    325:                       DO 150 I = M,1,-1
                    326:                           TEMP = ALPHA*B(I,J)
                    327:                           DO 140 K = I + 1,M
                    328:                               TEMP = TEMP - A(K,I)*B(K,J)
                    329:   140                     CONTINUE
                    330:                           IF (NOUNIT) TEMP = TEMP/A(I,I)
                    331:                           B(I,J) = TEMP
                    332:   150                 CONTINUE
                    333:   160             CONTINUE
                    334:               END IF
                    335:           END IF
                    336:       ELSE
                    337:           IF (LSAME(TRANSA,'N')) THEN
                    338: *
                    339: *           Form  B := alpha*B*inv( A ).
                    340: *
                    341:               IF (UPPER) THEN
                    342:                   DO 210 J = 1,N
                    343:                       IF (ALPHA.NE.ONE) THEN
                    344:                           DO 170 I = 1,M
                    345:                               B(I,J) = ALPHA*B(I,J)
                    346:   170                     CONTINUE
                    347:                       END IF
                    348:                       DO 190 K = 1,J - 1
                    349:                           IF (A(K,J).NE.ZERO) THEN
                    350:                               DO 180 I = 1,M
                    351:                                   B(I,J) = B(I,J) - A(K,J)*B(I,K)
                    352:   180                         CONTINUE
                    353:                           END IF
                    354:   190                 CONTINUE
                    355:                       IF (NOUNIT) THEN
                    356:                           TEMP = ONE/A(J,J)
                    357:                           DO 200 I = 1,M
                    358:                               B(I,J) = TEMP*B(I,J)
                    359:   200                     CONTINUE
                    360:                       END IF
                    361:   210             CONTINUE
                    362:               ELSE
                    363:                   DO 260 J = N,1,-1
                    364:                       IF (ALPHA.NE.ONE) THEN
                    365:                           DO 220 I = 1,M
                    366:                               B(I,J) = ALPHA*B(I,J)
                    367:   220                     CONTINUE
                    368:                       END IF
                    369:                       DO 240 K = J + 1,N
                    370:                           IF (A(K,J).NE.ZERO) THEN
                    371:                               DO 230 I = 1,M
                    372:                                   B(I,J) = B(I,J) - A(K,J)*B(I,K)
                    373:   230                         CONTINUE
                    374:                           END IF
                    375:   240                 CONTINUE
                    376:                       IF (NOUNIT) THEN
                    377:                           TEMP = ONE/A(J,J)
                    378:                           DO 250 I = 1,M
                    379:                               B(I,J) = TEMP*B(I,J)
                    380:   250                     CONTINUE
                    381:                       END IF
                    382:   260             CONTINUE
                    383:               END IF
                    384:           ELSE
                    385: *
1.7       bertrand  386: *           Form  B := alpha*B*inv( A**T ).
1.1       bertrand  387: *
                    388:               IF (UPPER) THEN
                    389:                   DO 310 K = N,1,-1
                    390:                       IF (NOUNIT) THEN
                    391:                           TEMP = ONE/A(K,K)
                    392:                           DO 270 I = 1,M
                    393:                               B(I,K) = TEMP*B(I,K)
                    394:   270                     CONTINUE
                    395:                       END IF
                    396:                       DO 290 J = 1,K - 1
                    397:                           IF (A(J,K).NE.ZERO) THEN
                    398:                               TEMP = A(J,K)
                    399:                               DO 280 I = 1,M
                    400:                                   B(I,J) = B(I,J) - TEMP*B(I,K)
                    401:   280                         CONTINUE
                    402:                           END IF
                    403:   290                 CONTINUE
                    404:                       IF (ALPHA.NE.ONE) THEN
                    405:                           DO 300 I = 1,M
                    406:                               B(I,K) = ALPHA*B(I,K)
                    407:   300                     CONTINUE
                    408:                       END IF
                    409:   310             CONTINUE
                    410:               ELSE
                    411:                   DO 360 K = 1,N
                    412:                       IF (NOUNIT) THEN
                    413:                           TEMP = ONE/A(K,K)
                    414:                           DO 320 I = 1,M
                    415:                               B(I,K) = TEMP*B(I,K)
                    416:   320                     CONTINUE
                    417:                       END IF
                    418:                       DO 340 J = K + 1,N
                    419:                           IF (A(J,K).NE.ZERO) THEN
                    420:                               TEMP = A(J,K)
                    421:                               DO 330 I = 1,M
                    422:                                   B(I,J) = B(I,J) - TEMP*B(I,K)
                    423:   330                         CONTINUE
                    424:                           END IF
                    425:   340                 CONTINUE
                    426:                       IF (ALPHA.NE.ONE) THEN
                    427:                           DO 350 I = 1,M
                    428:                               B(I,K) = ALPHA*B(I,K)
                    429:   350                     CONTINUE
                    430:                       END IF
                    431:   360             CONTINUE
                    432:               END IF
                    433:           END IF
                    434:       END IF
                    435: *
                    436:       RETURN
                    437: *
1.16    ! bertrand  438: *     End of DTRSM
1.1       bertrand  439: *
                    440:       END

CVSweb interface <joel.bertrand@systella.fr>