123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580581582583584585586587588589590591592593594595596597598599600601602603604605606607608609610611612613614615616617618619620621622623624625626627628629630631632633634635636637638639640641642643644645646647 |
- #include "datatypes.h"
- int ctbmv_(char *uplo, char *trans, char *diag, integer *n,
- integer *k, complex *a, integer *lda, complex *x, integer *incx,
- ftnlen uplo_len, ftnlen trans_len, ftnlen diag_len)
- {
-
- integer a_dim1, a_offset, i__1, i__2, i__3, i__4, i__5;
- complex q__1, q__2, q__3;
-
- void r_cnjg(complex *, complex *);
-
- integer i__, j, l, ix, jx, kx, info;
- complex temp;
- extern logical lsame_(char *, char *, ftnlen, ftnlen);
- integer kplus1;
- extern int xerbla_(char *, integer *, ftnlen);
- logical noconj, nounit;
-
- a_dim1 = *lda;
- a_offset = 1 + a_dim1;
- a -= a_offset;
- --x;
-
- info = 0;
- if (! lsame_(uplo, "U", (ftnlen)1, (ftnlen)1) && ! lsame_(uplo, "L", (
- ftnlen)1, (ftnlen)1)) {
- info = 1;
- } else if (! lsame_(trans, "N", (ftnlen)1, (ftnlen)1) && ! lsame_(trans,
- "T", (ftnlen)1, (ftnlen)1) && ! lsame_(trans, "C", (ftnlen)1, (
- ftnlen)1)) {
- info = 2;
- } else if (! lsame_(diag, "U", (ftnlen)1, (ftnlen)1) && ! lsame_(diag,
- "N", (ftnlen)1, (ftnlen)1)) {
- info = 3;
- } else if (*n < 0) {
- info = 4;
- } else if (*k < 0) {
- info = 5;
- } else if (*lda < *k + 1) {
- info = 7;
- } else if (*incx == 0) {
- info = 9;
- }
- if (info != 0) {
- xerbla_("CTBMV ", &info, (ftnlen)6);
- return 0;
- }
- if (*n == 0) {
- return 0;
- }
- noconj = lsame_(trans, "T", (ftnlen)1, (ftnlen)1);
- nounit = lsame_(diag, "N", (ftnlen)1, (ftnlen)1);
- if (*incx <= 0) {
- kx = 1 - (*n - 1) * *incx;
- } else if (*incx != 1) {
- kx = 1;
- }
- if (lsame_(trans, "N", (ftnlen)1, (ftnlen)1)) {
- if (lsame_(uplo, "U", (ftnlen)1, (ftnlen)1)) {
- kplus1 = *k + 1;
- if (*incx == 1) {
- i__1 = *n;
- for (j = 1; j <= i__1; ++j) {
- i__2 = j;
- if (x[i__2].r != 0.f || x[i__2].i != 0.f) {
- i__2 = j;
- temp.r = x[i__2].r, temp.i = x[i__2].i;
- l = kplus1 - j;
- i__2 = 1, i__3 = j - *k;
- i__4 = j - 1;
- for (i__ = max(i__2,i__3); i__ <= i__4; ++i__) {
- i__2 = i__;
- i__3 = i__;
- i__5 = l + i__ + j * a_dim1;
- q__2.r = temp.r * a[i__5].r - temp.i * a[i__5].i,
- q__2.i = temp.r * a[i__5].i + temp.i * a[
- i__5].r;
- q__1.r = x[i__3].r + q__2.r, q__1.i = x[i__3].i +
- q__2.i;
- x[i__2].r = q__1.r, x[i__2].i = q__1.i;
- }
- if (nounit) {
- i__4 = j;
- i__2 = j;
- i__3 = kplus1 + j * a_dim1;
- q__1.r = x[i__2].r * a[i__3].r - x[i__2].i * a[
- i__3].i, q__1.i = x[i__2].r * a[i__3].i +
- x[i__2].i * a[i__3].r;
- x[i__4].r = q__1.r, x[i__4].i = q__1.i;
- }
- }
- }
- } else {
- jx = kx;
- i__1 = *n;
- for (j = 1; j <= i__1; ++j) {
- i__4 = jx;
- if (x[i__4].r != 0.f || x[i__4].i != 0.f) {
- i__4 = jx;
- temp.r = x[i__4].r, temp.i = x[i__4].i;
- ix = kx;
- l = kplus1 - j;
- i__4 = 1, i__2 = j - *k;
- i__3 = j - 1;
- for (i__ = max(i__4,i__2); i__ <= i__3; ++i__) {
- i__4 = ix;
- i__2 = ix;
- i__5 = l + i__ + j * a_dim1;
- q__2.r = temp.r * a[i__5].r - temp.i * a[i__5].i,
- q__2.i = temp.r * a[i__5].i + temp.i * a[
- i__5].r;
- q__1.r = x[i__2].r + q__2.r, q__1.i = x[i__2].i +
- q__2.i;
- x[i__4].r = q__1.r, x[i__4].i = q__1.i;
- ix += *incx;
- }
- if (nounit) {
- i__3 = jx;
- i__4 = jx;
- i__2 = kplus1 + j * a_dim1;
- q__1.r = x[i__4].r * a[i__2].r - x[i__4].i * a[
- i__2].i, q__1.i = x[i__4].r * a[i__2].i +
- x[i__4].i * a[i__2].r;
- x[i__3].r = q__1.r, x[i__3].i = q__1.i;
- }
- }
- jx += *incx;
- if (j > *k) {
- kx += *incx;
- }
- }
- }
- } else {
- if (*incx == 1) {
- for (j = *n; j >= 1; --j) {
- i__1 = j;
- if (x[i__1].r != 0.f || x[i__1].i != 0.f) {
- i__1 = j;
- temp.r = x[i__1].r, temp.i = x[i__1].i;
- l = 1 - j;
- i__1 = *n, i__3 = j + *k;
- i__4 = j + 1;
- for (i__ = min(i__1,i__3); i__ >= i__4; --i__) {
- i__1 = i__;
- i__3 = i__;
- i__2 = l + i__ + j * a_dim1;
- q__2.r = temp.r * a[i__2].r - temp.i * a[i__2].i,
- q__2.i = temp.r * a[i__2].i + temp.i * a[
- i__2].r;
- q__1.r = x[i__3].r + q__2.r, q__1.i = x[i__3].i +
- q__2.i;
- x[i__1].r = q__1.r, x[i__1].i = q__1.i;
- }
- if (nounit) {
- i__4 = j;
- i__1 = j;
- i__3 = j * a_dim1 + 1;
- q__1.r = x[i__1].r * a[i__3].r - x[i__1].i * a[
- i__3].i, q__1.i = x[i__1].r * a[i__3].i +
- x[i__1].i * a[i__3].r;
- x[i__4].r = q__1.r, x[i__4].i = q__1.i;
- }
- }
- }
- } else {
- kx += (*n - 1) * *incx;
- jx = kx;
- for (j = *n; j >= 1; --j) {
- i__4 = jx;
- if (x[i__4].r != 0.f || x[i__4].i != 0.f) {
- i__4 = jx;
- temp.r = x[i__4].r, temp.i = x[i__4].i;
- ix = kx;
- l = 1 - j;
- i__4 = *n, i__1 = j + *k;
- i__3 = j + 1;
- for (i__ = min(i__4,i__1); i__ >= i__3; --i__) {
- i__4 = ix;
- i__1 = ix;
- i__2 = l + i__ + j * a_dim1;
- q__2.r = temp.r * a[i__2].r - temp.i * a[i__2].i,
- q__2.i = temp.r * a[i__2].i + temp.i * a[
- i__2].r;
- q__1.r = x[i__1].r + q__2.r, q__1.i = x[i__1].i +
- q__2.i;
- x[i__4].r = q__1.r, x[i__4].i = q__1.i;
- ix -= *incx;
- }
- if (nounit) {
- i__3 = jx;
- i__4 = jx;
- i__1 = j * a_dim1 + 1;
- q__1.r = x[i__4].r * a[i__1].r - x[i__4].i * a[
- i__1].i, q__1.i = x[i__4].r * a[i__1].i +
- x[i__4].i * a[i__1].r;
- x[i__3].r = q__1.r, x[i__3].i = q__1.i;
- }
- }
- jx -= *incx;
- if (*n - j >= *k) {
- kx -= *incx;
- }
- }
- }
- }
- } else {
- if (lsame_(uplo, "U", (ftnlen)1, (ftnlen)1)) {
- kplus1 = *k + 1;
- if (*incx == 1) {
- for (j = *n; j >= 1; --j) {
- i__3 = j;
- temp.r = x[i__3].r, temp.i = x[i__3].i;
- l = kplus1 - j;
- if (noconj) {
- if (nounit) {
- i__3 = kplus1 + j * a_dim1;
- q__1.r = temp.r * a[i__3].r - temp.i * a[i__3].i,
- q__1.i = temp.r * a[i__3].i + temp.i * a[
- i__3].r;
- temp.r = q__1.r, temp.i = q__1.i;
- }
- i__4 = 1, i__1 = j - *k;
- i__3 = max(i__4,i__1);
- for (i__ = j - 1; i__ >= i__3; --i__) {
- i__4 = l + i__ + j * a_dim1;
- i__1 = i__;
- q__2.r = a[i__4].r * x[i__1].r - a[i__4].i * x[
- i__1].i, q__2.i = a[i__4].r * x[i__1].i +
- a[i__4].i * x[i__1].r;
- q__1.r = temp.r + q__2.r, q__1.i = temp.i +
- q__2.i;
- temp.r = q__1.r, temp.i = q__1.i;
- }
- } else {
- if (nounit) {
- r_cnjg(&q__2, &a[kplus1 + j * a_dim1]);
- q__1.r = temp.r * q__2.r - temp.i * q__2.i,
- q__1.i = temp.r * q__2.i + temp.i *
- q__2.r;
- temp.r = q__1.r, temp.i = q__1.i;
- }
- i__4 = 1, i__1 = j - *k;
- i__3 = max(i__4,i__1);
- for (i__ = j - 1; i__ >= i__3; --i__) {
- r_cnjg(&q__3, &a[l + i__ + j * a_dim1]);
- i__4 = i__;
- q__2.r = q__3.r * x[i__4].r - q__3.i * x[i__4].i,
- q__2.i = q__3.r * x[i__4].i + q__3.i * x[
- i__4].r;
- q__1.r = temp.r + q__2.r, q__1.i = temp.i +
- q__2.i;
- temp.r = q__1.r, temp.i = q__1.i;
- }
- }
- i__3 = j;
- x[i__3].r = temp.r, x[i__3].i = temp.i;
- }
- } else {
- kx += (*n - 1) * *incx;
- jx = kx;
- for (j = *n; j >= 1; --j) {
- i__3 = jx;
- temp.r = x[i__3].r, temp.i = x[i__3].i;
- kx -= *incx;
- ix = kx;
- l = kplus1 - j;
- if (noconj) {
- if (nounit) {
- i__3 = kplus1 + j * a_dim1;
- q__1.r = temp.r * a[i__3].r - temp.i * a[i__3].i,
- q__1.i = temp.r * a[i__3].i + temp.i * a[
- i__3].r;
- temp.r = q__1.r, temp.i = q__1.i;
- }
- i__4 = 1, i__1 = j - *k;
- i__3 = max(i__4,i__1);
- for (i__ = j - 1; i__ >= i__3; --i__) {
- i__4 = l + i__ + j * a_dim1;
- i__1 = ix;
- q__2.r = a[i__4].r * x[i__1].r - a[i__4].i * x[
- i__1].i, q__2.i = a[i__4].r * x[i__1].i +
- a[i__4].i * x[i__1].r;
- q__1.r = temp.r + q__2.r, q__1.i = temp.i +
- q__2.i;
- temp.r = q__1.r, temp.i = q__1.i;
- ix -= *incx;
- }
- } else {
- if (nounit) {
- r_cnjg(&q__2, &a[kplus1 + j * a_dim1]);
- q__1.r = temp.r * q__2.r - temp.i * q__2.i,
- q__1.i = temp.r * q__2.i + temp.i *
- q__2.r;
- temp.r = q__1.r, temp.i = q__1.i;
- }
- i__4 = 1, i__1 = j - *k;
- i__3 = max(i__4,i__1);
- for (i__ = j - 1; i__ >= i__3; --i__) {
- r_cnjg(&q__3, &a[l + i__ + j * a_dim1]);
- i__4 = ix;
- q__2.r = q__3.r * x[i__4].r - q__3.i * x[i__4].i,
- q__2.i = q__3.r * x[i__4].i + q__3.i * x[
- i__4].r;
- q__1.r = temp.r + q__2.r, q__1.i = temp.i +
- q__2.i;
- temp.r = q__1.r, temp.i = q__1.i;
- ix -= *incx;
- }
- }
- i__3 = jx;
- x[i__3].r = temp.r, x[i__3].i = temp.i;
- jx -= *incx;
- }
- }
- } else {
- if (*incx == 1) {
- i__3 = *n;
- for (j = 1; j <= i__3; ++j) {
- i__4 = j;
- temp.r = x[i__4].r, temp.i = x[i__4].i;
- l = 1 - j;
- if (noconj) {
- if (nounit) {
- i__4 = j * a_dim1 + 1;
- q__1.r = temp.r * a[i__4].r - temp.i * a[i__4].i,
- q__1.i = temp.r * a[i__4].i + temp.i * a[
- i__4].r;
- temp.r = q__1.r, temp.i = q__1.i;
- }
- i__1 = *n, i__2 = j + *k;
- i__4 = min(i__1,i__2);
- for (i__ = j + 1; i__ <= i__4; ++i__) {
- i__1 = l + i__ + j * a_dim1;
- i__2 = i__;
- q__2.r = a[i__1].r * x[i__2].r - a[i__1].i * x[
- i__2].i, q__2.i = a[i__1].r * x[i__2].i +
- a[i__1].i * x[i__2].r;
- q__1.r = temp.r + q__2.r, q__1.i = temp.i +
- q__2.i;
- temp.r = q__1.r, temp.i = q__1.i;
- }
- } else {
- if (nounit) {
- r_cnjg(&q__2, &a[j * a_dim1 + 1]);
- q__1.r = temp.r * q__2.r - temp.i * q__2.i,
- q__1.i = temp.r * q__2.i + temp.i *
- q__2.r;
- temp.r = q__1.r, temp.i = q__1.i;
- }
- i__1 = *n, i__2 = j + *k;
- i__4 = min(i__1,i__2);
- for (i__ = j + 1; i__ <= i__4; ++i__) {
- r_cnjg(&q__3, &a[l + i__ + j * a_dim1]);
- i__1 = i__;
- q__2.r = q__3.r * x[i__1].r - q__3.i * x[i__1].i,
- q__2.i = q__3.r * x[i__1].i + q__3.i * x[
- i__1].r;
- q__1.r = temp.r + q__2.r, q__1.i = temp.i +
- q__2.i;
- temp.r = q__1.r, temp.i = q__1.i;
- }
- }
- i__4 = j;
- x[i__4].r = temp.r, x[i__4].i = temp.i;
- }
- } else {
- jx = kx;
- i__3 = *n;
- for (j = 1; j <= i__3; ++j) {
- i__4 = jx;
- temp.r = x[i__4].r, temp.i = x[i__4].i;
- kx += *incx;
- ix = kx;
- l = 1 - j;
- if (noconj) {
- if (nounit) {
- i__4 = j * a_dim1 + 1;
- q__1.r = temp.r * a[i__4].r - temp.i * a[i__4].i,
- q__1.i = temp.r * a[i__4].i + temp.i * a[
- i__4].r;
- temp.r = q__1.r, temp.i = q__1.i;
- }
- i__1 = *n, i__2 = j + *k;
- i__4 = min(i__1,i__2);
- for (i__ = j + 1; i__ <= i__4; ++i__) {
- i__1 = l + i__ + j * a_dim1;
- i__2 = ix;
- q__2.r = a[i__1].r * x[i__2].r - a[i__1].i * x[
- i__2].i, q__2.i = a[i__1].r * x[i__2].i +
- a[i__1].i * x[i__2].r;
- q__1.r = temp.r + q__2.r, q__1.i = temp.i +
- q__2.i;
- temp.r = q__1.r, temp.i = q__1.i;
- ix += *incx;
- }
- } else {
- if (nounit) {
- r_cnjg(&q__2, &a[j * a_dim1 + 1]);
- q__1.r = temp.r * q__2.r - temp.i * q__2.i,
- q__1.i = temp.r * q__2.i + temp.i *
- q__2.r;
- temp.r = q__1.r, temp.i = q__1.i;
- }
- i__1 = *n, i__2 = j + *k;
- i__4 = min(i__1,i__2);
- for (i__ = j + 1; i__ <= i__4; ++i__) {
- r_cnjg(&q__3, &a[l + i__ + j * a_dim1]);
- i__1 = ix;
- q__2.r = q__3.r * x[i__1].r - q__3.i * x[i__1].i,
- q__2.i = q__3.r * x[i__1].i + q__3.i * x[
- i__1].r;
- q__1.r = temp.r + q__2.r, q__1.i = temp.i +
- q__2.i;
- temp.r = q__1.r, temp.i = q__1.i;
- ix += *incx;
- }
- }
- i__4 = jx;
- x[i__4].r = temp.r, x[i__4].i = temp.i;
- jx += *incx;
- }
- }
- }
- }
- return 0;
- }
|