{-# LANGUAGE BangPatterns #-}
{-# LANGUAGE NoImplicitPrelude #-}
{-# LANGUAGE ScopedTypeVariables #-}
module Numeric.Floating.IEEE.Internal.Rounding.Rational where
import           Control.Exception (assert)
import           Data.Functor.Product
import           Data.Ratio
import           GHC.Float (expt)
import           Math.NumberTheory.Logarithms (integerLog2', integerLogBase')
import           MyPrelude
import           Numeric.Floating.IEEE.Internal.Base
import           Numeric.Floating.IEEE.Internal.Rounding.Common

default ()

-- |
-- Conversion from a rational number to floating-point value, with each rounding attributes.
fromRationalTiesToEven, fromRationalTiesToAway, fromRationalTowardPositive, fromRationalTowardNegative, fromRationalTowardZero :: RealFloat a => Rational -> a
fromRationalTiesToEven :: Rational -> a
fromRationalTiesToEven = RoundTiesToEven a -> a
forall a. RoundTiesToEven a -> a
roundTiesToEven (RoundTiesToEven a -> a)
-> (Rational -> RoundTiesToEven a) -> Rational -> a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Rational -> RoundTiesToEven a
forall a (f :: * -> *).
(RealFloat a, RoundingStrategy f) =>
Rational -> f a
fromRationalR
fromRationalTiesToAway :: Rational -> a
fromRationalTiesToAway = RoundTiesToAway a -> a
forall a. RoundTiesToAway a -> a
roundTiesToAway (RoundTiesToAway a -> a)
-> (Rational -> RoundTiesToAway a) -> Rational -> a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Rational -> RoundTiesToAway a
forall a (f :: * -> *).
(RealFloat a, RoundingStrategy f) =>
Rational -> f a
fromRationalR
fromRationalTowardPositive :: Rational -> a
fromRationalTowardPositive = RoundTowardPositive a -> a
forall a. RoundTowardPositive a -> a
roundTowardPositive (RoundTowardPositive a -> a)
-> (Rational -> RoundTowardPositive a) -> Rational -> a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Rational -> RoundTowardPositive a
forall a (f :: * -> *).
(RealFloat a, RoundingStrategy f) =>
Rational -> f a
fromRationalR
fromRationalTowardNegative :: Rational -> a
fromRationalTowardNegative = RoundTowardNegative a -> a
forall a. RoundTowardNegative a -> a
roundTowardNegative (RoundTowardNegative a -> a)
-> (Rational -> RoundTowardNegative a) -> Rational -> a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Rational -> RoundTowardNegative a
forall a (f :: * -> *).
(RealFloat a, RoundingStrategy f) =>
Rational -> f a
fromRationalR
fromRationalTowardZero :: Rational -> a
fromRationalTowardZero = RoundTowardZero a -> a
forall a. RoundTowardZero a -> a
roundTowardZero (RoundTowardZero a -> a)
-> (Rational -> RoundTowardZero a) -> Rational -> a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Rational -> RoundTowardZero a
forall a (f :: * -> *).
(RealFloat a, RoundingStrategy f) =>
Rational -> f a
fromRationalR
{-# INLINE fromRationalTiesToEven #-}
{-# INLINE fromRationalTiesToAway #-}
{-# INLINE fromRationalTowardPositive #-}
{-# INLINE fromRationalTowardNegative #-}
{-# INLINE fromRationalTowardZero #-}

fromRationalR :: (RealFloat a, RoundingStrategy f) => Rational -> f a
fromRationalR :: Rational -> f a
fromRationalR Rational
x = Integer -> Integer -> f a
forall a (f :: * -> *).
(RealFloat a, RoundingStrategy f) =>
Integer -> Integer -> f a
fromRatioR (Rational -> Integer
forall a. Ratio a -> a
numerator Rational
x) (Rational -> Integer
forall a. Ratio a -> a
denominator Rational
x)
{-# INLINE fromRationalR #-}

fromRatioR :: (RealFloat a, RoundingStrategy f)
           => Integer -- ^ numerator
           -> Integer -- ^ denominator
           -> f a
fromRatioR :: Integer -> Integer -> f a
fromRatioR Integer
0 !Integer
_ = a -> f a
forall (f :: * -> *) a. RoundingStrategy f => a -> f a
exact a
0
fromRatioR Integer
n Integer
0 | Integer
n Integer -> Integer -> Bool
forall a. Ord a => a -> a -> Bool
> Integer
0 = a -> f a
forall (f :: * -> *) a. RoundingStrategy f => a -> f a
exact (a
1 a -> a -> a
forall a. Fractional a => a -> a -> a
/ a
0) -- positive infinity
               | Bool
otherwise = a -> f a
forall (f :: * -> *) a. RoundingStrategy f => a -> f a
exact (- a
1 a -> a -> a
forall a. Fractional a => a -> a -> a
/ a
0) -- negative infinity
fromRatioR Integer
n Integer
d | Integer
d Integer -> Integer -> Bool
forall a. Ord a => a -> a -> Bool
< Integer
0 = [Char] -> f a
forall a. HasCallStack => [Char] -> a
error [Char]
"fromRatio: negative denominator"
               | Integer
n Integer -> Integer -> Bool
forall a. Ord a => a -> a -> Bool
< Integer
0 = a -> a
forall a. Num a => a -> a
negate (a -> a) -> f a -> f a
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
<$> Bool -> Integer -> Integer -> f a
forall (f :: * -> *) a.
(RealFloat a, RoundingStrategy f) =>
Bool -> Integer -> Integer -> f a
fromPositiveRatioR Bool
True (- Integer
n) Integer
d
               | Bool
otherwise = Bool -> Integer -> Integer -> f a
forall (f :: * -> *) a.
(RealFloat a, RoundingStrategy f) =>
Bool -> Integer -> Integer -> f a
fromPositiveRatioR Bool
False Integer
n Integer
d
{-# INLINE fromRatioR #-}

fromPositiveRatioR :: forall f a. (RealFloat a, RoundingStrategy f)
                   => Bool -- ^ True if the result will be negated
                   -> Integer -- ^ numerator (> 0)
                   -> Integer -- ^ denominator (> 0)
                   -> f a
fromPositiveRatioR :: Bool -> Integer -> Integer -> f a
fromPositiveRatioR !Bool
neg !Integer
n !Integer
d = Bool -> f a -> f a
forall a. HasCallStack => Bool -> a -> a
assert (Integer
n Integer -> Integer -> Bool
forall a. Ord a => a -> a -> Bool
> Integer
0 Bool -> Bool -> Bool
&& Integer
d Integer -> Integer -> Bool
forall a. Ord a => a -> a -> Bool
> Integer
0) f a
result
  where
    result :: f a
result = let e0 :: Int
                 e0 :: Int
e0 = if Integer
base Integer -> Integer -> Bool
forall a. Eq a => a -> a -> Bool
== Integer
2 then
                        Integer -> Int
integerLog2' Integer
n Int -> Int -> Int
forall a. Num a => a -> a -> a
- Integer -> Int
integerLog2' Integer
d Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int
fDigits
                      else
                        Integer -> Integer -> Int
integerLogBase' Integer
base Integer
n Int -> Int -> Int
forall a. Num a => a -> a -> a
- Integer -> Integer -> Int
integerLogBase' Integer
base Integer
d Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int
fDigits
                 q0, r0, d0 :: Integer
                 (!Integer
d0, (!Integer
q0, !Integer
r0)) =
                   if Int
e0 Int -> Int -> Bool
forall a. Ord a => a -> a -> Bool
>= Int
0 then
                     -- n = q0 * (d * base^e0) + r0, 0 <= r0 < d * base^e0
                     let d_ :: Integer
d_ = Integer -> Integer -> Int -> Integer
multiplyByExpt Integer
d Integer
base Int
e0
                     in (Integer
d_, Integer
n Integer -> Integer -> (Integer, Integer)
forall a. Integral a => a -> a -> (a, a)
`quotRem` Integer
d_)
                   else
                     -- n * base^(-e0) = q0 * d + r0, 0 <= r0 < d
                     (Integer
d, (Integer -> Integer -> Int -> Integer
multiplyByExpt Integer
n Integer
base (-Int
e0)) Integer -> Integer -> (Integer, Integer)
forall a. Integral a => a -> a -> (a, a)
`quotRem` Integer
d)
                 -- Invariant: n / d * base^^(-e0) = q0 + r0 / d0
                 !()
_ = Bool -> () -> ()
forall a. HasCallStack => Bool -> a -> a
assert (Integer
n Integer -> Integer -> Rational
forall a. Integral a => a -> a -> Ratio a
% Integer
d Rational -> Rational -> Rational
forall a. Num a => a -> a -> a
* Integer -> Rational
forall a. Num a => Integer -> a
fromInteger Integer
baseRational -> Int -> Rational
forall a b. (Fractional a, Integral b) => a -> b -> a
^^(-Int
e0) Rational -> Rational -> Bool
forall a. Eq a => a -> a -> Bool
== Integer -> Rational
forall a. Num a => Integer -> a
fromInteger Integer
q0 Rational -> Rational -> Rational
forall a. Num a => a -> a -> a
+ Integer
r0 Integer -> Integer -> Rational
forall a. Integral a => a -> a -> Ratio a
% Integer
d0) ()
                 !()
_ = Bool -> () -> ()
forall a. HasCallStack => Bool -> a -> a
assert (Integer
baseInteger -> Int -> Integer
forall a b. (Num a, Integral b) => a -> b -> a
^(Int
fDigitsInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
1) Integer -> Integer -> Bool
forall a. Ord a => a -> a -> Bool
<= Integer
q0 Bool -> Bool -> Bool
&& Integer
q0 Integer -> Integer -> Bool
forall a. Ord a => a -> a -> Bool
< Integer
baseInteger -> Int -> Integer
forall a b. (Num a, Integral b) => a -> b -> a
^(Int
fDigitsInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
1)) ()

                 q, r, d' :: Integer
                 e :: Int
                 (!Integer
q, !Integer
r, !Integer
d', !Int
e) =
                   if Integer
q0 Integer -> Integer -> Bool
forall a. Ord a => a -> a -> Bool
< Integer -> Int -> Integer
expt Integer
base Int
fDigits then
                     -- base^(fDigits-1) <= q0 < base^fDigits
                     (Integer
q0, Integer
r0, Integer
d0, Int
e0)
                   else
                     -- base^fDigits <= q0 < base^(fDigits+1)
                     let (Integer
q', Integer
r') = Integer
q0 Integer -> Integer -> (Integer, Integer)
forall a. Integral a => a -> a -> (a, a)
`quotRem` Integer
base
                     in (Integer
q', Integer
r' Integer -> Integer -> Integer
forall a. Num a => a -> a -> a
* Integer
d0 Integer -> Integer -> Integer
forall a. Num a => a -> a -> a
+ Integer
r0, Integer
base Integer -> Integer -> Integer
forall a. Num a => a -> a -> a
* Integer
d0, Int
e0 Int -> Int -> Int
forall a. Num a => a -> a -> a
+ Int
1)
                 -- Invariant: n / d * 2^^(-e) = q + r / d', base^(fDigits-1) <= q < base^fDigits, 0 <= r < d'
                 !()
_ = Bool -> () -> ()
forall a. HasCallStack => Bool -> a -> a
assert (Integer
n Integer -> Integer -> Rational
forall a. Integral a => a -> a -> Ratio a
% Integer
d Rational -> Rational -> Rational
forall a. Num a => a -> a -> a
* Integer -> Rational
forall a. Num a => Integer -> a
fromInteger Integer
baseRational -> Int -> Rational
forall a b. (Fractional a, Integral b) => a -> b -> a
^^(-Int
e) Rational -> Rational -> Bool
forall a. Eq a => a -> a -> Bool
== Integer -> Rational
forall a. Num a => Integer -> a
fromInteger Integer
q Rational -> Rational -> Rational
forall a. Num a => a -> a -> a
+ Integer
r Integer -> Integer -> Rational
forall a. Integral a => a -> a -> Ratio a
% Integer
d') ()
                 -- base^(e+fDigits-1) <= q * base^^e <= n/d < (q+1) * base^^e <= base^(e+fDigits)
                 -- In particular, base^(fDigits-1) <= q < base^fDigits
             in if Int
expMin Int -> Int -> Bool
forall a. Ord a => a -> a -> Bool
<= Int
e Int -> Int -> Int
forall a. Num a => a -> a -> a
+ Int
fDigits Bool -> Bool -> Bool
&& Int
e Int -> Int -> Int
forall a. Num a => a -> a -> a
+ Int
fDigits Int -> Int -> Bool
forall a. Ord a => a -> a -> Bool
<= Int
expMax then
                  -- normal: base^^(expMin-1) <= n/d < base^expMax
                  let towardzero_or_exact :: a
towardzero_or_exact = Integer -> Int -> a
forall a. RealFloat a => Integer -> Int -> a
encodeFloat Integer
q Int
e
                      awayfromzero :: a
awayfromzero = Integer -> Int -> a
forall a. RealFloat a => Integer -> Int -> a
encodeFloat (Integer
q Integer -> Integer -> Integer
forall a. Num a => a -> a -> a
+ Integer
1) Int
e -- may be infinity
                      parity :: Int
parity = Integer -> Int
forall a. Num a => Integer -> a
fromInteger Integer
q :: Int
                  in Bool -> Ordering -> Bool -> Int -> a -> a -> f a
forall (f :: * -> *) a.
RoundingStrategy f =>
Bool -> Ordering -> Bool -> Int -> a -> a -> f a
doRound
                       (Integer
r Integer -> Integer -> Bool
forall a. Eq a => a -> a -> Bool
== Integer
0)
                       (Integer -> Integer -> Ordering
forall a. Ord a => a -> a -> Ordering
compare (Integer
base Integer -> Integer -> Integer
forall a. Num a => a -> a -> a
* Integer
r) Integer
d')
                       Bool
neg
                       Int
parity
                       a
towardzero_or_exact
                       a
awayfromzero
                else
                  if Int
expMax Int -> Int -> Bool
forall a. Ord a => a -> a -> Bool
< Int
e Int -> Int -> Int
forall a. Num a => a -> a -> a
+ Int
fDigits then
                    -- overflow
                    let inf :: a
inf = a
1 a -> a -> a
forall a. Fractional a => a -> a -> a
/ a
0
                    in Ordering -> Bool -> Int -> a -> a -> f a
forall (f :: * -> *) a.
RoundingStrategy f =>
Ordering -> Bool -> Int -> a -> a -> f a
inexact Ordering
GT Bool
neg Int
1 a
forall a. RealFloat a => a
maxFinite a
inf
                  else
                    -- subnormal: 0 < n/d < base^^(expMin-1)
                    -- e + fDigits < expMin
                    let (Integer
q', Integer
r') = Integer -> Integer -> Int -> (Integer, Integer)
quotRemByExpt Integer
q Integer
base (Int
expMin Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int
fDigits Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int
e)
                        !()
_ = Bool -> () -> ()
forall a. HasCallStack => Bool -> a -> a
assert (Integer
q Integer -> Integer -> Bool
forall a. Eq a => a -> a -> Bool
== Integer
q' Integer -> Integer -> Integer
forall a. Num a => a -> a -> a
* Integer
baseInteger -> Int -> Integer
forall a b. (Num a, Integral b) => a -> b -> a
^(Int
expMinInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
fDigitsInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
e) Integer -> Integer -> Integer
forall a. Num a => a -> a -> a
+ Integer
r' Bool -> Bool -> Bool
&& Integer
0 Integer -> Integer -> Bool
forall a. Ord a => a -> a -> Bool
<= Integer
r' Bool -> Bool -> Bool
&& Integer
r' Integer -> Integer -> Bool
forall a. Ord a => a -> a -> Bool
< Integer
baseInteger -> Int -> Integer
forall a b. (Num a, Integral b) => a -> b -> a
^(Int
expMinInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
fDigitsInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
e)) ()
                        -- q = q' * base^(expMin-fDigits-e) + r', 0 <= r' < base^(expMin-fDigits-e)
                        -- n / d * base^^(-e) = q' * base^(expMin-fDigits-e) + r' + r / d'
                        -- n / d = q' * base^^(expMin - fDigits) + (r' + r / d') * base^^e
                        !()
_ = Bool -> () -> ()
forall a. HasCallStack => Bool -> a -> a
assert (Integer
n Integer -> Integer -> Rational
forall a. Integral a => a -> a -> Ratio a
% Integer
d Rational -> Rational -> Bool
forall a. Eq a => a -> a -> Bool
== Integer -> Rational
forall a. Num a => Integer -> a
fromInteger Integer
q' Rational -> Rational -> Rational
forall a. Num a => a -> a -> a
* Integer -> Rational
forall a. Num a => Integer -> a
fromInteger Integer
baseRational -> Int -> Rational
forall a b. (Fractional a, Integral b) => a -> b -> a
^^(Int
expMin Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int
fDigits) Rational -> Rational -> Rational
forall a. Num a => a -> a -> a
+ (Integer -> Rational
forall a. Num a => Integer -> a
fromInteger Integer
r' Rational -> Rational -> Rational
forall a. Num a => a -> a -> a
+ Integer
r Integer -> Integer -> Rational
forall a. Integral a => a -> a -> Ratio a
% Integer
d') Rational -> Rational -> Rational
forall a. Num a => a -> a -> a
* Integer -> Rational
forall a. Num a => Integer -> a
fromInteger Integer
baseRational -> Int -> Rational
forall a b. (Fractional a, Integral b) => a -> b -> a
^^Int
e) ()
                        -- rounding direction: (r' + r / d') * base^^e vs. base^^(expMin-fDigits-1)
                        towardzero :: a
towardzero = Integer -> Int -> a
forall a. RealFloat a => Integer -> Int -> a
encodeFloat Integer
q' (Int
expMin Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int
fDigits)
                        awayfromzero :: a
awayfromzero = Integer -> Int -> a
forall a. RealFloat a => Integer -> Int -> a
encodeFloat (Integer
q' Integer -> Integer -> Integer
forall a. Num a => a -> a -> a
+ Integer
1) (Int
expMin Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int
fDigits)
                        parity :: Int
parity = Integer -> Int
forall a. Num a => Integer -> a
fromInteger Integer
q' :: Int
                    in Bool -> Ordering -> Bool -> Int -> a -> a -> f a
forall (f :: * -> *) a.
RoundingStrategy f =>
Bool -> Ordering -> Bool -> Int -> a -> a -> f a
doRound
                         (Integer
r Integer -> Integer -> Bool
forall a. Eq a => a -> a -> Bool
== Integer
0 Bool -> Bool -> Bool
&& Integer
r' Integer -> Integer -> Bool
forall a. Eq a => a -> a -> Bool
== Integer
0)
                         (Integer -> Integer -> Integer -> Int -> Ordering
compareWithExpt Integer
base Integer
q Integer
r' (Int
expMin Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int
fDigits Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int
e Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int
1) Ordering -> Ordering -> Ordering
forall a. Semigroup a => a -> a -> a
<> if Integer
r Integer -> Integer -> Bool
forall a. Eq a => a -> a -> Bool
== Integer
0 then Ordering
EQ else Ordering
GT)
                         -- (compare r' (expt base (expMin - fDigits - e - 1)) <> if r == 0 then EQ else GT)
                         Bool
neg
                         Int
parity
                         a
towardzero
                         a
awayfromzero

    !base :: Integer
base = a -> Integer
forall a. RealFloat a => a -> Integer
floatRadix (a
forall a. HasCallStack => a
undefined :: a)
    !fDigits :: Int
fDigits = a -> Int
forall a. RealFloat a => a -> Int
floatDigits (a
forall a. HasCallStack => a
undefined :: a) -- 53 for Double
    (!Int
expMin, !Int
expMax) = a -> (Int, Int)
forall a. RealFloat a => a -> (Int, Int)
floatRange (a
forall a. HasCallStack => a
undefined :: a) -- (-1021, 1024) for Double
{-# INLINABLE [0] fromPositiveRatioR #-}
{-# SPECIALIZE
  fromPositiveRatioR :: RealFloat a => Bool -> Integer -> Integer -> RoundTiesToEven a
                      , RealFloat a => Bool -> Integer -> Integer -> RoundTiesToAway a
                      , RealFloat a => Bool -> Integer -> Integer -> RoundTowardPositive a
                      , RealFloat a => Bool -> Integer -> Integer -> RoundTowardZero a
                      , RealFloat a => Bool -> Integer -> Integer -> Product RoundTowardNegative RoundTowardPositive a
                      , RoundingStrategy f => Bool -> Integer -> Integer -> f Double
                      , RoundingStrategy f => Bool -> Integer -> Integer -> f Float
                      , Bool -> Integer -> Integer -> RoundTiesToEven Double
                      , Bool -> Integer -> Integer -> RoundTiesToAway Double
                      , Bool -> Integer -> Integer -> RoundTowardPositive Double
                      , Bool -> Integer -> Integer -> RoundTowardZero Double
                      , Bool -> Integer -> Integer -> RoundTiesToEven Float
                      , Bool -> Integer -> Integer -> RoundTiesToAway Float
                      , Bool -> Integer -> Integer -> RoundTowardPositive Float
                      , Bool -> Integer -> Integer -> RoundTowardZero Float
                      , Bool -> Integer -> Integer -> Product RoundTowardNegative RoundTowardPositive Double
                      , Bool -> Integer -> Integer -> Product RoundTowardNegative RoundTowardPositive Float
  #-}
{-# RULES
"fromPositiveRatioR/RoundTowardNegative"
  fromPositiveRatioR = \neg x y -> RoundTowardNegative (roundTowardPositive (fromPositiveRatioR (not neg) x y))
  #-}