319 SUBROUTINE dlarrd( RANGE, ORDER, N, VL, VU, IL, IU, GERS,
320 $ reltol, d, e, e2, pivmin, nsplit, isplit,
321 $ m, w, werr, wl, wu, iblock, indexw,
322 $ work, iwork, info )
330 CHARACTER order, range
331 INTEGER il, info, iu, m, n, nsplit
332 DOUBLE PRECISION pivmin, reltol, vl, vu, wl, wu
335 INTEGER iblock( * ), indexw( * ),
336 $ isplit( * ), iwork( * )
337 DOUBLE PRECISION d( * ), e( * ), e2( * ),
338 $ gers( * ), w( * ), werr( * ), work( * )
344 DOUBLE PRECISION zero, one, two, half, fudge
345 parameter( zero = 0.0d0, one = 1.0d0,
346 $ two = 2.0d0, half = one/two,
348 INTEGER allrng, valrng, indrng
349 parameter( allrng = 1, valrng = 2, indrng = 3 )
352 LOGICAL ncnvrg, toofew
353 INTEGER i, ib, ibegin, idiscl, idiscu, ie, iend, iinfo,
354 $ im, in, ioff, iout, irange, itmax, itmp1,
355 $ itmp2, iw, iwoff,
j, jblk, jdisc, je, jee, nb,
357 DOUBLE PRECISION atoli, eps, gl, gu, rtoli, tmp1, tmp2,
358 $ tnorm, uflow, wkill, wlu, wul
374 INTRINSIC abs, int, log, max, min
382 IF(
lsame( range,
'A' ) )
THEN
384 ELSE IF(
lsame( range,
'V' ) )
THEN
386 ELSE IF(
lsame( range,
'I' ) )
THEN
394 IF( irange.LE.0 )
THEN
396 ELSE IF( .NOT.(
lsame(order,
'B').OR.
lsame(order,
'E')) )
THEN
398 ELSE IF( n.LT.0 )
THEN
400 ELSE IF( irange.EQ.valrng )
THEN
403 ELSE IF( irange.EQ.indrng .AND.
404 $ ( il.LT.1 .OR. il.GT.max( 1, n ) ) )
THEN
406 ELSE IF( irange.EQ.indrng .AND.
407 $ ( iu.LT.min( n, il ) .OR. iu.GT.n ) )
THEN
425 IF( irange.EQ.indrng .AND. il.EQ.1 .AND. iu.EQ.n ) irange = 1
435 IF( (irange.EQ.allrng).OR.
436 $ ((irange.EQ.valrng).AND.(d(1).GT.vl).AND.(d(1).LE.vu)).OR.
437 $ ((irange.EQ.indrng).AND.(il.EQ.1).AND.(iu.EQ.1)) )
THEN
450 nb =
ilaenv( 1,
'DSTEBZ',
' ', n, -1, -1, -1 )
457 gl = min( gl, gers( 2*i - 1))
458 gu = max( gu, gers(2*i) )
461 tnorm = max( abs( gl ), abs( gu ) )
462 gl = gl - fudge*tnorm*eps*n - fudge*two*pivmin
463 gu = gu + fudge*tnorm*eps*n + fudge*two*pivmin
476 atoli = fudge*two*uflow + fudge*two*pivmin
478 IF( irange.EQ.indrng )
THEN
483 itmax = int( ( log( tnorm+pivmin )-log( pivmin ) ) /
498 CALL
dlaebz( 3, itmax, n, 2, 2, nb, atoli, rtoli, pivmin,
499 $ d, e, e2, iwork( 5 ), work( n+1 ), work( n+5 ), iout,
500 $ iwork, w, iblock, iinfo )
501 IF( iinfo .NE. 0 )
THEN
506 IF( iwork( 6 ).EQ.iu )
THEN
523 IF( nwl.LT.0 .OR. nwl.GE.n .OR. nwu.LT.1 .OR. nwu.GT.n )
THEN
528 ELSEIF( irange.EQ.valrng )
THEN
532 ELSEIF( irange.EQ.allrng )
THEN
548 DO 70 jblk = 1, nsplit
551 iend = isplit( jblk )
556 IF( wl.GE.d( ibegin )-pivmin )
558 IF( wu.GE.d( ibegin )-pivmin )
560 IF( irange.EQ.allrng .OR.
561 $ ( wl.LT.d( ibegin )-pivmin
562 $ .AND. wu.GE. d( ibegin )-pivmin ) )
THEN
626 DO 40
j = ibegin, iend
627 gl = min( gl, gers( 2*
j - 1))
628 gu = max( gu, gers(2*
j) )
636 gl = gl - fudge*tnorm*eps*in - fudge*pivmin
637 gu = gu + fudge*tnorm*eps*in + fudge*pivmin
639 IF( irange.GT.1 )
THEN
656 CALL
dlaebz( 1, 0, in, in, 1, nb, atoli, rtoli, pivmin,
657 $ d( ibegin ), e( ibegin ), e2( ibegin ),
658 $ idumma, work( n+1 ), work( n+2*in+1 ), im,
659 $ iwork, w( m+1 ), iblock( m+1 ), iinfo )
660 IF( iinfo .NE. 0 )
THEN
665 nwl = nwl + iwork( 1 )
666 nwu = nwu + iwork( in+1 )
667 iwoff = m - iwork( 1 )
670 itmax = int( ( log( gu-gl+pivmin )-log( pivmin ) ) /
672 CALL
dlaebz( 2, itmax, in, in, 1, nb, atoli, rtoli, pivmin,
673 $ d( ibegin ), e( ibegin ), e2( ibegin ),
674 $ idumma, work( n+1 ), work( n+2*in+1 ), iout,
675 $ iwork, w( m+1 ), iblock( m+1 ), iinfo )
676 IF( iinfo .NE. 0 )
THEN
686 tmp1 = half*( work(
j+n )+work(
j+in+n ) )
688 tmp2 = half*abs( work(
j+n )-work(
j+in+n ) )
689 IF(
j.GT.iout-iinfo )
THEN
696 DO 50 je = iwork(
j ) + 1 + iwoff,
697 $ iwork(
j+in ) + iwoff
700 indexw( je ) = je - iwoff
711 IF( irange.EQ.indrng )
THEN
712 idiscl = il - 1 - nwl
715 IF( idiscl.GT.0 )
THEN
720 IF( w( je ).LE.wlu .AND. idiscl.GT.0 )
THEN
725 werr( im ) = werr( je )
726 indexw( im ) = indexw( je )
727 iblock( im ) = iblock( je )
732 IF( idiscu.GT.0 )
THEN
737 IF( w( je ).GE.wul .AND. idiscu.GT.0 )
THEN
742 werr( im ) = werr( je )
743 indexw( im ) = indexw( je )
744 iblock( im ) = iblock( je )
751 werr( jee ) = werr( je )
752 indexw( jee ) = indexw( je )
753 iblock( jee ) = iblock( je )
758 IF( idiscl.GT.0 .OR. idiscu.GT.0 )
THEN
765 IF( idiscl.GT.0 )
THEN
767 DO 100 jdisc = 1, idiscl
770 IF( iblock( je ).NE.0 .AND.
771 $ ( w( je ).LT.wkill .OR. iw.EQ.0 ) )
THEN
779 IF( idiscu.GT.0 )
THEN
781 DO 120 jdisc = 1, idiscu
784 IF( iblock( je ).NE.0 .AND.
785 $ ( w( je ).GE.wkill .OR. iw.EQ.0 ) )
THEN
796 IF( iblock( je ).NE.0 )
THEN
799 werr( im ) = werr( je )
800 indexw( im ) = indexw( je )
801 iblock( im ) = iblock( je )
806 IF( idiscl.LT.0 .OR. idiscu.LT.0 )
THEN
811 IF(( irange.EQ.allrng .AND. m.NE.n ).OR.
812 $ ( irange.EQ.indrng .AND. m.NE.iu-il+1 ) )
THEN
820 IF(
lsame(order,
'E') .AND. nsplit.GT.1 )
THEN
825 IF( w(
j ).LT.tmp1 )
THEN
835 werr( ie ) = werr( je )
836 iblock( ie ) = iblock( je )
837 indexw( ie ) = indexw( je )
LOGICAL function lsame(CA, CB)
LSAME
subroutine dlarrd(RANGE, ORDER, N, VL, VU, IL, IU, GERS, RELTOL, D, E, E2, PIVMIN, NSPLIT, ISPLIT, M, W, WERR, WL, WU, IBLOCK, INDEXW, WORK, IWORK, INFO)
DLARRD computes the eigenvalues of a symmetric tridiagonal matrix to suitable accuracy.
subroutine dlaebz(IJOB, NITMAX, N, MMAX, MINP, NBMIN, ABSTOL, RELTOL, PIVMIN, D, E, E2, NVAL, AB, C, MOUT, NAB, WORK, IWORK, INFO)
DLAEBZ computes the number of eigenvalues of a real symmetric tridiagonal matrix which are less than ...
INTEGER function ilaenv(ISPEC, NAME, OPTS, N1, N2, N3, N4)
set ue cd $ADTTMP cat<< EOF > tmp f Program LinearEquations Implicit none Real j
DOUBLE PRECISION function dlamch(CMACH)
DLAMCH