{-# LANGUAGE BangPatterns #-}
{-# LANGUAGE DeriveFunctor #-}
{-# LANGUAGE NoImplicitPrelude #-}
module Numeric.Floating.IEEE.Internal.Rounding.Common where
import           Control.Exception (assert)
import           Data.Bits
import           Data.Functor.Product
import           Data.Int
import           GHC.Float (expt)
import           Math.NumberTheory.Logarithms (integerLog2')
import           MyPrelude
import           Numeric.Floating.IEEE.Internal.IntegerInternals

default ()

class Functor f => RoundingStrategy f where
  exact :: a -> f a
  inexact :: Ordering -- ^ LT -> toward-zero is the nearest, EQ -> midpoint, GT -> away-from-zero is the nearest
          -> Bool -- ^ negative (True -> negative, False -> positive)
          -> Int -- ^ parity (even -> toward-zero is even, odd -> toward-zero is odd)
          -> a -- ^ toward zero
          -> a -- ^ away from zero
          -> f a
  doRound :: Bool -- ^ exactness; if True, the Ordering must be LT
          -> Ordering -- ^ LT -> toward-zero is the nearest, EQ -> midpoint, GT -> away-from-zero is the nearest
          -> Bool -- ^ negative (True -> negative, False -> positive)
          -> Int -- ^ parity (even -> toward-zero is even, odd -> toward-zero is odd)
          -> a -- ^ toward zero
          -> a -- ^ away from zero
          -> f a
  exact a
x = Bool -> Ordering -> Bool -> Int -> a -> a -> f a
forall (f :: * -> *) a.
RoundingStrategy f =>
Bool -> Ordering -> Bool -> Int -> a -> a -> f a
doRound Bool
True Ordering
LT Bool
False Int
0 a
x a
x
  inexact Ordering
o Bool
neg Int
parity a
zero a
away = Bool -> Ordering -> Bool -> Int -> a -> a -> f a
forall (f :: * -> *) a.
RoundingStrategy f =>
Bool -> Ordering -> Bool -> Int -> a -> a -> f a
doRound Bool
False Ordering
o Bool
neg Int
parity a
zero a
away

newtype RoundTiesToEven a = RoundTiesToEven { RoundTiesToEven a -> a
roundTiesToEven :: a }
  deriving (a -> RoundTiesToEven b -> RoundTiesToEven a
(a -> b) -> RoundTiesToEven a -> RoundTiesToEven b
(forall a b. (a -> b) -> RoundTiesToEven a -> RoundTiesToEven b)
-> (forall a b. a -> RoundTiesToEven b -> RoundTiesToEven a)
-> Functor RoundTiesToEven
forall a b. a -> RoundTiesToEven b -> RoundTiesToEven a
forall a b. (a -> b) -> RoundTiesToEven a -> RoundTiesToEven b
forall (f :: * -> *).
(forall a b. (a -> b) -> f a -> f b)
-> (forall a b. a -> f b -> f a) -> Functor f
<$ :: a -> RoundTiesToEven b -> RoundTiesToEven a
$c<$ :: forall a b. a -> RoundTiesToEven b -> RoundTiesToEven a
fmap :: (a -> b) -> RoundTiesToEven a -> RoundTiesToEven b
$cfmap :: forall a b. (a -> b) -> RoundTiesToEven a -> RoundTiesToEven b
Functor)

instance RoundingStrategy RoundTiesToEven where
  exact :: a -> RoundTiesToEven a
exact = a -> RoundTiesToEven a
forall a. a -> RoundTiesToEven a
RoundTiesToEven
  inexact :: Ordering -> Bool -> Int -> a -> a -> RoundTiesToEven a
inexact Ordering
o Bool
_neg Int
parity a
zero a
away = a -> RoundTiesToEven a
forall a. a -> RoundTiesToEven a
RoundTiesToEven (a -> RoundTiesToEven a) -> a -> RoundTiesToEven a
forall a b. (a -> b) -> a -> b
$ case Ordering
o of
                                                        Ordering
LT -> a
zero
                                                        Ordering
EQ | Int -> Bool
forall a. Integral a => a -> Bool
even Int
parity -> a
zero
                                                           | Bool
otherwise -> a
away
                                                        Ordering
GT -> a
away
  doRound :: Bool -> Ordering -> Bool -> Int -> a -> a -> RoundTiesToEven a
doRound Bool
_ex Ordering
o Bool
_neg Int
parity a
zero a
away = a -> RoundTiesToEven a
forall a. a -> RoundTiesToEven a
RoundTiesToEven (a -> RoundTiesToEven a) -> a -> RoundTiesToEven a
forall a b. (a -> b) -> a -> b
$ case Ordering
o of
    Ordering
LT -> a
zero
    Ordering
EQ | Int -> Bool
forall a. Integral a => a -> Bool
even Int
parity -> a
zero
       | Bool
otherwise -> a
away
    Ordering
GT -> a
away
  {-# INLINE exact #-}
  {-# INLINE inexact #-}
  {-# INLINE doRound #-}

newtype RoundTiesToAway a = RoundTiesToAway { RoundTiesToAway a -> a
roundTiesToAway :: a }
  deriving (a -> RoundTiesToAway b -> RoundTiesToAway a
(a -> b) -> RoundTiesToAway a -> RoundTiesToAway b
(forall a b. (a -> b) -> RoundTiesToAway a -> RoundTiesToAway b)
-> (forall a b. a -> RoundTiesToAway b -> RoundTiesToAway a)
-> Functor RoundTiesToAway
forall a b. a -> RoundTiesToAway b -> RoundTiesToAway a
forall a b. (a -> b) -> RoundTiesToAway a -> RoundTiesToAway b
forall (f :: * -> *).
(forall a b. (a -> b) -> f a -> f b)
-> (forall a b. a -> f b -> f a) -> Functor f
<$ :: a -> RoundTiesToAway b -> RoundTiesToAway a
$c<$ :: forall a b. a -> RoundTiesToAway b -> RoundTiesToAway a
fmap :: (a -> b) -> RoundTiesToAway a -> RoundTiesToAway b
$cfmap :: forall a b. (a -> b) -> RoundTiesToAway a -> RoundTiesToAway b
Functor)

instance RoundingStrategy RoundTiesToAway where
  exact :: a -> RoundTiesToAway a
exact = a -> RoundTiesToAway a
forall a. a -> RoundTiesToAway a
RoundTiesToAway
  inexact :: Ordering -> Bool -> Int -> a -> a -> RoundTiesToAway a
inexact Ordering
o Bool
_neg Int
_parity a
zero a
away = a -> RoundTiesToAway a
forall a. a -> RoundTiesToAway a
RoundTiesToAway (a -> RoundTiesToAway a) -> a -> RoundTiesToAway a
forall a b. (a -> b) -> a -> b
$ case Ordering
o of
                                                         Ordering
LT -> a
zero
                                                         Ordering
EQ -> a
away
                                                         Ordering
GT -> a
away
  doRound :: Bool -> Ordering -> Bool -> Int -> a -> a -> RoundTiesToAway a
doRound Bool
_ex Ordering
o Bool
_neg Int
_parity a
zero a
away = a -> RoundTiesToAway a
forall a. a -> RoundTiesToAway a
RoundTiesToAway (a -> RoundTiesToAway a) -> a -> RoundTiesToAway a
forall a b. (a -> b) -> a -> b
$ case Ordering
o of
    Ordering
LT -> a
zero
    Ordering
EQ -> a
away
    Ordering
GT -> a
away
  {-# INLINE exact #-}
  {-# INLINE inexact #-}
  {-# INLINE doRound #-}

newtype RoundTowardPositive a = RoundTowardPositive { RoundTowardPositive a -> a
roundTowardPositive :: a }
  deriving (a -> RoundTowardPositive b -> RoundTowardPositive a
(a -> b) -> RoundTowardPositive a -> RoundTowardPositive b
(forall a b.
 (a -> b) -> RoundTowardPositive a -> RoundTowardPositive b)
-> (forall a b.
    a -> RoundTowardPositive b -> RoundTowardPositive a)
-> Functor RoundTowardPositive
forall a b. a -> RoundTowardPositive b -> RoundTowardPositive a
forall a b.
(a -> b) -> RoundTowardPositive a -> RoundTowardPositive b
forall (f :: * -> *).
(forall a b. (a -> b) -> f a -> f b)
-> (forall a b. a -> f b -> f a) -> Functor f
<$ :: a -> RoundTowardPositive b -> RoundTowardPositive a
$c<$ :: forall a b. a -> RoundTowardPositive b -> RoundTowardPositive a
fmap :: (a -> b) -> RoundTowardPositive a -> RoundTowardPositive b
$cfmap :: forall a b.
(a -> b) -> RoundTowardPositive a -> RoundTowardPositive b
Functor)

instance RoundingStrategy RoundTowardPositive where
  exact :: a -> RoundTowardPositive a
exact = a -> RoundTowardPositive a
forall a. a -> RoundTowardPositive a
RoundTowardPositive
  inexact :: Ordering -> Bool -> Int -> a -> a -> RoundTowardPositive a
inexact Ordering
_o Bool
neg Int
_parity a
zero a
away | Bool
neg = a -> RoundTowardPositive a
forall a. a -> RoundTowardPositive a
RoundTowardPositive a
zero
                                   | Bool
otherwise = a -> RoundTowardPositive a
forall a. a -> RoundTowardPositive a
RoundTowardPositive a
away
  doRound :: Bool -> Ordering -> Bool -> Int -> a -> a -> RoundTowardPositive a
doRound Bool
ex Ordering
_o Bool
neg Int
_parity a
zero a
away | Bool
ex Bool -> Bool -> Bool
|| Bool
neg = a -> RoundTowardPositive a
forall a. a -> RoundTowardPositive a
RoundTowardPositive a
zero
                                      | Bool
otherwise = a -> RoundTowardPositive a
forall a. a -> RoundTowardPositive a
RoundTowardPositive a
away
  {-# INLINE exact #-}
  {-# INLINE inexact #-}
  {-# INLINE doRound #-}

newtype RoundTowardNegative a = RoundTowardNegative { RoundTowardNegative a -> a
roundTowardNegative :: a }
  deriving (a -> RoundTowardNegative b -> RoundTowardNegative a
(a -> b) -> RoundTowardNegative a -> RoundTowardNegative b
(forall a b.
 (a -> b) -> RoundTowardNegative a -> RoundTowardNegative b)
-> (forall a b.
    a -> RoundTowardNegative b -> RoundTowardNegative a)
-> Functor RoundTowardNegative
forall a b. a -> RoundTowardNegative b -> RoundTowardNegative a
forall a b.
(a -> b) -> RoundTowardNegative a -> RoundTowardNegative b
forall (f :: * -> *).
(forall a b. (a -> b) -> f a -> f b)
-> (forall a b. a -> f b -> f a) -> Functor f
<$ :: a -> RoundTowardNegative b -> RoundTowardNegative a
$c<$ :: forall a b. a -> RoundTowardNegative b -> RoundTowardNegative a
fmap :: (a -> b) -> RoundTowardNegative a -> RoundTowardNegative b
$cfmap :: forall a b.
(a -> b) -> RoundTowardNegative a -> RoundTowardNegative b
Functor)

instance RoundingStrategy RoundTowardNegative where
  exact :: a -> RoundTowardNegative a
exact = a -> RoundTowardNegative a
forall a. a -> RoundTowardNegative a
RoundTowardNegative
  inexact :: Ordering -> Bool -> Int -> a -> a -> RoundTowardNegative a
inexact Ordering
_o Bool
neg Int
_parity a
zero a
away | Bool
neg = a -> RoundTowardNegative a
forall a. a -> RoundTowardNegative a
RoundTowardNegative a
away
                                   | Bool
otherwise = a -> RoundTowardNegative a
forall a. a -> RoundTowardNegative a
RoundTowardNegative a
zero
  doRound :: Bool -> Ordering -> Bool -> Int -> a -> a -> RoundTowardNegative a
doRound Bool
ex Ordering
_o Bool
neg Int
_parity a
zero a
away | Bool -> Bool
not Bool
ex Bool -> Bool -> Bool
&& Bool
neg = a -> RoundTowardNegative a
forall a. a -> RoundTowardNegative a
RoundTowardNegative a
away
                                      | Bool
otherwise = a -> RoundTowardNegative a
forall a. a -> RoundTowardNegative a
RoundTowardNegative a
zero
  {-# INLINE exact #-}
  {-# INLINE inexact #-}
  {-# INLINE doRound #-}

newtype RoundTowardZero a = RoundTowardZero { RoundTowardZero a -> a
roundTowardZero :: a }
  deriving (a -> RoundTowardZero b -> RoundTowardZero a
(a -> b) -> RoundTowardZero a -> RoundTowardZero b
(forall a b. (a -> b) -> RoundTowardZero a -> RoundTowardZero b)
-> (forall a b. a -> RoundTowardZero b -> RoundTowardZero a)
-> Functor RoundTowardZero
forall a b. a -> RoundTowardZero b -> RoundTowardZero a
forall a b. (a -> b) -> RoundTowardZero a -> RoundTowardZero b
forall (f :: * -> *).
(forall a b. (a -> b) -> f a -> f b)
-> (forall a b. a -> f b -> f a) -> Functor f
<$ :: a -> RoundTowardZero b -> RoundTowardZero a
$c<$ :: forall a b. a -> RoundTowardZero b -> RoundTowardZero a
fmap :: (a -> b) -> RoundTowardZero a -> RoundTowardZero b
$cfmap :: forall a b. (a -> b) -> RoundTowardZero a -> RoundTowardZero b
Functor)

instance RoundingStrategy RoundTowardZero where
  exact :: a -> RoundTowardZero a
exact = a -> RoundTowardZero a
forall a. a -> RoundTowardZero a
RoundTowardZero
  inexact :: Ordering -> Bool -> Int -> a -> a -> RoundTowardZero a
inexact Ordering
_o Bool
_neg Int
_parity a
zero a
_away = a -> RoundTowardZero a
forall a. a -> RoundTowardZero a
RoundTowardZero a
zero
  doRound :: Bool -> Ordering -> Bool -> Int -> a -> a -> RoundTowardZero a
doRound Bool
_ex Ordering
_o Bool
_neg Int
_parity a
zero a
_away = a -> RoundTowardZero a
forall a. a -> RoundTowardZero a
RoundTowardZero a
zero
  {-# INLINE exact #-}
  {-# INLINE inexact #-}
  {-# INLINE doRound #-}

instance (RoundingStrategy f, RoundingStrategy g) => RoundingStrategy (Product f g) where
  exact :: a -> Product f g a
exact a
x = f a -> g a -> Product f g a
forall k (f :: k -> *) (g :: k -> *) (a :: k).
f a -> g a -> Product f g a
Pair (a -> f a
forall (f :: * -> *) a. RoundingStrategy f => a -> f a
exact a
x) (a -> g a
forall (f :: * -> *) a. RoundingStrategy f => a -> f a
exact a
x)
  inexact :: Ordering -> Bool -> Int -> a -> a -> Product f g a
inexact Ordering
o Bool
neg Int
parity a
zero a
away = f a -> g a -> Product f g a
forall k (f :: k -> *) (g :: k -> *) (a :: k).
f a -> g a -> Product f g a
Pair (Ordering -> Bool -> Int -> a -> a -> f a
forall (f :: * -> *) a.
RoundingStrategy f =>
Ordering -> Bool -> Int -> a -> a -> f a
inexact Ordering
o Bool
neg Int
parity a
zero a
away) (Ordering -> Bool -> Int -> a -> a -> g a
forall (f :: * -> *) a.
RoundingStrategy f =>
Ordering -> Bool -> Int -> a -> a -> f a
inexact Ordering
o Bool
neg Int
parity a
zero a
away)
  doRound :: Bool -> Ordering -> Bool -> Int -> a -> a -> Product f g a
doRound Bool
ex Ordering
o Bool
neg Int
parity a
zero a
away = f a -> g a -> Product f g a
forall k (f :: k -> *) (g :: k -> *) (a :: k).
f a -> g a -> Product f g a
Pair (Bool -> Ordering -> Bool -> Int -> a -> a -> f a
forall (f :: * -> *) a.
RoundingStrategy f =>
Bool -> Ordering -> Bool -> Int -> a -> a -> f a
doRound Bool
ex Ordering
o Bool
neg Int
parity a
zero a
away) (Bool -> Ordering -> Bool -> Int -> a -> a -> g a
forall (f :: * -> *) a.
RoundingStrategy f =>
Bool -> Ordering -> Bool -> Int -> a -> a -> f a
doRound Bool
ex Ordering
o Bool
neg Int
parity a
zero a
away)
  {-# INLINE exact #-}
  {-# INLINE inexact #-}
  {-# INLINE doRound #-}

{-
from GHC.Float:
expt :: Integer -> Int -> Integer
expt base n = base ^ n
-}

quotRemByExpt :: Integer -- ^ the dividend @x@
              -> Integer -- ^ base
              -> Int -- ^ the exponent @e@ (must be non-negative)
              -> (Integer, Integer) -- ^ @x \`'quotRem'\` (base ^ e)@
quotRemByExpt :: Integer -> Integer -> Int -> (Integer, Integer)
quotRemByExpt Integer
x Integer
2 Int
n    = Bool -> (Integer, Integer) -> (Integer, Integer)
forall a. (?callStack::CallStack) => Bool -> a -> a
assert (Int
n Int -> Int -> Bool
forall a. Ord a => a -> a -> Bool
>= Int
0) (Integer
x Integer -> Int -> Integer
`unsafeShiftRInteger` Int
n, Integer
x Integer -> Integer -> Integer
forall a. Bits a => a -> a -> a
.&. (Int -> Integer
forall a. Bits a => Int -> a
bit Int
n Integer -> Integer -> Integer
forall a. Num a => a -> a -> a
- Integer
1))
quotRemByExpt Integer
x Integer
base Int
n = Integer
x Integer -> Integer -> (Integer, Integer)
forall a. Integral a => a -> a -> (a, a)
`quotRem` Integer -> Int -> Integer
expt Integer
base Int
n
{-# INLINE quotRemByExpt #-}

multiplyByExpt :: Integer -- ^ the multiplicand @x@
               -> Integer -- ^ base
               -> Int -- ^ the exponent @e@ (must be non-negative)
               -> Integer -- ^ @x * base ^ e@
multiplyByExpt :: Integer -> Integer -> Int -> Integer
multiplyByExpt Integer
x Integer
2 Int
n    = Bool -> Integer -> Integer
forall a. (?callStack::CallStack) => Bool -> a -> a
assert (Int
n Int -> Int -> Bool
forall a. Ord a => a -> a -> Bool
>= Int
0) (Integer
x Integer -> Int -> Integer
`unsafeShiftLInteger` Int
n)
multiplyByExpt Integer
x Integer
base Int
n = Integer
x Integer -> Integer -> Integer
forall a. Num a => a -> a -> a
* Integer -> Int -> Integer
expt Integer
base Int
n
{-# INLINE multiplyByExpt #-}

isDivisibleByExpt :: Integer -- ^ the dividend @x@
                  -> Integer -- ^ the base
                  -> Int -- ^ the exponent @e@ (must be non-negative)
                  -> Integer -- ^ the remainder @r@ (must be @x \`'rem'\` (base ^ e)@)
                  -> Bool -- ^ @r == 0@
isDivisibleByExpt :: Integer -> Integer -> Int -> Integer -> Bool
isDivisibleByExpt Integer
x Integer
2 Int
e Integer
r = Bool -> Bool -> Bool
forall a. (?callStack::CallStack) => Bool -> a -> a
assert (Integer
r Integer -> Integer -> Bool
forall a. Eq a => a -> a -> Bool
== Integer
x Integer -> Integer -> Integer
forall a. Integral a => a -> a -> a
`rem` (Integer
2 Integer -> Int -> Integer
forall a b. (Num a, Integral b) => a -> b -> a
^ Int
e)) (Bool -> Bool) -> Bool -> Bool
forall a b. (a -> b) -> a -> b
$ Integer
x Integer -> Integer -> Bool
forall a. Eq a => a -> a -> Bool
== Integer
0 Bool -> Bool -> Bool
|| Integer -> Int
Numeric.Floating.IEEE.Internal.IntegerInternals.countTrailingZerosInteger Integer
x Int -> Int -> Bool
forall a. Ord a => a -> a -> Bool
>= Int
e
isDivisibleByExpt Integer
x Integer
base Int
e Integer
r = Bool -> Bool -> Bool
forall a. (?callStack::CallStack) => Bool -> a -> a
assert (Integer
r Integer -> Integer -> Bool
forall a. Eq a => a -> a -> Bool
== Integer
x Integer -> Integer -> Integer
forall a. Integral a => a -> a -> a
`rem` (Integer
base Integer -> Int -> Integer
forall a b. (Num a, Integral b) => a -> b -> a
^ Int
e)) (Integer
r Integer -> Integer -> Bool
forall a. Eq a => a -> a -> Bool
== Integer
0)
{-# INLINE isDivisibleByExpt #-}

-- |
-- Assumption: @n >= 0@, @e >= 0@, and @r == n \`'rem'\` base^(e+1)@
--
-- Returns @compare r (base^e)@.
compareWithExpt :: Integer -- ^ base
                -> Integer -- ^ the number @n@ (must be non-negative)
                -> Integer -- ^ the remainder @r@ (must be @n \`'rem'\' base^(e+1)@)
                -> Int -- ^ the exponent @e@ (must be non-negative)
                -> Ordering
compareWithExpt :: Integer -> Integer -> Integer -> Int -> Ordering
compareWithExpt Integer
2 Integer
n Integer
r Int
e = Bool -> Ordering -> Ordering
forall a. (?callStack::CallStack) => Bool -> a -> a
assert (Integer
r Integer -> Integer -> Bool
forall a. Eq a => a -> a -> Bool
== Integer
n Integer -> Integer -> Integer
forall a. Integral a => a -> a -> a
`rem` Integer -> Int -> Integer
expt Integer
2 (Int
eInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
1)) (Ordering -> Ordering) -> Ordering -> Ordering
forall a b. (a -> b) -> a -> b
$
  if Integer
n Integer -> Integer -> Bool
forall a. Eq a => a -> a -> Bool
== Integer
0 Bool -> Bool -> Bool
|| Integer -> Int
integerLog2' Integer
n Int -> Int -> Bool
forall a. Ord a => a -> a -> Bool
< Int
e then
    -- If integerLog2 n < e (i.e. n < 2^e), it is trivial
    Ordering
LT
  else
    -- In this branch, n > 0 && integerLog2' n >= e
    let result :: Ordering
result = Integer -> Int -> Ordering
Numeric.Floating.IEEE.Internal.IntegerInternals.roundingMode Integer
n Int
e
        !()
_ = Bool -> () -> ()
forall a. (?callStack::CallStack) => Bool -> a -> a
assert (Ordering
result Ordering -> Ordering -> Bool
forall a. Eq a => a -> a -> Bool
== Integer -> Integer -> Ordering
forall a. Ord a => a -> a -> Ordering
compare Integer
r (Integer -> Int -> Integer
expt Integer
2 Int
e)) ()
    in Ordering
result
compareWithExpt Integer
base Integer
n Integer
r Int
e = Bool -> Ordering -> Ordering
forall a. (?callStack::CallStack) => Bool -> a -> a
assert (Integer
r Integer -> Integer -> Bool
forall a. Eq a => a -> a -> Bool
== Integer
n Integer -> Integer -> Integer
forall a. Integral a => a -> a -> a
`rem` Integer -> Int -> Integer
expt Integer
base (Int
eInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
1)) (Ordering -> Ordering) -> Ordering -> Ordering
forall a b. (a -> b) -> a -> b
$ Integer -> Integer -> Ordering
forall a. Ord a => a -> a -> Ordering
compare Integer
r (Integer -> Int -> Integer
expt Integer
base Int
e)
{-# INLINE compareWithExpt #-}