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

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>