BLAS/SRC/dsymv.f

   1 *> \brief \b DSYMV
   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 DSYMV(UPLO,N,ALPHA,A,LDA,X,INCX,BETA,Y,INCY)
  12 *
  13 *       .. Scalar Arguments ..
  14 *       DOUBLE PRECISION ALPHA,BETA
  15 *       INTEGER INCX,INCY,LDA,N
  16 *       CHARACTER UPLO
  17 *       ..
  18 *       .. Array Arguments ..
  19 *       DOUBLE PRECISION A(LDA,*),X(*),Y(*)
  20 *       ..
  21 *
  22 *
  23 *> \par Purpose:
  24 *  =============
  25 *>
  26 *> \verbatim
  27 *>
  28 *> DSYMV  performs the matrix-vector  operation
  29 *>
  30 *>    y := alpha*A*x + beta*y,
  31 *>
  32 *> where alpha and beta are scalars, x and y are n element vectors and
  33 *> A is an n by n symmetric matrix.
  34 *> \endverbatim
  35 *
  36 *  Arguments:
  37 *  ==========
  38 *
  39 *> \param[in] UPLO
  40 *> \verbatim
  41 *>          UPLO is CHARACTER*1
  42 *>           On entry, UPLO specifies whether the upper or lower
  43 *>           triangular part of the array A is to be referenced as
  44 *>           follows:
  45 *>
  46 *>              UPLO = 'U' or 'u'   Only the upper triangular part of A
  47 *>                                  is to be referenced.
  48 *>
  49 *>              UPLO = 'L' or 'l'   Only the lower triangular part of A
  50 *>                                  is to be referenced.
  51 *> \endverbatim
  52 *>
  53 *> \param[in] N
  54 *> \verbatim
  55 *>          N is INTEGER
  56 *>           On entry, N specifies the order of the matrix A.
  57 *>           N must be at least zero.
  58 *> \endverbatim
  59 *>
  60 *> \param[in] ALPHA
  61 *> \verbatim
  62 *>          ALPHA is DOUBLE PRECISION.
  63 *>           On entry, ALPHA specifies the scalar alpha.
  64 *> \endverbatim
  65 *>
  66 *> \param[in] A
  67 *> \verbatim
  68 *>          A is DOUBLE PRECISION array of DIMENSION ( LDA, n ).
  69 *>           Before entry with  UPLO = 'U' or 'u', the leading n by n
  70 *>           upper triangular part of the array A must contain the upper
  71 *>           triangular part of the symmetric matrix and the strictly
  72 *>           lower triangular part of A is not referenced.
  73 *>           Before entry with UPLO = 'L' or 'l', the leading n by n
  74 *>           lower triangular part of the array A must contain the lower
  75 *>           triangular part of the symmetric matrix and the strictly
  76 *>           upper triangular part of A is not referenced.
  77 *> \endverbatim
  78 *>
  79 *> \param[in] LDA
  80 *> \verbatim
  81 *>          LDA is INTEGER
  82 *>           On entry, LDA specifies the first dimension of A as declared
  83 *>           in the calling (sub) program. LDA must be at least
  84 *>           max( 1, n ).
  85 *> \endverbatim
  86 *>
  87 *> \param[in] X
  88 *> \verbatim
  89 *>          X is DOUBLE PRECISION array of dimension at least
  90 *>           ( 1 + ( n - 1 )*abs( INCX ) ).
  91 *>           Before entry, the incremented array X must contain the n
  92 *>           element vector x.
  93 *> \endverbatim
  94 *>
  95 *> \param[in] INCX
  96 *> \verbatim
  97 *>          INCX is INTEGER
  98 *>           On entry, INCX specifies the increment for the elements of
  99 *>           X. INCX must not be zero.
 100 *> \endverbatim
 101 *>
 102 *> \param[in] BETA
 103 *> \verbatim
 104 *>          BETA is DOUBLE PRECISION.
 105 *>           On entry, BETA specifies the scalar beta. When BETA is
 106 *>           supplied as zero then Y need not be set on input.
 107 *> \endverbatim
 108 *>
 109 *> \param[in,out] Y
 110 *> \verbatim
 111 *>          Y is DOUBLE PRECISION array of dimension at least
 112 *>           ( 1 + ( n - 1 )*abs( INCY ) ).
 113 *>           Before entry, the incremented array Y must contain the n
 114 *>           element vector y. On exit, Y is overwritten by the updated
 115 *>           vector y.
 116 *> \endverbatim
 117 *>
 118 *> \param[in] INCY
 119 *> \verbatim
 120 *>          INCY is INTEGER
 121 *>           On entry, INCY specifies the increment for the elements of
 122 *>           Y. INCY must not be zero.
 123 *> \endverbatim
 124 *
 125 *  Authors:
 126 *  ========
 127 *
 128 *> \author Univ. of Tennessee
 129 *> \author Univ. of California Berkeley
 130 *> \author Univ. of Colorado Denver
 131 *> \author NAG Ltd.
 132 *
 133 *> \date November 2011
 134 *
 135 *> \ingroup double_blas_level2
 136 *
 137 *> \par Further Details:
 138 *  =====================
 139 *>
 140 *> \verbatim
 141 *>
 142 *>  Level 2 Blas routine.
 143 *>  The vector and matrix arguments are not referenced when N = 0, or M = 0
 144 *>
 145 *>  -- Written on 22-October-1986.
 146 *>     Jack Dongarra, Argonne National Lab.
 147 *>     Jeremy Du Croz, Nag Central Office.
 148 *>     Sven Hammarling, Nag Central Office.
 149 *>     Richard Hanson, Sandia National Labs.
 150 *> \endverbatim
 151 *>
 152 *  =====================================================================
 153       SUBROUTINE DSYMV(UPLO,N,ALPHA,A,LDA,X,INCX,BETA,Y,INCY)
 154 *
 155 *  -- Reference BLAS level2 routine (version 3.4.0) --
 156 *  -- Reference BLAS is a software package provided by Univ. of Tennessee,    --
 157 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
 158 *     November 2011
 159 *
 160 *     .. Scalar Arguments ..
 161       DOUBLE PRECISION ALPHA,BETA
 162       INTEGER INCX,INCY,LDA,N
 163       CHARACTER UPLO
 164 *     ..
 165 *     .. Array Arguments ..
 166       DOUBLE PRECISION A(LDA,*),X(*),Y(*)
 167 *     ..
 168 *
 169 *  =====================================================================
 170 *
 171 *     .. Parameters ..
 172       DOUBLE PRECISION ONE,ZERO
 173       PARAMETER (ONE=1.0D+0,ZERO=0.0D+0)
 174 *     ..
 175 *     .. Local Scalars ..
 176       DOUBLE PRECISION TEMP1,TEMP2
 177       INTEGER I,INFO,IX,IY,J,JX,JY,KX,KY
 178 *     ..
 179 *     .. External Functions ..
 180       LOGICAL LSAME
 181       EXTERNAL LSAME
 182 *     ..
 183 *     .. External Subroutines ..
 184       EXTERNAL XERBLA
 185 *     ..
 186 *     .. Intrinsic Functions ..
 187       INTRINSIC MAX
 188 *     ..
 189 *
 190 *     Test the input parameters.
 191 *
 192       INFO = 0
 193       IF (.NOT.LSAME(UPLO,'U') .AND. .NOT.LSAME(UPLO,'L')) THEN
 194           INFO = 1
 195       ELSE IF (N.LT.0) THEN
 196           INFO = 2
 197       ELSE IF (LDA.LT.MAX(1,N)) THEN
 198           INFO = 5
 199       ELSE IF (INCX.EQ.0) THEN
 200           INFO = 7
 201       ELSE IF (INCY.EQ.0) THEN
 202           INFO = 10
 203       END IF
 204       IF (INFO.NE.0) THEN
 205           CALL XERBLA('DSYMV ',INFO)
 206           RETURN
 207       END IF
 208 *
 209 *     Quick return if possible.
 210 *
 211       IF ((N.EQ.0) .OR. ((ALPHA.EQ.ZERO).AND. (BETA.EQ.ONE))) RETURN
 212 *
 213 *     Set up the start points in  X  and  Y.
 214 *
 215       IF (INCX.GT.0) THEN
 216           KX = 1
 217       ELSE
 218           KX = 1 - (N-1)*INCX
 219       END IF
 220       IF (INCY.GT.0) THEN
 221           KY = 1
 222       ELSE
 223           KY = 1 - (N-1)*INCY
 224       END IF
 225 *
 226 *     Start the operations. In this version the elements of A are
 227 *     accessed sequentially with one pass through the triangular part
 228 *     of A.
 229 *
 230 *     First form  y := beta*y.
 231 *
 232       IF (BETA.NE.ONE) THEN
 233           IF (INCY.EQ.1) THEN
 234               IF (BETA.EQ.ZERO) THEN
 235                   DO 10 I = 1,N
 236                       Y(I) = ZERO
 237    10             CONTINUE
 238               ELSE
 239                   DO 20 I = 1,N
 240                       Y(I) = BETA*Y(I)
 241    20             CONTINUE
 242               END IF
 243           ELSE
 244               IY = KY
 245               IF (BETA.EQ.ZERO) THEN
 246                   DO 30 I = 1,N
 247                       Y(IY) = ZERO
 248                       IY = IY + INCY
 249    30             CONTINUE
 250               ELSE
 251                   DO 40 I = 1,N
 252                       Y(IY) = BETA*Y(IY)
 253                       IY = IY + INCY
 254    40             CONTINUE
 255               END IF
 256           END IF
 257       END IF
 258       IF (ALPHA.EQ.ZERO) RETURN
 259       IF (LSAME(UPLO,'U')) THEN
 260 *
 261 *        Form  y  when A is stored in upper triangle.
 262 *
 263           IF ((INCX.EQ.1) .AND. (INCY.EQ.1)) THEN
 264               DO 60 J = 1,N
 265                   TEMP1 = ALPHA*X(J)
 266                   TEMP2 = ZERO
 267                   DO 50 I = 1,J - 1
 268                       Y(I) = Y(I) + TEMP1*A(I,J)
 269                       TEMP2 = TEMP2 + A(I,J)*X(I)
 270    50             CONTINUE
 271                   Y(J) = Y(J) + TEMP1*A(J,J) + ALPHA*TEMP2
 272    60         CONTINUE
 273           ELSE
 274               JX = KX
 275               JY = KY
 276               DO 80 J = 1,N
 277                   TEMP1 = ALPHA*X(JX)
 278                   TEMP2 = ZERO
 279                   IX = KX
 280                   IY = KY
 281                   DO 70 I = 1,J - 1
 282                       Y(IY) = Y(IY) + TEMP1*A(I,J)
 283                       TEMP2 = TEMP2 + A(I,J)*X(IX)
 284                       IX = IX + INCX
 285                       IY = IY + INCY
 286    70             CONTINUE
 287                   Y(JY) = Y(JY) + TEMP1*A(J,J) + ALPHA*TEMP2
 288                   JX = JX + INCX
 289                   JY = JY + INCY
 290    80         CONTINUE
 291           END IF
 292       ELSE
 293 *
 294 *        Form  y  when A is stored in lower triangle.
 295 *
 296           IF ((INCX.EQ.1) .AND. (INCY.EQ.1)) THEN
 297               DO 100 J = 1,N
 298                   TEMP1 = ALPHA*X(J)
 299                   TEMP2 = ZERO
 300                   Y(J) = Y(J) + TEMP1*A(J,J)
 301                   DO 90 I = J + 1,N
 302                       Y(I) = Y(I) + TEMP1*A(I,J)
 303                       TEMP2 = TEMP2 + A(I,J)*X(I)
 304    90             CONTINUE
 305                   Y(J) = Y(J) + ALPHA*TEMP2
 306   100         CONTINUE
 307           ELSE
 308               JX = KX
 309               JY = KY
 310               DO 120 J = 1,N
 311                   TEMP1 = ALPHA*X(JX)
 312                   TEMP2 = ZERO
 313                   Y(JY) = Y(JY) + TEMP1*A(J,J)
 314                   IX = JX
 315                   IY = JY
 316                   DO 110 I = J + 1,N
 317                       IX = IX + INCX
 318                       IY = IY + INCY
 319                       Y(IY) = Y(IY) + TEMP1*A(I,J)
 320                       TEMP2 = TEMP2 + A(I,J)*X(IX)
 321   110             CONTINUE
 322                   Y(JY) = Y(JY) + ALPHA*TEMP2
 323                   JX = JX + INCX
 324                   JY = JY + INCY
 325   120         CONTINUE
 326           END IF
 327       END IF
 328 *
 329       RETURN
 330 *
 331 *     End of DSYMV .
 332 *
 333       END