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