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