File:  [local] / rpl / lapack / blas / dtrsm.f
Revision 1.16: download - view: text, annotated - select for diffs - revision graph
Mon Aug 7 08:38:44 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 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, 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, 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: *> \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: *  =====================================================================
  180:       SUBROUTINE DTRSM(SIDE,UPLO,TRANSA,DIAG,M,N,ALPHA,A,LDA,B,LDB)
  181: *
  182: *  -- Reference BLAS level3 routine --
  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: *
  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: *
  310: *           Form  B := alpha*inv( A**T )*B.
  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: *
  386: *           Form  B := alpha*B*inv( A**T ).
  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: *
  438: *     End of DTRSM
  439: *
  440:       END

CVSweb interface <joel.bertrand@systella.fr>