LAPACK 3.12.1
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches
slasq4.f
Go to the documentation of this file.
1*> \brief \b SLASQ4 computes an approximation to the smallest eigenvalue using values of d from the previous transform. Used by sbdsqr.
2*
3* =========== DOCUMENTATION ===========
4*
5* Online html documentation available at
6* http://www.netlib.org/lapack/explore-html/
7*
8*> Download SLASQ4 + dependencies
9*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/slasq4.f">
10*> [TGZ]</a>
11*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/slasq4.f">
12*> [ZIP]</a>
13*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/slasq4.f">
14*> [TXT]</a>
15*
16* Definition:
17* ===========
18*
19* SUBROUTINE SLASQ4( I0, N0, Z, PP, N0IN, DMIN, DMIN1, DMIN2, DN,
20* DN1, DN2, TAU, TTYPE, G )
21*
22* .. Scalar Arguments ..
23* INTEGER I0, N0, N0IN, PP, TTYPE
24* REAL DMIN, DMIN1, DMIN2, DN, DN1, DN2, G, TAU
25* ..
26* .. Array Arguments ..
27* REAL Z( * )
28* ..
29*
30*
31*> \par Purpose:
32* =============
33*>
34*> \verbatim
35*>
36*> SLASQ4 computes an approximation TAU to the smallest eigenvalue
37*> using values of d from the previous transform.
38*> \endverbatim
39*
40* Arguments:
41* ==========
42*
43*> \param[in] I0
44*> \verbatim
45*> I0 is INTEGER
46*> First index.
47*> \endverbatim
48*>
49*> \param[in] N0
50*> \verbatim
51*> N0 is INTEGER
52*> Last index.
53*> \endverbatim
54*>
55*> \param[in] Z
56*> \verbatim
57*> Z is REAL array, dimension ( 4*N0 )
58*> Z holds the qd array.
59*> \endverbatim
60*>
61*> \param[in] PP
62*> \verbatim
63*> PP is INTEGER
64*> PP=0 for ping, PP=1 for pong.
65*> \endverbatim
66*>
67*> \param[in] N0IN
68*> \verbatim
69*> N0IN is INTEGER
70*> The value of N0 at start of EIGTEST.
71*> \endverbatim
72*>
73*> \param[in] DMIN
74*> \verbatim
75*> DMIN is REAL
76*> Minimum value of d.
77*> \endverbatim
78*>
79*> \param[in] DMIN1
80*> \verbatim
81*> DMIN1 is REAL
82*> Minimum value of d, excluding D( N0 ).
83*> \endverbatim
84*>
85*> \param[in] DMIN2
86*> \verbatim
87*> DMIN2 is REAL
88*> Minimum value of d, excluding D( N0 ) and D( N0-1 ).
89*> \endverbatim
90*>
91*> \param[in] DN
92*> \verbatim
93*> DN is REAL
94*> d(N)
95*> \endverbatim
96*>
97*> \param[in] DN1
98*> \verbatim
99*> DN1 is REAL
100*> d(N-1)
101*> \endverbatim
102*>
103*> \param[in] DN2
104*> \verbatim
105*> DN2 is REAL
106*> d(N-2)
107*> \endverbatim
108*>
109*> \param[out] TAU
110*> \verbatim
111*> TAU is REAL
112*> This is the shift.
113*> \endverbatim
114*>
115*> \param[out] TTYPE
116*> \verbatim
117*> TTYPE is INTEGER
118*> Shift type.
119*> \endverbatim
120*>
121*> \param[in,out] G
122*> \verbatim
123*> G is REAL
124*> G is passed as an argument in order to save its value between
125*> calls to SLASQ4.
126*> \endverbatim
127*
128* Authors:
129* ========
130*
131*> \author Univ. of Tennessee
132*> \author Univ. of California Berkeley
133*> \author Univ. of Colorado Denver
134*> \author NAG Ltd.
135*
136*> \ingroup lasq4
137*
138*> \par Further Details:
139* =====================
140*>
141*> \verbatim
142*>
143*> CNST1 = 9/16
144*> \endverbatim
145*>
146* =====================================================================
147 SUBROUTINE slasq4( I0, N0, Z, PP, N0IN, DMIN, DMIN1, DMIN2, DN,
148 $ DN1, DN2, TAU, TTYPE, G )
149*
150* -- LAPACK computational routine --
151* -- LAPACK is a software package provided by Univ. of Tennessee, --
152* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
153*
154* .. Scalar Arguments ..
155 INTEGER I0, N0, N0IN, PP, TTYPE
156 REAL DMIN, DMIN1, DMIN2, DN, DN1, DN2, G, TAU
157* ..
158* .. Array Arguments ..
159 REAL Z( * )
160* ..
161*
162* =====================================================================
163*
164* .. Parameters ..
165 REAL CNST1, CNST2, CNST3
166 parameter( cnst1 = 0.5630e0, cnst2 = 1.010e0,
167 $ cnst3 = 1.050e0 )
168 REAL QURTR, THIRD, HALF, ZERO, ONE, TWO, HUNDRD
169 parameter( qurtr = 0.250e0, third = 0.3330e0,
170 $ half = 0.50e0, zero = 0.0e0, one = 1.0e0,
171 $ two = 2.0e0, hundrd = 100.0e0 )
172* ..
173* .. Local Scalars ..
174 INTEGER I4, NN, NP
175 REAL A2, B1, B2, GAM, GAP1, GAP2, S
176* ..
177* .. Intrinsic Functions ..
178 INTRINSIC max, min, sqrt
179* ..
180* .. Executable Statements ..
181*
182* A negative DMIN forces the shift to take that absolute value
183* TTYPE records the type of shift.
184*
185 IF( dmin.LE.zero ) THEN
186 tau = -dmin
187 ttype = -1
188 RETURN
189 END IF
190*
191 nn = 4*n0 + pp
192 IF( n0in.EQ.n0 ) THEN
193*
194* No eigenvalues deflated.
195*
196 IF( dmin.EQ.dn .OR. dmin.EQ.dn1 ) THEN
197*
198 b1 = sqrt( z( nn-3 ) )*sqrt( z( nn-5 ) )
199 b2 = sqrt( z( nn-7 ) )*sqrt( z( nn-9 ) )
200 a2 = z( nn-7 ) + z( nn-5 )
201*
202* Cases 2 and 3.
203*
204 IF( dmin.EQ.dn .AND. dmin1.EQ.dn1 ) THEN
205 gap2 = dmin2 - a2 - dmin2*qurtr
206 IF( gap2.GT.zero .AND. gap2.GT.b2 ) THEN
207 gap1 = a2 - dn - ( b2 / gap2 )*b2
208 ELSE
209 gap1 = a2 - dn - ( b1+b2 )
210 END IF
211 IF( gap1.GT.zero .AND. gap1.GT.b1 ) THEN
212 s = max( dn-( b1 / gap1 )*b1, half*dmin )
213 ttype = -2
214 ELSE
215 s = zero
216 IF( dn.GT.b1 )
217 $ s = dn - b1
218 IF( a2.GT.( b1+b2 ) )
219 $ s = min( s, a2-( b1+b2 ) )
220 s = max( s, third*dmin )
221 ttype = -3
222 END IF
223 ELSE
224*
225* Case 4.
226*
227 ttype = -4
228 s = qurtr*dmin
229 IF( dmin.EQ.dn ) THEN
230 gam = dn
231 a2 = zero
232 IF( z( nn-5 ) .GT. z( nn-7 ) )
233 $ RETURN
234 b2 = z( nn-5 ) / z( nn-7 )
235 np = nn - 9
236 ELSE
237 np = nn - 2*pp
238 gam = dn1
239 IF( z( np-4 ) .GT. z( np-2 ) )
240 $ RETURN
241 a2 = z( np-4 ) / z( np-2 )
242 IF( z( nn-9 ) .GT. z( nn-11 ) )
243 $ RETURN
244 b2 = z( nn-9 ) / z( nn-11 )
245 np = nn - 13
246 END IF
247*
248* Approximate contribution to norm squared from I < NN-1.
249*
250 a2 = a2 + b2
251 DO 10 i4 = np, 4*i0 - 1 + pp, -4
252 IF( b2.EQ.zero )
253 $ GO TO 20
254 b1 = b2
255 IF( z( i4 ) .GT. z( i4-2 ) )
256 $ RETURN
257 b2 = b2*( z( i4 ) / z( i4-2 ) )
258 a2 = a2 + b2
259 IF( hundrd*max( b2, b1 ).LT.a2 .OR. cnst1.LT.a2 )
260 $ GO TO 20
261 10 CONTINUE
262 20 CONTINUE
263 a2 = cnst3*a2
264*
265* Rayleigh quotient residual bound.
266*
267 IF( a2.LT.cnst1 )
268 $ s = gam*( one-sqrt( a2 ) ) / ( one+a2 )
269 END IF
270 ELSE IF( dmin.EQ.dn2 ) THEN
271*
272* Case 5.
273*
274 ttype = -5
275 s = qurtr*dmin
276*
277* Compute contribution to norm squared from I > NN-2.
278*
279 np = nn - 2*pp
280 b1 = z( np-2 )
281 b2 = z( np-6 )
282 gam = dn2
283 IF( z( np-8 ).GT.b2 .OR. z( np-4 ).GT.b1 )
284 $ RETURN
285 a2 = ( z( np-8 ) / b2 )*( one+z( np-4 ) / b1 )
286*
287* Approximate contribution to norm squared from I < NN-2.
288*
289 IF( n0-i0.GT.2 ) THEN
290 b2 = z( nn-13 ) / z( nn-15 )
291 a2 = a2 + b2
292 DO 30 i4 = nn - 17, 4*i0 - 1 + pp, -4
293 IF( b2.EQ.zero )
294 $ GO TO 40
295 b1 = b2
296 IF( z( i4 ) .GT. z( i4-2 ) )
297 $ RETURN
298 b2 = b2*( z( i4 ) / z( i4-2 ) )
299 a2 = a2 + b2
300 IF( hundrd*max( b2, b1 ).LT.a2 .OR. cnst1.LT.a2 )
301 $ GO TO 40
302 30 CONTINUE
303 40 CONTINUE
304 a2 = cnst3*a2
305 END IF
306*
307 IF( a2.LT.cnst1 )
308 $ s = gam*( one-sqrt( a2 ) ) / ( one+a2 )
309 ELSE
310*
311* Case 6, no information to guide us.
312*
313 IF( ttype.EQ.-6 ) THEN
314 g = g + third*( one-g )
315 ELSE IF( ttype.EQ.-18 ) THEN
316 g = qurtr*third
317 ELSE
318 g = qurtr
319 END IF
320 s = g*dmin
321 ttype = -6
322 END IF
323*
324 ELSE IF( n0in.EQ.( n0+1 ) ) THEN
325*
326* One eigenvalue just deflated. Use DMIN1, DN1 for DMIN and DN.
327*
328 IF( dmin1.EQ.dn1 .AND. dmin2.EQ.dn2 ) THEN
329*
330* Cases 7 and 8.
331*
332 ttype = -7
333 s = third*dmin1
334 IF( z( nn-5 ).GT.z( nn-7 ) )
335 $ RETURN
336 b1 = z( nn-5 ) / z( nn-7 )
337 b2 = b1
338 IF( b2.EQ.zero )
339 $ GO TO 60
340 DO 50 i4 = 4*n0 - 9 + pp, 4*i0 - 1 + pp, -4
341 a2 = b1
342 IF( z( i4 ).GT.z( i4-2 ) )
343 $ RETURN
344 b1 = b1*( z( i4 ) / z( i4-2 ) )
345 b2 = b2 + b1
346 IF( hundrd*max( b1, a2 ).LT.b2 )
347 $ GO TO 60
348 50 CONTINUE
349 60 CONTINUE
350 b2 = sqrt( cnst3*b2 )
351 a2 = dmin1 / ( one+b2**2 )
352 gap2 = half*dmin2 - a2
353 IF( gap2.GT.zero .AND. gap2.GT.b2*a2 ) THEN
354 s = max( s, a2*( one-cnst2*a2*( b2 / gap2 )*b2 ) )
355 ELSE
356 s = max( s, a2*( one-cnst2*b2 ) )
357 ttype = -8
358 END IF
359 ELSE
360*
361* Case 9.
362*
363 s = qurtr*dmin1
364 IF( dmin1.EQ.dn1 )
365 $ s = half*dmin1
366 ttype = -9
367 END IF
368*
369 ELSE IF( n0in.EQ.( n0+2 ) ) THEN
370*
371* Two eigenvalues deflated. Use DMIN2, DN2 for DMIN and DN.
372*
373* Cases 10 and 11.
374*
375 IF( dmin2.EQ.dn2 .AND. two*z( nn-5 ).LT.z( nn-7 ) ) THEN
376 ttype = -10
377 s = third*dmin2
378 IF( z( nn-5 ).GT.z( nn-7 ) )
379 $ RETURN
380 b1 = z( nn-5 ) / z( nn-7 )
381 b2 = b1
382 IF( b2.EQ.zero )
383 $ GO TO 80
384 DO 70 i4 = 4*n0 - 9 + pp, 4*i0 - 1 + pp, -4
385 IF( z( i4 ).GT.z( i4-2 ) )
386 $ RETURN
387 b1 = b1*( z( i4 ) / z( i4-2 ) )
388 b2 = b2 + b1
389 IF( hundrd*b1.LT.b2 )
390 $ GO TO 80
391 70 CONTINUE
392 80 CONTINUE
393 b2 = sqrt( cnst3*b2 )
394 a2 = dmin2 / ( one+b2**2 )
395 gap2 = z( nn-7 ) + z( nn-9 ) -
396 $ sqrt( z( nn-11 ) )*sqrt( z( nn-9 ) ) - a2
397 IF( gap2.GT.zero .AND. gap2.GT.b2*a2 ) THEN
398 s = max( s, a2*( one-cnst2*a2*( b2 / gap2 )*b2 ) )
399 ELSE
400 s = max( s, a2*( one-cnst2*b2 ) )
401 END IF
402 ELSE
403 s = qurtr*dmin2
404 ttype = -11
405 END IF
406 ELSE IF( n0in.GT.( n0+2 ) ) THEN
407*
408* Case 12, more than two eigenvalues deflated. No information.
409*
410 s = zero
411 ttype = -12
412 END IF
413*
414 tau = s
415 RETURN
416*
417* End of SLASQ4
418*
419 END
subroutine slasq4(i0, n0, z, pp, n0in, dmin, dmin1, dmin2, dn, dn1, dn2, tau, ttype, g)
SLASQ4 computes an approximation to the smallest eigenvalue using values of d from the previous trans...
Definition slasq4.f:149