File:  [local] / rpl / lapack / blas / ztpsv.f
Revision 1.16: download - view: text, annotated - select for diffs - revision graph
Mon Aug 7 08:38:46 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 ZTPSV
    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 ZTPSV(UPLO,TRANS,DIAG,N,AP,X,INCX)
   12: *
   13: *       .. Scalar Arguments ..
   14: *       INTEGER INCX,N
   15: *       CHARACTER DIAG,TRANS,UPLO
   16: *       ..
   17: *       .. Array Arguments ..
   18: *       COMPLEX*16 AP(*),X(*)
   19: *       ..
   20: *
   21: *
   22: *> \par Purpose:
   23: *  =============
   24: *>
   25: *> \verbatim
   26: *>
   27: *> ZTPSV  solves one of the systems of equations
   28: *>
   29: *>    A*x = b,   or   A**T*x = b,   or   A**H*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, supplied in packed form.
   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**H*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] AP
   85: *> \verbatim
   86: *>          AP is COMPLEX*16 array, dimension at least
   87: *>           ( ( n*( n + 1 ) )/2 ).
   88: *>           Before entry with  UPLO = 'U' or 'u', the array AP must
   89: *>           contain the upper triangular matrix packed sequentially,
   90: *>           column by column, so that AP( 1 ) contains a( 1, 1 ),
   91: *>           AP( 2 ) and AP( 3 ) contain a( 1, 2 ) and a( 2, 2 )
   92: *>           respectively, and so on.
   93: *>           Before entry with UPLO = 'L' or 'l', the array AP must
   94: *>           contain the lower triangular matrix packed sequentially,
   95: *>           column by column, so that AP( 1 ) contains a( 1, 1 ),
   96: *>           AP( 2 ) and AP( 3 ) contain a( 2, 1 ) and a( 3, 1 )
   97: *>           respectively, and so on.
   98: *>           Note that when  DIAG = 'U' or 'u', the diagonal elements of
   99: *>           A are not referenced, but are assumed to be unity.
  100: *> \endverbatim
  101: *>
  102: *> \param[in,out] X
  103: *> \verbatim
  104: *>          X is COMPLEX*16 array, dimension at least
  105: *>           ( 1 + ( n - 1 )*abs( INCX ) ).
  106: *>           Before entry, the incremented array X must contain the n
  107: *>           element right-hand side vector b. On exit, X is overwritten
  108: *>           with the solution vector x.
  109: *> \endverbatim
  110: *>
  111: *> \param[in] INCX
  112: *> \verbatim
  113: *>          INCX is INTEGER
  114: *>           On entry, INCX specifies the increment for the elements of
  115: *>           X. INCX must not be zero.
  116: *> \endverbatim
  117: *
  118: *  Authors:
  119: *  ========
  120: *
  121: *> \author Univ. of Tennessee
  122: *> \author Univ. of California Berkeley
  123: *> \author Univ. of Colorado Denver
  124: *> \author NAG Ltd.
  125: *
  126: *> \ingroup complex16_blas_level2
  127: *
  128: *> \par Further Details:
  129: *  =====================
  130: *>
  131: *> \verbatim
  132: *>
  133: *>  Level 2 Blas routine.
  134: *>
  135: *>  -- Written on 22-October-1986.
  136: *>     Jack Dongarra, Argonne National Lab.
  137: *>     Jeremy Du Croz, Nag Central Office.
  138: *>     Sven Hammarling, Nag Central Office.
  139: *>     Richard Hanson, Sandia National Labs.
  140: *> \endverbatim
  141: *>
  142: *  =====================================================================
  143:       SUBROUTINE ZTPSV(UPLO,TRANS,DIAG,N,AP,X,INCX)
  144: *
  145: *  -- Reference BLAS level2 routine --
  146: *  -- Reference BLAS is a software package provided by Univ. of Tennessee,    --
  147: *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
  148: *
  149: *     .. Scalar Arguments ..
  150:       INTEGER INCX,N
  151:       CHARACTER DIAG,TRANS,UPLO
  152: *     ..
  153: *     .. Array Arguments ..
  154:       COMPLEX*16 AP(*),X(*)
  155: *     ..
  156: *
  157: *  =====================================================================
  158: *
  159: *     .. Parameters ..
  160:       COMPLEX*16 ZERO
  161:       PARAMETER (ZERO= (0.0D+0,0.0D+0))
  162: *     ..
  163: *     .. Local Scalars ..
  164:       COMPLEX*16 TEMP
  165:       INTEGER I,INFO,IX,J,JX,K,KK,KX
  166:       LOGICAL NOCONJ,NOUNIT
  167: *     ..
  168: *     .. External Functions ..
  169:       LOGICAL LSAME
  170:       EXTERNAL LSAME
  171: *     ..
  172: *     .. External Subroutines ..
  173:       EXTERNAL XERBLA
  174: *     ..
  175: *     .. Intrinsic Functions ..
  176:       INTRINSIC DCONJG
  177: *     ..
  178: *
  179: *     Test the input parameters.
  180: *
  181:       INFO = 0
  182:       IF (.NOT.LSAME(UPLO,'U') .AND. .NOT.LSAME(UPLO,'L')) THEN
  183:           INFO = 1
  184:       ELSE IF (.NOT.LSAME(TRANS,'N') .AND. .NOT.LSAME(TRANS,'T') .AND.
  185:      +         .NOT.LSAME(TRANS,'C')) THEN
  186:           INFO = 2
  187:       ELSE IF (.NOT.LSAME(DIAG,'U') .AND. .NOT.LSAME(DIAG,'N')) THEN
  188:           INFO = 3
  189:       ELSE IF (N.LT.0) THEN
  190:           INFO = 4
  191:       ELSE IF (INCX.EQ.0) THEN
  192:           INFO = 7
  193:       END IF
  194:       IF (INFO.NE.0) THEN
  195:           CALL XERBLA('ZTPSV ',INFO)
  196:           RETURN
  197:       END IF
  198: *
  199: *     Quick return if possible.
  200: *
  201:       IF (N.EQ.0) RETURN
  202: *
  203:       NOCONJ = LSAME(TRANS,'T')
  204:       NOUNIT = LSAME(DIAG,'N')
  205: *
  206: *     Set up the start point in X if the increment is not unity. This
  207: *     will be  ( N - 1 )*INCX  too small for descending loops.
  208: *
  209:       IF (INCX.LE.0) THEN
  210:           KX = 1 - (N-1)*INCX
  211:       ELSE IF (INCX.NE.1) THEN
  212:           KX = 1
  213:       END IF
  214: *
  215: *     Start the operations. In this version the elements of AP are
  216: *     accessed sequentially with one pass through AP.
  217: *
  218:       IF (LSAME(TRANS,'N')) THEN
  219: *
  220: *        Form  x := inv( A )*x.
  221: *
  222:           IF (LSAME(UPLO,'U')) THEN
  223:               KK = (N* (N+1))/2
  224:               IF (INCX.EQ.1) THEN
  225:                   DO 20 J = N,1,-1
  226:                       IF (X(J).NE.ZERO) THEN
  227:                           IF (NOUNIT) X(J) = X(J)/AP(KK)
  228:                           TEMP = X(J)
  229:                           K = KK - 1
  230:                           DO 10 I = J - 1,1,-1
  231:                               X(I) = X(I) - TEMP*AP(K)
  232:                               K = K - 1
  233:    10                     CONTINUE
  234:                       END IF
  235:                       KK = KK - J
  236:    20             CONTINUE
  237:               ELSE
  238:                   JX = KX + (N-1)*INCX
  239:                   DO 40 J = N,1,-1
  240:                       IF (X(JX).NE.ZERO) THEN
  241:                           IF (NOUNIT) X(JX) = X(JX)/AP(KK)
  242:                           TEMP = X(JX)
  243:                           IX = JX
  244:                           DO 30 K = KK - 1,KK - J + 1,-1
  245:                               IX = IX - INCX
  246:                               X(IX) = X(IX) - TEMP*AP(K)
  247:    30                     CONTINUE
  248:                       END IF
  249:                       JX = JX - INCX
  250:                       KK = KK - J
  251:    40             CONTINUE
  252:               END IF
  253:           ELSE
  254:               KK = 1
  255:               IF (INCX.EQ.1) THEN
  256:                   DO 60 J = 1,N
  257:                       IF (X(J).NE.ZERO) THEN
  258:                           IF (NOUNIT) X(J) = X(J)/AP(KK)
  259:                           TEMP = X(J)
  260:                           K = KK + 1
  261:                           DO 50 I = J + 1,N
  262:                               X(I) = X(I) - TEMP*AP(K)
  263:                               K = K + 1
  264:    50                     CONTINUE
  265:                       END IF
  266:                       KK = KK + (N-J+1)
  267:    60             CONTINUE
  268:               ELSE
  269:                   JX = KX
  270:                   DO 80 J = 1,N
  271:                       IF (X(JX).NE.ZERO) THEN
  272:                           IF (NOUNIT) X(JX) = X(JX)/AP(KK)
  273:                           TEMP = X(JX)
  274:                           IX = JX
  275:                           DO 70 K = KK + 1,KK + N - J
  276:                               IX = IX + INCX
  277:                               X(IX) = X(IX) - TEMP*AP(K)
  278:    70                     CONTINUE
  279:                       END IF
  280:                       JX = JX + INCX
  281:                       KK = KK + (N-J+1)
  282:    80             CONTINUE
  283:               END IF
  284:           END IF
  285:       ELSE
  286: *
  287: *        Form  x := inv( A**T )*x  or  x := inv( A**H )*x.
  288: *
  289:           IF (LSAME(UPLO,'U')) THEN
  290:               KK = 1
  291:               IF (INCX.EQ.1) THEN
  292:                   DO 110 J = 1,N
  293:                       TEMP = X(J)
  294:                       K = KK
  295:                       IF (NOCONJ) THEN
  296:                           DO 90 I = 1,J - 1
  297:                               TEMP = TEMP - AP(K)*X(I)
  298:                               K = K + 1
  299:    90                     CONTINUE
  300:                           IF (NOUNIT) TEMP = TEMP/AP(KK+J-1)
  301:                       ELSE
  302:                           DO 100 I = 1,J - 1
  303:                               TEMP = TEMP - DCONJG(AP(K))*X(I)
  304:                               K = K + 1
  305:   100                     CONTINUE
  306:                           IF (NOUNIT) TEMP = TEMP/DCONJG(AP(KK+J-1))
  307:                       END IF
  308:                       X(J) = TEMP
  309:                       KK = KK + J
  310:   110             CONTINUE
  311:               ELSE
  312:                   JX = KX
  313:                   DO 140 J = 1,N
  314:                       TEMP = X(JX)
  315:                       IX = KX
  316:                       IF (NOCONJ) THEN
  317:                           DO 120 K = KK,KK + J - 2
  318:                               TEMP = TEMP - AP(K)*X(IX)
  319:                               IX = IX + INCX
  320:   120                     CONTINUE
  321:                           IF (NOUNIT) TEMP = TEMP/AP(KK+J-1)
  322:                       ELSE
  323:                           DO 130 K = KK,KK + J - 2
  324:                               TEMP = TEMP - DCONJG(AP(K))*X(IX)
  325:                               IX = IX + INCX
  326:   130                     CONTINUE
  327:                           IF (NOUNIT) TEMP = TEMP/DCONJG(AP(KK+J-1))
  328:                       END IF
  329:                       X(JX) = TEMP
  330:                       JX = JX + INCX
  331:                       KK = KK + J
  332:   140             CONTINUE
  333:               END IF
  334:           ELSE
  335:               KK = (N* (N+1))/2
  336:               IF (INCX.EQ.1) THEN
  337:                   DO 170 J = N,1,-1
  338:                       TEMP = X(J)
  339:                       K = KK
  340:                       IF (NOCONJ) THEN
  341:                           DO 150 I = N,J + 1,-1
  342:                               TEMP = TEMP - AP(K)*X(I)
  343:                               K = K - 1
  344:   150                     CONTINUE
  345:                           IF (NOUNIT) TEMP = TEMP/AP(KK-N+J)
  346:                       ELSE
  347:                           DO 160 I = N,J + 1,-1
  348:                               TEMP = TEMP - DCONJG(AP(K))*X(I)
  349:                               K = K - 1
  350:   160                     CONTINUE
  351:                           IF (NOUNIT) TEMP = TEMP/DCONJG(AP(KK-N+J))
  352:                       END IF
  353:                       X(J) = TEMP
  354:                       KK = KK - (N-J+1)
  355:   170             CONTINUE
  356:               ELSE
  357:                   KX = KX + (N-1)*INCX
  358:                   JX = KX
  359:                   DO 200 J = N,1,-1
  360:                       TEMP = X(JX)
  361:                       IX = KX
  362:                       IF (NOCONJ) THEN
  363:                           DO 180 K = KK,KK - (N- (J+1)),-1
  364:                               TEMP = TEMP - AP(K)*X(IX)
  365:                               IX = IX - INCX
  366:   180                     CONTINUE
  367:                           IF (NOUNIT) TEMP = TEMP/AP(KK-N+J)
  368:                       ELSE
  369:                           DO 190 K = KK,KK - (N- (J+1)),-1
  370:                               TEMP = TEMP - DCONJG(AP(K))*X(IX)
  371:                               IX = IX - INCX
  372:   190                     CONTINUE
  373:                           IF (NOUNIT) TEMP = TEMP/DCONJG(AP(KK-N+J))
  374:                       END IF
  375:                       X(JX) = TEMP
  376:                       JX = JX - INCX
  377:                       KK = KK - (N-J+1)
  378:   200             CONTINUE
  379:               END IF
  380:           END IF
  381:       END IF
  382: *
  383:       RETURN
  384: *
  385: *     End of ZTPSV
  386: *
  387:       END

CVSweb interface <joel.bertrand@systella.fr>