1 /* dgemv.f -- translated by f2c (version 19991025). You must link the resulting object file with the libraries: -lf2c -lm (in that order) */
2 #include "FLA_f2c.h"
3 /* Subroutine */
dgemv_(char * trans,integer * m,integer * n,doublereal * alpha,doublereal * a,integer * lda,doublereal * x,integer * incx,doublereal * beta,doublereal * y,integer * incy)4 int dgemv_(char *trans, integer *m, integer *n, doublereal * alpha, doublereal *a, integer *lda, doublereal *x, integer *incx, doublereal *beta, doublereal *y, integer *incy)
5 {
6 /* System generated locals */
7 integer a_dim1, a_offset, i__1, i__2;
8 /* Local variables */
9 integer info;
10 doublereal temp;
11 integer lenx, leny, i__, j;
12 extern logical lsame_(char *, char *);
13 integer ix, iy, jx, jy, kx, ky;
14 extern /* Subroutine */
15 int xerbla_(char *, integer *);
16 /* .. Scalar Arguments .. */
17 /* .. Array Arguments .. */
18 /* .. */
19 /* Purpose */
20 /* ======= */
21 /* DGEMV performs one of the matrix-vector operations */
22 /* y := alpha*A*x + beta*y, or y := alpha*A'*x + beta*y, */
23 /* where alpha and beta are scalars, x and y are vectors and A is an */
24 /* m by n matrix. */
25 /* Parameters */
26 /* ========== */
27 /* TRANS - CHARACTER*1. */
28 /* On entry, TRANS specifies the operation to be performed as */
29 /* follows: */
30 /* TRANS = 'N' or 'n' y := alpha*A*x + beta*y. */
31 /* TRANS = 'T' or 't' y := alpha*A'*x + beta*y. */
32 /* TRANS = 'C' or 'c' y := alpha*A'*x + beta*y. */
33 /* Unchanged on exit. */
34 /* M - INTEGER. */
35 /* On entry, M specifies the number of rows of the matrix A. */
36 /* M must be at least zero. */
37 /* Unchanged on exit. */
38 /* N - INTEGER. */
39 /* On entry, N specifies the number of columns of the matrix A. */
40 /* N must be at least zero. */
41 /* Unchanged on exit. */
42 /* ALPHA - DOUBLE PRECISION. */
43 /* On entry, ALPHA specifies the scalar alpha. */
44 /* Unchanged on exit. */
45 /* A - DOUBLE PRECISION array of DIMENSION ( LDA, n ). */
46 /* Before entry, the leading m by n part of the array A must */
47 /* contain the matrix of coefficients. */
48 /* Unchanged on exit. */
49 /* LDA - INTEGER. */
50 /* On entry, LDA specifies the first dimension of A as declared */
51 /* in the calling (sub) program. LDA must be at least */
52 /* max( 1, m ). */
53 /* Unchanged on exit. */
54 /* X - DOUBLE PRECISION array of DIMENSION at least */
55 /* ( 1 + ( n - 1 )*f2c_abs( INCX ) ) when TRANS = 'N' or 'n' */
56 /* and at least */
57 /* ( 1 + ( m - 1 )*f2c_abs( INCX ) ) otherwise. */
58 /* Before entry, the incremented array X must contain the */
59 /* vector x. */
60 /* Unchanged on exit. */
61 /* INCX - INTEGER. */
62 /* On entry, INCX specifies the increment for the elements of */
63 /* X. INCX must not be zero. */
64 /* Unchanged on exit. */
65 /* BETA - DOUBLE PRECISION. */
66 /* On entry, BETA specifies the scalar beta. When BETA is */
67 /* supplied as zero then Y need not be set on input. */
68 /* Unchanged on exit. */
69 /* Y - DOUBLE PRECISION array of DIMENSION at least */
70 /* ( 1 + ( m - 1 )*f2c_abs( INCY ) ) when TRANS = 'N' or 'n' */
71 /* and at least */
72 /* ( 1 + ( n - 1 )*f2c_abs( INCY ) ) otherwise. */
73 /* Before entry with BETA non-zero, the incremented array Y */
74 /* must contain the vector y. On exit, Y is overwritten by the */
75 /* updated vector y. */
76 /* INCY - INTEGER. */
77 /* On entry, INCY specifies the increment for the elements of */
78 /* Y. INCY must not be zero. */
79 /* Unchanged on exit. */
80 /* Level 2 Blas routine. */
81 /* -- Written on 22-October-1986. */
82 /* Jack Dongarra, Argonne National Lab. */
83 /* Jeremy Du Croz, Nag Central Office. */
84 /* Sven Hammarling, Nag Central Office. */
85 /* Richard Hanson, Sandia National Labs. */
86 /* .. Parameters .. */
87 /* .. Local Scalars .. */
88 /* .. External Functions .. */
89 /* .. External Subroutines .. */
90 /* .. Intrinsic Functions .. */
91 /* .. */
92 /* .. Executable Statements .. */
93 /* Test the input parameters. */
94 /* Parameter adjustments */
95 a_dim1 = *lda;
96 a_offset = 1 + a_dim1 * 1;
97 a -= a_offset;
98 --x;
99 --y;
100 /* Function Body */
101 info = 0;
102 if (! lsame_(trans, "N") && ! lsame_(trans, "T") && ! lsame_(trans, "C") )
103 {
104 info = 1;
105 }
106 else if (*m < 0)
107 {
108 info = 2;
109 }
110 else if (*n < 0)
111 {
112 info = 3;
113 }
114 else if (*lda < max(1,*m))
115 {
116 info = 6;
117 }
118 else if (*incx == 0)
119 {
120 info = 8;
121 }
122 else if (*incy == 0)
123 {
124 info = 11;
125 }
126 if (info != 0)
127 {
128 xerbla_("DGEMV ", &info);
129 return 0;
130 }
131 /* Quick return if possible. */
132 if (*m == 0 || *n == 0 || *alpha == 0. && *beta == 1.)
133 {
134 return 0;
135 }
136 /* Set LENX and LENY, the lengths of the vectors x and y, and set */
137 /* up the start points in X and Y. */
138 if (lsame_(trans, "N"))
139 {
140 lenx = *n;
141 leny = *m;
142 }
143 else
144 {
145 lenx = *m;
146 leny = *n;
147 }
148 if (*incx > 0)
149 {
150 kx = 1;
151 }
152 else
153 {
154 kx = 1 - (lenx - 1) * *incx;
155 }
156 if (*incy > 0)
157 {
158 ky = 1;
159 }
160 else
161 {
162 ky = 1 - (leny - 1) * *incy;
163 }
164 /* Start the operations. In this version the elements of A are */
165 /* accessed sequentially with one pass through A. */
166 /* First form y := beta*y. */
167 if (*beta != 1.)
168 {
169 if (*incy == 1)
170 {
171 if (*beta == 0.)
172 {
173 i__1 = leny;
174 for (i__ = 1;
175 i__ <= i__1;
176 ++i__)
177 {
178 y[i__] = 0.;
179 /* L10: */
180 }
181 }
182 else
183 {
184 i__1 = leny;
185 for (i__ = 1;
186 i__ <= i__1;
187 ++i__)
188 {
189 y[i__] = *beta * y[i__];
190 /* L20: */
191 }
192 }
193 }
194 else
195 {
196 iy = ky;
197 if (*beta == 0.)
198 {
199 i__1 = leny;
200 for (i__ = 1;
201 i__ <= i__1;
202 ++i__)
203 {
204 y[iy] = 0.;
205 iy += *incy;
206 /* L30: */
207 }
208 }
209 else
210 {
211 i__1 = leny;
212 for (i__ = 1;
213 i__ <= i__1;
214 ++i__)
215 {
216 y[iy] = *beta * y[iy];
217 iy += *incy;
218 /* L40: */
219 }
220 }
221 }
222 }
223 if (*alpha == 0.)
224 {
225 return 0;
226 }
227 if (lsame_(trans, "N"))
228 {
229 /* Form y := alpha*A*x + y. */
230 jx = kx;
231 if (*incy == 1)
232 {
233 i__1 = *n;
234 for (j = 1;
235 j <= i__1;
236 ++j)
237 {
238 if (x[jx] != 0.)
239 {
240 temp = *alpha * x[jx];
241 i__2 = *m;
242 for (i__ = 1;
243 i__ <= i__2;
244 ++i__)
245 {
246 y[i__] += temp * a[i__ + j * a_dim1];
247 /* L50: */
248 }
249 }
250 jx += *incx;
251 /* L60: */
252 }
253 }
254 else
255 {
256 i__1 = *n;
257 for (j = 1;
258 j <= i__1;
259 ++j)
260 {
261 if (x[jx] != 0.)
262 {
263 temp = *alpha * x[jx];
264 iy = ky;
265 i__2 = *m;
266 for (i__ = 1;
267 i__ <= i__2;
268 ++i__)
269 {
270 y[iy] += temp * a[i__ + j * a_dim1];
271 iy += *incy;
272 /* L70: */
273 }
274 }
275 jx += *incx;
276 /* L80: */
277 }
278 }
279 }
280 else
281 {
282 /* Form y := alpha*A'*x + y. */
283 jy = ky;
284 if (*incx == 1)
285 {
286 i__1 = *n;
287 for (j = 1;
288 j <= i__1;
289 ++j)
290 {
291 temp = 0.;
292 i__2 = *m;
293 for (i__ = 1;
294 i__ <= i__2;
295 ++i__)
296 {
297 temp += a[i__ + j * a_dim1] * x[i__];
298 /* L90: */
299 }
300 y[jy] += *alpha * temp;
301 jy += *incy;
302 /* L100: */
303 }
304 }
305 else
306 {
307 i__1 = *n;
308 for (j = 1;
309 j <= i__1;
310 ++j)
311 {
312 temp = 0.;
313 ix = kx;
314 i__2 = *m;
315 for (i__ = 1;
316 i__ <= i__2;
317 ++i__)
318 {
319 temp += a[i__ + j * a_dim1] * x[ix];
320 ix += *incx;
321 /* L110: */
322 }
323 y[jy] += *alpha * temp;
324 jy += *incy;
325 /* L120: */
326 }
327 }
328 }
329 return 0;
330 /* End of DGEMV . */
331 }
332 /* dgemv_ */
333
334