3 /* Subroutine */ int slasq5_(integer *i0, integer *n0, real *z__, integer *pp,
4 real *tau, real *dmin__, real *dmin1, real *dmin2, real *dn, real *
5 dnm1, real *dnm2, logical *ieee)
7 /* System generated locals */
17 /* -- LAPACK auxiliary routine (version 3.1) -- */
18 /* Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd.. */
21 /* .. Scalar Arguments .. */
23 /* .. Array Arguments .. */
29 /* SLASQ5 computes one dqds transform in ping-pong form, one */
30 /* version for IEEE machines another for non IEEE machines. */
35 /* I0 (input) INTEGER */
38 /* N0 (input) INTEGER */
41 /* Z (input) REAL array, dimension ( 4*N ) */
42 /* Z holds the qd array. EMIN is stored in Z(4*N0) to avoid */
43 /* an extra argument. */
45 /* PP (input) INTEGER */
46 /* PP=0 for ping, PP=1 for pong. */
48 /* TAU (input) REAL */
49 /* This is the shift. */
51 /* DMIN (output) REAL */
52 /* Minimum value of d. */
54 /* DMIN1 (output) REAL */
55 /* Minimum value of d, excluding D( N0 ). */
57 /* DMIN2 (output) REAL */
58 /* Minimum value of d, excluding D( N0 ) and D( N0-1 ). */
60 /* DN (output) REAL */
61 /* d(N0), the last value of d. */
63 /* DNM1 (output) REAL */
66 /* DNM2 (output) REAL */
69 /* IEEE (input) LOGICAL */
70 /* Flag for IEEE or non IEEE arithmetic. */
72 /* ===================================================================== */
76 /* .. Local Scalars .. */
78 /* .. Intrinsic Functions .. */
80 /* .. Executable Statements .. */
82 /* Parameter adjustments */
86 if (*n0 - *i0 - 1 <= 0) {
90 j4 = (*i0 << 2) + *pp - 3;
98 /* Code for IEEE arithmetic. */
102 for (j4 = *i0 << 2; j4 <= i__1; j4 += 4) {
103 z__[j4 - 2] = d__ + z__[j4 - 1];
104 temp = z__[j4 + 1] / z__[j4 - 2];
105 d__ = d__ * temp - *tau;
106 *dmin__ = dmin(*dmin__,d__);
107 z__[j4] = z__[j4 - 1] * temp;
110 emin = dmin(r__1,emin);
115 for (j4 = *i0 << 2; j4 <= i__1; j4 += 4) {
116 z__[j4 - 3] = d__ + z__[j4];
117 temp = z__[j4 + 2] / z__[j4 - 3];
118 d__ = d__ * temp - *tau;
119 *dmin__ = dmin(*dmin__,d__);
120 z__[j4 - 1] = z__[j4] * temp;
123 emin = dmin(r__1,emin);
128 /* Unroll last two steps. */
132 j4 = (*n0 - 2 << 2) - *pp;
133 j4p2 = j4 + (*pp << 1) - 1;
134 z__[j4 - 2] = *dnm2 + z__[j4p2];
135 z__[j4] = z__[j4p2 + 2] * (z__[j4p2] / z__[j4 - 2]);
136 *dnm1 = z__[j4p2 + 2] * (*dnm2 / z__[j4 - 2]) - *tau;
137 *dmin__ = dmin(*dmin__,*dnm1);
141 j4p2 = j4 + (*pp << 1) - 1;
142 z__[j4 - 2] = *dnm1 + z__[j4p2];
143 z__[j4] = z__[j4p2 + 2] * (z__[j4p2] / z__[j4 - 2]);
144 *dn = z__[j4p2 + 2] * (*dnm1 / z__[j4 - 2]) - *tau;
145 *dmin__ = dmin(*dmin__,*dn);
149 /* Code for non IEEE arithmetic. */
153 for (j4 = *i0 << 2; j4 <= i__1; j4 += 4) {
154 z__[j4 - 2] = d__ + z__[j4 - 1];
158 z__[j4] = z__[j4 + 1] * (z__[j4 - 1] / z__[j4 - 2]);
159 d__ = z__[j4 + 1] * (d__ / z__[j4 - 2]) - *tau;
161 *dmin__ = dmin(*dmin__,d__);
163 r__1 = emin, r__2 = z__[j4];
164 emin = dmin(r__1,r__2);
169 for (j4 = *i0 << 2; j4 <= i__1; j4 += 4) {
170 z__[j4 - 3] = d__ + z__[j4];
174 z__[j4 - 1] = z__[j4 + 2] * (z__[j4] / z__[j4 - 3]);
175 d__ = z__[j4 + 2] * (d__ / z__[j4 - 3]) - *tau;
177 *dmin__ = dmin(*dmin__,d__);
179 r__1 = emin, r__2 = z__[j4 - 1];
180 emin = dmin(r__1,r__2);
185 /* Unroll last two steps. */
189 j4 = (*n0 - 2 << 2) - *pp;
190 j4p2 = j4 + (*pp << 1) - 1;
191 z__[j4 - 2] = *dnm2 + z__[j4p2];
195 z__[j4] = z__[j4p2 + 2] * (z__[j4p2] / z__[j4 - 2]);
196 *dnm1 = z__[j4p2 + 2] * (*dnm2 / z__[j4 - 2]) - *tau;
198 *dmin__ = dmin(*dmin__,*dnm1);
202 j4p2 = j4 + (*pp << 1) - 1;
203 z__[j4 - 2] = *dnm1 + z__[j4p2];
207 z__[j4] = z__[j4p2 + 2] * (z__[j4p2] / z__[j4 - 2]);
208 *dn = z__[j4p2 + 2] * (*dnm1 / z__[j4 - 2]) - *tau;
210 *dmin__ = dmin(*dmin__,*dn);
215 z__[(*n0 << 2) - *pp] = emin;