A RetroSearch Logo

Home - News ( United States | United Kingdom | Italy | Germany ) - Football scores

Search Query:

Showing content from http://hackage.haskell.org/packages/archive/base/4.5.0.0/doc/html/src/GHC-Float.html below:

GHC/Float.lhs

\begin{code}





















#include "ieee-flpt.h"


module GHC.Float( module GHC.Float, Float(..), Double(..), Float#, Double#
                , double2Int, int2Double, float2Int, int2Float )
    where

import Data.Maybe

import Data.Bits
import GHC.Base
import GHC.List
import GHC.Enum
import GHC.Show
import GHC.Num
import GHC.Real
import GHC.Arr
import GHC.Float.RealFracMethods
import GHC.Float.ConversionUtils
import GHC.Integer.Logarithms ( integerLogBase# )
import GHC.Integer.Logarithms.Internals

infixr 8  **
\end{code} %********************************************************* %* * \subsection{Standard numeric classes} %* * %********************************************************* \begin{code}





class  (Fractional a) => Floating a  where
    pi                  :: a
    exp, log, sqrt      :: a -> a
    (**), logBase       :: a -> a -> a
    sin, cos, tan       :: a -> a
    asin, acos, atan    :: a -> a
    sinh, cosh, tanh    :: a -> a
    asinh, acosh, atanh :: a -> a

    
    
    
    
    
    x ** y              =  exp (log x * y)
    logBase x y         =  log y / log x
    sqrt x              =  x ** 0.5
    tan  x              =  sin  x / cos  x
    tanh x              =  sinh x / cosh x






class  (RealFrac a, Floating a) => RealFloat a  where
    
    
    floatRadix          :: a -> Integer
    
    
    floatDigits         :: a -> Int
    
    
    floatRange          :: a -> (Int,Int)
    
    
    
    
    
    
    
    
    
    
    
    decodeFloat         :: a -> (Integer,Int)
    
    
    
    
    
    
    
    encodeFloat         :: Integer -> Int -> a
    
    
    
    
    
    
    
    exponent            :: a -> Int
    
    
    
    
    significand         :: a -> a
    
    scaleFloat          :: Int -> a -> a
    
    isNaN               :: a -> Bool
    
    isInfinite          :: a -> Bool
    
    
    isDenormalized      :: a -> Bool
    
    isNegativeZero      :: a -> Bool
    
    isIEEE              :: a -> Bool
    
    
    
    
    
    
    
    
    
    atan2               :: a -> a -> a


    exponent x          =  if m == 0 then 0 else n + floatDigits x
                           where (m,n) = decodeFloat x

    significand x       =  encodeFloat m (negate (floatDigits x))
                           where (m,_) = decodeFloat x

    scaleFloat 0 x      =  x
    scaleFloat k x
      | isFix           =  x
      | otherwise       =  encodeFloat m (n + clamp b k)
                           where (m,n) = decodeFloat x
                                 (l,h) = floatRange x
                                 d     = floatDigits x
                                 b     = h  l + 4*d
                                 
                                 
                                 
                                 
                                 
                                 
                                 
                                 
                                 isFix = x == 0 || isNaN x || isInfinite x

    atan2 y x
      | x > 0            =  atan (y/x)
      | x == 0 && y > 0  =  pi/2
      | x <  0 && y > 0  =  pi + atan (y/x)
      |(x <= 0 && y < 0)            ||
       (x <  0 && isNegativeZero y) ||
       (isNegativeZero x && isNegativeZero y)
                         = atan2 (y) x
      | y == 0 && (x < 0 || isNegativeZero x)
                          =  pi    
      | x==0 && y==0      =  y     
      | otherwise         =  x + y 
\end{code} %********************************************************* %* * \subsection{Type @Float@} %* * %********************************************************* \begin{code}
instance  Num Float  where
    (+)         x y     =  plusFloat x y
    ()         x y     =  minusFloat x y
    negate      x       =  negateFloat x
    (*)         x y     =  timesFloat x y
    abs x | x >= 0.0    =  x
          | otherwise   =  negateFloat x
    signum x | x == 0.0  = 0
             | x > 0.0   = 1
             | otherwise = negate 1

    
    fromInteger i = F# (floatFromInteger i)

instance  Real Float  where
    toRational (F# x#)  =
        case decodeFloat_Int# x# of
          (# m#, e# #)
            | e# >=# 0#                                 ->
                    (smallInteger m# `shiftLInteger` e#) :% 1
            | (int2Word# m# `and#` 1##) `eqWord#` 0##   ->
                    case elimZerosInt# m# (negateInt# e#) of
                      (# n, d# #) -> n :% shiftLInteger 1 d#
            | otherwise                                 ->
                    smallInteger m# :% shiftLInteger 1 (negateInt# e#)

instance  Fractional Float  where
    (/) x y             =  divideFloat x y
    fromRational (n:%0)
        | n == 0        = 0/0
        | n < 0         = (1)/0
        | otherwise     = 1/0
    fromRational (n:%d)
        | n == 0        = encodeFloat 0 0
        | n < 0         = (fromRat'' minEx mantDigs (n) d)
        | otherwise     = fromRat'' minEx mantDigs n d
          where
            minEx       = FLT_MIN_EXP
            mantDigs    = FLT_MANT_DIG
    recip x             =  1.0 / x



instance  RealFrac Float  where

        
    
    
    


#if FLT_RADIX != 2
#error FLT_RADIX must be 2
#endif
    properFraction (F# x#)
      = case decodeFloat_Int# x# of
        (# m#, n# #) ->
            let m = I# m#
                n = I# n#
            in
            if n >= 0
            then (fromIntegral m * (2 ^ n), 0.0)
            else let i = if m >= 0 then                m `shiftR` negate n
                                   else negate (negate m `shiftR` negate n)
                     f = m  (i `shiftL` negate n)
                 in (fromIntegral i, encodeFloat (fromIntegral f) n)

    truncate x  = case properFraction x of
                     (n,_) -> n

    round x     = case properFraction x of
                     (n,r) -> let
                                m         = if r < 0.0 then n  1 else n + 1
                                half_down = abs r  0.5
                              in
                              case (compare half_down 0.0) of
                                LT -> n
                                EQ -> if even n then n else m
                                GT -> m

    ceiling x   = case properFraction x of
                    (n,r) -> if r > 0.0 then n + 1 else n

    floor x     = case properFraction x of
                    (n,r) -> if r < 0.0 then n  1 else n

instance  Floating Float  where
    pi                  =  3.141592653589793238
    exp x               =  expFloat x
    log x               =  logFloat x
    sqrt x              =  sqrtFloat x
    sin x               =  sinFloat x
    cos x               =  cosFloat x
    tan x               =  tanFloat x
    asin x              =  asinFloat x
    acos x              =  acosFloat x
    atan x              =  atanFloat x
    sinh x              =  sinhFloat x
    cosh x              =  coshFloat x
    tanh x              =  tanhFloat x
    (**) x y            =  powerFloat x y
    logBase x y         =  log y / log x

    asinh x = log (x + sqrt (1.0+x*x))
    acosh x = log (x + (x+1.0) * sqrt ((x1.0)/(x+1.0)))
    atanh x = 0.5 * log ((1.0+x) / (1.0x))

instance  RealFloat Float  where
    floatRadix _        =  FLT_RADIX        
    floatDigits _       =  FLT_MANT_DIG     
    floatRange _        =  (FLT_MIN_EXP, FLT_MAX_EXP) 

    decodeFloat (F# f#) = case decodeFloat_Int# f# of
                          (# i, e #) -> (smallInteger i, I# e)

    encodeFloat i (I# e) = F# (encodeFloatInteger i e)

    exponent x          = case decodeFloat x of
                            (m,n) -> if m == 0 then 0 else n + floatDigits x

    significand x       = case decodeFloat x of
                            (m,_) -> encodeFloat m (negate (floatDigits x))

    scaleFloat 0 x      = x
    scaleFloat k x
      | isFix           = x
      | otherwise       = case decodeFloat x of
                            (m,n) -> encodeFloat m (n + clamp bf k)
                        where bf = FLT_MAX_EXP  (FLT_MIN_EXP) + 4*FLT_MANT_DIG
                              isFix = x == 0 || isFloatFinite x == 0

    isNaN x          = 0 /= isFloatNaN x
    isInfinite x     = 0 /= isFloatInfinite x
    isDenormalized x = 0 /= isFloatDenormalized x
    isNegativeZero x = 0 /= isFloatNegativeZero x
    isIEEE _         = True

instance  Show Float  where
    showsPrec   x = showSignedFloat showFloat x
    showList = showList__ (showsPrec 0)
\end{code} %********************************************************* %* * \subsection{Type @Double@} %* * %********************************************************* \begin{code}
instance  Num Double  where
    (+)         x y     =  plusDouble x y
    ()         x y     =  minusDouble x y
    negate      x       =  negateDouble x
    (*)         x y     =  timesDouble x y
    abs x | x >= 0.0    =  x
          | otherwise   =  negateDouble x
    signum x | x == 0.0  = 0
             | x > 0.0   = 1
             | otherwise = negate 1

    
    fromInteger i = D# (doubleFromInteger i)


instance  Real Double  where
    toRational (D# x#)  =
        case decodeDoubleInteger x# of
          (# m, e# #)
            | e# >=# 0#                                     ->
                shiftLInteger m e# :% 1
            | (integerToWord m `and#` 1##) `eqWord#` 0##    ->
                case elimZerosInteger m (negateInt# e#) of
                    (# n, d# #) ->  n :% shiftLInteger 1 d#
            | otherwise                                     ->
                m :% shiftLInteger 1 (negateInt# e#)

instance  Fractional Double  where
    (/) x y             =  divideDouble x y
    fromRational (n:%0)
        | n == 0        = 0/0
        | n < 0         = (1)/0
        | otherwise     = 1/0
    fromRational (n:%d)
        | n == 0        = encodeFloat 0 0
        | n < 0         = (fromRat'' minEx mantDigs (n) d)
        | otherwise     = fromRat'' minEx mantDigs n d
          where
            minEx       = DBL_MIN_EXP
            mantDigs    = DBL_MANT_DIG
    recip x             =  1.0 / x

instance  Floating Double  where
    pi                  =  3.141592653589793238
    exp x               =  expDouble x
    log x               =  logDouble x
    sqrt x              =  sqrtDouble x
    sin  x              =  sinDouble x
    cos  x              =  cosDouble x
    tan  x              =  tanDouble x
    asin x              =  asinDouble x
    acos x              =  acosDouble x
    atan x              =  atanDouble x
    sinh x              =  sinhDouble x
    cosh x              =  coshDouble x
    tanh x              =  tanhDouble x
    (**) x y            =  powerDouble x y
    logBase x y         =  log y / log x

    asinh x = log (x + sqrt (1.0+x*x))
    acosh x = log (x + (x+1.0) * sqrt ((x1.0)/(x+1.0)))
    atanh x = 0.5 * log ((1.0+x) / (1.0x))



instance  RealFrac Double  where

        
    
    
    

    properFraction x
      = case (decodeFloat x)      of { (m,n) ->
        if n >= 0 then
            (fromInteger m * 2 ^ n, 0.0)
        else
            case (quotRem m (2^(negate n))) of { (w,r) ->
            (fromInteger w, encodeFloat r n)
            }
        }

    truncate x  = case properFraction x of
                     (n,_) -> n

    round x     = case properFraction x of
                     (n,r) -> let
                                m         = if r < 0.0 then n  1 else n + 1
                                half_down = abs r  0.5
                              in
                              case (compare half_down 0.0) of
                                LT -> n
                                EQ -> if even n then n else m
                                GT -> m

    ceiling x   = case properFraction x of
                    (n,r) -> if r > 0.0 then n + 1 else n

    floor x     = case properFraction x of
                    (n,r) -> if r < 0.0 then n  1 else n

instance  RealFloat Double  where
    floatRadix _        =  FLT_RADIX        
    floatDigits _       =  DBL_MANT_DIG     
    floatRange _        =  (DBL_MIN_EXP, DBL_MAX_EXP) 

    decodeFloat (D# x#)
      = case decodeDoubleInteger x#   of
          (# i, j #) -> (i, I# j)

    encodeFloat i (I# j) = D# (encodeDoubleInteger i j)

    exponent x          = case decodeFloat x of
                            (m,n) -> if m == 0 then 0 else n + floatDigits x

    significand x       = case decodeFloat x of
                            (m,_) -> encodeFloat m (negate (floatDigits x))

    scaleFloat 0 x      = x
    scaleFloat k x
      | isFix           = x
      | otherwise       = case decodeFloat x of
                            (m,n) -> encodeFloat m (n + clamp bd k)
                        where bd = DBL_MAX_EXP  (DBL_MIN_EXP) + 4*DBL_MANT_DIG
                              isFix = x == 0 || isDoubleFinite x == 0

    isNaN x             = 0 /= isDoubleNaN x
    isInfinite x        = 0 /= isDoubleInfinite x
    isDenormalized x    = 0 /= isDoubleDenormalized x
    isNegativeZero x    = 0 /= isDoubleNegativeZero x
    isIEEE _            = True

instance  Show Double  where
    showsPrec   x = showSignedFloat showFloat x
    showList = showList__ (showsPrec 0)
\end{code} %********************************************************* %* * \subsection{@Enum@ instances} %* * %********************************************************* The @Enum@ instances for Floats and Doubles are slightly unusual. The @toEnum@ function truncates numbers to Int. The definitions of @enumFrom@ and @enumFromThen@ allow floats to be used in arithmetic series: [0,0.1 .. 1.0]. However, roundoff errors make these somewhat dubious. This example may have either 10 or 11 elements, depending on how 0.1 is represented. NOTE: The instances for Float and Double do not make use of the default methods for @enumFromTo@ and @enumFromThenTo@, as these rely on there being a `non-lossy' conversion to and from Ints. Instead we make use of the 1.2 default methods (back in the days when Enum had Ord as a superclass) for these (@numericEnumFromTo@ and @numericEnumFromThenTo@ below.) \begin{code}
instance  Enum Float  where
    succ x         = x + 1
    pred x         = x  1
    toEnum         = int2Float
    fromEnum       = fromInteger . truncate   
    enumFrom       = numericEnumFrom
    enumFromTo     = numericEnumFromTo
    enumFromThen   = numericEnumFromThen
    enumFromThenTo = numericEnumFromThenTo

instance  Enum Double  where
    succ x         = x + 1
    pred x         = x  1
    toEnum         =  int2Double
    fromEnum       =  fromInteger . truncate   
    enumFrom       =  numericEnumFrom
    enumFromTo     =  numericEnumFromTo
    enumFromThen   =  numericEnumFromThen
    enumFromThenTo =  numericEnumFromThenTo
\end{code} %********************************************************* %* * \subsection{Printing floating point} %* * %********************************************************* \begin{code}



showFloat :: (RealFloat a) => a -> ShowS
showFloat x  =  showString (formatRealFloat FFGeneric Nothing x)



data FFFormat = FFExponent | FFFixed | FFGeneric

formatRealFloat :: (RealFloat a) => FFFormat -> Maybe Int -> a -> String
formatRealFloat fmt decs x
   | isNaN x                   = "NaN"
   | isInfinite x              = if x < 0 then "-Infinity" else "Infinity"
   | x < 0 || isNegativeZero x = '-':doFmt fmt (floatToDigits (toInteger base) (x))
   | otherwise                 = doFmt fmt (floatToDigits (toInteger base) x)
 where
  base = 10

  doFmt format (is, e) =
    let ds = map intToDigit is in
    case format of
     FFGeneric ->
      doFmt (if e < 0 || e > 7 then FFExponent else FFFixed)
            (is,e)
     FFExponent ->
      case decs of
       Nothing ->
        let show_e' = show (e1) in
        case ds of
          "0"     -> "0.0e0"
          [d]     -> d : ".0e" ++ show_e'
          (d:ds') -> d : '.' : ds' ++ "e" ++ show_e'
          []      -> error "formatRealFloat/doFmt/FFExponent: []"
       Just dec ->
        let dec' = max dec 1 in
        case is of
         [0] -> '0' :'.' : take dec' (repeat '0') ++ "e0"
         _ ->
          let
           (ei,is') = roundTo base (dec'+1) is
           (d:ds') = map intToDigit (if ei > 0 then init is' else is')
          in
          d:'.':ds' ++ 'e':show (e1+ei)
     FFFixed ->
      let
       mk0 ls = case ls of { "" -> "0" ; _ -> ls}
      in
      case decs of
       Nothing
          | e <= 0    -> "0." ++ replicate (e) '0' ++ ds
          | otherwise ->
             let
                f 0 s    rs  = mk0 (reverse s) ++ '.':mk0 rs
                f n s    ""  = f (n1) ('0':s) ""
                f n s (r:rs) = f (n1) (r:s) rs
             in
                f e "" ds
       Just dec ->
        let dec' = max dec 0 in
        if e >= 0 then
         let
          (ei,is') = roundTo base (dec' + e) is
          (ls,rs)  = splitAt (e+ei) (map intToDigit is')
         in
         mk0 ls ++ (if null rs then "" else '.':rs)
        else
         let
          (ei,is') = roundTo base dec' (replicate (e) 0 ++ is)
          d:ds' = map intToDigit (if ei > 0 then is' else 0:is')
         in
         d : (if null ds' then "" else '.':ds')


roundTo :: Int -> Int -> [Int] -> (Int,[Int])
roundTo base d is =
  case f d is of
    x@(0,_) -> x
    (1,xs)  -> (1, 1:xs)
    _       -> error "roundTo: bad Value"
 where
  b2 = base `div` 2

  f n []     = (0, replicate n 0)
  f 0 (x:_)  = (if x >= b2 then 1 else 0, [])
  f n (i:xs)
     | i' == base = (1,0:ds)
     | otherwise  = (0,i':ds)
      where
       (c,ds) = f (n1) xs
       i'     = c + i



















floatToDigits :: (RealFloat a) => Integer -> a -> ([Int], Int)
floatToDigits _ 0 = ([0], 0)
floatToDigits base x =
 let
  (f0, e0) = decodeFloat x
  (minExp0, _) = floatRange x
  p = floatDigits x
  b = floatRadix x
  minExp = minExp0  p 
  
  
  (f, e) =
   let n = minExp  e0 in
   if n > 0 then (f0 `quot` (expt b n), e0+n) else (f0, e0)
  (r, s, mUp, mDn) =
   if e >= 0 then
    let be = expt b e in
    if f == expt b (p1) then
      (f*be*b*2, 2*b, be*b, be)     
    else
      (f*be*2, 2, be, be)
   else
    if e > minExp && f == expt b (p1) then
      (f*b*2, expt b (e+1)*2, b, 1)
    else
      (f*2, expt b (e)*2, 1, 1)
  k :: Int
  k =
   let
    k0 :: Int
    k0 =
     if b == 2 && base == 10 then
        
        
        
        
        
        
        
        
        
        
        
        
        
        
        
        
        
        
        let lx = p  1 + e0
            k1 = (lx * 8651) `quot` 28738
        in if lx >= 0 then k1 + 1 else k1
     else
	
        
        ceiling ((log (fromInteger (f+1) :: Float) +
                 fromIntegral e * log (fromInteger b)) /
                   log (fromInteger base))


    fixup n =
      if n >= 0 then
        if r + mUp <= expt base n * s then n else fixup (n+1)
      else
        if expt base (n) * (r + mUp) <= s then n else fixup (n+1)
   in
   fixup k0

  gen ds rn sN mUpN mDnN =
   let
    (dn, rn') = (rn * base) `quotRem` sN
    mUpN' = mUpN * base
    mDnN' = mDnN * base
   in
   case (rn' < mDnN', rn' + mUpN' > sN) of
    (True,  False) -> dn : ds
    (False, True)  -> dn+1 : ds
    (True,  True)  -> if rn' * 2 < sN then dn : ds else dn+1 : ds
    (False, False) -> gen (dn:ds) rn' sN mUpN' mDnN'

  rds =
   if k >= 0 then
      gen [] r (s * expt base k) mUp mDn
   else
     let bk = expt base (k) in
     gen [] (r * bk) s (mUp * bk) (mDn * bk)
 in
 (map fromIntegral (reverse rds), k)

\end{code} %********************************************************* %* * \subsection{Converting from a Rational to a RealFloat %* * %********************************************************* [In response to a request for documentation of how fromRational works, Joe Fasel writes:] A quite reasonable request! This code was added to the Prelude just before the 1.2 release, when Lennart, working with an early version of hbi, noticed that (read . show) was not the identity for floating-point numbers. (There was a one-bit error about half the time.) The original version of the conversion function was in fact simply a floating-point divide, as you suggest above. The new version is, I grant you, somewhat denser. Unfortunately, Joe's code doesn't work! Here's an example: main = putStr (shows (1.82173691287639817263897126389712638972163e-300::Double) "\n") This program prints 0.0000000000000000 instead of 1.8217369128763981e-300 Here's Joe's code: \begin{pseudocode} fromRat :: (RealFloat a) => Rational -> a fromRat x = x' where x' = f e -- If the exponent of the nearest floating-point number to x -- is e, then the significand is the integer nearest xb^(-e), -- where b is the floating-point radix. We start with a good -- guess for e, and if it is correct, the exponent of the -- floating-point number we construct will again be e. If -- not, one more iteration is needed. f e = if e' == e then y else f e' where y = encodeFloat (round (x * (1 % b)^^e)) e (_,e') = decodeFloat y b = floatRadix x' -- We obtain a trial exponent by doing a floating-point -- division of x's numerator by its denominator. The -- result of this division may not itself be the ultimate -- result, because of an accumulation of three rounding -- errors. (s,e) = decodeFloat (fromInteger (numerator x) `asTypeOf` x' / fromInteger (denominator x)) \end{pseudocode} Now, here's Lennart's code (which works) \begin{code}


fromRat :: (RealFloat a) => Rational -> a


fromRat (n :% 0) | n > 0     =  1/0        
                 | n < 0     = 1/0        
                 | otherwise =  0/0        

fromRat (n :% d) | n > 0     = fromRat' (n :% d)
                 | n < 0     =  fromRat' ((n) :% d)
                 | otherwise = encodeFloat 0 0             









fromRat' :: (RealFloat a) => Rational -> a

fromRat' x = r
  where b = floatRadix r
        p = floatDigits r
        (minExp0, _) = floatRange r
        minExp = minExp0  p            
        xMax   = toRational (expt b p)
        p0 = (integerLogBase b (numerator x)  integerLogBase b (denominator x)  p) `max` minExp
        
        
        f = if p0 < 0 then 1 :% expt b (p0) else expt b p0 :% 1
        x0 = x / f
        
        
        (x', p') = if x0 >= xMax then (x0 / toRational b, p0+1) else (x0, p0)
        r = encodeFloat (round x') p'


minExpt, maxExpt :: Int
minExpt = 0
maxExpt = 1100

expt :: Integer -> Int -> Integer
expt base n =
    if base == 2 && n >= minExpt && n <= maxExpt then
        expts!n
    else
        if base == 10 && n <= maxExpt10 then
            expts10!n
        else
            base^n

expts :: Array Int Integer
expts = array (minExpt,maxExpt) [(n,2^n) | n <- [minExpt .. maxExpt]]

maxExpt10 :: Int
maxExpt10 = 324

expts10 :: Array Int Integer
expts10 = array (minExpt,maxExpt10) [(n,10^n) | n <- [minExpt .. maxExpt10]]







integerLogBase :: Integer -> Integer -> Int
integerLogBase b i
   | i < b     = 0
   | b == 2    = I# (integerLog2# i)
   | otherwise = I# (integerLogBase# b i)

\end{code} Unfortunately, the old conversion code was awfully slow due to a) a slow integer logarithm b) repeated calculation of gcd's For the case of Rational's coming from a Float or Double via toRational, we can exploit the fact that the denominator is a power of two, which for these brings a huge speedup since we need only shift and add instead of division. The below is an adaption of fromRat' for the conversion to Float or Double exploiting the known floatRadix and avoiding divisions as much as possible. \begin{code}

fromRat'' :: RealFloat a => Int -> Int -> Integer -> Integer -> a

fromRat'' minEx@(I# me#) mantDigs@(I# md#) n d =
    case integerLog2IsPowerOf2# d of
      (# ld#, pw# #)
        | pw# ==# 0# ->
          case integerLog2# n of
            ln# | ln# >=# (ld# +# me# -# 1#) ->
                  
                  
                  if ln# <# md#
                    then encodeFloat n (I# (negateInt# ld#))
                    else let n'  = n `shiftR` (I# (ln# +# 1# -# md#))
                             n'' = case roundingMode# n (ln# -# md#) of
                                    0# -> n'
                                    2# -> n' + 1
                                    _  -> case fromInteger n' .&. (1 :: Int) of
                                            0 -> n'
                                            _ -> n' + 1
                         in encodeFloat n'' (I# (ln# -# ld# +# 1# -# md#))
                | otherwise ->
                  
                  
                  
                  case ld# +# (me# -# md#) of
                    ld'# | ld'# <=# 0#  -> 
                           encodeFloat n (I# ((me# -# md#) -# ld'#))
                         | ld'# <=# ln#  ->
                           let n' = n `shiftR` (I# ld'#)
                           in case roundingMode# n (ld'# -# 1#) of
                                0# -> encodeFloat n' (minEx  mantDigs)
                                1# -> if fromInteger n' .&. (1 :: Int) == 0
                                        then encodeFloat n' (minExmantDigs)
                                        else encodeFloat (n' + 1) (minExmantDigs)
                                _  -> encodeFloat (n' + 1) (minExmantDigs)
                         | ld'# ># (ln# +# 1#)  -> encodeFloat 0 0 
                         | otherwise ->  
                           case integerLog2IsPowerOf2# n of
                            (# _, 0# #) -> encodeFloat 0 0  
                            (# _, _ #)  -> encodeFloat 1 (minEx  mantDigs)
        | otherwise ->
          let ln = I# (integerLog2# n)
              ld = I# ld#
              
              p0 = max minEx (ln  ld)
              (n', d')
                | p0 < mantDigs = (n `shiftL` (mantDigs  p0), d)
                | p0 == mantDigs = (n, d)
                | otherwise     = (n, d `shiftL` (p0  mantDigs))
              
              
              
              scale p a b
                | (b `shiftL` mantDigs) <= a = (p+1, a, b `shiftL` 1)
                | otherwise = (p, a, b)
              (p', n'', d'') = scale (p0mantDigs) n' d'
              
              rdq = case n'' `quotRem` d'' of
                     (q,r) -> case compare (r `shiftL` 1) d'' of
                                LT -> q
                                EQ -> if fromInteger q .&. (1 :: Int) == 0
                                        then q else q+1
                                GT -> q+1
          in  encodeFloat rdq p'
\end{code} %********************************************************* %* * \subsection{Floating point numeric primops} %* * %********************************************************* Definitions of the boxed PrimOps; these will be used in the case of partial applications, etc. \begin{code}
plusFloat, minusFloat, timesFloat, divideFloat :: Float -> Float -> Float
plusFloat   (F# x) (F# y) = F# (plusFloat# x y)
minusFloat  (F# x) (F# y) = F# (minusFloat# x y)
timesFloat  (F# x) (F# y) = F# (timesFloat# x y)
divideFloat (F# x) (F# y) = F# (divideFloat# x y)

negateFloat :: Float -> Float
negateFloat (F# x)        = F# (negateFloat# x)

gtFloat, geFloat, eqFloat, neFloat, ltFloat, leFloat :: Float -> Float -> Bool
gtFloat     (F# x) (F# y) = gtFloat# x y
geFloat     (F# x) (F# y) = geFloat# x y
eqFloat     (F# x) (F# y) = eqFloat# x y
neFloat     (F# x) (F# y) = neFloat# x y
ltFloat     (F# x) (F# y) = ltFloat# x y
leFloat     (F# x) (F# y) = leFloat# x y

expFloat, logFloat, sqrtFloat :: Float -> Float
sinFloat, cosFloat, tanFloat  :: Float -> Float
asinFloat, acosFloat, atanFloat  :: Float -> Float
sinhFloat, coshFloat, tanhFloat  :: Float -> Float
expFloat    (F# x) = F# (expFloat# x)
logFloat    (F# x) = F# (logFloat# x)
sqrtFloat   (F# x) = F# (sqrtFloat# x)
sinFloat    (F# x) = F# (sinFloat# x)
cosFloat    (F# x) = F# (cosFloat# x)
tanFloat    (F# x) = F# (tanFloat# x)
asinFloat   (F# x) = F# (asinFloat# x)
acosFloat   (F# x) = F# (acosFloat# x)
atanFloat   (F# x) = F# (atanFloat# x)
sinhFloat   (F# x) = F# (sinhFloat# x)
coshFloat   (F# x) = F# (coshFloat# x)
tanhFloat   (F# x) = F# (tanhFloat# x)

powerFloat :: Float -> Float -> Float
powerFloat  (F# x) (F# y) = F# (powerFloat# x y)




plusDouble, minusDouble, timesDouble, divideDouble :: Double -> Double -> Double
plusDouble   (D# x) (D# y) = D# (x +## y)
minusDouble  (D# x) (D# y) = D# (x -## y)
timesDouble  (D# x) (D# y) = D# (x *## y)
divideDouble (D# x) (D# y) = D# (x /## y)

negateDouble :: Double -> Double
negateDouble (D# x)        = D# (negateDouble# x)

gtDouble, geDouble, eqDouble, neDouble, leDouble, ltDouble :: Double -> Double -> Bool
gtDouble    (D# x) (D# y) = x >## y
geDouble    (D# x) (D# y) = x >=## y
eqDouble    (D# x) (D# y) = x ==## y
neDouble    (D# x) (D# y) = x /=## y
ltDouble    (D# x) (D# y) = x <## y
leDouble    (D# x) (D# y) = x <=## y

double2Float :: Double -> Float
double2Float (D# x) = F# (double2Float# x)

float2Double :: Float -> Double
float2Double (F# x) = D# (float2Double# x)

expDouble, logDouble, sqrtDouble :: Double -> Double
sinDouble, cosDouble, tanDouble  :: Double -> Double
asinDouble, acosDouble, atanDouble  :: Double -> Double
sinhDouble, coshDouble, tanhDouble  :: Double -> Double
expDouble    (D# x) = D# (expDouble# x)
logDouble    (D# x) = D# (logDouble# x)
sqrtDouble   (D# x) = D# (sqrtDouble# x)
sinDouble    (D# x) = D# (sinDouble# x)
cosDouble    (D# x) = D# (cosDouble# x)
tanDouble    (D# x) = D# (tanDouble# x)
asinDouble   (D# x) = D# (asinDouble# x)
acosDouble   (D# x) = D# (acosDouble# x)
atanDouble   (D# x) = D# (atanDouble# x)
sinhDouble   (D# x) = D# (sinhDouble# x)
coshDouble   (D# x) = D# (coshDouble# x)
tanhDouble   (D# x) = D# (tanhDouble# x)

powerDouble :: Double -> Double -> Double
powerDouble  (D# x) (D# y) = D# (x **## y)
\end{code} \begin{code}
foreign import ccall unsafe "isFloatNaN" isFloatNaN :: Float -> Int
foreign import ccall unsafe "isFloatInfinite" isFloatInfinite :: Float -> Int
foreign import ccall unsafe "isFloatDenormalized" isFloatDenormalized :: Float -> Int
foreign import ccall unsafe "isFloatNegativeZero" isFloatNegativeZero :: Float -> Int
foreign import ccall unsafe "isFloatFinite" isFloatFinite :: Float -> Int

foreign import ccall unsafe "isDoubleNaN" isDoubleNaN :: Double -> Int
foreign import ccall unsafe "isDoubleInfinite" isDoubleInfinite :: Double -> Int
foreign import ccall unsafe "isDoubleDenormalized" isDoubleDenormalized :: Double -> Int
foreign import ccall unsafe "isDoubleNegativeZero" isDoubleNegativeZero :: Double -> Int
foreign import ccall unsafe "isDoubleFinite" isDoubleFinite :: Double -> Int
\end{code} %********************************************************* %* * \subsection{Coercion rules} %* * %********************************************************* \begin{code}

\end{code} Note [realToFrac int-to-float] ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ Don found that the RULES for realToFrac/Int->Double and simliarly Float made a huge difference to some stream-fusion programs. Here's an example import Data.Array.Vector n = 40000000 main = do let c = replicateU n (2::Double) a = mapU realToFrac (enumFromToU 0 (n-1) ) :: UArr Double print (sumU (zipWithU (*) c a)) Without the RULE we get this loop body: case $wtoRational sc_sY4 of ww_aM7 { (# ww1_aM9, ww2_aMa #) -> case $wfromRat ww1_aM9 ww2_aMa of tpl_X1P { D# ipv_sW3 -> Main.$s$wfold (+# sc_sY4 1) (+# wild_X1i 1) (+## sc2_sY6 (*## 2.0 ipv_sW3)) And with the rule: Main.$s$wfold (+# sc_sXT 1) (+# wild_X1h 1) (+## sc2_sXV (*## 2.0 (int2Double# sc_sXT))) The running time of the program goes from 120 seconds to 0.198 seconds with the native backend, and 0.143 seconds with the C backend. A few more details in Trac #2251, and the patch message "Add RULES for realToFrac from Int". %********************************************************* %* * \subsection{Utils} %* * %********************************************************* \begin{code}
showSignedFloat :: (RealFloat a)
  => (a -> ShowS)       
  -> Int                
  -> a                  
  -> ShowS
showSignedFloat showPos p x
   | x < 0 || isNegativeZero x
       = showParen (p > 6) (showChar '-' . showPos (x))
   | otherwise = showPos x
\end{code} We need to prevent over/underflow of the exponent in encodeFloat when called from scaleFloat, hence we clamp the scaling parameter. We must have a large enough range to cover the maximum difference of exponents returned by decodeFloat. \begin{code}
clamp :: Int -> Int -> Int
clamp bd k = max (bd) (min bd k)
\end{code}

RetroSearch is an open source project built by @garambo | Open a GitHub Issue

Search and Browse the WWW like it's 1997 | Search results from DuckDuckGo

HTML: 3.2 | Encoding: UTF-8 | Version: 0.7.4