LAPACK  3.6.1 LAPACK: Linear Algebra PACKage
 subroutine dlasq4 ( integer I0, integer N0, double precision, dimension( * ) Z, integer PP, integer N0IN, double precision DMIN, double precision DMIN1, double precision DMIN2, double precision DN, double precision DN1, double precision DN2, double precision TAU, integer TTYPE, double precision G )

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

Purpose:
``` DLASQ4 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 DOUBLE PRECISION 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 DOUBLE PRECISION Minimum value of d.``` [in] DMIN1 ``` DMIN1 is DOUBLE PRECISION Minimum value of d, excluding D( N0 ).``` [in] DMIN2 ``` DMIN2 is DOUBLE PRECISION Minimum value of d, excluding D( N0 ) and D( N0-1 ).``` [in] DN ``` DN is DOUBLE PRECISION d(N)``` [in] DN1 ``` DN1 is DOUBLE PRECISION d(N-1)``` [in] DN2 ``` DN2 is DOUBLE PRECISION d(N-2)``` [out] TAU ``` TAU is DOUBLE PRECISION This is the shift.``` [out] TTYPE ``` TTYPE is INTEGER Shift type.``` [in,out] G ``` G is DOUBLE PRECISION G is passed as an argument in order to save its value between calls to DLASQ4.```
Date
June 2016
Further Details:
`  CNST1 = 9/16`

Definition at line 153 of file dlasq4.f.

153 *
154 * -- LAPACK computational routine (version 3.6.1) --
155 * -- LAPACK is a software package provided by Univ. of Tennessee, --
156 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
157 * June 2016
158 *
159 * .. Scalar Arguments ..
160  INTEGER i0, n0, n0in, pp, ttype
161  DOUBLE PRECISION dmin, dmin1, dmin2, dn, dn1, dn2, g, tau
162 * ..
163 * .. Array Arguments ..
164  DOUBLE PRECISION z( * )
165 * ..
166 *
167 * =====================================================================
168 *
169 * .. Parameters ..
170  DOUBLE PRECISION cnst1, cnst2, cnst3
171  parameter ( cnst1 = 0.5630d0, cnst2 = 1.010d0,
172  \$ cnst3 = 1.050d0 )
173  DOUBLE PRECISION qurtr, third, half, zero, one, two, hundrd
174  parameter ( qurtr = 0.250d0, third = 0.3330d0,
175  \$ half = 0.50d0, zero = 0.0d0, one = 1.0d0,
176  \$ two = 2.0d0, hundrd = 100.0d0 )
177 * ..
178 * .. Local Scalars ..
179  INTEGER i4, nn, np
180  DOUBLE PRECISION a2, b1, b2, gam, gap1, gap2, s
181 * ..
182 * .. Intrinsic Functions ..
183  INTRINSIC max, min, sqrt
184 * ..
185 * .. Executable Statements ..
186 *
187 * A negative DMIN forces the shift to take that absolute value
188 * TTYPE records the type of shift.
189 *
190  IF( dmin.LE.zero ) THEN
191  tau = -dmin
192  ttype = -1
193  RETURN
194  END IF
195 *
196  nn = 4*n0 + pp
197  IF( n0in.EQ.n0 ) THEN
198 *
199 * No eigenvalues deflated.
200 *
201  IF( dmin.EQ.dn .OR. dmin.EQ.dn1 ) THEN
202 *
203  b1 = sqrt( z( nn-3 ) )*sqrt( z( nn-5 ) )
204  b2 = sqrt( z( nn-7 ) )*sqrt( z( nn-9 ) )
205  a2 = z( nn-7 ) + z( nn-5 )
206 *
207 * Cases 2 and 3.
208 *
209  IF( dmin.EQ.dn .AND. dmin1.EQ.dn1 ) THEN
210  gap2 = dmin2 - a2 - dmin2*qurtr
211  IF( gap2.GT.zero .AND. gap2.GT.b2 ) THEN
212  gap1 = a2 - dn - ( b2 / gap2 )*b2
213  ELSE
214  gap1 = a2 - dn - ( b1+b2 )
215  END IF
216  IF( gap1.GT.zero .AND. gap1.GT.b1 ) THEN
217  s = max( dn-( b1 / gap1 )*b1, half*dmin )
218  ttype = -2
219  ELSE
220  s = zero
221  IF( dn.GT.b1 )
222  \$ s = dn - b1
223  IF( a2.GT.( b1+b2 ) )
224  \$ s = min( s, a2-( b1+b2 ) )
225  s = max( s, third*dmin )
226  ttype = -3
227  END IF
228  ELSE
229 *
230 * Case 4.
231 *
232  ttype = -4
233  s = qurtr*dmin
234  IF( dmin.EQ.dn ) THEN
235  gam = dn
236  a2 = zero
237  IF( z( nn-5 ) .GT. z( nn-7 ) )
238  \$ RETURN
239  b2 = z( nn-5 ) / z( nn-7 )
240  np = nn - 9
241  ELSE
242  np = nn - 2*pp
243  b2 = z( np-2 )
244  gam = dn1
245  IF( z( np-4 ) .GT. z( np-2 ) )
246  \$ RETURN
247  a2 = z( np-4 ) / z( np-2 )
248  IF( z( nn-9 ) .GT. z( nn-11 ) )
249  \$ RETURN
250  b2 = z( nn-9 ) / z( nn-11 )
251  np = nn - 13
252  END IF
253 *
254 * Approximate contribution to norm squared from I < NN-1.
255 *
256  a2 = a2 + b2
257  DO 10 i4 = np, 4*i0 - 1 + pp, -4
258  IF( b2.EQ.zero )
259  \$ GO TO 20
260  b1 = b2
261  IF( z( i4 ) .GT. z( i4-2 ) )
262  \$ RETURN
263  b2 = b2*( z( i4 ) / z( i4-2 ) )
264  a2 = a2 + b2
265  IF( hundrd*max( b2, b1 ).LT.a2 .OR. cnst1.LT.a2 )
266  \$ GO TO 20
267  10 CONTINUE
268  20 CONTINUE
269  a2 = cnst3*a2
270 *
271 * Rayleigh quotient residual bound.
272 *
273  IF( a2.LT.cnst1 )
274  \$ s = gam*( one-sqrt( a2 ) ) / ( one+a2 )
275  END IF
276  ELSE IF( dmin.EQ.dn2 ) THEN
277 *
278 * Case 5.
279 *
280  ttype = -5
281  s = qurtr*dmin
282 *
283 * Compute contribution to norm squared from I > NN-2.
284 *
285  np = nn - 2*pp
286  b1 = z( np-2 )
287  b2 = z( np-6 )
288  gam = dn2
289  IF( z( np-8 ).GT.b2 .OR. z( np-4 ).GT.b1 )
290  \$ RETURN
291  a2 = ( z( np-8 ) / b2 )*( one+z( np-4 ) / b1 )
292 *
293 * Approximate contribution to norm squared from I < NN-2.
294 *
295  IF( n0-i0.GT.2 ) THEN
296  b2 = z( nn-13 ) / z( nn-15 )
297  a2 = a2 + b2
298  DO 30 i4 = nn - 17, 4*i0 - 1 + pp, -4
299  IF( b2.EQ.zero )
300  \$ GO TO 40
301  b1 = b2
302  IF( z( i4 ) .GT. z( i4-2 ) )
303  \$ RETURN
304  b2 = b2*( z( i4 ) / z( i4-2 ) )
305  a2 = a2 + b2
306  IF( hundrd*max( b2, b1 ).LT.a2 .OR. cnst1.LT.a2 )
307  \$ GO TO 40
308  30 CONTINUE
309  40 CONTINUE
310  a2 = cnst3*a2
311  END IF
312 *
313  IF( a2.LT.cnst1 )
314  \$ s = gam*( one-sqrt( a2 ) ) / ( one+a2 )
315  ELSE
316 *
317 * Case 6, no information to guide us.
318 *
319  IF( ttype.EQ.-6 ) THEN
320  g = g + third*( one-g )
321  ELSE IF( ttype.EQ.-18 ) THEN
322  g = qurtr*third
323  ELSE
324  g = qurtr
325  END IF
326  s = g*dmin
327  ttype = -6
328  END IF
329 *
330  ELSE IF( n0in.EQ.( n0+1 ) ) THEN
331 *
332 * One eigenvalue just deflated. Use DMIN1, DN1 for DMIN and DN.
333 *
334  IF( dmin1.EQ.dn1 .AND. dmin2.EQ.dn2 ) THEN
335 *
336 * Cases 7 and 8.
337 *
338  ttype = -7
339  s = third*dmin1
340  IF( z( nn-5 ).GT.z( nn-7 ) )
341  \$ RETURN
342  b1 = z( nn-5 ) / z( nn-7 )
343  b2 = b1
344  IF( b2.EQ.zero )
345  \$ GO TO 60
346  DO 50 i4 = 4*n0 - 9 + pp, 4*i0 - 1 + pp, -4
347  a2 = b1
348  IF( z( i4 ).GT.z( i4-2 ) )
349  \$ RETURN
350  b1 = b1*( z( i4 ) / z( i4-2 ) )
351  b2 = b2 + b1
352  IF( hundrd*max( b1, a2 ).LT.b2 )
353  \$ GO TO 60
354  50 CONTINUE
355  60 CONTINUE
356  b2 = sqrt( cnst3*b2 )
357  a2 = dmin1 / ( one+b2**2 )
358  gap2 = half*dmin2 - a2
359  IF( gap2.GT.zero .AND. gap2.GT.b2*a2 ) THEN
360  s = max( s, a2*( one-cnst2*a2*( b2 / gap2 )*b2 ) )
361  ELSE
362  s = max( s, a2*( one-cnst2*b2 ) )
363  ttype = -8
364  END IF
365  ELSE
366 *
367 * Case 9.
368 *
369  s = qurtr*dmin1
370  IF( dmin1.EQ.dn1 )
371  \$ s = half*dmin1
372  ttype = -9
373  END IF
374 *
375  ELSE IF( n0in.EQ.( n0+2 ) ) THEN
376 *
377 * Two eigenvalues deflated. Use DMIN2, DN2 for DMIN and DN.
378 *
379 * Cases 10 and 11.
380 *
381  IF( dmin2.EQ.dn2 .AND. two*z( nn-5 ).LT.z( nn-7 ) ) THEN
382  ttype = -10
383  s = third*dmin2
384  IF( z( nn-5 ).GT.z( nn-7 ) )
385  \$ RETURN
386  b1 = z( nn-5 ) / z( nn-7 )
387  b2 = b1
388  IF( b2.EQ.zero )
389  \$ GO TO 80
390  DO 70 i4 = 4*n0 - 9 + pp, 4*i0 - 1 + pp, -4
391  IF( z( i4 ).GT.z( i4-2 ) )
392  \$ RETURN
393  b1 = b1*( z( i4 ) / z( i4-2 ) )
394  b2 = b2 + b1
395  IF( hundrd*b1.LT.b2 )
396  \$ GO TO 80
397  70 CONTINUE
398  80 CONTINUE
399  b2 = sqrt( cnst3*b2 )
400  a2 = dmin2 / ( one+b2**2 )
401  gap2 = z( nn-7 ) + z( nn-9 ) -
402  \$ sqrt( z( nn-11 ) )*sqrt( z( nn-9 ) ) - a2
403  IF( gap2.GT.zero .AND. gap2.GT.b2*a2 ) THEN
404  s = max( s, a2*( one-cnst2*a2*( b2 / gap2 )*b2 ) )
405  ELSE
406  s = max( s, a2*( one-cnst2*b2 ) )
407  END IF
408  ELSE
409  s = qurtr*dmin2
410  ttype = -11
411  END IF
412  ELSE IF( n0in.GT.( n0+2 ) ) THEN
413 *
414 * Case 12, more than two eigenvalues deflated. No information.
415 *
416  s = zero
417  ttype = -12
418  END IF
419 *
420  tau = s
421  RETURN
422 *
423 * End of DLASQ4
424 *

Here is the caller graph for this function: