LAPACK 3.12.1
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches

◆ slasq4()

subroutine slasq4 ( integer i0,
integer n0,
real, dimension( * ) z,
integer pp,
integer n0in,
real dmin,
real dmin1,
real dmin2,
real dn,
real dn1,
real dn2,
real tau,
integer ttype,
real g )

SLASQ4 computes an approximation to the smallest eigenvalue using values of d from the previous transform. Used by sbdsqr.

Download SLASQ4 + dependencies [TGZ] [ZIP] [TXT]

Purpose:
!>
!> SLASQ4 computes an approximation TAU to the smallest eigenvalue
!> using values of d from the previous transform.
!>
Parameters
[in]I0
!>          I0 is INTEGER
!>        First index.
!>
[in]N0
!>          N0 is INTEGER
!>        Last index.
!>
[in]Z
!>          Z is REAL array, dimension ( 4*N0 )
!>        Z holds the qd array.
!>
[in]PP
!>          PP is INTEGER
!>        PP=0 for ping, PP=1 for pong.
!>
[in]N0IN
!>          N0IN is INTEGER
!>        The value of N0 at start of EIGTEST.
!>
[in]DMIN
!>          DMIN is REAL
!>        Minimum value of d.
!>
[in]DMIN1
!>          DMIN1 is REAL
!>        Minimum value of d, excluding D( N0 ).
!>
[in]DMIN2
!>          DMIN2 is REAL
!>        Minimum value of d, excluding D( N0 ) and D( N0-1 ).
!>
[in]DN
!>          DN is REAL
!>        d(N)
!>
[in]DN1
!>          DN1 is REAL
!>        d(N-1)
!>
[in]DN2
!>          DN2 is REAL
!>        d(N-2)
!>
[out]TAU
!>          TAU is REAL
!>        This is the shift.
!>
[out]TTYPE
!>          TTYPE is INTEGER
!>        Shift type.
!>
[in,out]G
!>          G is REAL
!>        G is passed as an argument in order to save its value between
!>        calls to SLASQ4.
!>
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Further Details:
!>
!>  CNST1 = 9/16
!>

Definition at line 147 of file slasq4.f.

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*
Here is the caller graph for this function: