37 #ifndef TEMPLATE_LAPACK_TGEVC_HEADER 38 #define TEMPLATE_LAPACK_TGEVC_HEADER 262 integer a_dim1, a_offset, b_dim1, b_offset, vl_dim1, vl_offset, vr_dim1,
263 vr_offset, i__1, i__2, i__3, i__4, i__5;
264 Treal d__1, d__2, d__3, d__4, d__5, d__6;
267 Treal dmin__, temp, suma[4] , sumb[4]
269 Treal cim2a, cim2b, cre2a, cre2b, temp2, bdiag[2];
286 Treal acoefa, bcoefa, cimaga, cimagb;
289 Treal bcoefi, ascale, bscale, creala;
294 Treal salfar, safmin;
295 Treal xscale, bignum;
301 #define suma_ref(a_1,a_2) suma[(a_2)*2 + a_1 - 3] 302 #define sumb_ref(a_1,a_2) sumb[(a_2)*2 + a_1 - 3] 303 #define a_ref(a_1,a_2) a[(a_2)*a_dim1 + a_1] 304 #define b_ref(a_1,a_2) b[(a_2)*b_dim1 + a_1] 305 #define vl_ref(a_1,a_2) vl[(a_2)*vl_dim1 + a_1] 306 #define vr_ref(a_1,a_2) vr[(a_2)*vr_dim1 + a_1] 307 #define sum_ref(a_1,a_2) sum[(a_2)*2 + a_1 - 3] 312 a_offset = 1 + a_dim1 * 1;
315 b_offset = 1 + b_dim1 * 1;
318 vl_offset = 1 + vl_dim1 * 1;
321 vr_offset = 1 + vr_dim1 * 1;
365 }
else if (ihwmny < 0) {
386 for (j = 1; j <= i__1; ++j) {
392 if (
a_ref(j + 1, j) != 0.) {
397 if (select[j] || select[j + 1]) {
417 for (j = 1; j <= i__1; ++j) {
418 if (
a_ref(j + 1, j) != 0.) {
424 if (
a_ref(j + 2, j + 1) != 0.) {
436 }
else if ( ( compl_AAAA && *ldvl < *n ) || *ldvl < 1) {
438 }
else if ( ( compr && *ldvr < *n ) || *ldvr < 1) {
440 }
else if (*mm < im) {
462 small = safmin * *n / ulp;
464 bignum = 1. / (safmin * *n);
480 for (j = 2; j <= i__1; ++j) {
483 if (
a_ref(j, j - 1) == 0.) {
489 for (i__ = 1; i__ <= i__2; ++i__) {
495 work[*n + j] = temp2;
499 for (i__ = iend + 1; i__ <= i__2; ++i__) {
509 ascale = 1. /
maxMACRO(anorm,safmin);
510 bscale = 1. /
maxMACRO(bnorm,safmin);
521 for (je = 1; je <= i__1; ++je) {
534 if (
a_ref(je + 1, je) != 0.) {
542 ilcomp = select[je] || select[je + 1];
554 if ((d__1 =
a_ref(je, je),
absMACRO(d__1)) <= safmin && (d__2 =
561 for (jr = 1; jr <= i__2; ++jr) {
573 for (jr = 1; jr <= i__2; ++jr) {
574 work[(*n << 1) + jr] = 0.;
587 d__3 = (d__1 =
a_ref(je, je),
absMACRO(d__1)) * ascale, d__4 = (
591 salfar = temp *
a_ref(je, je) * ascale;
592 sbeta = temp *
b_ref(je, je) * bscale;
593 acoef = sbeta * ascale;
594 bcoefr = salfar * bscale;
615 d__1 = scale, d__2 = 1. / (safmin *
maxMACRO(d__3,d__4));
618 acoef = ascale * (scale * sbeta);
620 acoef = scale * acoef;
623 bcoefr = bscale * (scale * salfar);
625 bcoefr = scale * bcoefr;
633 work[(*n << 1) + je] = 1.;
639 d__1 = safmin * 100.;
641 acoef, &temp, &bcoefr, &temp2, &bcoefi);
653 if (acoefa * ulp < safmin && acoefa >= safmin) {
654 scale = safmin / ulp / acoefa;
656 if (bcoefa * ulp < safmin && bcoefa >= safmin) {
658 d__1 = scale, d__2 = safmin / ulp / bcoefa;
661 if (safmin * acoefa > ascale) {
662 scale = ascale / (safmin * acoefa);
664 if (safmin * bcoefa > bscale) {
666 d__1 = scale, d__2 = bscale / (safmin * bcoefa);
670 acoef = scale * acoef;
672 bcoefr = scale * bcoefr;
673 bcoefi = scale * bcoefi;
679 temp = acoef *
a_ref(je + 1, je);
680 temp2r = acoef *
a_ref(je, je) - bcoefr *
b_ref(je, je);
681 temp2i = -bcoefi *
b_ref(je, je);
683 work[(*n << 1) + je] = 1.;
684 work[*n * 3 + je] = 0.;
685 work[(*n << 1) + je + 1] = -temp2r / temp;
686 work[*n * 3 + je + 1] = -temp2i / temp;
688 work[(*n << 1) + je + 1] = 1.;
689 work[*n * 3 + je + 1] = 0.;
690 temp = acoef *
a_ref(je, je + 1);
691 work[(*n << 1) + je] = (bcoefr *
b_ref(je + 1, je + 1) -
692 acoef *
a_ref(je + 1, je + 1)) / temp;
693 work[*n * 3 + je] = bcoefi *
b_ref(je + 1, je + 1) / temp;
696 d__5 = (d__1 = work[(*n << 1) + je],
absMACRO(d__1)) + (d__2 =
697 work[*n * 3 + je],
absMACRO(d__2)), d__6 = (d__3 = work[(*
698 n << 1) + je + 1],
absMACRO(d__3)) + (d__4 = work[*n * 3 +
704 d__1 = ulp * acoefa * anorm, d__2 = ulp * bcoefa * bnorm, d__1 =
717 for (j = je + nw; j <= i__2; ++j) {
724 bdiag[0] =
b_ref(j, j);
726 if (
a_ref(j + 1, j) != 0.) {
728 bdiag[1] =
b_ref(j + 1, j + 1);
737 d__1 = work[j], d__2 = work[*n + j], d__1 =
maxMACRO(d__1,d__2),
738 d__2 = acoefa * work[j] + bcoefa * work[*n + j];
742 d__1 = temp, d__2 = work[j + 1], d__1 =
maxMACRO(d__1,d__2),
743 d__2 = work[*n + j + 1], d__1 =
maxMACRO(d__1,d__2),
744 d__2 = acoefa * work[j + 1] + bcoefa * work[*n +
748 if (temp > bignum * xscale) {
750 for (jw = 0; jw <= i__3; ++jw) {
752 for (jr = je; jr <= i__4; ++jr) {
753 work[(jw + 2) * *n + jr] = xscale * work[(jw + 2)
794 for (jw = 1; jw <= i__3; ++jw) {
809 for (ja = 1; ja <= i__4; ++ja) {
814 for (jr = je; jr <= i__5; ++jr) {
816 + ja - 1) * work[(jw + 1) * *n + jr];
818 + ja - 1) * work[(jw + 1) * *n + jr];
839 for (ja = 1; ja <= i__3; ++ja) {
857 bdiag, &bdiag[1], sum, &c__2, &bcoefr, &bcoefi, &
858 work[(*n << 1) + j], n, &scale, &temp, &iinfo);
861 for (jw = 0; jw <= i__3; ++jw) {
863 for (jr = je; jr <= i__4; ++jr) {
864 work[(jw + 2) * *n + jr] = scale * work[(jw + 2) *
883 for (jw = 0; jw <= i__2; ++jw) {
886 (jw + 2) * *n + je], &c__1, &c_b37, &work[(jw + 4)
904 for (j = ibeg; j <= i__2; ++j) {
913 for (j = ibeg; j <= i__2; ++j) {
925 for (jw = 0; jw <= i__2; ++jw) {
927 for (jr = ibeg; jr <= i__3; ++jr) {
935 ieig = ieig + nw - 1;
950 for (je = *n; je >= 1; --je) {
966 if (
a_ref(je, je - 1) != 0.) {
974 ilcomp = select[je] || select[je - 1];
986 if ((d__1 =
a_ref(je, je),
absMACRO(d__1)) <= safmin && (d__2 =
993 for (jr = 1; jr <= i__1; ++jr) {
1005 for (jw = 0; jw <= i__1; ++jw) {
1007 for (jr = 1; jr <= i__2; ++jr) {
1008 work[(jw + 2) * *n + jr] = 0.;
1023 d__3 = (d__1 =
a_ref(je, je),
absMACRO(d__1)) * ascale, d__4 = (
1027 salfar = temp *
a_ref(je, je) * ascale;
1028 sbeta = temp *
b_ref(je, je) * bscale;
1029 acoef = sbeta * ascale;
1030 bcoefr = salfar * bscale;
1051 d__1 = scale, d__2 = 1. / (safmin *
maxMACRO(d__3,d__4));
1054 acoef = ascale * (scale * sbeta);
1056 acoef = scale * acoef;
1059 bcoefr = bscale * (scale * salfar);
1061 bcoefr = scale * bcoefr;
1069 work[(*n << 1) + je] = 1.;
1076 for (jr = 1; jr <= i__1; ++jr) {
1077 work[(*n << 1) + jr] = bcoefr *
b_ref(jr, je) - acoef *
1085 d__1 = safmin * 100.;
1087 ldb, &d__1, &acoef, &temp, &bcoefr, &temp2, &bcoefi);
1098 if (acoefa * ulp < safmin && acoefa >= safmin) {
1099 scale = safmin / ulp / acoefa;
1101 if (bcoefa * ulp < safmin && bcoefa >= safmin) {
1103 d__1 = scale, d__2 = safmin / ulp / bcoefa;
1106 if (safmin * acoefa > ascale) {
1107 scale = ascale / (safmin * acoefa);
1109 if (safmin * bcoefa > bscale) {
1111 d__1 = scale, d__2 = bscale / (safmin * bcoefa);
1115 acoef = scale * acoef;
1117 bcoefr = scale * bcoefr;
1118 bcoefi = scale * bcoefi;
1125 temp = acoef *
a_ref(je, je - 1);
1126 temp2r = acoef *
a_ref(je, je) - bcoefr *
b_ref(je, je);
1127 temp2i = -bcoefi *
b_ref(je, je);
1129 work[(*n << 1) + je] = 1.;
1130 work[*n * 3 + je] = 0.;
1131 work[(*n << 1) + je - 1] = -temp2r / temp;
1132 work[*n * 3 + je - 1] = -temp2i / temp;
1134 work[(*n << 1) + je - 1] = 1.;
1135 work[*n * 3 + je - 1] = 0.;
1136 temp = acoef *
a_ref(je - 1, je);
1137 work[(*n << 1) + je] = (bcoefr *
b_ref(je - 1, je - 1) -
1138 acoef *
a_ref(je - 1, je - 1)) / temp;
1139 work[*n * 3 + je] = bcoefi *
b_ref(je - 1, je - 1) / temp;
1143 d__5 = (d__1 = work[(*n << 1) + je],
absMACRO(d__1)) + (d__2 =
1144 work[*n * 3 + je],
absMACRO(d__2)), d__6 = (d__3 = work[(*
1145 n << 1) + je - 1],
absMACRO(d__3)) + (d__4 = work[*n * 3 +
1152 creala = acoef * work[(*n << 1) + je - 1];
1153 cimaga = acoef * work[*n * 3 + je - 1];
1154 crealb = bcoefr * work[(*n << 1) + je - 1] - bcoefi * work[*n
1156 cimagb = bcoefi * work[(*n << 1) + je - 1] + bcoefr * work[*n
1158 cre2a = acoef * work[(*n << 1) + je];
1159 cim2a = acoef * work[*n * 3 + je];
1160 cre2b = bcoefr * work[(*n << 1) + je] - bcoefi * work[*n * 3
1162 cim2b = bcoefi * work[(*n << 1) + je] + bcoefr * work[*n * 3
1165 for (jr = 1; jr <= i__1; ++jr) {
1166 work[(*n << 1) + jr] = -creala *
a_ref(jr, je - 1) +
1167 crealb *
b_ref(jr, je - 1) - cre2a *
a_ref(jr, je)
1168 + cre2b *
b_ref(jr, je);
1169 work[*n * 3 + jr] = -cimaga *
a_ref(jr, je - 1) + cimagb *
1171 cim2b *
b_ref(jr, je);
1177 d__1 = ulp * acoefa * anorm, d__2 = ulp * bcoefa * bnorm, d__1 =
1184 for (j = je - nw; j >= 1; --j) {
1189 if (! il2by2 && j > 1) {
1190 if (
a_ref(j, j - 1) != 0.) {
1195 bdiag[0] =
b_ref(j, j);
1198 bdiag[1] =
b_ref(j + 1, j + 1);
1206 lda, bdiag, &bdiag[1], &work[(*n << 1) + j], n, &
1207 bcoefr, &bcoefi, sum, &c__2, &scale, &temp, &iinfo);
1211 for (jw = 0; jw <= i__1; ++jw) {
1213 for (jr = 1; jr <= i__2; ++jr) {
1214 work[(jw + 2) * *n + jr] = scale * work[(jw + 2) *
1222 d__1 = scale * xmax;
1226 for (jw = 1; jw <= i__1; ++jw) {
1228 for (ja = 1; ja <= i__2; ++ja) {
1229 work[(jw + 1) * *n + j + ja - 1] =
sum_ref(ja, jw);
1242 temp = acoefa * work[j] + bcoefa * work[*n + j];
1245 d__1 = temp, d__2 = acoefa * work[j + 1] + bcoefa *
1252 if (temp > bignum * xscale) {
1255 for (jw = 0; jw <= i__1; ++jw) {
1257 for (jr = 1; jr <= i__2; ++jr) {
1258 work[(jw + 2) * *n + jr] = xscale * work[(jw
1273 for (ja = 1; ja <= i__1; ++ja) {
1275 creala = acoef * work[(*n << 1) + j + ja - 1];
1276 cimaga = acoef * work[*n * 3 + j + ja - 1];
1277 crealb = bcoefr * work[(*n << 1) + j + ja - 1] -
1278 bcoefi * work[*n * 3 + j + ja - 1];
1279 cimagb = bcoefi * work[(*n << 1) + j + ja - 1] +
1280 bcoefr * work[*n * 3 + j + ja - 1];
1282 for (jr = 1; jr <= i__2; ++jr) {
1283 work[(*n << 1) + jr] = work[(*n << 1) + jr] -
1284 creala *
a_ref(jr, j + ja - 1) +
1285 crealb *
b_ref(jr, j + ja - 1);
1286 work[*n * 3 + jr] = work[*n * 3 + jr] -
1287 cimaga *
a_ref(jr, j + ja - 1) +
1288 cimagb *
b_ref(jr, j + ja - 1);
1292 creala = acoef * work[(*n << 1) + j + ja - 1];
1293 crealb = bcoefr * work[(*n << 1) + j + ja - 1];
1295 for (jr = 1; jr <= i__2; ++jr) {
1296 work[(*n << 1) + jr] = work[(*n << 1) + jr] -
1297 creala *
a_ref(jr, j + ja - 1) +
1298 crealb *
b_ref(jr, j + ja - 1);
1318 for (jw = 0; jw <= i__1; ++jw) {
1320 for (jr = 1; jr <= i__2; ++jr) {
1321 work[(jw + 4) * *n + jr] = work[(jw + 2) * *n + 1] *
1331 for (jc = 2; jc <= i__2; ++jc) {
1333 for (jr = 1; jr <= i__3; ++jr) {
1334 work[(jw + 4) * *n + jr] += work[(jw + 2) * *n +
1344 for (jw = 0; jw <= i__1; ++jw) {
1346 for (jr = 1; jr <= i__2; ++jr) {
1347 vr_ref(jr, ieig + jw) = work[(jw + 4) * *n + jr];
1356 for (jw = 0; jw <= i__1; ++jw) {
1358 for (jr = 1; jr <= i__2; ++jr) {
1359 vr_ref(jr, ieig + jw) = work[(jw + 2) * *n + jr];
1373 for (j = 1; j <= i__1; ++j) {
1382 for (j = 1; j <= i__1; ++j) {
1390 if (xmax > safmin) {
1393 for (jw = 0; jw <= i__1; ++jw) {
1395 for (jr = 1; jr <= i__2; ++jr) {
#define absMACRO(x)
Definition: template_blas_common.h:47
int template_lapack_lacpy(const char *uplo, const integer *m, const integer *n, const Treal *a, const integer *lda, Treal *b, const integer *ldb)
Definition: template_lapack_lacpy.h:42
int integer
Definition: template_blas_common.h:40
#define maxMACRO(a, b)
Definition: template_blas_common.h:45
int template_lapack_labad(Treal *small, Treal *large)
Definition: template_lapack_labad.h:42
#define sum_ref(a_1, a_2)
#define minMACRO(a, b)
Definition: template_blas_common.h:46
int template_lapack_lag2(const Treal *a, const integer *lda, const Treal *b, const integer *ldb, const Treal *safmin, Treal *scale1, Treal *scale2, Treal *wr1, Treal *wr2, Treal *wi)
Definition: template_lapack_lag2.h:42
int template_blas_erbla(const char *srname, integer *info)
Definition: template_blas_common.cc:146
#define suma_ref(a_1, a_2)
Treal template_lapack_lamch(const char *cmach, Treal dummyReal)
Definition: template_lapack_lamch.h:202
int template_blas_gemv(const char *trans, const integer *m, const integer *n, const Treal *alpha, const Treal *a, const integer *lda, const Treal *x, const integer *incx, const Treal *beta, Treal *y, const integer *incy)
Definition: template_blas_gemv.h:43
#define sumb_ref(a_1, a_2)
bool logical
Definition: template_blas_common.h:41
side
Definition: Matrix.h:75
#define TRUE_
Definition: template_lapack_common.h:42
int template_lapack_laln2(const logical *ltrans, const integer *na, const integer *nw, const Treal *smin, const Treal *ca, const Treal *a, const integer *lda, const Treal *d1, const Treal *d2, const Treal *b, const integer *ldb, const Treal *wr, const Treal *wi, Treal *x, const integer *ldx, Treal *scale, Treal *xnorm, integer *info)
Definition: template_lapack_laln2.h:42
#define FALSE_
Definition: template_lapack_common.h:43
logical template_blas_lsame(const char *ca, const char *cb)
Definition: template_blas_common.cc:46
int template_lapack_tgevc(const char *side, const char *howmny, const logical *select, const integer *n, const Treal *a, const integer *lda, const Treal *b, const integer *ldb, Treal *vl, const integer *ldvl, Treal *vr, const integer *ldvr, const integer *mm, integer *m, Treal *work, integer *info)
Definition: template_lapack_tgevc.h:46