BLAS/SRC/dtrmv.f

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