File:  [local] / rpl / lapack / blas / dtrsv.f
Revision 1.11: download - view: text, annotated - select for diffs - revision graph
Mon Jan 27 09:28:13 2014 UTC (10 years, 3 months ago) by bertrand
Branches: MAIN
CVS tags: rpl-4_1_24, rpl-4_1_23, rpl-4_1_22, rpl-4_1_21, rpl-4_1_20, rpl-4_1_19, rpl-4_1_18, rpl-4_1_17, HEAD
Cohérence.

    1: *> \brief \b DTRSV
    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 DTRSV(UPLO,TRANS,DIAG,N,A,LDA,X,INCX)
   12:    13: *       .. Scalar Arguments ..
   14: *       INTEGER INCX,LDA,N
   15: *       CHARACTER DIAG,TRANS,UPLO
   16: *       ..
   17: *       .. Array Arguments ..
   18: *       DOUBLE PRECISION A(LDA,*),X(*)
   19: *       ..
   20: *  
   21: *
   22: *> \par Purpose:
   23: *  =============
   24: *>
   25: *> \verbatim
   26: *>
   27: *> DTRSV  solves one of the systems of equations
   28: *>
   29: *>    A*x = b,   or   A**T*x = b,
   30: *>
   31: *> where b and x are n element vectors and A is an n by n unit, or
   32: *> non-unit, upper or lower triangular matrix.
   33: *>
   34: *> No test for singularity or near-singularity is included in this
   35: *> routine. Such tests must be performed before calling this routine.
   36: *> \endverbatim
   37: *
   38: *  Arguments:
   39: *  ==========
   40: *
   41: *> \param[in] UPLO
   42: *> \verbatim
   43: *>          UPLO is CHARACTER*1
   44: *>           On entry, UPLO specifies whether the matrix is an upper or
   45: *>           lower triangular matrix as follows:
   46: *>
   47: *>              UPLO = 'U' or 'u'   A is an upper triangular matrix.
   48: *>
   49: *>              UPLO = 'L' or 'l'   A is a lower triangular matrix.
   50: *> \endverbatim
   51: *>
   52: *> \param[in] TRANS
   53: *> \verbatim
   54: *>          TRANS is CHARACTER*1
   55: *>           On entry, TRANS specifies the equations to be solved as
   56: *>           follows:
   57: *>
   58: *>              TRANS = 'N' or 'n'   A*x = b.
   59: *>
   60: *>              TRANS = 'T' or 't'   A**T*x = b.
   61: *>
   62: *>              TRANS = 'C' or 'c'   A**T*x = b.
   63: *> \endverbatim
   64: *>
   65: *> \param[in] DIAG
   66: *> \verbatim
   67: *>          DIAG is CHARACTER*1
   68: *>           On entry, DIAG specifies whether or not A is unit
   69: *>           triangular as follows:
   70: *>
   71: *>              DIAG = 'U' or 'u'   A is assumed to be unit triangular.
   72: *>
   73: *>              DIAG = 'N' or 'n'   A is not assumed to be unit
   74: *>                                  triangular.
   75: *> \endverbatim
   76: *>
   77: *> \param[in] N
   78: *> \verbatim
   79: *>          N is INTEGER
   80: *>           On entry, N specifies the order of the matrix A.
   81: *>           N must be at least zero.
   82: *> \endverbatim
   83: *>
   84: *> \param[in] A
   85: *> \verbatim
   86: *>          A is DOUBLE PRECISION array of DIMENSION ( LDA, n ).
   87: *>           Before entry with  UPLO = 'U' or 'u', the leading n by n
   88: *>           upper triangular part of the array A must contain the upper
   89: *>           triangular matrix and the strictly lower triangular part of
   90: *>           A is not referenced.
   91: *>           Before entry with UPLO = 'L' or 'l', the leading n by n
   92: *>           lower triangular part of the array A must contain the lower
   93: *>           triangular matrix and the strictly upper triangular part of
   94: *>           A is not referenced.
   95: *>           Note that when  DIAG = 'U' or 'u', the diagonal elements of
   96: *>           A are not referenced either, but are assumed to be unity.
   97: *> \endverbatim
   98: *>
   99: *> \param[in] LDA
  100: *> \verbatim
  101: *>          LDA is INTEGER
  102: *>           On entry, LDA specifies the first dimension of A as declared
  103: *>           in the calling (sub) program. LDA must be at least
  104: *>           max( 1, n ).
  105: *> \endverbatim
  106: *>
  107: *> \param[in,out] X
  108: *> \verbatim
  109: *>          X is DOUBLE PRECISION array of dimension at least
  110: *>           ( 1 + ( n - 1 )*abs( INCX ) ).
  111: *>           Before entry, the incremented array X must contain the n
  112: *>           element right-hand side vector b. On exit, X is overwritten
  113: *>           with the solution vector x.
  114: *> \endverbatim
  115: *>
  116: *> \param[in] INCX
  117: *> \verbatim
  118: *>          INCX is INTEGER
  119: *>           On entry, INCX specifies the increment for the elements of
  120: *>           X. INCX must not be zero.
  121: *>
  122: *>  Level 2 Blas routine.
  123: *>
  124: *>  -- Written on 22-October-1986.
  125: *>     Jack Dongarra, Argonne National Lab.
  126: *>     Jeremy Du Croz, Nag Central Office.
  127: *>     Sven Hammarling, Nag Central Office.
  128: *>     Richard Hanson, Sandia National Labs.
  129: *> \endverbatim
  130: *
  131: *  Authors:
  132: *  ========
  133: *
  134: *> \author Univ. of Tennessee 
  135: *> \author Univ. of California Berkeley 
  136: *> \author Univ. of Colorado Denver 
  137: *> \author NAG Ltd. 
  138: *
  139: *> \date November 2011
  140: *
  141: *> \ingroup double_blas_level1
  142: *
  143: *  =====================================================================
  144:       SUBROUTINE DTRSV(UPLO,TRANS,DIAG,N,A,LDA,X,INCX)
  145: *
  146: *  -- Reference BLAS level1 routine (version 3.4.0) --
  147: *  -- Reference BLAS is a software package provided by Univ. of Tennessee,    --
  148: *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
  149: *     November 2011
  150: *
  151: *     .. Scalar Arguments ..
  152:       INTEGER INCX,LDA,N
  153:       CHARACTER DIAG,TRANS,UPLO
  154: *     ..
  155: *     .. Array Arguments ..
  156:       DOUBLE PRECISION A(LDA,*),X(*)
  157: *     ..
  158: *
  159: *  =====================================================================
  160: *
  161: *     .. Parameters ..
  162:       DOUBLE PRECISION ZERO
  163:       PARAMETER (ZERO=0.0D+0)
  164: *     ..
  165: *     .. Local Scalars ..
  166:       DOUBLE PRECISION TEMP
  167:       INTEGER I,INFO,IX,J,JX,KX
  168:       LOGICAL NOUNIT
  169: *     ..
  170: *     .. External Functions ..
  171:       LOGICAL LSAME
  172:       EXTERNAL LSAME
  173: *     ..
  174: *     .. External Subroutines ..
  175:       EXTERNAL XERBLA
  176: *     ..
  177: *     .. Intrinsic Functions ..
  178:       INTRINSIC MAX
  179: *     ..
  180: *
  181: *     Test the input parameters.
  182: *
  183:       INFO = 0
  184:       IF (.NOT.LSAME(UPLO,'U') .AND. .NOT.LSAME(UPLO,'L')) THEN
  185:           INFO = 1
  186:       ELSE IF (.NOT.LSAME(TRANS,'N') .AND. .NOT.LSAME(TRANS,'T') .AND.
  187:      +         .NOT.LSAME(TRANS,'C')) THEN
  188:           INFO = 2
  189:       ELSE IF (.NOT.LSAME(DIAG,'U') .AND. .NOT.LSAME(DIAG,'N')) THEN
  190:           INFO = 3
  191:       ELSE IF (N.LT.0) THEN
  192:           INFO = 4
  193:       ELSE IF (LDA.LT.MAX(1,N)) THEN
  194:           INFO = 6
  195:       ELSE IF (INCX.EQ.0) THEN
  196:           INFO = 8
  197:       END IF
  198:       IF (INFO.NE.0) THEN
  199:           CALL XERBLA('DTRSV ',INFO)
  200:           RETURN
  201:       END IF
  202: *
  203: *     Quick return if possible.
  204: *
  205:       IF (N.EQ.0) RETURN
  206: *
  207:       NOUNIT = LSAME(DIAG,'N')
  208: *
  209: *     Set up the start point in X if the increment is not unity. This
  210: *     will be  ( N - 1 )*INCX  too small for descending loops.
  211: *
  212:       IF (INCX.LE.0) THEN
  213:           KX = 1 - (N-1)*INCX
  214:       ELSE IF (INCX.NE.1) THEN
  215:           KX = 1
  216:       END IF
  217: *
  218: *     Start the operations. In this version the elements of A are
  219: *     accessed sequentially with one pass through A.
  220: *
  221:       IF (LSAME(TRANS,'N')) THEN
  222: *
  223: *        Form  x := inv( A )*x.
  224: *
  225:           IF (LSAME(UPLO,'U')) THEN
  226:               IF (INCX.EQ.1) THEN
  227:                   DO 20 J = N,1,-1
  228:                       IF (X(J).NE.ZERO) THEN
  229:                           IF (NOUNIT) X(J) = X(J)/A(J,J)
  230:                           TEMP = X(J)
  231:                           DO 10 I = J - 1,1,-1
  232:                               X(I) = X(I) - TEMP*A(I,J)
  233:    10                     CONTINUE
  234:                       END IF
  235:    20             CONTINUE
  236:               ELSE
  237:                   JX = KX + (N-1)*INCX
  238:                   DO 40 J = N,1,-1
  239:                       IF (X(JX).NE.ZERO) THEN
  240:                           IF (NOUNIT) X(JX) = X(JX)/A(J,J)
  241:                           TEMP = X(JX)
  242:                           IX = JX
  243:                           DO 30 I = J - 1,1,-1
  244:                               IX = IX - INCX
  245:                               X(IX) = X(IX) - TEMP*A(I,J)
  246:    30                     CONTINUE
  247:                       END IF
  248:                       JX = JX - INCX
  249:    40             CONTINUE
  250:               END IF
  251:           ELSE
  252:               IF (INCX.EQ.1) THEN
  253:                   DO 60 J = 1,N
  254:                       IF (X(J).NE.ZERO) THEN
  255:                           IF (NOUNIT) X(J) = X(J)/A(J,J)
  256:                           TEMP = X(J)
  257:                           DO 50 I = J + 1,N
  258:                               X(I) = X(I) - TEMP*A(I,J)
  259:    50                     CONTINUE
  260:                       END IF
  261:    60             CONTINUE
  262:               ELSE
  263:                   JX = KX
  264:                   DO 80 J = 1,N
  265:                       IF (X(JX).NE.ZERO) THEN
  266:                           IF (NOUNIT) X(JX) = X(JX)/A(J,J)
  267:                           TEMP = X(JX)
  268:                           IX = JX
  269:                           DO 70 I = J + 1,N
  270:                               IX = IX + INCX
  271:                               X(IX) = X(IX) - TEMP*A(I,J)
  272:    70                     CONTINUE
  273:                       END IF
  274:                       JX = JX + INCX
  275:    80             CONTINUE
  276:               END IF
  277:           END IF
  278:       ELSE
  279: *
  280: *        Form  x := inv( A**T )*x.
  281: *
  282:           IF (LSAME(UPLO,'U')) THEN
  283:               IF (INCX.EQ.1) THEN
  284:                   DO 100 J = 1,N
  285:                       TEMP = X(J)
  286:                       DO 90 I = 1,J - 1
  287:                           TEMP = TEMP - A(I,J)*X(I)
  288:    90                 CONTINUE
  289:                       IF (NOUNIT) TEMP = TEMP/A(J,J)
  290:                       X(J) = TEMP
  291:   100             CONTINUE
  292:               ELSE
  293:                   JX = KX
  294:                   DO 120 J = 1,N
  295:                       TEMP = X(JX)
  296:                       IX = KX
  297:                       DO 110 I = 1,J - 1
  298:                           TEMP = TEMP - A(I,J)*X(IX)
  299:                           IX = IX + INCX
  300:   110                 CONTINUE
  301:                       IF (NOUNIT) TEMP = TEMP/A(J,J)
  302:                       X(JX) = TEMP
  303:                       JX = JX + INCX
  304:   120             CONTINUE
  305:               END IF
  306:           ELSE
  307:               IF (INCX.EQ.1) THEN
  308:                   DO 140 J = N,1,-1
  309:                       TEMP = X(J)
  310:                       DO 130 I = N,J + 1,-1
  311:                           TEMP = TEMP - A(I,J)*X(I)
  312:   130                 CONTINUE
  313:                       IF (NOUNIT) TEMP = TEMP/A(J,J)
  314:                       X(J) = TEMP
  315:   140             CONTINUE
  316:               ELSE
  317:                   KX = KX + (N-1)*INCX
  318:                   JX = KX
  319:                   DO 160 J = N,1,-1
  320:                       TEMP = X(JX)
  321:                       IX = KX
  322:                       DO 150 I = N,J + 1,-1
  323:                           TEMP = TEMP - A(I,J)*X(IX)
  324:                           IX = IX - INCX
  325:   150                 CONTINUE
  326:                       IF (NOUNIT) TEMP = TEMP/A(J,J)
  327:                       X(JX) = TEMP
  328:                       JX = JX - INCX
  329:   160             CONTINUE
  330:               END IF
  331:           END IF
  332:       END IF
  333: *
  334:       RETURN
  335: *
  336: *     End of DTRSV .
  337: *
  338:       END

CVSweb interface <joel.bertrand@systella.fr>