Annotation of rpl/lapack/blas/dtpmv.f, revision 1.12

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

CVSweb interface <joel.bertrand@systella.fr>