00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013 #include "f2c.h"
00014 #include "blaswrap.h"
00015
00016 int zhpmv_(char *uplo, integer *n, doublecomplex *alpha,
00017 doublecomplex *ap, doublecomplex *x, integer *incx, doublecomplex *
00018 beta, doublecomplex *y, integer *incy)
00019 {
00020
00021 integer i__1, i__2, i__3, i__4, i__5;
00022 doublereal d__1;
00023 doublecomplex z__1, z__2, z__3, z__4;
00024
00025
00026 void d_cnjg(doublecomplex *, doublecomplex *);
00027
00028
00029 integer i__, j, k, kk, ix, iy, jx, jy, kx, ky, info;
00030 doublecomplex temp1, temp2;
00031 extern logical lsame_(char *, char *);
00032 extern int xerbla_(char *, integer *);
00033
00034
00035
00036
00037
00038
00039
00040
00041
00042
00043
00044
00045
00046
00047
00048
00049
00050
00051
00052
00053
00054
00055
00056
00057
00058
00059
00060
00061
00062
00063
00064
00065
00066
00067
00068
00069
00070
00071
00072
00073
00074
00075
00076
00077
00078
00079
00080
00081
00082
00083
00084
00085
00086
00087
00088
00089
00090
00091
00092
00093
00094
00095
00096
00097
00098
00099
00100
00101
00102
00103
00104
00105
00106
00107
00108
00109
00110
00111
00112
00113
00114
00115
00116
00117
00118
00119
00120
00121
00122
00123
00124
00125
00126
00127
00128
00129
00130
00131
00132
00133
00134
00135
00136
00137
00138
00139
00140
00141 --y;
00142 --x;
00143 --ap;
00144
00145
00146 info = 0;
00147 if (! lsame_(uplo, "U") && ! lsame_(uplo, "L")) {
00148 info = 1;
00149 } else if (*n < 0) {
00150 info = 2;
00151 } else if (*incx == 0) {
00152 info = 6;
00153 } else if (*incy == 0) {
00154 info = 9;
00155 }
00156 if (info != 0) {
00157 xerbla_("ZHPMV ", &info);
00158 return 0;
00159 }
00160
00161
00162
00163 if (*n == 0 || alpha->r == 0. && alpha->i == 0. && (beta->r == 1. &&
00164 beta->i == 0.)) {
00165 return 0;
00166 }
00167
00168
00169
00170 if (*incx > 0) {
00171 kx = 1;
00172 } else {
00173 kx = 1 - (*n - 1) * *incx;
00174 }
00175 if (*incy > 0) {
00176 ky = 1;
00177 } else {
00178 ky = 1 - (*n - 1) * *incy;
00179 }
00180
00181
00182
00183
00184
00185
00186 if (beta->r != 1. || beta->i != 0.) {
00187 if (*incy == 1) {
00188 if (beta->r == 0. && beta->i == 0.) {
00189 i__1 = *n;
00190 for (i__ = 1; i__ <= i__1; ++i__) {
00191 i__2 = i__;
00192 y[i__2].r = 0., y[i__2].i = 0.;
00193
00194 }
00195 } else {
00196 i__1 = *n;
00197 for (i__ = 1; i__ <= i__1; ++i__) {
00198 i__2 = i__;
00199 i__3 = i__;
00200 z__1.r = beta->r * y[i__3].r - beta->i * y[i__3].i,
00201 z__1.i = beta->r * y[i__3].i + beta->i * y[i__3]
00202 .r;
00203 y[i__2].r = z__1.r, y[i__2].i = z__1.i;
00204
00205 }
00206 }
00207 } else {
00208 iy = ky;
00209 if (beta->r == 0. && beta->i == 0.) {
00210 i__1 = *n;
00211 for (i__ = 1; i__ <= i__1; ++i__) {
00212 i__2 = iy;
00213 y[i__2].r = 0., y[i__2].i = 0.;
00214 iy += *incy;
00215
00216 }
00217 } else {
00218 i__1 = *n;
00219 for (i__ = 1; i__ <= i__1; ++i__) {
00220 i__2 = iy;
00221 i__3 = iy;
00222 z__1.r = beta->r * y[i__3].r - beta->i * y[i__3].i,
00223 z__1.i = beta->r * y[i__3].i + beta->i * y[i__3]
00224 .r;
00225 y[i__2].r = z__1.r, y[i__2].i = z__1.i;
00226 iy += *incy;
00227
00228 }
00229 }
00230 }
00231 }
00232 if (alpha->r == 0. && alpha->i == 0.) {
00233 return 0;
00234 }
00235 kk = 1;
00236 if (lsame_(uplo, "U")) {
00237
00238
00239
00240 if (*incx == 1 && *incy == 1) {
00241 i__1 = *n;
00242 for (j = 1; j <= i__1; ++j) {
00243 i__2 = j;
00244 z__1.r = alpha->r * x[i__2].r - alpha->i * x[i__2].i, z__1.i =
00245 alpha->r * x[i__2].i + alpha->i * x[i__2].r;
00246 temp1.r = z__1.r, temp1.i = z__1.i;
00247 temp2.r = 0., temp2.i = 0.;
00248 k = kk;
00249 i__2 = j - 1;
00250 for (i__ = 1; i__ <= i__2; ++i__) {
00251 i__3 = i__;
00252 i__4 = i__;
00253 i__5 = k;
00254 z__2.r = temp1.r * ap[i__5].r - temp1.i * ap[i__5].i,
00255 z__2.i = temp1.r * ap[i__5].i + temp1.i * ap[i__5]
00256 .r;
00257 z__1.r = y[i__4].r + z__2.r, z__1.i = y[i__4].i + z__2.i;
00258 y[i__3].r = z__1.r, y[i__3].i = z__1.i;
00259 d_cnjg(&z__3, &ap[k]);
00260 i__3 = i__;
00261 z__2.r = z__3.r * x[i__3].r - z__3.i * x[i__3].i, z__2.i =
00262 z__3.r * x[i__3].i + z__3.i * x[i__3].r;
00263 z__1.r = temp2.r + z__2.r, z__1.i = temp2.i + z__2.i;
00264 temp2.r = z__1.r, temp2.i = z__1.i;
00265 ++k;
00266
00267 }
00268 i__2 = j;
00269 i__3 = j;
00270 i__4 = kk + j - 1;
00271 d__1 = ap[i__4].r;
00272 z__3.r = d__1 * temp1.r, z__3.i = d__1 * temp1.i;
00273 z__2.r = y[i__3].r + z__3.r, z__2.i = y[i__3].i + z__3.i;
00274 z__4.r = alpha->r * temp2.r - alpha->i * temp2.i, z__4.i =
00275 alpha->r * temp2.i + alpha->i * temp2.r;
00276 z__1.r = z__2.r + z__4.r, z__1.i = z__2.i + z__4.i;
00277 y[i__2].r = z__1.r, y[i__2].i = z__1.i;
00278 kk += j;
00279
00280 }
00281 } else {
00282 jx = kx;
00283 jy = ky;
00284 i__1 = *n;
00285 for (j = 1; j <= i__1; ++j) {
00286 i__2 = jx;
00287 z__1.r = alpha->r * x[i__2].r - alpha->i * x[i__2].i, z__1.i =
00288 alpha->r * x[i__2].i + alpha->i * x[i__2].r;
00289 temp1.r = z__1.r, temp1.i = z__1.i;
00290 temp2.r = 0., temp2.i = 0.;
00291 ix = kx;
00292 iy = ky;
00293 i__2 = kk + j - 2;
00294 for (k = kk; k <= i__2; ++k) {
00295 i__3 = iy;
00296 i__4 = iy;
00297 i__5 = k;
00298 z__2.r = temp1.r * ap[i__5].r - temp1.i * ap[i__5].i,
00299 z__2.i = temp1.r * ap[i__5].i + temp1.i * ap[i__5]
00300 .r;
00301 z__1.r = y[i__4].r + z__2.r, z__1.i = y[i__4].i + z__2.i;
00302 y[i__3].r = z__1.r, y[i__3].i = z__1.i;
00303 d_cnjg(&z__3, &ap[k]);
00304 i__3 = ix;
00305 z__2.r = z__3.r * x[i__3].r - z__3.i * x[i__3].i, z__2.i =
00306 z__3.r * x[i__3].i + z__3.i * x[i__3].r;
00307 z__1.r = temp2.r + z__2.r, z__1.i = temp2.i + z__2.i;
00308 temp2.r = z__1.r, temp2.i = z__1.i;
00309 ix += *incx;
00310 iy += *incy;
00311
00312 }
00313 i__2 = jy;
00314 i__3 = jy;
00315 i__4 = kk + j - 1;
00316 d__1 = ap[i__4].r;
00317 z__3.r = d__1 * temp1.r, z__3.i = d__1 * temp1.i;
00318 z__2.r = y[i__3].r + z__3.r, z__2.i = y[i__3].i + z__3.i;
00319 z__4.r = alpha->r * temp2.r - alpha->i * temp2.i, z__4.i =
00320 alpha->r * temp2.i + alpha->i * temp2.r;
00321 z__1.r = z__2.r + z__4.r, z__1.i = z__2.i + z__4.i;
00322 y[i__2].r = z__1.r, y[i__2].i = z__1.i;
00323 jx += *incx;
00324 jy += *incy;
00325 kk += j;
00326
00327 }
00328 }
00329 } else {
00330
00331
00332
00333 if (*incx == 1 && *incy == 1) {
00334 i__1 = *n;
00335 for (j = 1; j <= i__1; ++j) {
00336 i__2 = j;
00337 z__1.r = alpha->r * x[i__2].r - alpha->i * x[i__2].i, z__1.i =
00338 alpha->r * x[i__2].i + alpha->i * x[i__2].r;
00339 temp1.r = z__1.r, temp1.i = z__1.i;
00340 temp2.r = 0., temp2.i = 0.;
00341 i__2 = j;
00342 i__3 = j;
00343 i__4 = kk;
00344 d__1 = ap[i__4].r;
00345 z__2.r = d__1 * temp1.r, z__2.i = d__1 * temp1.i;
00346 z__1.r = y[i__3].r + z__2.r, z__1.i = y[i__3].i + z__2.i;
00347 y[i__2].r = z__1.r, y[i__2].i = z__1.i;
00348 k = kk + 1;
00349 i__2 = *n;
00350 for (i__ = j + 1; i__ <= i__2; ++i__) {
00351 i__3 = i__;
00352 i__4 = i__;
00353 i__5 = k;
00354 z__2.r = temp1.r * ap[i__5].r - temp1.i * ap[i__5].i,
00355 z__2.i = temp1.r * ap[i__5].i + temp1.i * ap[i__5]
00356 .r;
00357 z__1.r = y[i__4].r + z__2.r, z__1.i = y[i__4].i + z__2.i;
00358 y[i__3].r = z__1.r, y[i__3].i = z__1.i;
00359 d_cnjg(&z__3, &ap[k]);
00360 i__3 = i__;
00361 z__2.r = z__3.r * x[i__3].r - z__3.i * x[i__3].i, z__2.i =
00362 z__3.r * x[i__3].i + z__3.i * x[i__3].r;
00363 z__1.r = temp2.r + z__2.r, z__1.i = temp2.i + z__2.i;
00364 temp2.r = z__1.r, temp2.i = z__1.i;
00365 ++k;
00366
00367 }
00368 i__2 = j;
00369 i__3 = j;
00370 z__2.r = alpha->r * temp2.r - alpha->i * temp2.i, z__2.i =
00371 alpha->r * temp2.i + alpha->i * temp2.r;
00372 z__1.r = y[i__3].r + z__2.r, z__1.i = y[i__3].i + z__2.i;
00373 y[i__2].r = z__1.r, y[i__2].i = z__1.i;
00374 kk += *n - j + 1;
00375
00376 }
00377 } else {
00378 jx = kx;
00379 jy = ky;
00380 i__1 = *n;
00381 for (j = 1; j <= i__1; ++j) {
00382 i__2 = jx;
00383 z__1.r = alpha->r * x[i__2].r - alpha->i * x[i__2].i, z__1.i =
00384 alpha->r * x[i__2].i + alpha->i * x[i__2].r;
00385 temp1.r = z__1.r, temp1.i = z__1.i;
00386 temp2.r = 0., temp2.i = 0.;
00387 i__2 = jy;
00388 i__3 = jy;
00389 i__4 = kk;
00390 d__1 = ap[i__4].r;
00391 z__2.r = d__1 * temp1.r, z__2.i = d__1 * temp1.i;
00392 z__1.r = y[i__3].r + z__2.r, z__1.i = y[i__3].i + z__2.i;
00393 y[i__2].r = z__1.r, y[i__2].i = z__1.i;
00394 ix = jx;
00395 iy = jy;
00396 i__2 = kk + *n - j;
00397 for (k = kk + 1; k <= i__2; ++k) {
00398 ix += *incx;
00399 iy += *incy;
00400 i__3 = iy;
00401 i__4 = iy;
00402 i__5 = k;
00403 z__2.r = temp1.r * ap[i__5].r - temp1.i * ap[i__5].i,
00404 z__2.i = temp1.r * ap[i__5].i + temp1.i * ap[i__5]
00405 .r;
00406 z__1.r = y[i__4].r + z__2.r, z__1.i = y[i__4].i + z__2.i;
00407 y[i__3].r = z__1.r, y[i__3].i = z__1.i;
00408 d_cnjg(&z__3, &ap[k]);
00409 i__3 = ix;
00410 z__2.r = z__3.r * x[i__3].r - z__3.i * x[i__3].i, z__2.i =
00411 z__3.r * x[i__3].i + z__3.i * x[i__3].r;
00412 z__1.r = temp2.r + z__2.r, z__1.i = temp2.i + z__2.i;
00413 temp2.r = z__1.r, temp2.i = z__1.i;
00414
00415 }
00416 i__2 = jy;
00417 i__3 = jy;
00418 z__2.r = alpha->r * temp2.r - alpha->i * temp2.i, z__2.i =
00419 alpha->r * temp2.i + alpha->i * temp2.r;
00420 z__1.r = y[i__3].r + z__2.r, z__1.i = y[i__3].i + z__2.i;
00421 y[i__2].r = z__1.r, y[i__2].i = z__1.i;
00422 jx += *incx;
00423 jy += *incy;
00424 kk += *n - j + 1;
00425
00426 }
00427 }
00428 }
00429
00430 return 0;
00431
00432
00433
00434 }