1 /*
2
3 BLIS
4 An object-based framework for developing high-performance BLAS-like
5 libraries.
6
7 Copyright (C) 2014, The University of Texas at Austin
8
9 Redistribution and use in source and binary forms, with or without
10 modification, are permitted provided that the following conditions are
11 met:
12 - Redistributions of source code must retain the above copyright
13 notice, this list of conditions and the following disclaimer.
14 - Redistributions in binary form must reproduce the above copyright
15 notice, this list of conditions and the following disclaimer in the
16 documentation and/or other materials provided with the distribution.
17 - Neither the name(s) of the copyright holder(s) nor the names of its
18 contributors may be used to endorse or promote products derived
19 from this software without specific prior written permission.
20
21 THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
22 "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
23 LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
24 A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT
25 HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
26 SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
27 LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
28 DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
29 THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
30 (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
31 OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
32
33 */
34
35 #include "blis.h"
36
37 #ifdef BLIS_ENABLE_BLAS
38
39 /* cgbmv.f -- translated by f2c (version 19991025).
40 You must link the resulting object file with the libraries:
41 -lf2c -lm (in that order)
42 */
43
PASTEF77(c,gbmv)44 /* Subroutine */ int PASTEF77(c,gbmv)(const bla_character *trans, const bla_integer *m, const bla_integer *n, const bla_integer *kl, const bla_integer *ku, const bla_scomplex *alpha, const bla_scomplex *a, const bla_integer *lda, const bla_scomplex *x, const bla_integer *incx, const bla_scomplex *beta, bla_scomplex *y, const bla_integer *incy)
45 {
46 /* System generated locals */
47 bla_integer a_dim1, a_offset, i__1, i__2, i__3, i__4, i__5, i__6;
48 bla_scomplex q__1, q__2, q__3;
49
50 /* Builtin functions */
51 //void bla_r_cnjg(bla_scomplex *, bla_scomplex *);
52
53 /* Local variables */
54 bla_integer info;
55 bla_scomplex temp;
56 bla_integer lenx, leny, i__, j, k;
57 //extern bla_logical PASTEF770(lsame)(bla_character *, bla_character *, ftnlen, ftnlen);
58 bla_integer ix, iy, jx, jy, kx, ky;
59 //extern /* Subroutine */ int PASTEF770(xerbla)(bla_character *, bla_integer *, ftnlen);
60 bla_logical noconj;
61 bla_integer kup1;
62
63 /* .. Scalar Arguments .. */
64 /* .. Array Arguments .. */
65 /* .. */
66
67 /* Purpose */
68 /* ======= */
69
70 /* CGBMV performs one of the matrix-vector operations */
71
72 /* y := alpha*A*x + beta*y, or y := alpha*A'*x + beta*y, or */
73
74 /* y := alpha*conjg( A' )*x + beta*y, */
75
76 /* where alpha and beta are scalars, x and y are vectors and A is an */
77 /* m by n band matrix, with kl sub-diagonals and ku super-diagonals. */
78
79 /* Parameters */
80 /* ========== */
81
82 /* TRANS - CHARACTER*1. */
83 /* On entry, TRANS specifies the operation to be performed as */
84 /* follows: */
85
86 /* TRANS = 'N' or 'n' y := alpha*A*x + beta*y. */
87
88 /* TRANS = 'T' or 't' y := alpha*A'*x + beta*y. */
89
90 /* TRANS = 'C' or 'c' y := alpha*conjg( A' )*x + beta*y. */
91
92 /* Unchanged on exit. */
93
94 /* M - INTEGER. */
95 /* On entry, M specifies the number of rows of the matrix A. */
96 /* M must be at least zero. */
97 /* Unchanged on exit. */
98
99 /* N - INTEGER. */
100 /* On entry, N specifies the number of columns of the matrix A. */
101 /* N must be at least zero. */
102 /* Unchanged on exit. */
103
104 /* KL - INTEGER. */
105 /* On entry, KL specifies the number of sub-diagonals of the */
106 /* matrix A. KL must satisfy 0 .le. KL. */
107 /* Unchanged on exit. */
108
109 /* KU - INTEGER. */
110 /* On entry, KU specifies the number of super-diagonals of the */
111 /* matrix A. KU must satisfy 0 .le. KU. */
112 /* Unchanged on exit. */
113
114 /* ALPHA - COMPLEX . */
115 /* On entry, ALPHA specifies the scalar alpha. */
116 /* Unchanged on exit. */
117
118 /* A - COMPLEX array of DIMENSION ( LDA, n ). */
119 /* Before entry, the leading ( kl + ku + 1 ) by n part of the */
120 /* array A must contain the matrix of coefficients, supplied */
121 /* column by column, with the leading diagonal of the matrix in */
122 /* row ( ku + 1 ) of the array, the first super-diagonal */
123 /* starting at position 2 in row ku, the first sub-diagonal */
124 /* starting at position 1 in row ( ku + 2 ), and so on. */
125 /* Elements in the array A that do not correspond to elements */
126 /* in the band matrix (such as the top left ku by ku triangle) */
127 /* are not referenced. */
128 /* The following program segment will transfer a band matrix */
129 /* from conventional full matrix storage to band storage: */
130
131 /* DO 20, J = 1, N */
132 /* K = KU + 1 - J */
133 /* DO 10, I = MAX( 1, J - KU ), MIN( M, J + KL ) */
134 /* A( K + I, J ) = matrix( I, J ) */
135 /* 10 CONTINUE */
136 /* 20 CONTINUE */
137
138 /* Unchanged on exit. */
139
140 /* LDA - INTEGER. */
141 /* On entry, LDA specifies the first dimension of A as declared */
142 /* in the calling (sub) program. LDA must be at least */
143 /* ( kl + ku + 1 ). */
144 /* Unchanged on exit. */
145
146 /* X - COMPLEX array of DIMENSION at least */
147 /* ( 1 + ( n - 1 )*abs( INCX ) ) when TRANS = 'N' or 'n' */
148 /* and at least */
149 /* ( 1 + ( m - 1 )*abs( INCX ) ) otherwise. */
150 /* Before entry, the incremented array X must contain the */
151 /* vector x. */
152 /* Unchanged on exit. */
153
154 /* INCX - INTEGER. */
155 /* On entry, INCX specifies the increment for the elements of */
156 /* X. INCX must not be zero. */
157 /* Unchanged on exit. */
158
159 /* BETA - COMPLEX . */
160 /* On entry, BETA specifies the scalar beta. When BETA is */
161 /* supplied as zero then Y need not be set on input. */
162 /* Unchanged on exit. */
163
164 /* Y - COMPLEX array of DIMENSION at least */
165 /* ( 1 + ( m - 1 )*abs( INCY ) ) when TRANS = 'N' or 'n' */
166 /* and at least */
167 /* ( 1 + ( n - 1 )*abs( INCY ) ) otherwise. */
168 /* Before entry, the incremented array Y must contain the */
169 /* vector y. On exit, Y is overwritten by the updated vector y. */
170
171
172 /* INCY - INTEGER. */
173 /* On entry, INCY specifies the increment for the elements of */
174 /* Y. INCY must not be zero. */
175 /* Unchanged on exit. */
176
177
178 /* Level 2 Blas routine. */
179
180 /* -- Written on 22-October-1986. */
181 /* Jack Dongarra, Argonne National Lab. */
182 /* Jeremy Du Croz, Nag Central Office. */
183 /* Sven Hammarling, Nag Central Office. */
184 /* Richard Hanson, Sandia National Labs. */
185
186
187 /* .. Parameters .. */
188 /* .. Local Scalars .. */
189 /* .. External Functions .. */
190 /* .. External Subroutines .. */
191 /* .. Intrinsic Functions .. */
192 /* .. */
193 /* .. Executable Statements .. */
194
195 /* Test the input parameters. */
196
197 /* Parameter adjustments */
198 a_dim1 = *lda;
199 a_offset = 1 + a_dim1 * 1;
200 a -= a_offset;
201 --x;
202 --y;
203
204 /* Function Body */
205 info = 0;
206 if (! PASTEF770(lsame)(trans, "N", (ftnlen)1, (ftnlen)1) && ! PASTEF770(lsame)(trans, "T", (
207 ftnlen)1, (ftnlen)1) && ! PASTEF770(lsame)(trans, "C", (ftnlen)1, (ftnlen)1)
208 ) {
209 info = 1;
210 } else if (*m < 0) {
211 info = 2;
212 } else if (*n < 0) {
213 info = 3;
214 } else if (*kl < 0) {
215 info = 4;
216 } else if (*ku < 0) {
217 info = 5;
218 } else if (*lda < *kl + *ku + 1) {
219 info = 8;
220 } else if (*incx == 0) {
221 info = 10;
222 } else if (*incy == 0) {
223 info = 13;
224 }
225 if (info != 0) {
226 PASTEF770(xerbla)("CGBMV ", &info, (ftnlen)6);
227 return 0;
228 }
229
230 /* Quick return if possible. */
231
232 if (*m == 0 || *n == 0 || (bli_creal(*alpha) == 0.f && bli_cimag(*alpha) == 0.f && (bli_creal(*beta)
233 == 1.f && bli_cimag(*beta) == 0.f))) {
234 return 0;
235 }
236
237 noconj = PASTEF770(lsame)(trans, "T", (ftnlen)1, (ftnlen)1);
238
239 /* Set LENX and LENY, the lengths of the vectors x and y, and set */
240 /* up the start points in X and Y. */
241
242 if (PASTEF770(lsame)(trans, "N", (ftnlen)1, (ftnlen)1)) {
243 lenx = *n;
244 leny = *m;
245 } else {
246 lenx = *m;
247 leny = *n;
248 }
249 if (*incx > 0) {
250 kx = 1;
251 } else {
252 kx = 1 - (lenx - 1) * *incx;
253 }
254 if (*incy > 0) {
255 ky = 1;
256 } else {
257 ky = 1 - (leny - 1) * *incy;
258 }
259
260 /* Start the operations. In this version the elements of A are */
261 /* accessed sequentially with one pass through the band part of A. */
262
263 /* First form y := beta*y. */
264
265 if (bli_creal(*beta) != 1.f || bli_cimag(*beta) != 0.f) {
266 if (*incy == 1) {
267 if (bli_creal(*beta) == 0.f && bli_cimag(*beta) == 0.f) {
268 i__1 = leny;
269 for (i__ = 1; i__ <= i__1; ++i__) {
270 i__2 = i__;
271 bli_csets( (0.f), (0.f), y[i__2] );
272 /* L10: */
273 }
274 } else {
275 i__1 = leny;
276 for (i__ = 1; i__ <= i__1; ++i__) {
277 i__2 = i__;
278 i__3 = i__;
279 bli_csets( (bli_creal(*beta) * bli_creal(y[i__3]) - bli_cimag(*beta) * bli_cimag(y[i__3])), (bli_creal(*beta) * bli_cimag(y[i__3]) + bli_cimag(*beta) * bli_creal(y[i__3])), q__1 );
280 bli_csets( (bli_creal(q__1)), (bli_cimag(q__1)), y[i__2] );
281 /* L20: */
282 }
283 }
284 } else {
285 iy = ky;
286 if (bli_creal(*beta) == 0.f && bli_cimag(*beta) == 0.f) {
287 i__1 = leny;
288 for (i__ = 1; i__ <= i__1; ++i__) {
289 i__2 = iy;
290 bli_csets( (0.f), (0.f), y[i__2] );
291 iy += *incy;
292 /* L30: */
293 }
294 } else {
295 i__1 = leny;
296 for (i__ = 1; i__ <= i__1; ++i__) {
297 i__2 = iy;
298 i__3 = iy;
299 bli_csets( (bli_creal(*beta) * bli_creal(y[i__3]) - bli_cimag(*beta) * bli_cimag(y[i__3])), (bli_creal(*beta) * bli_cimag(y[i__3]) + bli_cimag(*beta) * bli_creal(y[i__3])), q__1 );
300 bli_csets( (bli_creal(q__1)), (bli_cimag(q__1)), y[i__2] );
301 iy += *incy;
302 /* L40: */
303 }
304 }
305 }
306 }
307 if (bli_creal(*alpha) == 0.f && bli_cimag(*alpha) == 0.f) {
308 return 0;
309 }
310 kup1 = *ku + 1;
311 if (PASTEF770(lsame)(trans, "N", (ftnlen)1, (ftnlen)1)) {
312
313 /* Form y := alpha*A*x + y. */
314
315 jx = kx;
316 if (*incy == 1) {
317 i__1 = *n;
318 for (j = 1; j <= i__1; ++j) {
319 i__2 = jx;
320 if (bli_creal(x[i__2]) != 0.f || bli_cimag(x[i__2]) != 0.f) {
321 i__2 = jx;
322 bli_csets( (bli_creal(*alpha) * bli_creal(x[i__2]) - bli_cimag(*alpha) * bli_cimag(x[i__2])), (bli_creal(*alpha) * bli_cimag(x[i__2]) + bli_cimag(*alpha) * bli_creal(x[i__2])), q__1 );
323 bli_csets( (bli_creal(q__1)), (bli_cimag(q__1)), temp );
324 k = kup1 - j;
325 /* Computing MAX */
326 i__2 = 1, i__3 = j - *ku;
327 /* Computing MIN */
328 i__5 = *m, i__6 = j + *kl;
329 i__4 = f2c_min(i__5,i__6);
330 for (i__ = f2c_max(i__2,i__3); i__ <= i__4; ++i__) {
331 i__2 = i__;
332 i__3 = i__;
333 i__5 = k + i__ + j * a_dim1;
334 bli_csets( (bli_creal(temp) * bli_creal(a[i__5]) - bli_cimag(temp) * bli_cimag(a[i__5])), (bli_creal(temp) * bli_cimag(a[i__5]) + bli_cimag(temp) * bli_creal(a[i__5])), q__2 );
335 bli_csets( (bli_creal(y[i__3]) + bli_creal(q__2)), (bli_cimag(y[i__3]) + bli_cimag(q__2)), q__1 );
336 bli_csets( (bli_creal(q__1)), (bli_cimag(q__1)), y[i__2] );
337 /* L50: */
338 }
339 }
340 jx += *incx;
341 /* L60: */
342 }
343 } else {
344 i__1 = *n;
345 for (j = 1; j <= i__1; ++j) {
346 i__4 = jx;
347 if (bli_creal(x[i__4]) != 0.f || bli_cimag(x[i__4]) != 0.f) {
348 i__4 = jx;
349 bli_csets( (bli_creal(*alpha) * bli_creal(x[i__4]) - bli_cimag(*alpha) * bli_cimag(x[i__4])), (bli_creal(*alpha) * bli_cimag(x[i__4]) + bli_cimag(*alpha) * bli_creal(x[i__4])), q__1 );
350 bli_csets( (bli_creal(q__1)), (bli_cimag(q__1)), temp );
351 iy = ky;
352 k = kup1 - j;
353 /* Computing MAX */
354 i__4 = 1, i__2 = j - *ku;
355 /* Computing MIN */
356 i__5 = *m, i__6 = j + *kl;
357 i__3 = f2c_min(i__5,i__6);
358 for (i__ = f2c_max(i__4,i__2); i__ <= i__3; ++i__) {
359 i__4 = iy;
360 i__2 = iy;
361 i__5 = k + i__ + j * a_dim1;
362 bli_csets( (bli_creal(temp) * bli_creal(a[i__5]) - bli_cimag(temp) * bli_cimag(a[i__5])), (bli_creal(temp) * bli_cimag(a[i__5]) + bli_cimag(temp) * bli_creal(a[i__5])), q__2 );
363 bli_csets( (bli_creal(y[i__2]) + bli_creal(q__2)), (bli_cimag(y[i__2]) + bli_cimag(q__2)), q__1 );
364 bli_csets( (bli_creal(q__1)), (bli_cimag(q__1)), y[i__4] );
365 iy += *incy;
366 /* L70: */
367 }
368 }
369 jx += *incx;
370 if (j > *ku) {
371 ky += *incy;
372 }
373 /* L80: */
374 }
375 }
376 } else {
377
378 /* Form y := alpha*A'*x + y or y := alpha*conjg( A' )*x + y. */
379
380 jy = ky;
381 if (*incx == 1) {
382 i__1 = *n;
383 for (j = 1; j <= i__1; ++j) {
384 bli_csets( (0.f), (0.f), temp );
385 k = kup1 - j;
386 if (noconj) {
387 /* Computing MAX */
388 i__3 = 1, i__4 = j - *ku;
389 /* Computing MIN */
390 i__5 = *m, i__6 = j + *kl;
391 i__2 = f2c_min(i__5,i__6);
392 for (i__ = f2c_max(i__3,i__4); i__ <= i__2; ++i__) {
393 i__3 = k + i__ + j * a_dim1;
394 i__4 = i__;
395 bli_csets( (bli_creal(a[i__3]) * bli_creal(x[i__4]) - bli_cimag(a[i__3]) * bli_cimag(x[i__4])), (bli_creal(a[i__3]) * bli_cimag(x[i__4]) + bli_cimag(a[i__3]) * bli_creal(x[i__4])), q__2 );
396 bli_csets( (bli_creal(temp) + bli_creal(q__2)), (bli_cimag(temp) + bli_cimag(q__2)), q__1 );
397 bli_csets( (bli_creal(q__1)), (bli_cimag(q__1)), temp );
398 /* L90: */
399 }
400 } else {
401 /* Computing MAX */
402 i__2 = 1, i__3 = j - *ku;
403 /* Computing MIN */
404 i__5 = *m, i__6 = j + *kl;
405 i__4 = f2c_min(i__5,i__6);
406 for (i__ = f2c_max(i__2,i__3); i__ <= i__4; ++i__) {
407 bla_r_cnjg(&q__3, &a[k + i__ + j * a_dim1]);
408 i__2 = i__;
409 bli_csets( (bli_creal(q__3) * bli_creal(x[i__2]) - bli_cimag(q__3) * bli_cimag(x[i__2])), (bli_creal(q__3) * bli_cimag(x[i__2]) + bli_cimag(q__3) * bli_creal(x[i__2])), q__2 );
410 bli_csets( (bli_creal(temp) + bli_creal(q__2)), (bli_cimag(temp) + bli_cimag(q__2)), q__1 );
411 bli_csets( (bli_creal(q__1)), (bli_cimag(q__1)), temp );
412 /* L100: */
413 }
414 }
415 i__4 = jy;
416 i__2 = jy;
417 bli_csets( (bli_creal(*alpha) * bli_creal(temp) - bli_cimag(*alpha) * bli_cimag(temp)), (bli_creal(*alpha) * bli_cimag(temp) + bli_cimag(*alpha) * bli_creal(temp)), q__2 );
418 bli_csets( (bli_creal(y[i__2]) + bli_creal(q__2)), (bli_cimag(y[i__2]) + bli_cimag(q__2)), q__1 );
419 bli_csets( (bli_creal(q__1)), (bli_cimag(q__1)), y[i__4] );
420 jy += *incy;
421 /* L110: */
422 }
423 } else {
424 i__1 = *n;
425 for (j = 1; j <= i__1; ++j) {
426 bli_csets( (0.f), (0.f), temp );
427 ix = kx;
428 k = kup1 - j;
429 if (noconj) {
430 /* Computing MAX */
431 i__4 = 1, i__2 = j - *ku;
432 /* Computing MIN */
433 i__5 = *m, i__6 = j + *kl;
434 i__3 = f2c_min(i__5,i__6);
435 for (i__ = f2c_max(i__4,i__2); i__ <= i__3; ++i__) {
436 i__4 = k + i__ + j * a_dim1;
437 i__2 = ix;
438 bli_csets( (bli_creal(a[i__4]) * bli_creal(x[i__2]) - bli_cimag(a[i__4]) * bli_cimag(x[i__2])), (bli_creal(a[i__4]) * bli_cimag(x[i__2]) + bli_cimag(a[i__4]) * bli_creal(x[i__2])), q__2 );
439 bli_csets( (bli_creal(temp) + bli_creal(q__2)), (bli_cimag(temp) + bli_cimag(q__2)), q__1 );
440 bli_csets( (bli_creal(q__1)), (bli_cimag(q__1)), temp );
441 ix += *incx;
442 /* L120: */
443 }
444 } else {
445 /* Computing MAX */
446 i__3 = 1, i__4 = j - *ku;
447 /* Computing MIN */
448 i__5 = *m, i__6 = j + *kl;
449 i__2 = f2c_min(i__5,i__6);
450 for (i__ = f2c_max(i__3,i__4); i__ <= i__2; ++i__) {
451 bla_r_cnjg(&q__3, &a[k + i__ + j * a_dim1]);
452 i__3 = ix;
453 bli_csets( (bli_creal(q__3) * bli_creal(x[i__3]) - bli_cimag(q__3) * bli_cimag(x[i__3])), (bli_creal(q__3) * bli_cimag(x[i__3]) + bli_cimag(q__3) * bli_creal(x[i__3])), q__2 );
454 bli_csets( (bli_creal(temp) + bli_creal(q__2)), (bli_cimag(temp) + bli_cimag(q__2)), q__1 );
455 bli_csets( (bli_creal(q__1)), (bli_cimag(q__1)), temp );
456 ix += *incx;
457 /* L130: */
458 }
459 }
460 i__2 = jy;
461 i__3 = jy;
462 bli_csets( (bli_creal(*alpha) * bli_creal(temp) - bli_cimag(*alpha) * bli_cimag(temp)), (bli_creal(*alpha) * bli_cimag(temp) + bli_cimag(*alpha) * bli_creal(temp)), q__2 );
463 bli_csets( (bli_creal(y[i__3]) + bli_creal(q__2)), (bli_cimag(y[i__3]) + bli_cimag(q__2)), q__1 );
464 bli_csets( (bli_creal(q__1)), (bli_cimag(q__1)), y[i__2] );
465 jy += *incy;
466 if (j > *ku) {
467 kx += *incx;
468 }
469 /* L140: */
470 }
471 }
472 }
473
474 return 0;
475
476 /* End of CGBMV . */
477
478 } /* cgbmv_ */
479
480 /* dgbmv.f -- translated by f2c (version 19991025).
481 You must link the resulting object file with the libraries:
482 -lf2c -lm (in that order)
483 */
484
PASTEF77(d,gbmv)485 /* Subroutine */ int PASTEF77(d,gbmv)(const bla_character *trans, const bla_integer *m, const bla_integer *n, const bla_integer *kl, const bla_integer *ku, const bla_double *alpha, const bla_double *a, const bla_integer *lda, const bla_double *x, const bla_integer *incx, const bla_double *beta, bla_double *y, const bla_integer *incy)
486 {
487 /* System generated locals */
488 bla_integer a_dim1, a_offset, i__1, i__2, i__3, i__4, i__5, i__6;
489
490 /* Local variables */
491 bla_integer info;
492 bla_double temp;
493 bla_integer lenx, leny, i__, j, k;
494 //extern bla_logical PASTEF770(lsame)(bla_character *, bla_character *, ftnlen, ftnlen);
495 bla_integer ix, iy, jx, jy, kx, ky;
496 //extern /* Subroutine */ int PASTEF770(xerbla)(bla_character *, bla_integer *, ftnlen);
497 bla_integer kup1;
498
499 /* .. Scalar Arguments .. */
500 /* .. Array Arguments .. */
501 /* .. */
502
503 /* Purpose */
504 /* ======= */
505
506 /* DGBMV performs one of the matrix-vector operations */
507
508 /* y := alpha*A*x + beta*y, or y := alpha*A'*x + beta*y, */
509
510 /* where alpha and beta are scalars, x and y are vectors and A is an */
511 /* m by n band matrix, with kl sub-diagonals and ku super-diagonals. */
512
513 /* Parameters */
514 /* ========== */
515
516 /* TRANS - CHARACTER*1. */
517 /* On entry, TRANS specifies the operation to be performed as */
518 /* follows: */
519
520 /* TRANS = 'N' or 'n' y := alpha*A*x + beta*y. */
521
522 /* TRANS = 'T' or 't' y := alpha*A'*x + beta*y. */
523
524 /* TRANS = 'C' or 'c' y := alpha*A'*x + beta*y. */
525
526 /* Unchanged on exit. */
527
528 /* M - INTEGER. */
529 /* On entry, M specifies the number of rows of the matrix A. */
530 /* M must be at least zero. */
531 /* Unchanged on exit. */
532
533 /* N - INTEGER. */
534 /* On entry, N specifies the number of columns of the matrix A. */
535 /* N must be at least zero. */
536 /* Unchanged on exit. */
537
538 /* KL - INTEGER. */
539 /* On entry, KL specifies the number of sub-diagonals of the */
540 /* matrix A. KL must satisfy 0 .le. KL. */
541 /* Unchanged on exit. */
542
543 /* KU - INTEGER. */
544 /* On entry, KU specifies the number of super-diagonals of the */
545 /* matrix A. KU must satisfy 0 .le. KU. */
546 /* Unchanged on exit. */
547
548 /* ALPHA - DOUBLE PRECISION. */
549 /* On entry, ALPHA specifies the scalar alpha. */
550 /* Unchanged on exit. */
551
552 /* A - DOUBLE PRECISION array of DIMENSION ( LDA, n ). */
553 /* Before entry, the leading ( kl + ku + 1 ) by n part of the */
554 /* array A must contain the matrix of coefficients, supplied */
555 /* column by column, with the leading diagonal of the matrix in */
556 /* row ( ku + 1 ) of the array, the first super-diagonal */
557 /* starting at position 2 in row ku, the first sub-diagonal */
558 /* starting at position 1 in row ( ku + 2 ), and so on. */
559 /* Elements in the array A that do not correspond to elements */
560 /* in the band matrix (such as the top left ku by ku triangle) */
561 /* are not referenced. */
562 /* The following program segment will transfer a band matrix */
563 /* from conventional full matrix storage to band storage: */
564
565 /* DO 20, J = 1, N */
566 /* K = KU + 1 - J */
567 /* DO 10, I = MAX( 1, J - KU ), MIN( M, J + KL ) */
568 /* A( K + I, J ) = matrix( I, J ) */
569 /* 10 CONTINUE */
570 /* 20 CONTINUE */
571
572 /* Unchanged on exit. */
573
574 /* LDA - INTEGER. */
575 /* On entry, LDA specifies the first dimension of A as declared */
576 /* in the calling (sub) program. LDA must be at least */
577 /* ( kl + ku + 1 ). */
578 /* Unchanged on exit. */
579
580 /* X - DOUBLE PRECISION array of DIMENSION at least */
581 /* ( 1 + ( n - 1 )*abs( INCX ) ) when TRANS = 'N' or 'n' */
582 /* and at least */
583 /* ( 1 + ( m - 1 )*abs( INCX ) ) otherwise. */
584 /* Before entry, the incremented array X must contain the */
585 /* vector x. */
586 /* Unchanged on exit. */
587
588 /* INCX - INTEGER. */
589 /* On entry, INCX specifies the increment for the elements of */
590 /* X. INCX must not be zero. */
591 /* Unchanged on exit. */
592
593 /* BETA - DOUBLE PRECISION. */
594 /* On entry, BETA specifies the scalar beta. When BETA is */
595 /* supplied as zero then Y need not be set on input. */
596 /* Unchanged on exit. */
597
598 /* Y - DOUBLE PRECISION array of DIMENSION at least */
599 /* ( 1 + ( m - 1 )*abs( INCY ) ) when TRANS = 'N' or 'n' */
600 /* and at least */
601 /* ( 1 + ( n - 1 )*abs( INCY ) ) otherwise. */
602 /* Before entry, the incremented array Y must contain the */
603 /* vector y. On exit, Y is overwritten by the updated vector y. */
604
605 /* INCY - INTEGER. */
606 /* On entry, INCY specifies the increment for the elements of */
607 /* Y. INCY must not be zero. */
608 /* Unchanged on exit. */
609
610
611 /* Level 2 Blas routine. */
612
613 /* -- Written on 22-October-1986. */
614 /* Jack Dongarra, Argonne National Lab. */
615 /* Jeremy Du Croz, Nag Central Office. */
616 /* Sven Hammarling, Nag Central Office. */
617 /* Richard Hanson, Sandia National Labs. */
618
619 /* .. Parameters .. */
620 /* .. Local Scalars .. */
621 /* .. External Functions .. */
622 /* .. External Subroutines .. */
623 /* .. Intrinsic Functions .. */
624 /* .. */
625 /* .. Executable Statements .. */
626
627 /* Test the input parameters. */
628
629 /* Parameter adjustments */
630 a_dim1 = *lda;
631 a_offset = 1 + a_dim1 * 1;
632 a -= a_offset;
633 --x;
634 --y;
635
636 /* Function Body */
637 info = 0;
638 if (! PASTEF770(lsame)(trans, "N", (ftnlen)1, (ftnlen)1) && ! PASTEF770(lsame)(trans, "T", (
639 ftnlen)1, (ftnlen)1) && ! PASTEF770(lsame)(trans, "C", (ftnlen)1, (ftnlen)1)
640 ) {
641 info = 1;
642 } else if (*m < 0) {
643 info = 2;
644 } else if (*n < 0) {
645 info = 3;
646 } else if (*kl < 0) {
647 info = 4;
648 } else if (*ku < 0) {
649 info = 5;
650 } else if (*lda < *kl + *ku + 1) {
651 info = 8;
652 } else if (*incx == 0) {
653 info = 10;
654 } else if (*incy == 0) {
655 info = 13;
656 }
657 if (info != 0) {
658 PASTEF770(xerbla)("DGBMV ", &info, (ftnlen)6);
659 return 0;
660 }
661
662 /* Quick return if possible. */
663
664 if (*m == 0 || *n == 0 || (*alpha == 0. && *beta == 1.)) {
665 return 0;
666 }
667
668 /* Set LENX and LENY, the lengths of the vectors x and y, and set */
669 /* up the start points in X and Y. */
670
671 if (PASTEF770(lsame)(trans, "N", (ftnlen)1, (ftnlen)1)) {
672 lenx = *n;
673 leny = *m;
674 } else {
675 lenx = *m;
676 leny = *n;
677 }
678 if (*incx > 0) {
679 kx = 1;
680 } else {
681 kx = 1 - (lenx - 1) * *incx;
682 }
683 if (*incy > 0) {
684 ky = 1;
685 } else {
686 ky = 1 - (leny - 1) * *incy;
687 }
688
689 /* Start the operations. In this version the elements of A are */
690 /* accessed sequentially with one pass through the band part of A. */
691
692 /* First form y := beta*y. */
693
694 if (*beta != 1.) {
695 if (*incy == 1) {
696 if (*beta == 0.) {
697 i__1 = leny;
698 for (i__ = 1; i__ <= i__1; ++i__) {
699 y[i__] = 0.;
700 /* L10: */
701 }
702 } else {
703 i__1 = leny;
704 for (i__ = 1; i__ <= i__1; ++i__) {
705 y[i__] = *beta * y[i__];
706 /* L20: */
707 }
708 }
709 } else {
710 iy = ky;
711 if (*beta == 0.) {
712 i__1 = leny;
713 for (i__ = 1; i__ <= i__1; ++i__) {
714 y[iy] = 0.;
715 iy += *incy;
716 /* L30: */
717 }
718 } else {
719 i__1 = leny;
720 for (i__ = 1; i__ <= i__1; ++i__) {
721 y[iy] = *beta * y[iy];
722 iy += *incy;
723 /* L40: */
724 }
725 }
726 }
727 }
728 if (*alpha == 0.) {
729 return 0;
730 }
731 kup1 = *ku + 1;
732 if (PASTEF770(lsame)(trans, "N", (ftnlen)1, (ftnlen)1)) {
733
734 /* Form y := alpha*A*x + y. */
735
736 jx = kx;
737 if (*incy == 1) {
738 i__1 = *n;
739 for (j = 1; j <= i__1; ++j) {
740 if (x[jx] != 0.) {
741 temp = *alpha * x[jx];
742 k = kup1 - j;
743 /* Computing MAX */
744 i__2 = 1, i__3 = j - *ku;
745 /* Computing MIN */
746 i__5 = *m, i__6 = j + *kl;
747 i__4 = f2c_min(i__5,i__6);
748 for (i__ = f2c_max(i__2,i__3); i__ <= i__4; ++i__) {
749 y[i__] += temp * a[k + i__ + j * a_dim1];
750 /* L50: */
751 }
752 }
753 jx += *incx;
754 /* L60: */
755 }
756 } else {
757 i__1 = *n;
758 for (j = 1; j <= i__1; ++j) {
759 if (x[jx] != 0.) {
760 temp = *alpha * x[jx];
761 iy = ky;
762 k = kup1 - j;
763 /* Computing MAX */
764 i__4 = 1, i__2 = j - *ku;
765 /* Computing MIN */
766 i__5 = *m, i__6 = j + *kl;
767 i__3 = f2c_min(i__5,i__6);
768 for (i__ = f2c_max(i__4,i__2); i__ <= i__3; ++i__) {
769 y[iy] += temp * a[k + i__ + j * a_dim1];
770 iy += *incy;
771 /* L70: */
772 }
773 }
774 jx += *incx;
775 if (j > *ku) {
776 ky += *incy;
777 }
778 /* L80: */
779 }
780 }
781 } else {
782
783 /* Form y := alpha*A'*x + y. */
784
785 jy = ky;
786 if (*incx == 1) {
787 i__1 = *n;
788 for (j = 1; j <= i__1; ++j) {
789 temp = 0.;
790 k = kup1 - j;
791 /* Computing MAX */
792 i__3 = 1, i__4 = j - *ku;
793 /* Computing MIN */
794 i__5 = *m, i__6 = j + *kl;
795 i__2 = f2c_min(i__5,i__6);
796 for (i__ = f2c_max(i__3,i__4); i__ <= i__2; ++i__) {
797 temp += a[k + i__ + j * a_dim1] * x[i__];
798 /* L90: */
799 }
800 y[jy] += *alpha * temp;
801 jy += *incy;
802 /* L100: */
803 }
804 } else {
805 i__1 = *n;
806 for (j = 1; j <= i__1; ++j) {
807 temp = 0.;
808 ix = kx;
809 k = kup1 - j;
810 /* Computing MAX */
811 i__2 = 1, i__3 = j - *ku;
812 /* Computing MIN */
813 i__5 = *m, i__6 = j + *kl;
814 i__4 = f2c_min(i__5,i__6);
815 for (i__ = f2c_max(i__2,i__3); i__ <= i__4; ++i__) {
816 temp += a[k + i__ + j * a_dim1] * x[ix];
817 ix += *incx;
818 /* L110: */
819 }
820 y[jy] += *alpha * temp;
821 jy += *incy;
822 if (j > *ku) {
823 kx += *incx;
824 }
825 /* L120: */
826 }
827 }
828 }
829
830 return 0;
831
832 /* End of DGBMV . */
833
834 } /* dgbmv_ */
835
836 /* sgbmv.f -- translated by f2c (version 19991025).
837 You must link the resulting object file with the libraries:
838 -lf2c -lm (in that order)
839 */
840
PASTEF77(s,gbmv)841 /* Subroutine */ int PASTEF77(s,gbmv)(const bla_character *trans, const bla_integer *m, const bla_integer *n, const bla_integer *kl, const bla_integer *ku, const bla_real *alpha, const bla_real *a, const bla_integer *lda, const bla_real *x, const bla_integer * incx, const bla_real *beta, bla_real *y, const bla_integer *incy)
842 {
843 /* System generated locals */
844 bla_integer a_dim1, a_offset, i__1, i__2, i__3, i__4, i__5, i__6;
845
846 /* Local variables */
847 bla_integer info;
848 bla_real temp;
849 bla_integer lenx, leny, i__, j, k;
850 //extern bla_logical PASTEF770(lsame)(bla_character *, bla_character *, ftnlen, ftnlen);
851 bla_integer ix, iy, jx, jy, kx, ky;
852 //extern /* Subroutine */ int PASTEF770(xerbla)(bla_character *, bla_integer *, ftnlen);
853 bla_integer kup1;
854
855 /* .. Scalar Arguments .. */
856 /* .. Array Arguments .. */
857 /* .. */
858
859 /* Purpose */
860 /* ======= */
861
862 /* SGBMV performs one of the matrix-vector operations */
863
864 /* y := alpha*A*x + beta*y, or y := alpha*A'*x + beta*y, */
865
866 /* where alpha and beta are scalars, x and y are vectors and A is an */
867 /* m by n band matrix, with kl sub-diagonals and ku super-diagonals. */
868
869 /* Parameters */
870 /* ========== */
871
872 /* TRANS - CHARACTER*1. */
873 /* On entry, TRANS specifies the operation to be performed as */
874 /* follows: */
875
876 /* TRANS = 'N' or 'n' y := alpha*A*x + beta*y. */
877
878 /* TRANS = 'T' or 't' y := alpha*A'*x + beta*y. */
879
880 /* TRANS = 'C' or 'c' y := alpha*A'*x + beta*y. */
881
882 /* Unchanged on exit. */
883
884 /* M - INTEGER. */
885 /* On entry, M specifies the number of rows of the matrix A. */
886 /* M must be at least zero. */
887 /* Unchanged on exit. */
888
889 /* N - INTEGER. */
890 /* On entry, N specifies the number of columns of the matrix A. */
891 /* N must be at least zero. */
892 /* Unchanged on exit. */
893
894 /* KL - INTEGER. */
895 /* On entry, KL specifies the number of sub-diagonals of the */
896 /* matrix A. KL must satisfy 0 .le. KL. */
897 /* Unchanged on exit. */
898
899 /* KU - INTEGER. */
900 /* On entry, KU specifies the number of super-diagonals of the */
901 /* matrix A. KU must satisfy 0 .le. KU. */
902 /* Unchanged on exit. */
903
904 /* ALPHA - REAL . */
905 /* On entry, ALPHA specifies the scalar alpha. */
906 /* Unchanged on exit. */
907
908 /* A - REAL array of DIMENSION ( LDA, n ). */
909 /* Before entry, the leading ( kl + ku + 1 ) by n part of the */
910 /* array A must contain the matrix of coefficients, supplied */
911 /* column by column, with the leading diagonal of the matrix in */
912 /* row ( ku + 1 ) of the array, the first super-diagonal */
913 /* starting at position 2 in row ku, the first sub-diagonal */
914 /* starting at position 1 in row ( ku + 2 ), and so on. */
915 /* Elements in the array A that do not correspond to elements */
916 /* in the band matrix (such as the top left ku by ku triangle) */
917 /* are not referenced. */
918 /* The following program segment will transfer a band matrix */
919 /* from conventional full matrix storage to band storage: */
920
921 /* DO 20, J = 1, N */
922 /* K = KU + 1 - J */
923 /* DO 10, I = MAX( 1, J - KU ), MIN( M, J + KL ) */
924 /* A( K + I, J ) = matrix( I, J ) */
925 /* 10 CONTINUE */
926 /* 20 CONTINUE */
927
928 /* Unchanged on exit. */
929
930 /* LDA - INTEGER. */
931 /* On entry, LDA specifies the first dimension of A as declared */
932 /* in the calling (sub) program. LDA must be at least */
933 /* ( kl + ku + 1 ). */
934 /* Unchanged on exit. */
935
936 /* X - REAL array of DIMENSION at least */
937 /* ( 1 + ( n - 1 )*abs( INCX ) ) when TRANS = 'N' or 'n' */
938 /* and at least */
939 /* ( 1 + ( m - 1 )*abs( INCX ) ) otherwise. */
940 /* Before entry, the incremented array X must contain the */
941 /* vector x. */
942 /* Unchanged on exit. */
943
944 /* INCX - INTEGER. */
945 /* On entry, INCX specifies the increment for the elements of */
946 /* X. INCX must not be zero. */
947 /* Unchanged on exit. */
948
949 /* BETA - REAL . */
950 /* On entry, BETA specifies the scalar beta. When BETA is */
951 /* supplied as zero then Y need not be set on input. */
952 /* Unchanged on exit. */
953
954 /* Y - REAL array of DIMENSION at least */
955 /* ( 1 + ( m - 1 )*abs( INCY ) ) when TRANS = 'N' or 'n' */
956 /* and at least */
957 /* ( 1 + ( n - 1 )*abs( INCY ) ) otherwise. */
958 /* Before entry, the incremented array Y must contain the */
959 /* vector y. On exit, Y is overwritten by the updated vector y. */
960
961 /* INCY - INTEGER. */
962 /* On entry, INCY specifies the increment for the elements of */
963 /* Y. INCY must not be zero. */
964 /* Unchanged on exit. */
965
966
967 /* Level 2 Blas routine. */
968
969 /* -- Written on 22-October-1986. */
970 /* Jack Dongarra, Argonne National Lab. */
971 /* Jeremy Du Croz, Nag Central Office. */
972 /* Sven Hammarling, Nag Central Office. */
973 /* Richard Hanson, Sandia National Labs. */
974
975 /* .. Parameters .. */
976 /* .. Local Scalars .. */
977 /* .. External Functions .. */
978 /* .. External Subroutines .. */
979 /* .. Intrinsic Functions .. */
980 /* .. */
981 /* .. Executable Statements .. */
982
983 /* Test the input parameters. */
984
985 /* Parameter adjustments */
986 a_dim1 = *lda;
987 a_offset = 1 + a_dim1 * 1;
988 a -= a_offset;
989 --x;
990 --y;
991
992 /* Function Body */
993 info = 0;
994 if (! PASTEF770(lsame)(trans, "N", (ftnlen)1, (ftnlen)1) && ! PASTEF770(lsame)(trans, "T", (
995 ftnlen)1, (ftnlen)1) && ! PASTEF770(lsame)(trans, "C", (ftnlen)1, (ftnlen)1)
996 ) {
997 info = 1;
998 } else if (*m < 0) {
999 info = 2;
1000 } else if (*n < 0) {
1001 info = 3;
1002 } else if (*kl < 0) {
1003 info = 4;
1004 } else if (*ku < 0) {
1005 info = 5;
1006 } else if (*lda < *kl + *ku + 1) {
1007 info = 8;
1008 } else if (*incx == 0) {
1009 info = 10;
1010 } else if (*incy == 0) {
1011 info = 13;
1012 }
1013 if (info != 0) {
1014 PASTEF770(xerbla)("SGBMV ", &info, (ftnlen)6);
1015 return 0;
1016 }
1017
1018 /* Quick return if possible. */
1019
1020 if (*m == 0 || *n == 0 || (*alpha == 0.f && *beta == 1.f)) {
1021 return 0;
1022 }
1023
1024 /* Set LENX and LENY, the lengths of the vectors x and y, and set */
1025 /* up the start points in X and Y. */
1026
1027 if (PASTEF770(lsame)(trans, "N", (ftnlen)1, (ftnlen)1)) {
1028 lenx = *n;
1029 leny = *m;
1030 } else {
1031 lenx = *m;
1032 leny = *n;
1033 }
1034 if (*incx > 0) {
1035 kx = 1;
1036 } else {
1037 kx = 1 - (lenx - 1) * *incx;
1038 }
1039 if (*incy > 0) {
1040 ky = 1;
1041 } else {
1042 ky = 1 - (leny - 1) * *incy;
1043 }
1044
1045 /* Start the operations. In this version the elements of A are */
1046 /* accessed sequentially with one pass through the band part of A. */
1047
1048 /* First form y := beta*y. */
1049
1050 if (*beta != 1.f) {
1051 if (*incy == 1) {
1052 if (*beta == 0.f) {
1053 i__1 = leny;
1054 for (i__ = 1; i__ <= i__1; ++i__) {
1055 y[i__] = 0.f;
1056 /* L10: */
1057 }
1058 } else {
1059 i__1 = leny;
1060 for (i__ = 1; i__ <= i__1; ++i__) {
1061 y[i__] = *beta * y[i__];
1062 /* L20: */
1063 }
1064 }
1065 } else {
1066 iy = ky;
1067 if (*beta == 0.f) {
1068 i__1 = leny;
1069 for (i__ = 1; i__ <= i__1; ++i__) {
1070 y[iy] = 0.f;
1071 iy += *incy;
1072 /* L30: */
1073 }
1074 } else {
1075 i__1 = leny;
1076 for (i__ = 1; i__ <= i__1; ++i__) {
1077 y[iy] = *beta * y[iy];
1078 iy += *incy;
1079 /* L40: */
1080 }
1081 }
1082 }
1083 }
1084 if (*alpha == 0.f) {
1085 return 0;
1086 }
1087 kup1 = *ku + 1;
1088 if (PASTEF770(lsame)(trans, "N", (ftnlen)1, (ftnlen)1)) {
1089
1090 /* Form y := alpha*A*x + y. */
1091
1092 jx = kx;
1093 if (*incy == 1) {
1094 i__1 = *n;
1095 for (j = 1; j <= i__1; ++j) {
1096 if (x[jx] != 0.f) {
1097 temp = *alpha * x[jx];
1098 k = kup1 - j;
1099 /* Computing MAX */
1100 i__2 = 1, i__3 = j - *ku;
1101 /* Computing MIN */
1102 i__5 = *m, i__6 = j + *kl;
1103 i__4 = f2c_min(i__5,i__6);
1104 for (i__ = f2c_max(i__2,i__3); i__ <= i__4; ++i__) {
1105 y[i__] += temp * a[k + i__ + j * a_dim1];
1106 /* L50: */
1107 }
1108 }
1109 jx += *incx;
1110 /* L60: */
1111 }
1112 } else {
1113 i__1 = *n;
1114 for (j = 1; j <= i__1; ++j) {
1115 if (x[jx] != 0.f) {
1116 temp = *alpha * x[jx];
1117 iy = ky;
1118 k = kup1 - j;
1119 /* Computing MAX */
1120 i__4 = 1, i__2 = j - *ku;
1121 /* Computing MIN */
1122 i__5 = *m, i__6 = j + *kl;
1123 i__3 = f2c_min(i__5,i__6);
1124 for (i__ = f2c_max(i__4,i__2); i__ <= i__3; ++i__) {
1125 y[iy] += temp * a[k + i__ + j * a_dim1];
1126 iy += *incy;
1127 /* L70: */
1128 }
1129 }
1130 jx += *incx;
1131 if (j > *ku) {
1132 ky += *incy;
1133 }
1134 /* L80: */
1135 }
1136 }
1137 } else {
1138
1139 /* Form y := alpha*A'*x + y. */
1140
1141 jy = ky;
1142 if (*incx == 1) {
1143 i__1 = *n;
1144 for (j = 1; j <= i__1; ++j) {
1145 temp = 0.f;
1146 k = kup1 - j;
1147 /* Computing MAX */
1148 i__3 = 1, i__4 = j - *ku;
1149 /* Computing MIN */
1150 i__5 = *m, i__6 = j + *kl;
1151 i__2 = f2c_min(i__5,i__6);
1152 for (i__ = f2c_max(i__3,i__4); i__ <= i__2; ++i__) {
1153 temp += a[k + i__ + j * a_dim1] * x[i__];
1154 /* L90: */
1155 }
1156 y[jy] += *alpha * temp;
1157 jy += *incy;
1158 /* L100: */
1159 }
1160 } else {
1161 i__1 = *n;
1162 for (j = 1; j <= i__1; ++j) {
1163 temp = 0.f;
1164 ix = kx;
1165 k = kup1 - j;
1166 /* Computing MAX */
1167 i__2 = 1, i__3 = j - *ku;
1168 /* Computing MIN */
1169 i__5 = *m, i__6 = j + *kl;
1170 i__4 = f2c_min(i__5,i__6);
1171 for (i__ = f2c_max(i__2,i__3); i__ <= i__4; ++i__) {
1172 temp += a[k + i__ + j * a_dim1] * x[ix];
1173 ix += *incx;
1174 /* L110: */
1175 }
1176 y[jy] += *alpha * temp;
1177 jy += *incy;
1178 if (j > *ku) {
1179 kx += *incx;
1180 }
1181 /* L120: */
1182 }
1183 }
1184 }
1185
1186 return 0;
1187
1188 /* End of SGBMV . */
1189
1190 } /* sgbmv_ */
1191
1192 /* zgbmv.f -- translated by f2c (version 19991025).
1193 You must link the resulting object file with the libraries:
1194 -lf2c -lm (in that order)
1195 */
1196
PASTEF77(z,gbmv)1197 /* Subroutine */ int PASTEF77(z,gbmv)(const bla_character *trans, const bla_integer *m, const bla_integer *n, const bla_integer *kl, const bla_integer *ku, const bla_dcomplex *alpha, const bla_dcomplex *a, const bla_integer *lda, const bla_dcomplex *x, const bla_integer *incx, const bla_dcomplex *beta, bla_dcomplex * y, const bla_integer *incy)
1198 {
1199 /* System generated locals */
1200 bla_integer a_dim1, a_offset, i__1, i__2, i__3, i__4, i__5, i__6;
1201 bla_dcomplex z__1, z__2, z__3;
1202
1203 /* Builtin functions */
1204 //void bla_d_cnjg(bla_dcomplex *, bla_dcomplex *);
1205
1206 /* Local variables */
1207 bla_integer info;
1208 bla_dcomplex temp;
1209 bla_integer lenx, leny, i__, j, k;
1210 //extern bla_logical PASTEF770(lsame)(bla_character *, bla_character *, ftnlen, ftnlen);
1211 bla_integer ix, iy, jx, jy, kx, ky;
1212 //extern /* Subroutine */ int PASTEF770(xerbla)(bla_character *, bla_integer *, ftnlen);
1213 bla_logical noconj;
1214 bla_integer kup1;
1215
1216 /* .. Scalar Arguments .. */
1217 /* .. Array Arguments .. */
1218 /* .. */
1219
1220 /* Purpose */
1221 /* ======= */
1222
1223 /* ZGBMV performs one of the matrix-vector operations */
1224
1225 /* y := alpha*A*x + beta*y, or y := alpha*A'*x + beta*y, or */
1226
1227 /* y := alpha*conjg( A' )*x + beta*y, */
1228
1229 /* where alpha and beta are scalars, x and y are vectors and A is an */
1230 /* m by n band matrix, with kl sub-diagonals and ku super-diagonals. */
1231
1232 /* Parameters */
1233 /* ========== */
1234
1235 /* TRANS - CHARACTER*1. */
1236 /* On entry, TRANS specifies the operation to be performed as */
1237 /* follows: */
1238
1239 /* TRANS = 'N' or 'n' y := alpha*A*x + beta*y. */
1240
1241 /* TRANS = 'T' or 't' y := alpha*A'*x + beta*y. */
1242
1243 /* TRANS = 'C' or 'c' y := alpha*conjg( A' )*x + beta*y. */
1244
1245 /* Unchanged on exit. */
1246
1247 /* M - INTEGER. */
1248 /* On entry, M specifies the number of rows of the matrix A. */
1249 /* M must be at least zero. */
1250 /* Unchanged on exit. */
1251
1252 /* N - INTEGER. */
1253 /* On entry, N specifies the number of columns of the matrix A. */
1254 /* N must be at least zero. */
1255 /* Unchanged on exit. */
1256
1257 /* KL - INTEGER. */
1258 /* On entry, KL specifies the number of sub-diagonals of the */
1259 /* matrix A. KL must satisfy 0 .le. KL. */
1260 /* Unchanged on exit. */
1261
1262 /* KU - INTEGER. */
1263 /* On entry, KU specifies the number of super-diagonals of the */
1264 /* matrix A. KU must satisfy 0 .le. KU. */
1265 /* Unchanged on exit. */
1266
1267 /* ALPHA - COMPLEX*16 . */
1268 /* On entry, ALPHA specifies the scalar alpha. */
1269 /* Unchanged on exit. */
1270
1271 /* A - COMPLEX*16 array of DIMENSION ( LDA, n ). */
1272 /* Before entry, the leading ( kl + ku + 1 ) by n part of the */
1273 /* array A must contain the matrix of coefficients, supplied */
1274 /* column by column, with the leading diagonal of the matrix in */
1275 /* row ( ku + 1 ) of the array, the first super-diagonal */
1276 /* starting at position 2 in row ku, the first sub-diagonal */
1277 /* starting at position 1 in row ( ku + 2 ), and so on. */
1278 /* Elements in the array A that do not correspond to elements */
1279 /* in the band matrix (such as the top left ku by ku triangle) */
1280 /* are not referenced. */
1281 /* The following program segment will transfer a band matrix */
1282 /* from conventional full matrix storage to band storage: */
1283
1284 /* DO 20, J = 1, N */
1285 /* K = KU + 1 - J */
1286 /* DO 10, I = MAX( 1, J - KU ), MIN( M, J + KL ) */
1287 /* A( K + I, J ) = matrix( I, J ) */
1288 /* 10 CONTINUE */
1289 /* 20 CONTINUE */
1290
1291 /* Unchanged on exit. */
1292
1293 /* LDA - INTEGER. */
1294 /* On entry, LDA specifies the first dimension of A as declared */
1295 /* in the calling (sub) program. LDA must be at least */
1296 /* ( kl + ku + 1 ). */
1297 /* Unchanged on exit. */
1298
1299 /* X - COMPLEX*16 array of DIMENSION at least */
1300 /* ( 1 + ( n - 1 )*abs( INCX ) ) when TRANS = 'N' or 'n' */
1301 /* and at least */
1302 /* ( 1 + ( m - 1 )*abs( INCX ) ) otherwise. */
1303 /* Before entry, the incremented array X must contain the */
1304 /* vector x. */
1305 /* Unchanged on exit. */
1306
1307 /* INCX - INTEGER. */
1308 /* On entry, INCX specifies the increment for the elements of */
1309 /* X. INCX must not be zero. */
1310 /* Unchanged on exit. */
1311
1312 /* BETA - COMPLEX*16 . */
1313 /* On entry, BETA specifies the scalar beta. When BETA is */
1314 /* supplied as zero then Y need not be set on input. */
1315 /* Unchanged on exit. */
1316
1317 /* Y - COMPLEX*16 array of DIMENSION at least */
1318 /* ( 1 + ( m - 1 )*abs( INCY ) ) when TRANS = 'N' or 'n' */
1319 /* and at least */
1320 /* ( 1 + ( n - 1 )*abs( INCY ) ) otherwise. */
1321 /* Before entry, the incremented array Y must contain the */
1322 /* vector y. On exit, Y is overwritten by the updated vector y. */
1323
1324
1325 /* INCY - INTEGER. */
1326 /* On entry, INCY specifies the increment for the elements of */
1327 /* Y. INCY must not be zero. */
1328 /* Unchanged on exit. */
1329
1330
1331 /* Level 2 Blas routine. */
1332
1333 /* -- Written on 22-October-1986. */
1334 /* Jack Dongarra, Argonne National Lab. */
1335 /* Jeremy Du Croz, Nag Central Office. */
1336 /* Sven Hammarling, Nag Central Office. */
1337 /* Richard Hanson, Sandia National Labs. */
1338
1339
1340 /* .. Parameters .. */
1341 /* .. Local Scalars .. */
1342 /* .. External Functions .. */
1343 /* .. External Subroutines .. */
1344 /* .. Intrinsic Functions .. */
1345 /* .. */
1346 /* .. Executable Statements .. */
1347
1348 /* Test the input parameters. */
1349
1350 /* Parameter adjustments */
1351 a_dim1 = *lda;
1352 a_offset = 1 + a_dim1 * 1;
1353 a -= a_offset;
1354 --x;
1355 --y;
1356
1357 /* Function Body */
1358 info = 0;
1359 if (! PASTEF770(lsame)(trans, "N", (ftnlen)1, (ftnlen)1) && ! PASTEF770(lsame)(trans, "T", (
1360 ftnlen)1, (ftnlen)1) && ! PASTEF770(lsame)(trans, "C", (ftnlen)1, (ftnlen)1)
1361 ) {
1362 info = 1;
1363 } else if (*m < 0) {
1364 info = 2;
1365 } else if (*n < 0) {
1366 info = 3;
1367 } else if (*kl < 0) {
1368 info = 4;
1369 } else if (*ku < 0) {
1370 info = 5;
1371 } else if (*lda < *kl + *ku + 1) {
1372 info = 8;
1373 } else if (*incx == 0) {
1374 info = 10;
1375 } else if (*incy == 0) {
1376 info = 13;
1377 }
1378 if (info != 0) {
1379 PASTEF770(xerbla)("ZGBMV ", &info, (ftnlen)6);
1380 return 0;
1381 }
1382
1383 /* Quick return if possible. */
1384
1385 if (*m == 0 || *n == 0 || (bli_zreal(*alpha) == 0. && bli_zimag(*alpha) == 0. && (bli_zreal(*beta) ==
1386 1. && bli_zimag(*beta) == 0.))) {
1387 return 0;
1388 }
1389
1390 noconj = PASTEF770(lsame)(trans, "T", (ftnlen)1, (ftnlen)1);
1391
1392 /* Set LENX and LENY, the lengths of the vectors x and y, and set */
1393 /* up the start points in X and Y. */
1394
1395 if (PASTEF770(lsame)(trans, "N", (ftnlen)1, (ftnlen)1)) {
1396 lenx = *n;
1397 leny = *m;
1398 } else {
1399 lenx = *m;
1400 leny = *n;
1401 }
1402 if (*incx > 0) {
1403 kx = 1;
1404 } else {
1405 kx = 1 - (lenx - 1) * *incx;
1406 }
1407 if (*incy > 0) {
1408 ky = 1;
1409 } else {
1410 ky = 1 - (leny - 1) * *incy;
1411 }
1412
1413 /* Start the operations. In this version the elements of A are */
1414 /* accessed sequentially with one pass through the band part of A. */
1415
1416 /* First form y := beta*y. */
1417
1418 if (bli_zreal(*beta) != 1. || bli_zimag(*beta) != 0.) {
1419 if (*incy == 1) {
1420 if (bli_zreal(*beta) == 0. && bli_zimag(*beta) == 0.) {
1421 i__1 = leny;
1422 for (i__ = 1; i__ <= i__1; ++i__) {
1423 i__2 = i__;
1424 bli_zsets( (0.), (0.), y[i__2] );
1425 /* L10: */
1426 }
1427 } else {
1428 i__1 = leny;
1429 for (i__ = 1; i__ <= i__1; ++i__) {
1430 i__2 = i__;
1431 i__3 = i__;
1432 bli_zsets( (bli_zreal(*beta) * bli_zreal(y[i__3]) - bli_zimag(*beta) * bli_zimag(y[i__3])), (bli_zreal(*beta) * bli_zimag(y[i__3]) + bli_zimag(*beta) * bli_zreal(y[i__3])), z__1 );
1433 bli_zsets( (bli_zreal(z__1)), (bli_zimag(z__1)), y[i__2] );
1434 /* L20: */
1435 }
1436 }
1437 } else {
1438 iy = ky;
1439 if (bli_zreal(*beta) == 0. && bli_zimag(*beta) == 0.) {
1440 i__1 = leny;
1441 for (i__ = 1; i__ <= i__1; ++i__) {
1442 i__2 = iy;
1443 bli_zsets( (0.), (0.), y[i__2] );
1444 iy += *incy;
1445 /* L30: */
1446 }
1447 } else {
1448 i__1 = leny;
1449 for (i__ = 1; i__ <= i__1; ++i__) {
1450 i__2 = iy;
1451 i__3 = iy;
1452 bli_zsets( (bli_zreal(*beta) * bli_zreal(y[i__3]) - bli_zimag(*beta) * bli_zimag(y[i__3])), (bli_zreal(*beta) * bli_zimag(y[i__3]) + bli_zimag(*beta) * bli_zreal(y[i__3])), z__1 );
1453 bli_zsets( (bli_zreal(z__1)), (bli_zimag(z__1)), y[i__2] );
1454 iy += *incy;
1455 /* L40: */
1456 }
1457 }
1458 }
1459 }
1460 if (bli_zreal(*alpha) == 0. && bli_zimag(*alpha) == 0.) {
1461 return 0;
1462 }
1463 kup1 = *ku + 1;
1464 if (PASTEF770(lsame)(trans, "N", (ftnlen)1, (ftnlen)1)) {
1465
1466 /* Form y := alpha*A*x + y. */
1467
1468 jx = kx;
1469 if (*incy == 1) {
1470 i__1 = *n;
1471 for (j = 1; j <= i__1; ++j) {
1472 i__2 = jx;
1473 if (bli_zreal(x[i__2]) != 0. || bli_zimag(x[i__2]) != 0.) {
1474 i__2 = jx;
1475 bli_zsets( (bli_zreal(*alpha) * bli_zreal(x[i__2]) - bli_zimag(*alpha) * bli_zimag(x[i__2])), (bli_zreal(*alpha) * bli_zimag(x[i__2]) + bli_zimag(*alpha) * bli_zreal(x[i__2])), z__1 );
1476 bli_zsets( (bli_zreal(z__1)), (bli_zimag(z__1)), temp );
1477 k = kup1 - j;
1478 /* Computing MAX */
1479 i__2 = 1, i__3 = j - *ku;
1480 /* Computing MIN */
1481 i__5 = *m, i__6 = j + *kl;
1482 i__4 = f2c_min(i__5,i__6);
1483 for (i__ = f2c_max(i__2,i__3); i__ <= i__4; ++i__) {
1484 i__2 = i__;
1485 i__3 = i__;
1486 i__5 = k + i__ + j * a_dim1;
1487 bli_zsets( (bli_zreal(temp) * bli_zreal(a[i__5]) - bli_zimag(temp) * bli_zimag(a[i__5])), (bli_zreal(temp) * bli_zimag(a[i__5]) + bli_zimag(temp) * bli_zreal(a[i__5])), z__2 );
1488 bli_zsets( (bli_zreal(y[i__3]) + bli_zreal(z__2)), (bli_zimag(y[i__3]) + bli_zimag(z__2)), z__1 );
1489 bli_zsets( (bli_zreal(z__1)), (bli_zimag(z__1)), y[i__2] );
1490 /* L50: */
1491 }
1492 }
1493 jx += *incx;
1494 /* L60: */
1495 }
1496 } else {
1497 i__1 = *n;
1498 for (j = 1; j <= i__1; ++j) {
1499 i__4 = jx;
1500 if (bli_zreal(x[i__4]) != 0. || bli_zimag(x[i__4]) != 0.) {
1501 i__4 = jx;
1502 bli_zsets( (bli_zreal(*alpha) * bli_zreal(x[i__4]) - bli_zimag(*alpha) * bli_zimag(x[i__4])), (bli_zreal(*alpha) * bli_zimag(x[i__4]) + bli_zimag(*alpha) * bli_zreal(x[i__4])), z__1 );
1503 bli_zsets( (bli_zreal(z__1)), (bli_zimag(z__1)), temp );
1504 iy = ky;
1505 k = kup1 - j;
1506 /* Computing MAX */
1507 i__4 = 1, i__2 = j - *ku;
1508 /* Computing MIN */
1509 i__5 = *m, i__6 = j + *kl;
1510 i__3 = f2c_min(i__5,i__6);
1511 for (i__ = f2c_max(i__4,i__2); i__ <= i__3; ++i__) {
1512 i__4 = iy;
1513 i__2 = iy;
1514 i__5 = k + i__ + j * a_dim1;
1515 bli_zsets( (bli_zreal(temp) * bli_zreal(a[i__5]) - bli_zimag(temp) * bli_zimag(a[i__5])), (bli_zreal(temp) * bli_zimag(a[i__5]) + bli_zimag(temp) * bli_zreal(a[i__5])), z__2 );
1516 bli_zsets( (bli_zreal(y[i__2]) + bli_zreal(z__2)), (bli_zimag(y[i__2]) + bli_zimag(z__2)), z__1 );
1517 bli_zsets( (bli_zreal(z__1)), (bli_zimag(z__1)), y[i__4] );
1518 iy += *incy;
1519 /* L70: */
1520 }
1521 }
1522 jx += *incx;
1523 if (j > *ku) {
1524 ky += *incy;
1525 }
1526 /* L80: */
1527 }
1528 }
1529 } else {
1530
1531 /* Form y := alpha*A'*x + y or y := alpha*conjg( A' )*x + y. */
1532
1533 jy = ky;
1534 if (*incx == 1) {
1535 i__1 = *n;
1536 for (j = 1; j <= i__1; ++j) {
1537 bli_zsets( (0.), (0.), temp );
1538 k = kup1 - j;
1539 if (noconj) {
1540 /* Computing MAX */
1541 i__3 = 1, i__4 = j - *ku;
1542 /* Computing MIN */
1543 i__5 = *m, i__6 = j + *kl;
1544 i__2 = f2c_min(i__5,i__6);
1545 for (i__ = f2c_max(i__3,i__4); i__ <= i__2; ++i__) {
1546 i__3 = k + i__ + j * a_dim1;
1547 i__4 = i__;
1548 bli_zsets( (bli_zreal(a[i__3]) * bli_zreal(x[i__4]) - bli_zimag(a[i__3]) * bli_zimag(x[i__4])), (bli_zreal(a[i__3]) * bli_zimag(x[i__4]) + bli_zimag(a[i__3]) * bli_zreal(x[i__4])), z__2 );
1549 bli_zsets( (bli_zreal(temp) + bli_zreal(z__2)), (bli_zimag(temp) + bli_zimag(z__2)), z__1 );
1550 bli_zsets( (bli_zreal(z__1)), (bli_zimag(z__1)), temp );
1551 /* L90: */
1552 }
1553 } else {
1554 /* Computing MAX */
1555 i__2 = 1, i__3 = j - *ku;
1556 /* Computing MIN */
1557 i__5 = *m, i__6 = j + *kl;
1558 i__4 = f2c_min(i__5,i__6);
1559 for (i__ = f2c_max(i__2,i__3); i__ <= i__4; ++i__) {
1560 bla_d_cnjg(&z__3, &a[k + i__ + j * a_dim1]);
1561 i__2 = i__;
1562 bli_zsets( (bli_zreal(z__3) * bli_zreal(x[i__2]) - bli_zimag(z__3) * bli_zimag(x[i__2])), (bli_zreal(z__3) * bli_zimag(x[i__2]) + bli_zimag(z__3) * bli_zreal(x[i__2])), z__2 );
1563 bli_zsets( (bli_zreal(temp) + bli_zreal(z__2)), (bli_zimag(temp) + bli_zimag(z__2)), z__1 );
1564 bli_zsets( (bli_zreal(z__1)), (bli_zimag(z__1)), temp );
1565 /* L100: */
1566 }
1567 }
1568 i__4 = jy;
1569 i__2 = jy;
1570 bli_zsets( (bli_zreal(*alpha) * bli_zreal(temp) - bli_zimag(*alpha) * bli_zimag(temp)), (bli_zreal(*alpha) * bli_zimag(temp) + bli_zimag(*alpha) * bli_zreal(temp)), z__2 );
1571 bli_zsets( (bli_zreal(y[i__2]) + bli_zreal(z__2)), (bli_zimag(y[i__2]) + bli_zimag(z__2)), z__1 );
1572 bli_zsets( (bli_zreal(z__1)), (bli_zimag(z__1)), y[i__4] );
1573 jy += *incy;
1574 /* L110: */
1575 }
1576 } else {
1577 i__1 = *n;
1578 for (j = 1; j <= i__1; ++j) {
1579 bli_zsets( (0.), (0.), temp );
1580 ix = kx;
1581 k = kup1 - j;
1582 if (noconj) {
1583 /* Computing MAX */
1584 i__4 = 1, i__2 = j - *ku;
1585 /* Computing MIN */
1586 i__5 = *m, i__6 = j + *kl;
1587 i__3 = f2c_min(i__5,i__6);
1588 for (i__ = f2c_max(i__4,i__2); i__ <= i__3; ++i__) {
1589 i__4 = k + i__ + j * a_dim1;
1590 i__2 = ix;
1591 bli_zsets( (bli_zreal(a[i__4]) * bli_zreal(x[i__2]) - bli_zimag(a[i__4]) * bli_zimag(x[i__2])), (bli_zreal(a[i__4]) * bli_zimag(x[i__2]) + bli_zimag(a[i__4]) * bli_zreal(x[i__2])), z__2 );
1592 bli_zsets( (bli_zreal(temp) + bli_zreal(z__2)), (bli_zimag(temp) + bli_zimag(z__2)), z__1 );
1593 bli_zsets( (bli_zreal(z__1)), (bli_zimag(z__1)), temp );
1594 ix += *incx;
1595 /* L120: */
1596 }
1597 } else {
1598 /* Computing MAX */
1599 i__3 = 1, i__4 = j - *ku;
1600 /* Computing MIN */
1601 i__5 = *m, i__6 = j + *kl;
1602 i__2 = f2c_min(i__5,i__6);
1603 for (i__ = f2c_max(i__3,i__4); i__ <= i__2; ++i__) {
1604 bla_d_cnjg(&z__3, &a[k + i__ + j * a_dim1]);
1605 i__3 = ix;
1606 bli_zsets( (bli_zreal(z__3) * bli_zreal(x[i__3]) - bli_zimag(z__3) * bli_zimag(x[i__3])), (bli_zreal(z__3) * bli_zimag(x[i__3]) + bli_zimag(z__3) * bli_zreal(x[i__3])), z__2 );
1607 bli_zsets( (bli_zreal(temp) + bli_zreal(z__2)), (bli_zimag(temp) + bli_zimag(z__2)), z__1 );
1608 bli_zsets( (bli_zreal(z__1)), (bli_zimag(z__1)), temp );
1609 ix += *incx;
1610 /* L130: */
1611 }
1612 }
1613 i__2 = jy;
1614 i__3 = jy;
1615 bli_zsets( (bli_zreal(*alpha) * bli_zreal(temp) - bli_zimag(*alpha) * bli_zimag(temp)), (bli_zreal(*alpha) * bli_zimag(temp) + bli_zimag(*alpha) * bli_zreal(temp)), z__2 );
1616 bli_zsets( (bli_zreal(y[i__3]) + bli_zreal(z__2)), (bli_zimag(y[i__3]) + bli_zimag(z__2)), z__1 );
1617 bli_zsets( (bli_zreal(z__1)), (bli_zimag(z__1)), y[i__2] );
1618 jy += *incy;
1619 if (j > *ku) {
1620 kx += *incx;
1621 }
1622 /* L140: */
1623 }
1624 }
1625 }
1626
1627 return 0;
1628
1629 /* End of ZGBMV . */
1630
1631 } /* zgbmv_ */
1632
1633 #endif
1634
1635