dlamc2()

subroutine dlamc2	(	integer	beta,
		integer	t,
		logical	rnd,
		double precision	eps,
		integer	emin,
		double precision	rmin,
		integer	emax,
		double precision	rmax
	)

Definition at line 326 of file tools.f.

*
*  -- LAPACK auxiliary routine (version 2.0) --
*     Univ. of Tennessee, Univ. of California Berkeley, NAG Ltd.,
*     Courant Institute, Argonne National Lab, and Rice University
*     October 31, 1992
*
*     .. Scalar Arguments ..
      LOGICAL            RND
      INTEGER            BETA, EMAX, EMIN, T
      DOUBLE PRECISION   EPS, RMAX, RMIN
*     ..
*
*  Purpose
*  =======
*
*  DLAMC2 determines the machine parameters specified in its argument
*  list.
*
*  Arguments
*  =========
*
*  BETA    (output) INTEGER
*          The base of the machine.
*
*  T       (output) INTEGER
*          The number of ( BETA ) digits in the mantissa.
*
*  RND     (output) LOGICAL
*          Specifies whether proper rounding  ( RND = .TRUE. )  or
*          chopping  ( RND = .FALSE. )  occurs in addition. This may not
*          be a reliable guide to the way in which the machine performs
*          its arithmetic.
*
*  EPS     (output) DOUBLE PRECISION
*          The smallest positive number such that
*
*             fl( 1.0 - EPS ) .LT. 1.0,
*
*          where fl denotes the computed value.
*
*  EMIN    (output) INTEGER
*          The minimum exponent before (gradual) underflow occurs.
*
*  RMIN    (output) DOUBLE PRECISION
*          The smallest normalized number for the machine, given by
*          BASE**( EMIN - 1 ), where  BASE  is the floating point value
*          of BETA.
*
*  EMAX    (output) INTEGER
*          The maximum exponent before overflow occurs.
*
*  RMAX    (output) DOUBLE PRECISION
*          The largest positive number for the machine, given by
*          BASE**EMAX * ( 1 - EPS ), where  BASE  is the floating point
*          value of BETA.
*
*  Further Details
*  ===============
*
*  The computation of  EPS  is based on a routine PARANOIA by
*  W. Kahan of the University of California at Berkeley.
*
* =====================================================================
*
*     .. Local Scalars ..
      LOGICAL            FIRST, IEEE, IWARN, LIEEE1, LRND
      INTEGER            GNMIN, GPMIN, I, LBETA, LEMAX, LEMIN, LT,
     $                   NGNMIN, NGPMIN
      DOUBLE PRECISION   A, B, C, HALF, LEPS, LRMAX, LRMIN, ONE, RBASE,
     $                   SIXTH, SMALL, THIRD, TWO, ZERO
*     ..
*     .. External Functions ..
      DOUBLE PRECISION   DLAMC3
      EXTERNAL           dlamc3
*     ..
*     .. External Subroutines ..
      EXTERNAL           dlamc1, dlamc4, dlamc5
*     ..
*     .. Intrinsic Functions ..
      INTRINSIC          abs, max, min
*     ..
*     .. Save statement ..
      SAVE               first, iwarn, lbeta, lemax, lemin, leps, lrmax,
     $                   lrmin, lt
*     ..
*     .. Data statements ..
      DATA               first / .true. / , iwarn / .false. /
*     ..
*     .. Executable Statements ..
*
      IF( first ) THEN
         first = .false.
         zero = 0
         one = 1
         two = 2
*
*        LBETA, LT, LRND, LEPS, LEMIN and LRMIN  are the local values of
*        BETA, T, RND, EPS, EMIN and RMIN.
*
*        Throughout this routine  we use the function  DLAMC3  to ensure
*        that relevant values are stored  and not held in registers,  or
*        are not affected by optimizers.
*
*        DLAMC1 returns the parameters  LBETA, LT, LRND and LIEEE1.
*
         CALL dlamc1( lbeta, lt, lrnd, lieee1 )
*
*        Start to find EPS.
*
         b = lbeta
         a = b**( -lt )
         leps = a
*
*        Try some tricks to see whether or not this is the correct  EPS.
*
         b = two / 3
         half = one / 2
         sixth = dlamc3( b, -half )
         third = dlamc3( sixth, sixth )
         b = dlamc3( third, -half )
         b = dlamc3( b, sixth )
         b = abs( b )
         IF( b.LT.leps )
     $      b = leps
*
         leps = 1
*
*+       WHILE( ( LEPS.GT.B ).AND.( B.GT.ZERO ) )LOOP
   10    CONTINUE
         IF( ( leps.GT.b ) .AND. ( b.GT.zero ) ) THEN
            leps = b
            c = dlamc3( half*leps, ( two**5 )*( leps**2 ) )
            c = dlamc3( half, -c )
            b = dlamc3( half, c )
            c = dlamc3( half, -b )
            b = dlamc3( half, c )
            GO TO 10
         END IF
*+       END WHILE
*
         IF( a.LT.leps )
     $      leps = a
*
*        Computation of EPS complete.
*
*        Now find  EMIN.  Let A = + or - 1, and + or - (1 + BASE**(-3)).
*        Keep dividing  A by BETA until (gradual) underflow occurs. This
*        is detected when we cannot recover the previous A.
*
         rbase = one / lbeta
         small = one
         DO 20 i = 1, 3
            small = dlamc3( small*rbase, zero )
   20    CONTINUE
         a = dlamc3( one, small )
         CALL dlamc4( ngpmin, one, lbeta )
         CALL dlamc4( ngnmin, -one, lbeta )
         CALL dlamc4( gpmin, a, lbeta )
         CALL dlamc4( gnmin, -a, lbeta )
         ieee = .false.
*
         IF( ( ngpmin.EQ.ngnmin ) .AND. ( gpmin.EQ.gnmin ) ) THEN
            IF( ngpmin.EQ.gpmin ) THEN
               lemin = ngpmin
*            ( Non twos-complement machines, no gradual underflow;
*              e.g.,  VAX )
            ELSE IF( ( gpmin-ngpmin ).EQ.3 ) THEN
               lemin = ngpmin - 1 + lt
               ieee = .true.
*            ( Non twos-complement machines, with gradual underflow;
*              e.g., IEEE standard followers )
            ELSE
               lemin = min( ngpmin, gpmin )
*            ( A guess; no known machine )
               iwarn = .true.
            END IF
*
         ELSE IF( ( ngpmin.EQ.gpmin ) .AND. ( ngnmin.EQ.gnmin ) ) THEN
            IF( abs( ngpmin-ngnmin ).EQ.1 ) THEN
               lemin = max( ngpmin, ngnmin )
*            ( Twos-complement machines, no gradual underflow;
*              e.g., CYBER 205 )
            ELSE
               lemin = min( ngpmin, ngnmin )
*            ( A guess; no known machine )
               iwarn = .true.
            END IF
*
         ELSE IF( ( abs( ngpmin-ngnmin ).EQ.1 ) .AND.
     $            ( gpmin.EQ.gnmin ) ) THEN
            IF( ( gpmin-min( ngpmin, ngnmin ) ).EQ.3 ) THEN
               lemin = max( ngpmin, ngnmin ) - 1 + lt
*            ( Twos-complement machines with gradual underflow;
*              no known machine )
            ELSE
               lemin = min( ngpmin, ngnmin )
*            ( A guess; no known machine )
               iwarn = .true.
            END IF
*
         ELSE
            lemin = min( ngpmin, ngnmin, gpmin, gnmin )
*         ( A guess; no known machine )
            iwarn = .true.
         END IF
***
* Comment out this if block if EMIN is ok
         IF( iwarn ) THEN
            first = .true.
            WRITE( 6, fmt = 9999 )lemin
         END IF
***
*
*        Assume IEEE arithmetic if we found denormalised  numbers above,
*        or if arithmetic seems to round in the  IEEE style,  determined
*        in routine DLAMC1. A true IEEE machine should have both  things
*        true; however, faulty machines may have one or the other.
*
         ieee = ieee .OR. lieee1
*
*        Compute  RMIN by successive division by  BETA. We could compute
*        RMIN as BASE**( EMIN - 1 ),  but some machines underflow during
*        this computation.
*
         lrmin = 1
         DO 30 i = 1, 1 - lemin
            lrmin = dlamc3( lrmin*rbase, zero )
   30    CONTINUE
*
*        Finally, call DLAMC5 to compute EMAX and RMAX.
*
         CALL dlamc5( lbeta, lt, lemin, ieee, lemax, lrmax )
      END IF
*
      beta = lbeta
      t = lt
      rnd = lrnd
      eps = leps
      emin = lemin
      rmin = lrmin
      emax = lemax
      rmax = lrmax
*
      RETURN
*
 9999 FORMAT( / / ' WARNING. The value EMIN may be incorrect:-',
     $      '  EMIN = ', i8, /
     $      ' If, after inspection, the value EMIN looks',
     $      ' acceptable please comment out ',
     $      / ' the IF block as marked within the code of routine',
     $      ' DLAMC2,', / ' otherwise supply EMIN explicitly.', / )
*
*     End of DLAMC2
*

Here is the call graph for this function:

Here is the caller graph for this function: