1 SUBROUTINE DLAED6( KNITER, ORGATI, RHO, D, Z, FINIT, TAU, INFO )
2 *
3 * -- LAPACK routine (version 3.2) --
4 * -- LAPACK is a software package provided by Univ. of Tennessee, --
5 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
6 * February 2007
7 *
8 * .. Scalar Arguments ..
9 LOGICAL ORGATI
10 INTEGER INFO, KNITER
11 DOUBLE PRECISION FINIT, RHO, TAU
12 * ..
13 * .. Array Arguments ..
14 DOUBLE PRECISION D( 3 ), Z( 3 )
15 * ..
16 *
17 * Purpose
18 * =======
19 *
20 * DLAED6 computes the positive or negative root (closest to the origin)
21 * of
22 * z(1) z(2) z(3)
23 * f(x) = rho + --------- + ---------- + ---------
24 * d(1)-x d(2)-x d(3)-x
25 *
26 * It is assumed that
27 *
28 * if ORGATI = .true. the root is between d(2) and d(3);
29 * otherwise it is between d(1) and d(2)
30 *
31 * This routine will be called by DLAED4 when necessary. In most cases,
32 * the root sought is the smallest in magnitude, though it might not be
33 * in some extremely rare situations.
34 *
35 * Arguments
36 * =========
37 *
38 * KNITER (input) INTEGER
39 * Refer to DLAED4 for its significance.
40 *
41 * ORGATI (input) LOGICAL
42 * If ORGATI is true, the needed root is between d(2) and
43 * d(3); otherwise it is between d(1) and d(2). See
44 * DLAED4 for further details.
45 *
46 * RHO (input) DOUBLE PRECISION
47 * Refer to the equation f(x) above.
48 *
49 * D (input) DOUBLE PRECISION array, dimension (3)
50 * D satisfies d(1) < d(2) < d(3).
51 *
52 * Z (input) DOUBLE PRECISION array, dimension (3)
53 * Each of the elements in z must be positive.
54 *
55 * FINIT (input) DOUBLE PRECISION
56 * The value of f at 0. It is more accurate than the one
57 * evaluated inside this routine (if someone wants to do
58 * so).
59 *
60 * TAU (output) DOUBLE PRECISION
61 * The root of the equation f(x).
62 *
63 * INFO (output) INTEGER
64 * = 0: successful exit
65 * > 0: if INFO = 1, failure to converge
66 *
67 * Further Details
68 * ===============
69 *
70 * 30/06/99: Based on contributions by
71 * Ren-Cang Li, Computer Science Division, University of California
72 * at Berkeley, USA
73 *
74 * 10/02/03: This version has a few statements commented out for thread
75 * safety (machine parameters are computed on each entry). SJH.
76 *
77 * 05/10/06: Modified from a new version of Ren-Cang Li, use
78 * Gragg-Thornton-Warner cubic convergent scheme for better stability.
79 *
80 * =====================================================================
81 *
82 * .. Parameters ..
83 INTEGER MAXIT
84 PARAMETER ( MAXIT = 40 )
85 DOUBLE PRECISION ZERO, ONE, TWO, THREE, FOUR, EIGHT
86 PARAMETER ( ZERO = 0.0D0, ONE = 1.0D0, TWO = 2.0D0,
87 $ THREE = 3.0D0, FOUR = 4.0D0, EIGHT = 8.0D0 )
88 * ..
89 * .. External Functions ..
90 DOUBLE PRECISION DLAMCH
91 EXTERNAL DLAMCH
92 * ..
93 * .. Local Arrays ..
94 DOUBLE PRECISION DSCALE( 3 ), ZSCALE( 3 )
95 * ..
96 * .. Local Scalars ..
97 LOGICAL SCALE
98 INTEGER I, ITER, NITER
99 DOUBLE PRECISION A, B, BASE, C, DDF, DF, EPS, ERRETM, ETA, F,
100 $ FC, SCLFAC, SCLINV, SMALL1, SMALL2, SMINV1,
101 $ SMINV2, TEMP, TEMP1, TEMP2, TEMP3, TEMP4,
102 $ LBD, UBD
103 * ..
104 * .. Intrinsic Functions ..
105 INTRINSIC ABS, INT, LOG, MAX, MIN, SQRT
106 * ..
107 * .. Executable Statements ..
108 *
109 INFO = 0
110 *
111 IF( ORGATI ) THEN
112 LBD = D(2)
113 UBD = D(3)
114 ELSE
115 LBD = D(1)
116 UBD = D(2)
117 END IF
118 IF( FINIT .LT. ZERO )THEN
119 LBD = ZERO
120 ELSE
121 UBD = ZERO
122 END IF
123 *
124 NITER = 1
125 TAU = ZERO
126 IF( KNITER.EQ.2 ) THEN
127 IF( ORGATI ) THEN
128 TEMP = ( D( 3 )-D( 2 ) ) / TWO
129 C = RHO + Z( 1 ) / ( ( D( 1 )-D( 2 ) )-TEMP )
130 A = C*( D( 2 )+D( 3 ) ) + Z( 2 ) + Z( 3 )
131 B = C*D( 2 )*D( 3 ) + Z( 2 )*D( 3 ) + Z( 3 )*D( 2 )
132 ELSE
133 TEMP = ( D( 1 )-D( 2 ) ) / TWO
134 C = RHO + Z( 3 ) / ( ( D( 3 )-D( 2 ) )-TEMP )
135 A = C*( D( 1 )+D( 2 ) ) + Z( 1 ) + Z( 2 )
136 B = C*D( 1 )*D( 2 ) + Z( 1 )*D( 2 ) + Z( 2 )*D( 1 )
137 END IF
138 TEMP = MAX( ABS( A ), ABS( B ), ABS( C ) )
139 A = A / TEMP
140 B = B / TEMP
141 C = C / TEMP
142 IF( C.EQ.ZERO ) THEN
143 TAU = B / A
144 ELSE IF( A.LE.ZERO ) THEN
145 TAU = ( A-SQRT( ABS( A*A-FOUR*B*C ) ) ) / ( TWO*C )
146 ELSE
147 TAU = TWO*B / ( A+SQRT( ABS( A*A-FOUR*B*C ) ) )
148 END IF
149 IF( TAU .LT. LBD .OR. TAU .GT. UBD )
150 $ TAU = ( LBD+UBD )/TWO
151 IF( D(1).EQ.TAU .OR. D(2).EQ.TAU .OR. D(3).EQ.TAU ) THEN
152 TAU = ZERO
153 ELSE
154 TEMP = FINIT + TAU*Z(1)/( D(1)*( D( 1 )-TAU ) ) +
155 $ TAU*Z(2)/( D(2)*( D( 2 )-TAU ) ) +
156 $ TAU*Z(3)/( D(3)*( D( 3 )-TAU ) )
157 IF( TEMP .LE. ZERO )THEN
158 LBD = TAU
159 ELSE
160 UBD = TAU
161 END IF
162 IF( ABS( FINIT ).LE.ABS( TEMP ) )
163 $ TAU = ZERO
164 END IF
165 END IF
166 *
167 * get machine parameters for possible scaling to avoid overflow
168 *
169 * modified by Sven: parameters SMALL1, SMINV1, SMALL2,
170 * SMINV2, EPS are not SAVEd anymore between one call to the
171 * others but recomputed at each call
172 *
173 EPS = DLAMCH( 'Epsilon' )
174 BASE = DLAMCH( 'Base' )
175 SMALL1 = BASE**( INT( LOG( DLAMCH( 'SafMin' ) ) / LOG( BASE ) /
176 $ THREE ) )
177 SMINV1 = ONE / SMALL1
178 SMALL2 = SMALL1*SMALL1
179 SMINV2 = SMINV1*SMINV1
180 *
181 * Determine if scaling of inputs necessary to avoid overflow
182 * when computing 1/TEMP**3
183 *
184 IF( ORGATI ) THEN
185 TEMP = MIN( ABS( D( 2 )-TAU ), ABS( D( 3 )-TAU ) )
186 ELSE
187 TEMP = MIN( ABS( D( 1 )-TAU ), ABS( D( 2 )-TAU ) )
188 END IF
189 SCALE = .FALSE.
190 IF( TEMP.LE.SMALL1 ) THEN
191 SCALE = .TRUE.
192 IF( TEMP.LE.SMALL2 ) THEN
193 *
194 * Scale up by power of radix nearest 1/SAFMIN**(2/3)
195 *
196 SCLFAC = SMINV2
197 SCLINV = SMALL2
198 ELSE
199 *
200 * Scale up by power of radix nearest 1/SAFMIN**(1/3)
201 *
202 SCLFAC = SMINV1
203 SCLINV = SMALL1
204 END IF
205 *
206 * Scaling up safe because D, Z, TAU scaled elsewhere to be O(1)
207 *
208 DO 10 I = 1, 3
209 DSCALE( I ) = D( I )*SCLFAC
210 ZSCALE( I ) = Z( I )*SCLFAC
211 10 CONTINUE
212 TAU = TAU*SCLFAC
213 LBD = LBD*SCLFAC
214 UBD = UBD*SCLFAC
215 ELSE
216 *
217 * Copy D and Z to DSCALE and ZSCALE
218 *
219 DO 20 I = 1, 3
220 DSCALE( I ) = D( I )
221 ZSCALE( I ) = Z( I )
222 20 CONTINUE
223 END IF
224 *
225 FC = ZERO
226 DF = ZERO
227 DDF = ZERO
228 DO 30 I = 1, 3
229 TEMP = ONE / ( DSCALE( I )-TAU )
230 TEMP1 = ZSCALE( I )*TEMP
231 TEMP2 = TEMP1*TEMP
232 TEMP3 = TEMP2*TEMP
233 FC = FC + TEMP1 / DSCALE( I )
234 DF = DF + TEMP2
235 DDF = DDF + TEMP3
236 30 CONTINUE
237 F = FINIT + TAU*FC
238 *
239 IF( ABS( F ).LE.ZERO )
240 $ GO TO 60
241 IF( F .LE. ZERO )THEN
242 LBD = TAU
243 ELSE
244 UBD = TAU
245 END IF
246 *
247 * Iteration begins -- Use Gragg-Thornton-Warner cubic convergent
248 * scheme
249 *
250 * It is not hard to see that
251 *
252 * 1) Iterations will go up monotonically
253 * if FINIT < 0;
254 *
255 * 2) Iterations will go down monotonically
256 * if FINIT > 0.
257 *
258 ITER = NITER + 1
259 *
260 DO 50 NITER = ITER, MAXIT
261 *
262 IF( ORGATI ) THEN
263 TEMP1 = DSCALE( 2 ) - TAU
264 TEMP2 = DSCALE( 3 ) - TAU
265 ELSE
266 TEMP1 = DSCALE( 1 ) - TAU
267 TEMP2 = DSCALE( 2 ) - TAU
268 END IF
269 A = ( TEMP1+TEMP2 )*F - TEMP1*TEMP2*DF
270 B = TEMP1*TEMP2*F
271 C = F - ( TEMP1+TEMP2 )*DF + TEMP1*TEMP2*DDF
272 TEMP = MAX( ABS( A ), ABS( B ), ABS( C ) )
273 A = A / TEMP
274 B = B / TEMP
275 C = C / TEMP
276 IF( C.EQ.ZERO ) THEN
277 ETA = B / A
278 ELSE IF( A.LE.ZERO ) THEN
279 ETA = ( A-SQRT( ABS( A*A-FOUR*B*C ) ) ) / ( TWO*C )
280 ELSE
281 ETA = TWO*B / ( A+SQRT( ABS( A*A-FOUR*B*C ) ) )
282 END IF
283 IF( F*ETA.GE.ZERO ) THEN
284 ETA = -F / DF
285 END IF
286 *
287 TAU = TAU + ETA
288 IF( TAU .LT. LBD .OR. TAU .GT. UBD )
289 $ TAU = ( LBD + UBD )/TWO
290 *
291 FC = ZERO
292 ERRETM = ZERO
293 DF = ZERO
294 DDF = ZERO
295 DO 40 I = 1, 3
296 TEMP = ONE / ( DSCALE( I )-TAU )
297 TEMP1 = ZSCALE( I )*TEMP
298 TEMP2 = TEMP1*TEMP
299 TEMP3 = TEMP2*TEMP
300 TEMP4 = TEMP1 / DSCALE( I )
301 FC = FC + TEMP4
302 ERRETM = ERRETM + ABS( TEMP4 )
303 DF = DF + TEMP2
304 DDF = DDF + TEMP3
305 40 CONTINUE
306 F = FINIT + TAU*FC
307 ERRETM = EIGHT*( ABS( FINIT )+ABS( TAU )*ERRETM ) +
308 $ ABS( TAU )*DF
309 IF( ABS( F ).LE.EPS*ERRETM )
310 $ GO TO 60
311 IF( F .LE. ZERO )THEN
312 LBD = TAU
313 ELSE
314 UBD = TAU
315 END IF
316 50 CONTINUE
317 INFO = 1
318 60 CONTINUE
319 *
320 * Undo scaling
321 *
322 IF( SCALE )
323 $ TAU = TAU*SCLINV
324 RETURN
325 *
326 * End of DLAED6
327 *
328 END
2 *
3 * -- LAPACK routine (version 3.2) --
4 * -- LAPACK is a software package provided by Univ. of Tennessee, --
5 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
6 * February 2007
7 *
8 * .. Scalar Arguments ..
9 LOGICAL ORGATI
10 INTEGER INFO, KNITER
11 DOUBLE PRECISION FINIT, RHO, TAU
12 * ..
13 * .. Array Arguments ..
14 DOUBLE PRECISION D( 3 ), Z( 3 )
15 * ..
16 *
17 * Purpose
18 * =======
19 *
20 * DLAED6 computes the positive or negative root (closest to the origin)
21 * of
22 * z(1) z(2) z(3)
23 * f(x) = rho + --------- + ---------- + ---------
24 * d(1)-x d(2)-x d(3)-x
25 *
26 * It is assumed that
27 *
28 * if ORGATI = .true. the root is between d(2) and d(3);
29 * otherwise it is between d(1) and d(2)
30 *
31 * This routine will be called by DLAED4 when necessary. In most cases,
32 * the root sought is the smallest in magnitude, though it might not be
33 * in some extremely rare situations.
34 *
35 * Arguments
36 * =========
37 *
38 * KNITER (input) INTEGER
39 * Refer to DLAED4 for its significance.
40 *
41 * ORGATI (input) LOGICAL
42 * If ORGATI is true, the needed root is between d(2) and
43 * d(3); otherwise it is between d(1) and d(2). See
44 * DLAED4 for further details.
45 *
46 * RHO (input) DOUBLE PRECISION
47 * Refer to the equation f(x) above.
48 *
49 * D (input) DOUBLE PRECISION array, dimension (3)
50 * D satisfies d(1) < d(2) < d(3).
51 *
52 * Z (input) DOUBLE PRECISION array, dimension (3)
53 * Each of the elements in z must be positive.
54 *
55 * FINIT (input) DOUBLE PRECISION
56 * The value of f at 0. It is more accurate than the one
57 * evaluated inside this routine (if someone wants to do
58 * so).
59 *
60 * TAU (output) DOUBLE PRECISION
61 * The root of the equation f(x).
62 *
63 * INFO (output) INTEGER
64 * = 0: successful exit
65 * > 0: if INFO = 1, failure to converge
66 *
67 * Further Details
68 * ===============
69 *
70 * 30/06/99: Based on contributions by
71 * Ren-Cang Li, Computer Science Division, University of California
72 * at Berkeley, USA
73 *
74 * 10/02/03: This version has a few statements commented out for thread
75 * safety (machine parameters are computed on each entry). SJH.
76 *
77 * 05/10/06: Modified from a new version of Ren-Cang Li, use
78 * Gragg-Thornton-Warner cubic convergent scheme for better stability.
79 *
80 * =====================================================================
81 *
82 * .. Parameters ..
83 INTEGER MAXIT
84 PARAMETER ( MAXIT = 40 )
85 DOUBLE PRECISION ZERO, ONE, TWO, THREE, FOUR, EIGHT
86 PARAMETER ( ZERO = 0.0D0, ONE = 1.0D0, TWO = 2.0D0,
87 $ THREE = 3.0D0, FOUR = 4.0D0, EIGHT = 8.0D0 )
88 * ..
89 * .. External Functions ..
90 DOUBLE PRECISION DLAMCH
91 EXTERNAL DLAMCH
92 * ..
93 * .. Local Arrays ..
94 DOUBLE PRECISION DSCALE( 3 ), ZSCALE( 3 )
95 * ..
96 * .. Local Scalars ..
97 LOGICAL SCALE
98 INTEGER I, ITER, NITER
99 DOUBLE PRECISION A, B, BASE, C, DDF, DF, EPS, ERRETM, ETA, F,
100 $ FC, SCLFAC, SCLINV, SMALL1, SMALL2, SMINV1,
101 $ SMINV2, TEMP, TEMP1, TEMP2, TEMP3, TEMP4,
102 $ LBD, UBD
103 * ..
104 * .. Intrinsic Functions ..
105 INTRINSIC ABS, INT, LOG, MAX, MIN, SQRT
106 * ..
107 * .. Executable Statements ..
108 *
109 INFO = 0
110 *
111 IF( ORGATI ) THEN
112 LBD = D(2)
113 UBD = D(3)
114 ELSE
115 LBD = D(1)
116 UBD = D(2)
117 END IF
118 IF( FINIT .LT. ZERO )THEN
119 LBD = ZERO
120 ELSE
121 UBD = ZERO
122 END IF
123 *
124 NITER = 1
125 TAU = ZERO
126 IF( KNITER.EQ.2 ) THEN
127 IF( ORGATI ) THEN
128 TEMP = ( D( 3 )-D( 2 ) ) / TWO
129 C = RHO + Z( 1 ) / ( ( D( 1 )-D( 2 ) )-TEMP )
130 A = C*( D( 2 )+D( 3 ) ) + Z( 2 ) + Z( 3 )
131 B = C*D( 2 )*D( 3 ) + Z( 2 )*D( 3 ) + Z( 3 )*D( 2 )
132 ELSE
133 TEMP = ( D( 1 )-D( 2 ) ) / TWO
134 C = RHO + Z( 3 ) / ( ( D( 3 )-D( 2 ) )-TEMP )
135 A = C*( D( 1 )+D( 2 ) ) + Z( 1 ) + Z( 2 )
136 B = C*D( 1 )*D( 2 ) + Z( 1 )*D( 2 ) + Z( 2 )*D( 1 )
137 END IF
138 TEMP = MAX( ABS( A ), ABS( B ), ABS( C ) )
139 A = A / TEMP
140 B = B / TEMP
141 C = C / TEMP
142 IF( C.EQ.ZERO ) THEN
143 TAU = B / A
144 ELSE IF( A.LE.ZERO ) THEN
145 TAU = ( A-SQRT( ABS( A*A-FOUR*B*C ) ) ) / ( TWO*C )
146 ELSE
147 TAU = TWO*B / ( A+SQRT( ABS( A*A-FOUR*B*C ) ) )
148 END IF
149 IF( TAU .LT. LBD .OR. TAU .GT. UBD )
150 $ TAU = ( LBD+UBD )/TWO
151 IF( D(1).EQ.TAU .OR. D(2).EQ.TAU .OR. D(3).EQ.TAU ) THEN
152 TAU = ZERO
153 ELSE
154 TEMP = FINIT + TAU*Z(1)/( D(1)*( D( 1 )-TAU ) ) +
155 $ TAU*Z(2)/( D(2)*( D( 2 )-TAU ) ) +
156 $ TAU*Z(3)/( D(3)*( D( 3 )-TAU ) )
157 IF( TEMP .LE. ZERO )THEN
158 LBD = TAU
159 ELSE
160 UBD = TAU
161 END IF
162 IF( ABS( FINIT ).LE.ABS( TEMP ) )
163 $ TAU = ZERO
164 END IF
165 END IF
166 *
167 * get machine parameters for possible scaling to avoid overflow
168 *
169 * modified by Sven: parameters SMALL1, SMINV1, SMALL2,
170 * SMINV2, EPS are not SAVEd anymore between one call to the
171 * others but recomputed at each call
172 *
173 EPS = DLAMCH( 'Epsilon' )
174 BASE = DLAMCH( 'Base' )
175 SMALL1 = BASE**( INT( LOG( DLAMCH( 'SafMin' ) ) / LOG( BASE ) /
176 $ THREE ) )
177 SMINV1 = ONE / SMALL1
178 SMALL2 = SMALL1*SMALL1
179 SMINV2 = SMINV1*SMINV1
180 *
181 * Determine if scaling of inputs necessary to avoid overflow
182 * when computing 1/TEMP**3
183 *
184 IF( ORGATI ) THEN
185 TEMP = MIN( ABS( D( 2 )-TAU ), ABS( D( 3 )-TAU ) )
186 ELSE
187 TEMP = MIN( ABS( D( 1 )-TAU ), ABS( D( 2 )-TAU ) )
188 END IF
189 SCALE = .FALSE.
190 IF( TEMP.LE.SMALL1 ) THEN
191 SCALE = .TRUE.
192 IF( TEMP.LE.SMALL2 ) THEN
193 *
194 * Scale up by power of radix nearest 1/SAFMIN**(2/3)
195 *
196 SCLFAC = SMINV2
197 SCLINV = SMALL2
198 ELSE
199 *
200 * Scale up by power of radix nearest 1/SAFMIN**(1/3)
201 *
202 SCLFAC = SMINV1
203 SCLINV = SMALL1
204 END IF
205 *
206 * Scaling up safe because D, Z, TAU scaled elsewhere to be O(1)
207 *
208 DO 10 I = 1, 3
209 DSCALE( I ) = D( I )*SCLFAC
210 ZSCALE( I ) = Z( I )*SCLFAC
211 10 CONTINUE
212 TAU = TAU*SCLFAC
213 LBD = LBD*SCLFAC
214 UBD = UBD*SCLFAC
215 ELSE
216 *
217 * Copy D and Z to DSCALE and ZSCALE
218 *
219 DO 20 I = 1, 3
220 DSCALE( I ) = D( I )
221 ZSCALE( I ) = Z( I )
222 20 CONTINUE
223 END IF
224 *
225 FC = ZERO
226 DF = ZERO
227 DDF = ZERO
228 DO 30 I = 1, 3
229 TEMP = ONE / ( DSCALE( I )-TAU )
230 TEMP1 = ZSCALE( I )*TEMP
231 TEMP2 = TEMP1*TEMP
232 TEMP3 = TEMP2*TEMP
233 FC = FC + TEMP1 / DSCALE( I )
234 DF = DF + TEMP2
235 DDF = DDF + TEMP3
236 30 CONTINUE
237 F = FINIT + TAU*FC
238 *
239 IF( ABS( F ).LE.ZERO )
240 $ GO TO 60
241 IF( F .LE. ZERO )THEN
242 LBD = TAU
243 ELSE
244 UBD = TAU
245 END IF
246 *
247 * Iteration begins -- Use Gragg-Thornton-Warner cubic convergent
248 * scheme
249 *
250 * It is not hard to see that
251 *
252 * 1) Iterations will go up monotonically
253 * if FINIT < 0;
254 *
255 * 2) Iterations will go down monotonically
256 * if FINIT > 0.
257 *
258 ITER = NITER + 1
259 *
260 DO 50 NITER = ITER, MAXIT
261 *
262 IF( ORGATI ) THEN
263 TEMP1 = DSCALE( 2 ) - TAU
264 TEMP2 = DSCALE( 3 ) - TAU
265 ELSE
266 TEMP1 = DSCALE( 1 ) - TAU
267 TEMP2 = DSCALE( 2 ) - TAU
268 END IF
269 A = ( TEMP1+TEMP2 )*F - TEMP1*TEMP2*DF
270 B = TEMP1*TEMP2*F
271 C = F - ( TEMP1+TEMP2 )*DF + TEMP1*TEMP2*DDF
272 TEMP = MAX( ABS( A ), ABS( B ), ABS( C ) )
273 A = A / TEMP
274 B = B / TEMP
275 C = C / TEMP
276 IF( C.EQ.ZERO ) THEN
277 ETA = B / A
278 ELSE IF( A.LE.ZERO ) THEN
279 ETA = ( A-SQRT( ABS( A*A-FOUR*B*C ) ) ) / ( TWO*C )
280 ELSE
281 ETA = TWO*B / ( A+SQRT( ABS( A*A-FOUR*B*C ) ) )
282 END IF
283 IF( F*ETA.GE.ZERO ) THEN
284 ETA = -F / DF
285 END IF
286 *
287 TAU = TAU + ETA
288 IF( TAU .LT. LBD .OR. TAU .GT. UBD )
289 $ TAU = ( LBD + UBD )/TWO
290 *
291 FC = ZERO
292 ERRETM = ZERO
293 DF = ZERO
294 DDF = ZERO
295 DO 40 I = 1, 3
296 TEMP = ONE / ( DSCALE( I )-TAU )
297 TEMP1 = ZSCALE( I )*TEMP
298 TEMP2 = TEMP1*TEMP
299 TEMP3 = TEMP2*TEMP
300 TEMP4 = TEMP1 / DSCALE( I )
301 FC = FC + TEMP4
302 ERRETM = ERRETM + ABS( TEMP4 )
303 DF = DF + TEMP2
304 DDF = DDF + TEMP3
305 40 CONTINUE
306 F = FINIT + TAU*FC
307 ERRETM = EIGHT*( ABS( FINIT )+ABS( TAU )*ERRETM ) +
308 $ ABS( TAU )*DF
309 IF( ABS( F ).LE.EPS*ERRETM )
310 $ GO TO 60
311 IF( F .LE. ZERO )THEN
312 LBD = TAU
313 ELSE
314 UBD = TAU
315 END IF
316 50 CONTINUE
317 INFO = 1
318 60 CONTINUE
319 *
320 * Undo scaling
321 *
322 IF( SCALE )
323 $ TAU = TAU*SCLINV
324 RETURN
325 *
326 * End of DLAED6
327 *
328 END