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