113 SUBROUTINE dlasq2( N, Z, INFO )
124 DOUBLE PRECISION Z( * )
130 DOUBLE PRECISION CBIAS
131 parameter ( cbias = 1.50d0 )
132 DOUBLE PRECISION ZERO, HALF, ONE, TWO, FOUR, HUNDRD
133 parameter ( zero = 0.0d0, half = 0.5d0, one = 1.0d0,
134 $ two = 2.0d0, four = 4.0d0, hundrd = 100.0d0 )
138 INTEGER I0, I1, I4, IINFO, IPN4, ITER, IWHILA, IWHILB,
139 $ k, kmin, n0, n1, nbig, ndiv, nfail, pp, splt,
141 DOUBLE PRECISION D, DEE, DEEMIN, DESIG, DMIN, DMIN1, DMIN2, DN,
142 $ dn1, dn2, e, emax, emin, eps, g, oldemn, qmax,
143 $ qmin, s, safmin, sigma, t, tau, temp, tol,
144 $ tol2, trace, zmax, tempe, tempq
151 DOUBLE PRECISION DLAMCH
152 EXTERNAL dlamch, ilaenv
155 INTRINSIC abs, dble, max, min, sqrt
163 eps = dlamch(
'Precision' )
164 safmin = dlamch(
'Safe minimum' )
170 CALL xerbla(
'DLASQ2', 1 )
172 ELSE IF( n.EQ.0 )
THEN
174 ELSE IF( n.EQ.1 )
THEN
178 IF( z( 1 ).LT.zero )
THEN
180 CALL xerbla(
'DLASQ2', 2 )
183 ELSE IF( n.EQ.2 )
THEN
187 IF( z( 2 ).LT.zero .OR. z( 3 ).LT.zero )
THEN
189 CALL xerbla(
'DLASQ2', 2 )
191 ELSE IF( z( 3 ).GT.z( 1 ) )
THEN
196 z( 5 ) = z( 1 ) + z( 2 ) + z( 3 )
197 IF( z( 2 ).GT.z( 3 )*tol2 )
THEN
198 t = half*( ( z( 1 )-z( 3 ) )+z( 2 ) )
199 s = z( 3 )*( z( 2 ) / t )
201 s = z( 3 )*( z( 2 ) / ( t*( one+sqrt( one+s / t ) ) ) )
203 s = z( 3 )*( z( 2 ) / ( t+sqrt( t )*sqrt( t+s ) ) )
205 t = z( 1 ) + ( s+z( 2 ) )
206 z( 3 ) = z( 3 )*( z( 1 ) / t )
210 z( 6 ) = z( 2 ) + z( 1 )
223 DO 10 k = 1, 2*( n-1 ), 2
224 IF( z( k ).LT.zero )
THEN
226 CALL xerbla(
'DLASQ2', 2 )
228 ELSE IF( z( k+1 ).LT.zero )
THEN
230 CALL xerbla(
'DLASQ2', 2 )
235 qmax = max( qmax, z( k ) )
236 emin = min( emin, z( k+1 ) )
237 zmax = max( qmax, zmax, z( k+1 ) )
239 IF( z( 2*n-1 ).LT.zero )
THEN
240 info = -( 200+2*n-1 )
241 CALL xerbla(
'DLASQ2', 2 )
245 qmax = max( qmax, z( 2*n-1 ) )
246 zmax = max( qmax, zmax )
254 CALL dlasrt(
'D', n, z, iinfo )
263 IF( trace.EQ.zero )
THEN
270 ieee = ilaenv( 10,
'DLASQ2',
'N', 1, 2, 3, 4 ).EQ.1 .AND.
271 $ ilaenv( 11,
'DLASQ2',
'N', 1, 2, 3, 4 ).EQ.1
279 z( 2*k-3 ) = z( k-1 )
287 IF( cbias*z( 4*i0-3 ).LT.z( 4*n0-3 ) )
THEN
289 DO 40 i4 = 4*i0, 2*( i0+n0-1 ), 4
291 z( i4-3 ) = z( ipn4-i4-3 )
292 z( ipn4-i4-3 ) = temp
294 z( i4-1 ) = z( ipn4-i4-5 )
295 z( ipn4-i4-5 ) = temp
306 DO 50 i4 = 4*( n0-1 ) + pp, 4*i0 + pp, -4
307 IF( z( i4-1 ).LE.tol2*d )
THEN
311 d = z( i4-3 )*( d / ( d+z( i4-1 ) ) )
317 emin = z( 4*i0+pp+1 )
319 DO 60 i4 = 4*i0 + pp, 4*( n0-1 ) + pp, 4
320 z( i4-2*pp-2 ) = d + z( i4-1 )
321 IF( z( i4-1 ).LE.tol2*d )
THEN
326 ELSE IF( safmin*z( i4+1 ).LT.z( i4-2*pp-2 ) .AND.
327 $ safmin*z( i4-2*pp-2 ).LT.z( i4+1 ) )
THEN
328 temp = z( i4+1 ) / z( i4-2*pp-2 )
329 z( i4-2*pp ) = z( i4-1 )*temp
332 z( i4-2*pp ) = z( i4+1 )*( z( i4-1 ) / z( i4-2*pp-2 ) )
333 d = z( i4+1 )*( d / z( i4-2*pp-2 ) )
335 emin = min( emin, z( i4-2*pp ) )
341 qmax = z( 4*i0-pp-2 )
342 DO 70 i4 = 4*i0 - pp + 2, 4*n0 - pp - 2, 4
343 qmax = max( qmax, z( i4 ) )
366 DO 160 iwhila = 1, n + 1
381 IF( sigma.LT.zero )
THEN
391 emin = abs( z( 4*n0-5 ) )
397 DO 90 i4 = 4*n0, 8, -4
398 IF( z( i4-5 ).LE.zero )
400 IF( qmin.GE.four*emax )
THEN
401 qmin = min( qmin, z( i4-3 ) )
402 emax = max( emax, z( i4-5 ) )
404 qmax = max( qmax, z( i4-7 )+z( i4-5 ) )
405 emin = min( emin, z( i4-5 ) )
413 IF( n0-i0.GT.1 )
THEN
417 DO 110 i4 = 4*i0+1, 4*n0-3, 4
418 dee = z( i4 )*( dee /( dee+z( i4-2 ) ) )
419 IF( dee.LE.deemin )
THEN
424 IF( (kmin-i0)*2.LT.n0-kmin .AND.
425 $ deemin.LE.half*z(4*n0-3) )
THEN
428 DO 120 i4 = 4*i0, 2*( i0+n0-1 ), 4
430 z( i4-3 ) = z( ipn4-i4-3 )
431 z( ipn4-i4-3 ) = temp
433 z( i4-2 ) = z( ipn4-i4-2 )
434 z( ipn4-i4-2 ) = temp
436 z( i4-1 ) = z( ipn4-i4-5 )
437 z( ipn4-i4-5 ) = temp
439 z( i4 ) = z( ipn4-i4-4 )
440 z( ipn4-i4-4 ) = temp
447 dmin = -max( zero, qmin-two*sqrt( qmin )*sqrt( emax ) )
455 nbig = 100*( n0-i0+1 )
456 DO 140 iwhilb = 1, nbig
462 CALL dlasq3( i0, n0, z, pp, dmin, sigma, desig, qmax, nfail,
463 $ iter, ndiv, ieee, ttype, dmin1, dmin2, dn, dn1,
470 IF( pp.EQ.0 .AND. n0-i0.GE.3 )
THEN
471 IF( z( 4*n0 ).LE.tol2*qmax .OR.
472 $ z( 4*n0-1 ).LE.tol2*sigma )
THEN
477 DO 130 i4 = 4*i0, 4*( n0-3 ), 4
478 IF( z( i4 ).LE.tol2*z( i4-3 ) .OR.
479 $ z( i4-1 ).LE.tol2*sigma )
THEN
486 qmax = max( qmax, z( i4+1 ) )
487 emin = min( emin, z( i4-1 ) )
488 oldemn = min( oldemn, z( i4 ) )
509 z( 4*i0-3 ) = z( 4*i0-3 ) + sigma
512 z( 4*k-5 ) = z( 4*k-5 ) * (tempq / z( 4*k-7 ))
514 z( 4*k-3 ) = z( 4*k-3 ) + sigma + tempe - z( 4*k-5 )
521 DO WHILE( ( i1.GE.2 ) .AND. ( z(4*i1-5).GE.zero ) )
529 z( 2*k-1 ) = z( 4*k-3 )
536 z( 2*k ) = z( 4*k-1 )
564 CALL dlasrt(
'D', n, z, iinfo )
575 z( 2*n+3 ) = dble( iter )
576 z( 2*n+4 ) = dble( ndiv ) / dble( n**2 )
577 z( 2*n+5 ) = hundrd*nfail / dble( iter )
subroutine dlasq2(N, Z, INFO)
DLASQ2 computes all the eigenvalues of the symmetric positive definite tridiagonal matrix associated ...
subroutine xerbla(SRNAME, INFO)
XERBLA
subroutine dlasq3(I0, N0, Z, PP, DMIN, SIGMA, DESIG, QMAX, NFAIL, ITER, NDIV, IEEE, TTYPE, DMIN1, DMIN2, DN, DN1, DN2, G, TAU)
DLASQ3 checks for deflation, computes a shift and calls dqds. Used by sbdsqr.
subroutine dlasrt(ID, N, D, INFO)
DLASRT sorts numbers in increasing or decreasing order.