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

1.8       bertrand    1: *> \brief \b ZTBSV
                      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 ZTBSV(UPLO,TRANS,DIAG,N,K,A,LDA,X,INCX)
1.13      bertrand   12: *
1.8       bertrand   13: *       .. Scalar Arguments ..
                     14: *       INTEGER INCX,K,LDA,N
                     15: *       CHARACTER DIAG,TRANS,UPLO
                     16: *       ..
                     17: *       .. Array Arguments ..
                     18: *       COMPLEX*16 A(LDA,*),X(*)
                     19: *       ..
1.13      bertrand   20: *
1.8       bertrand   21: *
                     22: *> \par Purpose:
                     23: *  =============
                     24: *>
                     25: *> \verbatim
                     26: *>
                     27: *> ZTBSV  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 band matrix, with ( k + 1 )
                     33: *> diagonals.
                     34: *>
                     35: *> No test for singularity or near-singularity is included in this
                     36: *> routine. Such tests must be performed before calling this routine.
                     37: *> \endverbatim
                     38: *
                     39: *  Arguments:
                     40: *  ==========
                     41: *
                     42: *> \param[in] UPLO
                     43: *> \verbatim
                     44: *>          UPLO is CHARACTER*1
                     45: *>           On entry, UPLO specifies whether the matrix is an upper or
                     46: *>           lower triangular matrix as follows:
                     47: *>
                     48: *>              UPLO = 'U' or 'u'   A is an upper triangular matrix.
                     49: *>
                     50: *>              UPLO = 'L' or 'l'   A is a lower triangular matrix.
                     51: *> \endverbatim
                     52: *>
                     53: *> \param[in] TRANS
                     54: *> \verbatim
                     55: *>          TRANS is CHARACTER*1
                     56: *>           On entry, TRANS specifies the equations to be solved as
                     57: *>           follows:
                     58: *>
                     59: *>              TRANS = 'N' or 'n'   A*x = b.
                     60: *>
                     61: *>              TRANS = 'T' or 't'   A**T*x = b.
                     62: *>
                     63: *>              TRANS = 'C' or 'c'   A**H*x = b.
                     64: *> \endverbatim
                     65: *>
                     66: *> \param[in] DIAG
                     67: *> \verbatim
                     68: *>          DIAG is CHARACTER*1
                     69: *>           On entry, DIAG specifies whether or not A is unit
                     70: *>           triangular as follows:
                     71: *>
                     72: *>              DIAG = 'U' or 'u'   A is assumed to be unit triangular.
                     73: *>
                     74: *>              DIAG = 'N' or 'n'   A is not assumed to be unit
                     75: *>                                  triangular.
                     76: *> \endverbatim
                     77: *>
                     78: *> \param[in] N
                     79: *> \verbatim
                     80: *>          N is INTEGER
                     81: *>           On entry, N specifies the order of the matrix A.
                     82: *>           N must be at least zero.
                     83: *> \endverbatim
                     84: *>
                     85: *> \param[in] K
                     86: *> \verbatim
                     87: *>          K is INTEGER
                     88: *>           On entry with UPLO = 'U' or 'u', K specifies the number of
                     89: *>           super-diagonals of the matrix A.
                     90: *>           On entry with UPLO = 'L' or 'l', K specifies the number of
                     91: *>           sub-diagonals of the matrix A.
                     92: *>           K must satisfy  0 .le. K.
                     93: *> \endverbatim
                     94: *>
                     95: *> \param[in] A
                     96: *> \verbatim
1.14      bertrand   97: *>          A is COMPLEX*16 array, dimension ( LDA, N )
1.8       bertrand   98: *>           Before entry with UPLO = 'U' or 'u', the leading ( k + 1 )
                     99: *>           by n part of the array A must contain the upper triangular
                    100: *>           band part of the matrix of coefficients, supplied column by
                    101: *>           column, with the leading diagonal of the matrix in row
                    102: *>           ( k + 1 ) of the array, the first super-diagonal starting at
                    103: *>           position 2 in row k, and so on. The top left k by k triangle
                    104: *>           of the array A is not referenced.
                    105: *>           The following program segment will transfer an upper
                    106: *>           triangular band matrix from conventional full matrix storage
                    107: *>           to band storage:
                    108: *>
                    109: *>                 DO 20, J = 1, N
                    110: *>                    M = K + 1 - J
                    111: *>                    DO 10, I = MAX( 1, J - K ), J
                    112: *>                       A( M + I, J ) = matrix( I, J )
                    113: *>              10    CONTINUE
                    114: *>              20 CONTINUE
                    115: *>
                    116: *>           Before entry with UPLO = 'L' or 'l', the leading ( k + 1 )
                    117: *>           by n part of the array A must contain the lower triangular
                    118: *>           band part of the matrix of coefficients, supplied column by
                    119: *>           column, with the leading diagonal of the matrix in row 1 of
                    120: *>           the array, the first sub-diagonal starting at position 1 in
                    121: *>           row 2, and so on. The bottom right k by k triangle of the
                    122: *>           array A is not referenced.
                    123: *>           The following program segment will transfer a lower
                    124: *>           triangular band matrix from conventional full matrix storage
                    125: *>           to band storage:
                    126: *>
                    127: *>                 DO 20, J = 1, N
                    128: *>                    M = 1 - J
                    129: *>                    DO 10, I = J, MIN( N, J + K )
                    130: *>                       A( M + I, J ) = matrix( I, J )
                    131: *>              10    CONTINUE
                    132: *>              20 CONTINUE
                    133: *>
                    134: *>           Note that when DIAG = 'U' or 'u' the elements of the array A
                    135: *>           corresponding to the diagonal elements of the matrix are not
                    136: *>           referenced, but are assumed to be unity.
                    137: *> \endverbatim
                    138: *>
                    139: *> \param[in] LDA
                    140: *> \verbatim
                    141: *>          LDA is INTEGER
                    142: *>           On entry, LDA specifies the first dimension of A as declared
                    143: *>           in the calling (sub) program. LDA must be at least
                    144: *>           ( k + 1 ).
                    145: *> \endverbatim
                    146: *>
                    147: *> \param[in,out] X
                    148: *> \verbatim
1.14      bertrand  149: *>          X is COMPLEX*16 array, dimension at least
1.8       bertrand  150: *>           ( 1 + ( n - 1 )*abs( INCX ) ).
                    151: *>           Before entry, the incremented array X must contain the n
                    152: *>           element right-hand side vector b. On exit, X is overwritten
                    153: *>           with the solution vector x.
                    154: *> \endverbatim
                    155: *>
                    156: *> \param[in] INCX
                    157: *> \verbatim
                    158: *>          INCX is INTEGER
                    159: *>           On entry, INCX specifies the increment for the elements of
                    160: *>           X. INCX must not be zero.
                    161: *> \endverbatim
                    162: *
                    163: *  Authors:
                    164: *  ========
                    165: *
1.13      bertrand  166: *> \author Univ. of Tennessee
                    167: *> \author Univ. of California Berkeley
                    168: *> \author Univ. of Colorado Denver
                    169: *> \author NAG Ltd.
1.8       bertrand  170: *
                    171: *> \ingroup complex16_blas_level2
                    172: *
                    173: *> \par Further Details:
                    174: *  =====================
                    175: *>
                    176: *> \verbatim
                    177: *>
                    178: *>  Level 2 Blas routine.
                    179: *>
                    180: *>  -- Written on 22-October-1986.
                    181: *>     Jack Dongarra, Argonne National Lab.
                    182: *>     Jeremy Du Croz, Nag Central Office.
                    183: *>     Sven Hammarling, Nag Central Office.
                    184: *>     Richard Hanson, Sandia National Labs.
                    185: *> \endverbatim
                    186: *>
                    187: *  =====================================================================
1.1       bertrand  188:       SUBROUTINE ZTBSV(UPLO,TRANS,DIAG,N,K,A,LDA,X,INCX)
1.8       bertrand  189: *
1.16    ! bertrand  190: *  -- Reference BLAS level2 routine --
1.8       bertrand  191: *  -- Reference BLAS is a software package provided by Univ. of Tennessee,    --
                    192: *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
                    193: *
1.1       bertrand  194: *     .. Scalar Arguments ..
                    195:       INTEGER INCX,K,LDA,N
                    196:       CHARACTER DIAG,TRANS,UPLO
                    197: *     ..
                    198: *     .. Array Arguments ..
1.8       bertrand  199:       COMPLEX*16 A(LDA,*),X(*)
1.1       bertrand  200: *     ..
                    201: *
                    202: *  =====================================================================
                    203: *
                    204: *     .. Parameters ..
1.8       bertrand  205:       COMPLEX*16 ZERO
1.1       bertrand  206:       PARAMETER (ZERO= (0.0D+0,0.0D+0))
                    207: *     ..
                    208: *     .. Local Scalars ..
1.8       bertrand  209:       COMPLEX*16 TEMP
1.1       bertrand  210:       INTEGER I,INFO,IX,J,JX,KPLUS1,KX,L
                    211:       LOGICAL NOCONJ,NOUNIT
                    212: *     ..
                    213: *     .. External Functions ..
                    214:       LOGICAL LSAME
                    215:       EXTERNAL LSAME
                    216: *     ..
                    217: *     .. External Subroutines ..
                    218:       EXTERNAL XERBLA
                    219: *     ..
                    220: *     .. Intrinsic Functions ..
                    221:       INTRINSIC DCONJG,MAX,MIN
                    222: *     ..
                    223: *
                    224: *     Test the input parameters.
                    225: *
                    226:       INFO = 0
                    227:       IF (.NOT.LSAME(UPLO,'U') .AND. .NOT.LSAME(UPLO,'L')) THEN
                    228:           INFO = 1
                    229:       ELSE IF (.NOT.LSAME(TRANS,'N') .AND. .NOT.LSAME(TRANS,'T') .AND.
                    230:      +         .NOT.LSAME(TRANS,'C')) THEN
                    231:           INFO = 2
                    232:       ELSE IF (.NOT.LSAME(DIAG,'U') .AND. .NOT.LSAME(DIAG,'N')) THEN
                    233:           INFO = 3
                    234:       ELSE IF (N.LT.0) THEN
                    235:           INFO = 4
                    236:       ELSE IF (K.LT.0) THEN
                    237:           INFO = 5
                    238:       ELSE IF (LDA.LT. (K+1)) THEN
                    239:           INFO = 7
                    240:       ELSE IF (INCX.EQ.0) THEN
                    241:           INFO = 9
                    242:       END IF
                    243:       IF (INFO.NE.0) THEN
                    244:           CALL XERBLA('ZTBSV ',INFO)
                    245:           RETURN
                    246:       END IF
                    247: *
                    248: *     Quick return if possible.
                    249: *
                    250:       IF (N.EQ.0) RETURN
                    251: *
                    252:       NOCONJ = LSAME(TRANS,'T')
                    253:       NOUNIT = LSAME(DIAG,'N')
                    254: *
                    255: *     Set up the start point in X if the increment is not unity. This
                    256: *     will be  ( N - 1 )*INCX  too small for descending loops.
                    257: *
                    258:       IF (INCX.LE.0) THEN
                    259:           KX = 1 - (N-1)*INCX
                    260:       ELSE IF (INCX.NE.1) THEN
                    261:           KX = 1
                    262:       END IF
                    263: *
                    264: *     Start the operations. In this version the elements of A are
                    265: *     accessed by sequentially with one pass through A.
                    266: *
                    267:       IF (LSAME(TRANS,'N')) THEN
                    268: *
                    269: *        Form  x := inv( A )*x.
                    270: *
                    271:           IF (LSAME(UPLO,'U')) THEN
                    272:               KPLUS1 = K + 1
                    273:               IF (INCX.EQ.1) THEN
                    274:                   DO 20 J = N,1,-1
                    275:                       IF (X(J).NE.ZERO) THEN
                    276:                           L = KPLUS1 - J
                    277:                           IF (NOUNIT) X(J) = X(J)/A(KPLUS1,J)
                    278:                           TEMP = X(J)
                    279:                           DO 10 I = J - 1,MAX(1,J-K),-1
                    280:                               X(I) = X(I) - TEMP*A(L+I,J)
                    281:    10                     CONTINUE
                    282:                       END IF
                    283:    20             CONTINUE
                    284:               ELSE
                    285:                   KX = KX + (N-1)*INCX
                    286:                   JX = KX
                    287:                   DO 40 J = N,1,-1
                    288:                       KX = KX - INCX
                    289:                       IF (X(JX).NE.ZERO) THEN
                    290:                           IX = KX
                    291:                           L = KPLUS1 - J
                    292:                           IF (NOUNIT) X(JX) = X(JX)/A(KPLUS1,J)
                    293:                           TEMP = X(JX)
                    294:                           DO 30 I = J - 1,MAX(1,J-K),-1
                    295:                               X(IX) = X(IX) - TEMP*A(L+I,J)
                    296:                               IX = IX - INCX
                    297:    30                     CONTINUE
                    298:                       END IF
                    299:                       JX = JX - INCX
                    300:    40             CONTINUE
                    301:               END IF
                    302:           ELSE
                    303:               IF (INCX.EQ.1) THEN
                    304:                   DO 60 J = 1,N
                    305:                       IF (X(J).NE.ZERO) THEN
                    306:                           L = 1 - J
                    307:                           IF (NOUNIT) X(J) = X(J)/A(1,J)
                    308:                           TEMP = X(J)
                    309:                           DO 50 I = J + 1,MIN(N,J+K)
                    310:                               X(I) = X(I) - TEMP*A(L+I,J)
                    311:    50                     CONTINUE
                    312:                       END IF
                    313:    60             CONTINUE
                    314:               ELSE
                    315:                   JX = KX
                    316:                   DO 80 J = 1,N
                    317:                       KX = KX + INCX
                    318:                       IF (X(JX).NE.ZERO) THEN
                    319:                           IX = KX
                    320:                           L = 1 - J
                    321:                           IF (NOUNIT) X(JX) = X(JX)/A(1,J)
                    322:                           TEMP = X(JX)
                    323:                           DO 70 I = J + 1,MIN(N,J+K)
                    324:                               X(IX) = X(IX) - TEMP*A(L+I,J)
                    325:                               IX = IX + INCX
                    326:    70                     CONTINUE
                    327:                       END IF
                    328:                       JX = JX + INCX
                    329:    80             CONTINUE
                    330:               END IF
                    331:           END IF
                    332:       ELSE
                    333: *
1.7       bertrand  334: *        Form  x := inv( A**T )*x  or  x := inv( A**H )*x.
1.1       bertrand  335: *
                    336:           IF (LSAME(UPLO,'U')) THEN
                    337:               KPLUS1 = K + 1
                    338:               IF (INCX.EQ.1) THEN
                    339:                   DO 110 J = 1,N
                    340:                       TEMP = X(J)
                    341:                       L = KPLUS1 - J
                    342:                       IF (NOCONJ) THEN
                    343:                           DO 90 I = MAX(1,J-K),J - 1
                    344:                               TEMP = TEMP - A(L+I,J)*X(I)
                    345:    90                     CONTINUE
                    346:                           IF (NOUNIT) TEMP = TEMP/A(KPLUS1,J)
                    347:                       ELSE
                    348:                           DO 100 I = MAX(1,J-K),J - 1
                    349:                               TEMP = TEMP - DCONJG(A(L+I,J))*X(I)
                    350:   100                     CONTINUE
                    351:                           IF (NOUNIT) TEMP = TEMP/DCONJG(A(KPLUS1,J))
                    352:                       END IF
                    353:                       X(J) = TEMP
                    354:   110             CONTINUE
                    355:               ELSE
                    356:                   JX = KX
                    357:                   DO 140 J = 1,N
                    358:                       TEMP = X(JX)
                    359:                       IX = KX
                    360:                       L = KPLUS1 - J
                    361:                       IF (NOCONJ) THEN
                    362:                           DO 120 I = MAX(1,J-K),J - 1
                    363:                               TEMP = TEMP - A(L+I,J)*X(IX)
                    364:                               IX = IX + INCX
                    365:   120                     CONTINUE
                    366:                           IF (NOUNIT) TEMP = TEMP/A(KPLUS1,J)
                    367:                       ELSE
                    368:                           DO 130 I = MAX(1,J-K),J - 1
                    369:                               TEMP = TEMP - DCONJG(A(L+I,J))*X(IX)
                    370:                               IX = IX + INCX
                    371:   130                     CONTINUE
                    372:                           IF (NOUNIT) TEMP = TEMP/DCONJG(A(KPLUS1,J))
                    373:                       END IF
                    374:                       X(JX) = TEMP
                    375:                       JX = JX + INCX
                    376:                       IF (J.GT.K) KX = KX + INCX
                    377:   140             CONTINUE
                    378:               END IF
                    379:           ELSE
                    380:               IF (INCX.EQ.1) THEN
                    381:                   DO 170 J = N,1,-1
                    382:                       TEMP = X(J)
                    383:                       L = 1 - J
                    384:                       IF (NOCONJ) THEN
                    385:                           DO 150 I = MIN(N,J+K),J + 1,-1
                    386:                               TEMP = TEMP - A(L+I,J)*X(I)
                    387:   150                     CONTINUE
                    388:                           IF (NOUNIT) TEMP = TEMP/A(1,J)
                    389:                       ELSE
                    390:                           DO 160 I = MIN(N,J+K),J + 1,-1
                    391:                               TEMP = TEMP - DCONJG(A(L+I,J))*X(I)
                    392:   160                     CONTINUE
                    393:                           IF (NOUNIT) TEMP = TEMP/DCONJG(A(1,J))
                    394:                       END IF
                    395:                       X(J) = TEMP
                    396:   170             CONTINUE
                    397:               ELSE
                    398:                   KX = KX + (N-1)*INCX
                    399:                   JX = KX
                    400:                   DO 200 J = N,1,-1
                    401:                       TEMP = X(JX)
                    402:                       IX = KX
                    403:                       L = 1 - J
                    404:                       IF (NOCONJ) THEN
                    405:                           DO 180 I = MIN(N,J+K),J + 1,-1
                    406:                               TEMP = TEMP - A(L+I,J)*X(IX)
                    407:                               IX = IX - INCX
                    408:   180                     CONTINUE
                    409:                           IF (NOUNIT) TEMP = TEMP/A(1,J)
                    410:                       ELSE
                    411:                           DO 190 I = MIN(N,J+K),J + 1,-1
                    412:                               TEMP = TEMP - DCONJG(A(L+I,J))*X(IX)
                    413:                               IX = IX - INCX
                    414:   190                     CONTINUE
                    415:                           IF (NOUNIT) TEMP = TEMP/DCONJG(A(1,J))
                    416:                       END IF
                    417:                       X(JX) = TEMP
                    418:                       JX = JX - INCX
                    419:                       IF ((N-J).GE.K) KX = KX - INCX
                    420:   200             CONTINUE
                    421:               END IF
                    422:           END IF
                    423:       END IF
                    424: *
                    425:       RETURN
                    426: *
1.16    ! bertrand  427: *     End of ZTBSV
1.1       bertrand  428: *
                    429:       END

CVSweb interface <joel.bertrand@systella.fr>