share/lapack/blas/fortran/zgemv.f

   1       SUBROUTINE ZGEMV ( TRANS, M, N, ALPHA, A, LDA, X, INCX,
   2      $                   BETA, Y, INCY )
   3 *     .. Scalar Arguments ..
   4       COMPLEX*16         ALPHA, BETA
   5       INTEGER            INCX, INCY, LDA, M, N
   6       CHARACTER*1        TRANS
   7 *     .. Array Arguments ..
   8       COMPLEX*16         A( LDA, * ), X( * ), Y( * )
   9 *     ..
  10 *
  11 *  Purpose
  12 *  =======
  13 *
  14 *  ZGEMV  performs one of the matrix-vector operations
  15 *
  16 *     y := alpha*A*x + beta*y,   or   y := alpha*A'*x + beta*y,   or
  17 *
  18 *     y := alpha*conjg( A' )*x + beta*y,
  19 *
  20 *  where alpha and beta are scalars, x and y are vectors and A is an
  21 *  m by n matrix.
  22 *
  23 *  Parameters
  24 *  ==========
  25 *
  26 *  TRANS  - CHARACTER*1.
  27 *           On entry, TRANS specifies the operation to be performed as
  28 *           follows:
  29 *
  30 *              TRANS = 'N' or 'n'   y := alpha*A*x + beta*y.
  31 *
  32 *              TRANS = 'T' or 't'   y := alpha*A'*x + beta*y.
  33 *
  34 *              TRANS = 'C' or 'c'   y := alpha*conjg( A' )*x + beta*y.
  35 *
  36 *           Unchanged on exit.
  37 *
  38 *  M      - INTEGER.
  39 *           On entry, M specifies the number of rows of the matrix A.
  40 *           M must be at least zero.
  41 *           Unchanged on exit.
  42 *
  43 *  N      - INTEGER.
  44 *           On entry, N specifies the number of columns of the matrix A.
  45 *           N must be at least zero.
  46 *           Unchanged on exit.
  47 *
  48 *  ALPHA  - COMPLEX*16      .
  49 *           On entry, ALPHA specifies the scalar alpha.
  50 *           Unchanged on exit.
  51 *
  52 *  A      - COMPLEX*16       array of DIMENSION ( LDA, n ).
  53 *           Before entry, the leading m by n part of the array A must
  54 *           contain the matrix of coefficients.
  55 *           Unchanged on exit.
  56 *
  57 *  LDA    - INTEGER.
  58 *           On entry, LDA specifies the first dimension of A as declared
  59 *           in the calling (sub) program. LDA must be at least
  60 *           max( 1, m ).
  61 *           Unchanged on exit.
  62 *
  63 *  X      - COMPLEX*16       array of DIMENSION at least
  64 *           ( 1 + ( n - 1 )*abs( INCX ) ) when TRANS = 'N' or 'n'
  65 *           and at least
  66 *           ( 1 + ( m - 1 )*abs( INCX ) ) otherwise.
  67 *           Before entry, the incremented array X must contain the
  68 *           vector x.
  69 *           Unchanged on exit.
  70 *
  71 *  INCX   - INTEGER.
  72 *           On entry, INCX specifies the increment for the elements of
  73 *           X. INCX must not be zero.
  74 *           Unchanged on exit.
  75 *
  76 *  BETA   - COMPLEX*16      .
  77 *           On entry, BETA specifies the scalar beta. When BETA is
  78 *           supplied as zero then Y need not be set on input.
  79 *           Unchanged on exit.
  80 *
  81 *  Y      - COMPLEX*16       array of DIMENSION at least
  82 *           ( 1 + ( m - 1 )*abs( INCY ) ) when TRANS = 'N' or 'n'
  83 *           and at least
  84 *           ( 1 + ( n - 1 )*abs( INCY ) ) otherwise.
  85 *           Before entry with BETA non-zero, the incremented array Y
  86 *           must contain the vector y. On exit, Y is overwritten by the
  87 *           updated vector y.
  88 *
  89 *  INCY   - INTEGER.
  90 *           On entry, INCY specifies the increment for the elements of
  91 *           Y. INCY must not be zero.
  92 *           Unchanged on exit.
  93 *
  94 *
  95 *  Level 2 Blas routine.
  96 *
  97 *  -- Written on 22-October-1986.
  98 *     Jack Dongarra, Argonne National Lab.
  99 *     Jeremy Du Croz, Nag Central Office.
 100 *     Sven Hammarling, Nag Central Office.
 101 *     Richard Hanson, Sandia National Labs.
 102 *
 103 *
 104 *     .. Parameters ..
 105       COMPLEX*16         ONE
 106       PARAMETER        ( ONE  = ( 1.0D+0, 0.0D+0 ) )
 107       COMPLEX*16         ZERO
 108       PARAMETER        ( ZERO = ( 0.0D+0, 0.0D+0 ) )
 109 *     .. Local Scalars ..
 110       COMPLEX*16         TEMP
 111       INTEGER            I, INFO, IX, IY, J, JX, JY, KX, KY, LENX, LENY
 112       LOGICAL            NOCONJ
 113 *     .. External Functions ..
 114       LOGICAL            LSAME
 115       EXTERNAL           LSAME
 116 *     .. External Subroutines ..
 117       EXTERNAL           XERBLA
 118 *     .. Intrinsic Functions ..
 119       INTRINSIC          DCONJG, MAX
 120 *     ..
 121 *     .. Executable Statements ..
 122 *
 123 *     Test the input parameters.
 124 *
 125       INFO = 0
 126       IF     ( .NOT.LSAME( TRANS, 'N' ).AND.
 127      $         .NOT.LSAME( TRANS, 'T' ).AND.
 128      $         .NOT.LSAME( TRANS, 'C' )      )THEN
 129          INFO = 1
 130       ELSE IF( M.LT.0 )THEN
 131          INFO = 2
 132       ELSE IF( N.LT.0 )THEN
 133          INFO = 3
 134       ELSE IF( LDA.LT.MAX( 1, M ) )THEN
 135          INFO = 6
 136       ELSE IF( INCX.EQ.0 )THEN
 137          INFO = 8
 138       ELSE IF( INCY.EQ.0 )THEN
 139          INFO = 11
 140       END IF
 141       IF( INFO.NE.0 )THEN
 142          CALL XERBLA( 'ZGEMV ', INFO )
 143          RETURN
 144       END IF
 145 *
 146 *     Quick return if possible.
 147 *
 148       IF( ( M.EQ.0 ).OR.( N.EQ.0 ).OR.
 149      $    ( ( ALPHA.EQ.ZERO ).AND.( BETA.EQ.ONE ) ) )
 150      $   RETURN
 151 *
 152       NOCONJ = LSAME( TRANS, 'T' )
 153 *
 154 *     Set  LENX  and  LENY, the lengths of the vectors x and y, and set
 155 *     up the start points in  X  and  Y.
 156 *
 157       IF( LSAME( TRANS, 'N' ) )THEN
 158          LENX = N
 159          LENY = M
 160       ELSE
 161          LENX = M
 162          LENY = N
 163       END IF
 164       IF( INCX.GT.0 )THEN
 165          KX = 1
 166       ELSE
 167          KX = 1 - ( LENX - 1 )*INCX
 168       END IF
 169       IF( INCY.GT.0 )THEN
 170          KY = 1
 171       ELSE
 172          KY = 1 - ( LENY - 1 )*INCY
 173       END IF
 174 *
 175 *     Start the operations. In this version the elements of A are
 176 *     accessed sequentially with one pass through A.
 177 *
 178 *     First form  y := beta*y.
 179 *
 180       IF( BETA.NE.ONE )THEN
 181          IF( INCY.EQ.1 )THEN
 182             IF( BETA.EQ.ZERO )THEN
 183                DO 10, I = 1, LENY
 184                   Y( I ) = ZERO
 185    10          CONTINUE
 186             ELSE
 187                DO 20, I = 1, LENY
 188                   Y( I ) = BETA*Y( I )
 189    20          CONTINUE
 190             END IF
 191          ELSE
 192             IY = KY
 193             IF( BETA.EQ.ZERO )THEN
 194                DO 30, I = 1, LENY
 195                   Y( IY ) = ZERO
 196                   IY      = IY   + INCY
 197    30          CONTINUE
 198             ELSE
 199                DO 40, I = 1, LENY
 200                   Y( IY ) = BETA*Y( IY )
 201                   IY      = IY           + INCY
 202    40          CONTINUE
 203             END IF
 204          END IF
 205       END IF
 206       IF( ALPHA.EQ.ZERO )
 207      $   RETURN
 208       IF( LSAME( TRANS, 'N' ) )THEN
 209 *
 210 *        Form  y := alpha*A*x + y.
 211 *
 212          JX = KX
 213          IF( INCY.EQ.1 )THEN
 214             DO 60, J = 1, N
 215                IF( X( JX ).NE.ZERO )THEN
 216                   TEMP = ALPHA*X( JX )
 217                   DO 50, I = 1, M
 218                      Y( I ) = Y( I ) + TEMP*A( I, J )
 219    50             CONTINUE
 220                END IF
 221                JX = JX + INCX
 222    60       CONTINUE
 223          ELSE
 224             DO 80, J = 1, N
 225                IF( X( JX ).NE.ZERO )THEN
 226                   TEMP = ALPHA*X( JX )
 227                   IY   = KY
 228                   DO 70, I = 1, M
 229                      Y( IY ) = Y( IY ) + TEMP*A( I, J )
 230                      IY      = IY      + INCY
 231    70             CONTINUE
 232                END IF
 233                JX = JX + INCX
 234    80       CONTINUE
 235          END IF
 236       ELSE
 237 *
 238 *        Form  y := alpha*A'*x + y  or  y := alpha*conjg( A' )*x + y.
 239 *
 240          JY = KY
 241          IF( INCX.EQ.1 )THEN
 242             DO 110, J = 1, N
 243                TEMP = ZERO
 244                IF( NOCONJ )THEN
 245                   DO 90, I = 1, M
 246                      TEMP = TEMP + A( I, J )*X( I )
 247    90             CONTINUE
 248                ELSE
 249                   DO 100, I = 1, M
 250                      TEMP = TEMP + DCONJG( A( I, J ) )*X( I )
 251   100             CONTINUE
 252                END IF
 253                Y( JY ) = Y( JY ) + ALPHA*TEMP
 254                JY      = JY      + INCY
 255   110       CONTINUE
 256          ELSE
 257             DO 140, J = 1, N
 258                TEMP = ZERO
 259                IX   = KX
 260                IF( NOCONJ )THEN
 261                   DO 120, I = 1, M
 262                      TEMP = TEMP + A( I, J )*X( IX )
 263                      IX   = IX   + INCX
 264   120             CONTINUE
 265                ELSE
 266                   DO 130, I = 1, M
 267                      TEMP = TEMP + DCONJG( A( I, J ) )*X( IX )
 268                      IX   = IX   + INCX
 269   130             CONTINUE
 270                END IF
 271                Y( JY ) = Y( JY ) + ALPHA*TEMP
 272                JY      = JY      + INCY
 273   140       CONTINUE
 274          END IF
 275       END IF
 276 *
 277       RETURN
 278 *
 279 *     End of ZGEMV .
 280 *
 281       END