35 #ifndef TEMPLATE_LAPACK_LARRE_HEADER
36 #define TEMPLATE_LAPACK_LARRE_HEADER
42 *e, Treal *e2, Treal *rtol1, Treal *rtol2, Treal *
45 Treal *gers, Treal *pivmin, Treal *work,
integer *
50 Treal d__1, d__2, d__3;
61 Treal eps, tau, tmp, rtl;
83 Treal isrght, bsrtol, dpivot;
312 if (irange == 1 || ( irange == 3 && d__[1] > *vl && d__[1] <= *vu ) ||
313 ( irange == 2 && *il == 1 && *iu == 1 ) ) {
338 for (i__ = 1; i__ <= i__1; ++i__) {
341 eabs = (d__1 = e[i__],
absMACRO(d__1));
346 gers[(i__ << 1) - 1] = d__[i__] - tmp1;
348 d__1 = gl, d__2 = gers[(i__ << 1) - 1];
350 gers[i__ * 2] = d__[i__] + tmp1;
352 d__1 =
gu, d__2 = gers[i__ * 2];
361 d__1 = 1., d__2 = d__3 * d__3;
362 *pivmin = safmin *
maxMACRO(d__1,d__2);
374 usedqd = irange == 1 && ! forceb;
375 if (irange == 1 && ! forceb) {
386 template_lapack_larrd(range,
"B", n, vl, vu, il, iu, &gers[1], &bsrtol, &d__[1], &e[
387 1], &e2[1], pivmin, nsplit, &isplit[1], &mm, &w[1], &werr[1],
388 vl, vu, &iblock[1], &indexw[1], &work[1], &iwork[1], &iinfo);
395 for (i__ = mm + 1; i__ <= i__1; ++i__) {
408 for (jblk = 1; jblk <= i__1; ++jblk) {
410 in = iend - ibegin + 1;
413 if (irange == 1 || ( irange == 3 && d__[ibegin] > *vl && d__[ibegin]
414 <= *vu ) || ( irange == 2 && iblock[wbegin] == jblk ) ) {
440 for (i__ = ibegin; i__ <= i__2; ++i__) {
442 d__1 = gers[(i__ << 1) - 1];
445 d__1 = gers[i__ * 2];
450 if (! (irange == 1 && ! forceb)) {
454 for (i__ = wbegin; i__ <= i__2; ++i__) {
455 if (iblock[i__] == jblk) {
471 usedqd = (Treal) mb > in * .5 && ! forceb;
472 wend = wbegin + mb - 1;
478 for (i__ = wbegin; i__ <= i__2; ++i__) {
480 d__1 = 0., d__2 = w[i__ + 1] - werr[i__ + 1] - (w[i__] +
486 d__1 = 0., d__2 = *vu - sigma - (w[wend] + werr[wend]);
489 indl = indexw[wbegin];
493 if ( ( irange == 1 && ! forceb ) || usedqd) {
497 rtl, &tmp, &tmp1, &iinfo);
503 d__2 = gl, d__3 = tmp - tmp1 - eps * 100. * (d__1 = tmp - tmp1,
507 rtl, &tmp, &tmp1, &iinfo);
513 d__2 =
gu, d__3 = tmp + tmp1 + eps * 100. * (d__1 = tmp + tmp1,
517 spdiam = isrght - isleft;
522 d__2 = gl, d__3 = w[wbegin] - werr[wbegin] - eps * 100. * (d__1 =
523 w[wbegin] - werr[wbegin],
absMACRO(d__1));
526 d__2 =
gu, d__3 = w[wend] + werr[wend] + eps * 100. * (d__1 = w[
527 wend] + werr[wend],
absMACRO(d__1));
538 if (irange == 1 && ! forceb) {
546 wend = wbegin + mb - 1;
548 s1 = isleft + spdiam * .25;
549 s2 = isrght - spdiam * .25;
555 s1 = isleft + spdiam * .25;
556 s2 = isrght - spdiam * .25;
559 s1 =
maxMACRO(isleft,*vl) + tmp * .25;
560 s2 =
minMACRO(isrght,*vu) - tmp * .25;
566 cnt, &cnt1, &cnt2, &iinfo);
571 }
else if (cnt1 - indl >= indu - cnt2) {
572 if (irange == 1 && ! forceb) {
585 if (irange == 1 && ! forceb) {
606 tau = spdiam * eps * *n + *pivmin * 2.;
609 clwdth = w[wend] + werr[wend] - w[wbegin] - werr[wbegin];
610 avgap = (d__1 = clwdth / (Treal) (wend - wbegin),
absMACRO(
617 d__1 = tau, d__2 = werr[wbegin];
621 d__1 = wgap[wend - 1];
624 d__1 = tau, d__2 = werr[wend];
632 for (idum = 1; idum <= 6; ++idum) {
636 dpivot = d__[ibegin] - sigma;
641 for (i__ = 1; i__ <= i__2; ++i__) {
642 work[(in << 1) + i__] = 1. / work[i__];
643 tmp = e[j] * work[(in << 1) + i__];
644 work[in + i__] = tmp;
645 dpivot = d__[j + 1] - sigma - tmp * e[j];
646 work[i__ + 1] = dpivot;
648 d__1 = dmax__, d__2 =
absMACRO(dpivot);
654 if (dmax__ > spdiam * 64.) {
659 if (usedqd && ! norep) {
663 for (i__ = 1; i__ <= i__2; ++i__) {
664 tmp = sgndef * work[i__];
678 sigma = gl - spdiam * 2. * eps * *n - *pivmin * 4.;
680 sigma = gu + spdiam * 2. * eps * *n + *pivmin * 4.;
683 sigma -= sgndef * tau;
711 for (i__ = 1; i__ <= 4; ++i__) {
715 i__2 = (in << 1) - 1;
718 for (i__ = 1; i__ <= i__2; ++i__) {
719 d__[ibegin + i__ - 1] *= eps * 8. * work[i__] + 1.;
720 e[ibegin + i__ - 1] *= eps * 8. * work[in + i__] + 1.;
723 d__[iend] *= eps * 4. * work[in] + 1.;
739 for (j = wbegin; j <= i__2; ++j) {
741 werr[j] += (d__1 = w[j],
absMACRO(d__1)) * eps;
747 for (i__ = ibegin; i__ <= i__2; ++i__) {
750 work[i__] = d__[i__] * (d__1 * d__1);
756 rtol2, &i__2, &w[wbegin], &wgap[wbegin], &werr[wbegin], &
757 work[(*n << 1) + 1], &iwork[1], pivmin, &spdiam, &in, &
766 d__1 = 0., d__2 = *vu - sigma - (w[wend] + werr[wend]);
769 for (i__ = indl; i__ <= i__2; ++i__) {
790 for (i__ = 1; i__ <= i__2; ++i__) {
791 work[(i__ << 1) - 1] = (d__1 = d__[j],
absMACRO(d__1));
792 work[i__ * 2] = e[j] * e[j] * work[(i__ << 1) - 1];
796 work[(in << 1) - 1] = (d__1 = d__[iend],
absMACRO(d__1));
808 for (i__ = 1; i__ <= i__2; ++i__) {
809 if (work[i__] < 0.) {
818 for (i__ = indl; i__ <= i__2; ++i__) {
820 w[*m] = work[in - i__ + 1];
827 for (i__ = indl; i__ <= i__2; ++i__) {
836 for (i__ = *m - mb + 1; i__ <= i__2; ++i__) {
838 werr[i__] = rtol * (d__1 = w[i__],
absMACRO(d__1));
842 for (i__ = *m - mb + 1; i__ <= i__2; ++i__) {
845 d__1 = 0., d__2 = w[i__ + 1] - werr[i__ + 1] - (w[i__] + werr[
851 d__1 = 0., d__2 = *vu - sigma - (w[*m] + werr[*m]);