-- | A collection of generic numerical and list manipulation functions.
module Goal.Core.Util
    ( -- * List Manipulation
      takeEvery
    , breakEvery
    , kFold
    , kFold'
    -- * Numeric
    , roundSD
    , toPi
    , circularDistance
    , integrate
    , logistic
    , logit
    , square
    , triangularNumber
    -- ** List Numerics
    , average
    , weightedAverage
    , circularAverage
    , weightedCircularAverage
    , range
    , discretizeFunction
    , logSumExp
    , logIntegralExp
    -- * Tracing
    , traceGiven
    -- * TypeNats
    , finiteInt
    , natValInt
    , Triangular
    -- ** Type Rationals
    , Rat
    , type (/)
    , ratVal
    ) where


--- Imports ---


-- Unqualified --

import Numeric
import Data.Ratio
import Data.Proxy
import Debug.Trace
import Data.Finite
import GHC.TypeNats

-- Qualified --

import qualified Numeric.GSL.Integration as I
import qualified Data.List as L

--- General Functions ---


-- | Takes every nth element, starting with the head of the list.
takeEvery :: Int -> [x] -> [x]
{-# INLINE takeEvery #-}
takeEvery :: Int -> [x] -> [x]
takeEvery Int
m = ((Int, x) -> x) -> [(Int, x)] -> [x]
forall a b. (a -> b) -> [a] -> [b]
map (Int, x) -> x
forall a b. (a, b) -> b
snd ([(Int, x)] -> [x]) -> ([x] -> [(Int, x)]) -> [x] -> [x]
forall b c a. (b -> c) -> (a -> b) -> a -> c
. ((Int, x) -> Bool) -> [(Int, x)] -> [(Int, x)]
forall a. (a -> Bool) -> [a] -> [a]
filter (\(Int
x,x
_) -> Int -> Int -> Int
forall a. Integral a => a -> a -> a
mod Int
x Int
m Int -> Int -> Bool
forall a. Eq a => a -> a -> Bool
== Int
0) ([(Int, x)] -> [(Int, x)])
-> ([x] -> [(Int, x)]) -> [x] -> [(Int, x)]
forall b c a. (b -> c) -> (a -> b) -> a -> c
. [Int] -> [x] -> [(Int, x)]
forall a b. [a] -> [b] -> [(a, b)]
zip [Int
0..]

-- | Break the list up into lists of length n.
breakEvery :: Int -> [x] -> [[x]]
{-# INLINE breakEvery #-}
breakEvery :: Int -> [x] -> [[x]]
breakEvery Int
_ [] = []
breakEvery Int
n [x]
xs = Int -> [x] -> [x]
forall a. Int -> [a] -> [a]
take Int
n [x]
xs [x] -> [[x]] -> [[x]]
forall a. a -> [a] -> [a]
: Int -> [x] -> [[x]]
forall x. Int -> [x] -> [[x]]
breakEvery Int
n (Int -> [x] -> [x]
forall a. Int -> [a] -> [a]
drop Int
n [x]
xs)

-- | Runs traceShow on the given element.
traceGiven :: Show a => a -> a
{-# INLINE traceGiven #-}
traceGiven :: a -> a
traceGiven a
a = a -> a -> a
forall a b. Show a => a -> b -> b
traceShow a
a a
a


--- Numeric ---

-- | Numerically integrates a 1-d function over an interval.
integrate
    :: Double -- ^ Error Tolerance
    -> (Double -> Double) -- ^ Function
    -> Double -- ^ Interval beginning
    -> Double -- ^ Interval end
    -> (Double,Double) -- ^ Integral
{-# INLINE integrate #-}
integrate :: Double
-> (Double -> Double) -> Double -> Double -> (Double, Double)
integrate Double
errbnd = Double
-> Int
-> (Double -> Double)
-> Double
-> Double
-> (Double, Double)
I.integrateQAGS Double
errbnd Int
10000

-- | Rounds the number to the specified significant digit.
roundSD :: RealFloat x => Int -> x -> x
{-# INLINE roundSD #-}
roundSD :: Int -> x -> x
roundSD Int
n x
x =
    let n' :: Int
        n' :: Int
n' = x -> Int
forall a b. (RealFrac a, Integral b) => a -> b
round (x -> Int) -> x -> Int
forall a b. (a -> b) -> a -> b
$ x
10x -> Int -> x
forall a b. (Num a, Integral b) => a -> b -> a
^Int
n x -> x -> x
forall a. Num a => a -> a -> a
* x
x
     in Int -> x
forall a b. (Integral a, Num b) => a -> b
fromIntegral Int
n'x -> x -> x
forall a. Fractional a => a -> a -> a
/x
10x -> Int -> x
forall a b. (Num a, Integral b) => a -> b -> a
^Int
n

-- | Value of a point on a circle, minus rotations.
toPi :: RealFloat x => x -> x
{-# INLINE toPi #-}
toPi :: x -> x
toPi x
x =
    let xpi :: x
xpi = x
x x -> x -> x
forall a. Fractional a => a -> a -> a
/ (x
2x -> x -> x
forall a. Num a => a -> a -> a
*x
forall a. Floating a => a
pi)
        f :: x
f = x
xpi x -> x -> x
forall a. Num a => a -> a -> a
- Int -> x
forall a b. (Integral a, Num b) => a -> b
fromIntegral (x -> Int
forall a b. (RealFrac a, Integral b) => a -> b
floor x
xpi :: Int)
     in x
2 x -> x -> x
forall a. Num a => a -> a -> a
* x
forall a. Floating a => a
pi x -> x -> x
forall a. Num a => a -> a -> a
* x
f

-- | Distance between two points on a circle, removing rotations.
circularDistance :: RealFloat x => x -> x -> x
{-# INLINE circularDistance #-}
circularDistance :: x -> x -> x
circularDistance x
x x
y =
    let x' :: x
x' = x -> x
forall x. RealFloat x => x -> x
toPi x
x
        y' :: x
y' = x -> x
forall x. RealFloat x => x -> x
toPi x
y
     in x -> x -> x
forall a. Ord a => a -> a -> a
min (x -> x
forall x. RealFloat x => x -> x
toPi (x -> x) -> x -> x
forall a b. (a -> b) -> a -> b
$ x
x' x -> x -> x
forall a. Num a => a -> a -> a
- x
y') (x -> x
forall x. RealFloat x => x -> x
toPi (x -> x) -> x -> x
forall a b. (a -> b) -> a -> b
$ x
y' x -> x -> x
forall a. Num a => a -> a -> a
- x
x')

-- | A standard sigmoid function.
logistic :: Floating x => x -> x
{-# INLINE logistic #-}
logistic :: x -> x
logistic x
x = x
1 x -> x -> x
forall a. Fractional a => a -> a -> a
/ (x
1 x -> x -> x
forall a. Num a => a -> a -> a
+ x -> x
forall a. Floating a => a -> a
exp (x -> x
forall a. Num a => a -> a
negate x
x))

-- | The inverse of the logistic.
logit :: Floating x => x -> x
{-# INLINE logit #-}
logit :: x -> x
logit x
x = x -> x
forall a. Floating a => a -> a
log (x -> x) -> x -> x
forall a b. (a -> b) -> a -> b
$ x
x x -> x -> x
forall a. Fractional a => a -> a -> a
/ (x
1 x -> x -> x
forall a. Num a => a -> a -> a
- x
x)

-- | The square of a number (for avoiding endless default values).
square :: Floating x => x -> x
{-# INLINE square #-}
square :: x -> x
square x
x = x
xx -> Int -> x
forall a b. (Num a, Integral b) => a -> b -> a
^(Int
2::Int)

-- | Triangular number.
triangularNumber :: Int -> Int
{-# INLINE triangularNumber #-}
triangularNumber :: Int -> Int
triangularNumber Int
n = (Int -> Int -> Int) -> Int -> Int -> Int
forall a b c. (a -> b -> c) -> b -> a -> c
flip Int -> Int -> Int
forall a. Integral a => a -> a -> a
div Int
2 (Int -> Int) -> Int -> Int
forall a b. (a -> b) -> a -> b
$ Int
n Int -> Int -> Int
forall a. Num a => a -> a -> a
* (Int
nInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
1)


-- Lists --

-- | Average value of a 'Traversable' of 'Fractional's.
average :: (Foldable f, Fractional x) => f x -> x
{-# INLINE average #-}
average :: f x -> x
average = (x -> x -> x) -> (x, x) -> x
forall a b c. (a -> b -> c) -> (a, b) -> c
uncurry x -> x -> x
forall a. Fractional a => a -> a -> a
(/) ((x, x) -> x) -> (f x -> (x, x)) -> f x -> x
forall b c a. (b -> c) -> (a -> b) -> a -> c
. ((x, x) -> x -> (x, x)) -> (x, x) -> f x -> (x, x)
forall (t :: Type -> Type) b a.
Foldable t =>
(b -> a -> b) -> b -> t a -> b
L.foldl' (\(x
s,x
c) x
e -> (x
ex -> x -> x
forall a. Num a => a -> a -> a
+x
s,x
cx -> x -> x
forall a. Num a => a -> a -> a
+x
1)) (x
0,x
0)

-- | Weighted Average given a 'Traversable' of (weight,value) pairs.
weightedAverage :: (Foldable f, Fractional x) => f (x,x) -> x
{-# INLINE weightedAverage #-}
weightedAverage :: f (x, x) -> x
weightedAverage = (x -> x -> x) -> (x, x) -> x
forall a b c. (a -> b -> c) -> (a, b) -> c
uncurry x -> x -> x
forall a. Fractional a => a -> a -> a
(/) ((x, x) -> x) -> (f (x, x) -> (x, x)) -> f (x, x) -> x
forall b c a. (b -> c) -> (a -> b) -> a -> c
. ((x, x) -> (x, x) -> (x, x)) -> (x, x) -> f (x, x) -> (x, x)
forall (t :: Type -> Type) b a.
Foldable t =>
(b -> a -> b) -> b -> t a -> b
L.foldl' (\(x
sm,x
nrm) (x
w,x
x) -> (x
sm x -> x -> x
forall a. Num a => a -> a -> a
+ x
wx -> x -> x
forall a. Num a => a -> a -> a
*x
x,x
nrm x -> x -> x
forall a. Num a => a -> a -> a
+ x
w)) (x
0,x
0)

-- | Circular average value of a 'Traversable' of radians.
circularAverage :: (Traversable f, RealFloat x) => f x -> x
{-# INLINE circularAverage #-}
circularAverage :: f x -> x
circularAverage f x
rds =
    let snmu :: x
snmu = f x -> x
forall (f :: Type -> Type) x.
(Foldable f, Fractional x) =>
f x -> x
average (f x -> x) -> f x -> x
forall a b. (a -> b) -> a -> b
$ x -> x
forall a. Floating a => a -> a
sin (x -> x) -> f x -> f x
forall (f :: Type -> Type) a b. Functor f => (a -> b) -> f a -> f b
<$> f x
rds
        csmu :: x
csmu = f x -> x
forall (f :: Type -> Type) x.
(Foldable f, Fractional x) =>
f x -> x
average (f x -> x) -> f x -> x
forall a b. (a -> b) -> a -> b
$ x -> x
forall a. Floating a => a -> a
cos (x -> x) -> f x -> f x
forall (f :: Type -> Type) a b. Functor f => (a -> b) -> f a -> f b
<$> f x
rds
     in x -> x -> x
forall a. RealFloat a => a -> a -> a
atan2 x
snmu x
csmu

-- | Returns k (training,validation) pairs. k should be greater than or equal to 2.
kFold :: Int -> [x] -> [([x],[x])]
{-# INLINE kFold #-}
kFold :: Int -> [x] -> [([x], [x])]
kFold Int
k [x]
xs =
    let nvls :: Int
nvls = Double -> Int
forall a b. (RealFrac a, Integral b) => a -> b
ceiling (Double -> Int) -> (Int -> Double) -> Int -> Int
forall b c a. (b -> c) -> (a -> b) -> a -> c
. (Double -> Double -> Double
forall a. Fractional a => a -> a -> a
/(Int -> Double
forall a b. (Integral a, Num b) => a -> b
fromIntegral Int
k :: Double)) (Double -> Double) -> (Int -> Double) -> Int -> Double
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Int -> Double
forall a b. (Integral a, Num b) => a -> b
fromIntegral (Int -> Int) -> Int -> Int
forall a b. (a -> b) -> a -> b
$ [x] -> Int
forall (t :: Type -> Type) a. Foldable t => t a -> Int
length [x]
xs
     in (([[x]], [[x]]) -> Maybe (([x], [x]), ([[x]], [[x]])))
-> ([[x]], [[x]]) -> [([x], [x])]
forall b a. (b -> Maybe (a, b)) -> b -> [a]
L.unfoldr ([[x]], [[x]]) -> Maybe (([x], [x]), ([[x]], [[x]]))
forall a. ([[a]], [[a]]) -> Maybe (([a], [a]), ([[a]], [[a]]))
unfoldFun ([], Int -> [x] -> [[x]]
forall x. Int -> [x] -> [[x]]
breakEvery Int
nvls [x]
xs)
    where unfoldFun :: ([[a]], [[a]]) -> Maybe (([a], [a]), ([[a]], [[a]]))
unfoldFun ([[a]]
_,[]) = Maybe (([a], [a]), ([[a]], [[a]]))
forall a. Maybe a
Nothing
          unfoldFun ([[a]]
hds,[a]
tl:[[a]]
tls) = (([a], [a]), ([[a]], [[a]])) -> Maybe (([a], [a]), ([[a]], [[a]]))
forall a. a -> Maybe a
Just (([[a]] -> [a]
forall (t :: Type -> Type) a. Foldable t => t [a] -> [a]
concat ([[a]] -> [a]) -> [[a]] -> [a]
forall a b. (a -> b) -> a -> b
$ [[a]]
hds [[a]] -> [[a]] -> [[a]]
forall a. [a] -> [a] -> [a]
++ [[a]]
tls,[a]
tl),([a]
tl[a] -> [[a]] -> [[a]]
forall a. a -> [a] -> [a]
:[[a]]
hds,[[a]]
tls))

-- | Returns k (training,test,validation) pairs for early stopping algorithms. k
-- should be greater than or equal to 3.
kFold' :: Int -> [x] -> [([x],[x],[x])]
{-# INLINE kFold' #-}
kFold' :: Int -> [x] -> [([x], [x], [x])]
kFold' Int
k [x]
xs =
    let nvls :: Int
nvls = Double -> Int
forall a b. (RealFrac a, Integral b) => a -> b
ceiling (Double -> Int) -> (Int -> Double) -> Int -> Int
forall b c a. (b -> c) -> (a -> b) -> a -> c
. (Double -> Double -> Double
forall a. Fractional a => a -> a -> a
/(Int -> Double
forall a b. (Integral a, Num b) => a -> b
fromIntegral Int
k :: Double)) (Double -> Double) -> (Int -> Double) -> Int -> Double
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Int -> Double
forall a b. (Integral a, Num b) => a -> b
fromIntegral (Int -> Int) -> Int -> Int
forall a b. (a -> b) -> a -> b
$ [x] -> Int
forall (t :: Type -> Type) a. Foldable t => t a -> Int
length [x]
xs
        brks :: [[x]]
brks = Int -> [x] -> [[x]]
forall x. Int -> [x] -> [[x]]
breakEvery Int
nvls [x]
xs
     in (([[x]], [[x]]) -> Maybe (([x], [x], [x]), ([[x]], [[x]])))
-> ([[x]], [[x]]) -> [([x], [x], [x])]
forall b a. (b -> Maybe (a, b)) -> b -> [a]
L.unfoldr ([[x]], [[x]]) -> Maybe (([x], [x], [x]), ([[x]], [[x]]))
forall a. ([[a]], [[a]]) -> Maybe (([a], [a], [a]), ([[a]], [[a]]))
unfoldFun ([], [[x]]
brks)
    where unfoldFun :: ([[a]], [[a]]) -> Maybe (([a], [a], [a]), ([[a]], [[a]]))
unfoldFun ([[a]]
hds,[a]
tl:[a]
tl':[[a]]
tls) = (([a], [a], [a]), ([[a]], [[a]]))
-> Maybe (([a], [a], [a]), ([[a]], [[a]]))
forall a. a -> Maybe a
Just (([[a]] -> [a]
forall (t :: Type -> Type) a. Foldable t => t [a] -> [a]
concat ([[a]] -> [a]) -> [[a]] -> [a]
forall a b. (a -> b) -> a -> b
$ [[a]]
hds [[a]] -> [[a]] -> [[a]]
forall a. [a] -> [a] -> [a]
++ [[a]]
tls,[a]
tl,[a]
tl'),([a]
tl[a] -> [[a]] -> [[a]]
forall a. a -> [a] -> [a]
:[[a]]
hds,[a]
tl'[a] -> [[a]] -> [[a]]
forall a. a -> [a] -> [a]
:[[a]]
tls))
          unfoldFun ([[a]]
hds,[a]
tl:[[a]]
tls) =
              let ([a]
tl0:[[a]]
hds') = [[a]] -> [[a]]
forall a. [a] -> [a]
reverse [[a]]
hds
               in (([a], [a], [a]), ([[a]], [[a]]))
-> Maybe (([a], [a], [a]), ([[a]], [[a]]))
forall a. a -> Maybe a
Just (([[a]] -> [a]
forall (t :: Type -> Type) a. Foldable t => t [a] -> [a]
concat ([[a]] -> [a]) -> [[a]] -> [a]
forall a b. (a -> b) -> a -> b
$ [[a]] -> [[a]]
forall a. [a] -> [a]
reverse [[a]]
hds' [[a]] -> [[a]] -> [[a]]
forall a. [a] -> [a] -> [a]
++ [[a]]
tls,[a]
tl,[a]
tl0),([a]
tl[a] -> [[a]] -> [[a]]
forall a. a -> [a] -> [a]
:[[a]]
hds,[[a]]
tls))
          unfoldFun ([[a]]
_,[]) = Maybe (([a], [a], [a]), ([[a]], [[a]]))
forall a. Maybe a
Nothing


-- | Weighted Circular average value of a 'Traversable' of radians.
weightedCircularAverage :: (Traversable f, RealFloat x) => f (x,x) -> x
{-# INLINE weightedCircularAverage #-}
weightedCircularAverage :: f (x, x) -> x
weightedCircularAverage f (x, x)
wxs =
    let snmu :: x
snmu = f (x, x) -> x
forall (f :: Type -> Type) x.
(Foldable f, Fractional x) =>
f (x, x) -> x
weightedAverage (f (x, x) -> x) -> f (x, x) -> x
forall a b. (a -> b) -> a -> b
$ (x, x) -> (x, x)
forall b a. Floating b => (a, b) -> (a, b)
sinPair ((x, x) -> (x, x)) -> f (x, x) -> f (x, x)
forall (f :: Type -> Type) a b. Functor f => (a -> b) -> f a -> f b
<$> f (x, x)
wxs
        csmu :: x
csmu = f (x, x) -> x
forall (f :: Type -> Type) x.
(Foldable f, Fractional x) =>
f (x, x) -> x
weightedAverage (f (x, x) -> x) -> f (x, x) -> x
forall a b. (a -> b) -> a -> b
$ (x, x) -> (x, x)
forall b a. Floating b => (a, b) -> (a, b)
cosPair ((x, x) -> (x, x)) -> f (x, x) -> f (x, x)
forall (f :: Type -> Type) a b. Functor f => (a -> b) -> f a -> f b
<$> f (x, x)
wxs
     in x -> x -> x
forall a. RealFloat a => a -> a -> a
atan2 x
snmu x
csmu
    where sinPair :: (a, b) -> (a, b)
sinPair (a
w,b
rd) = (a
w,b -> b
forall a. Floating a => a -> a
sin b
rd)
          cosPair :: (a, b) -> (a, b)
cosPair (a
w,b
rd) = (a
w,b -> b
forall a. Floating a => a -> a
cos b
rd)

-- | Returns n numbers which uniformly partitions the interval [mn,mx].
range
    :: RealFloat x => x -> x -> Int -> [x]
{-# INLINE range #-}
range :: x -> x -> Int -> [x]
range x
_ x
_ Int
0 = []
range x
mn x
mx Int
1 = [(x
mn x -> x -> x
forall a. Num a => a -> a -> a
+ x
mx) x -> x -> x
forall a. Fractional a => a -> a -> a
/ x
2]
range x
mn x
mx Int
n =
    [ x
x x -> x -> x
forall a. Num a => a -> a -> a
* x
mx x -> x -> x
forall a. Num a => a -> a -> a
+ (x
1 x -> x -> x
forall a. Num a => a -> a -> a
- x
x) x -> x -> x
forall a. Num a => a -> a -> a
* x
mn | x
x <- (x -> x -> x
forall a. Fractional a => a -> a -> a
/ (Int -> x
forall a b. (Integral a, Num b) => a -> b
fromIntegral Int
n x -> x -> x
forall a. Num a => a -> a -> a
- x
1)) (x -> x) -> (Int -> x) -> Int -> x
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Int -> x
forall a b. (Integral a, Num b) => a -> b
fromIntegral (Int -> x) -> [Int] -> [x]
forall (f :: Type -> Type) a b. Functor f => (a -> b) -> f a -> f b
<$> [Int
0 .. Int
nInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
1] ]

-- | Takes range information in the form of a minimum, maximum, and sample count,
-- a function to sample, and returns a list of pairs (x,f(x)) over the specified
-- range.
discretizeFunction :: Double -> Double -> Int -> (Double -> Double) -> [(Double,Double)]
{-# INLINE discretizeFunction #-}
discretizeFunction :: Double -> Double -> Int -> (Double -> Double) -> [(Double, Double)]
discretizeFunction Double
mn Double
mx Int
n Double -> Double
f =
    let rng :: [Double]
rng = Double -> Double -> Int -> [Double]
forall x. RealFloat x => x -> x -> Int -> [x]
range Double
mn Double
mx Int
n
    in [Double] -> [Double] -> [(Double, Double)]
forall a b. [a] -> [b] -> [(a, b)]
zip [Double]
rng ([Double] -> [(Double, Double)]) -> [Double] -> [(Double, Double)]
forall a b. (a -> b) -> a -> b
$ Double -> Double
f (Double -> Double) -> [Double] -> [Double]
forall (f :: Type -> Type) a b. Functor f => (a -> b) -> f a -> f b
<$> [Double]
rng

-- | Given a set of values, computes the "soft maximum" by way of taking the
-- exponential of every value, summing the results, and then taking the
-- logarithm. Incorporates some tricks to improve numerical stability.
logSumExp :: (Ord x, Floating x, Traversable f) => f x -> x
{-# INLINE logSumExp #-}
logSumExp :: f x -> x
logSumExp f x
xs =
    let mx :: x
mx = f x -> x
forall (t :: Type -> Type) a. (Foldable t, Ord a) => t a -> a
maximum f x
xs
     in (x -> x -> x
forall a. Num a => a -> a -> a
+ x
mx) (x -> x) -> (f x -> x) -> f x -> x
forall b c a. (b -> c) -> (a -> b) -> a -> c
. x -> x
forall a. Floating a => a -> a
log1p (x -> x) -> (f x -> x) -> f x -> x
forall b c a. (b -> c) -> (a -> b) -> a -> c
. x -> x -> x
forall a. Num a => a -> a -> a
subtract x
1 (x -> x) -> (f x -> x) -> f x -> x
forall b c a. (b -> c) -> (a -> b) -> a -> c
. f x -> x
forall (t :: Type -> Type) a. (Foldable t, Num a) => t a -> a
sum (f x -> x) -> f x -> x
forall a b. (a -> b) -> a -> b
$ x -> x
forall a. Floating a => a -> a
exp (x -> x) -> (x -> x) -> x -> x
forall b c a. (b -> c) -> (a -> b) -> a -> c
. x -> x -> x
forall a. Num a => a -> a -> a
subtract x
mx (x -> x) -> f x -> f x
forall (f :: Type -> Type) a b. Functor f => (a -> b) -> f a -> f b
<$> f x
xs

-- | Given a function, computes the "soft maximum" of the function by computing
-- the integral of the exponential of the function, and taking the logarithm of
-- the result. The maximum is first approximated on a given set of samples to
-- improve numerical stability. Pro tip: If you want to compute the normalizer
-- of a an exponential family probability density, provide the unnormalized
-- log-density to this function.
logIntegralExp
    :: Traversable f
    => Double -- ^ Error Tolerance
    -> (Double -> Double) -- ^ Function
    -> Double -- ^ Interval beginning
    -> Double -- ^ Interval end
    -> f Double -- ^ Samples (for approximating the max)
    -> Double -- ^ Log-Integral-Exp
{-# INLINE logIntegralExp #-}
logIntegralExp :: Double
-> (Double -> Double) -> Double -> Double -> f Double -> Double
logIntegralExp Double
err Double -> Double
f Double
mnbnd Double
mxbnd f Double
xsmps =
    let mx :: Double
mx = f Double -> Double
forall (t :: Type -> Type) a. (Foldable t, Ord a) => t a -> a
maximum (f Double -> Double) -> f Double -> Double
forall a b. (a -> b) -> a -> b
$ Double -> Double
f (Double -> Double) -> f Double -> f Double
forall (f :: Type -> Type) a b. Functor f => (a -> b) -> f a -> f b
<$> f Double
xsmps
        expf :: Double -> Double
expf Double
x = Double -> Double
forall a. Floating a => a -> a
exp (Double -> Double) -> Double -> Double
forall a b. (a -> b) -> a -> b
$ Double -> Double
f Double
x Double -> Double -> Double
forall a. Num a => a -> a -> a
- Double
mx
     in (Double -> Double -> Double
forall a. Num a => a -> a -> a
+ Double
mx) (Double -> Double)
-> ((Double, Double) -> Double) -> (Double, Double) -> Double
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Double -> Double
forall a. Floating a => a -> a
log1p (Double -> Double)
-> ((Double, Double) -> Double) -> (Double, Double) -> Double
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Double -> Double -> Double
forall a. Num a => a -> a -> a
subtract Double
1 (Double -> Double)
-> ((Double, Double) -> Double) -> (Double, Double) -> Double
forall b c a. (b -> c) -> (a -> b) -> a -> c
. (Double, Double) -> Double
forall a b. (a, b) -> a
fst ((Double, Double) -> Double) -> (Double, Double) -> Double
forall a b. (a -> b) -> a -> b
$ Double
-> (Double -> Double) -> Double -> Double -> (Double, Double)
integrate Double
err Double -> Double
expf Double
mnbnd Double
mxbnd


--- TypeLits ---


-- | Type-level triangular number.
type Triangular n = Div (n * (n + 1)) 2

-- | Type level rational numbers. This implementation does not currently permit negative numbers.
data Rat (n :: Nat) (d :: Nat)

-- | Infix 'Rat'.
type (/) n d = Rat n d

-- | Recover a rational value from a 'Proxy'.
ratVal :: (KnownNat n, KnownNat d) => Proxy (n / d) -> Rational
{-# INLINE ratVal #-}
ratVal :: Proxy (n / d) -> Rational
ratVal = Proxy n -> Proxy d -> Proxy (n / d) -> Rational
forall (n :: Nat) (d :: Nat).
(KnownNat n, KnownNat d) =>
Proxy n -> Proxy d -> Proxy (n / d) -> Rational
ratVal0 Proxy n
forall k (t :: k). Proxy t
Proxy Proxy d
forall k (t :: k). Proxy t
Proxy


-- | 'natVal and 'fromIntegral'.
natValInt :: KnownNat n => Proxy n -> Int
{-# INLINE natValInt #-}
natValInt :: Proxy n -> Int
natValInt = Natural -> Int
forall a b. (Integral a, Num b) => a -> b
fromIntegral (Natural -> Int) -> (Proxy n -> Natural) -> Proxy n -> Int
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Proxy n -> Natural
forall (n :: Nat) (proxy :: Nat -> Type).
KnownNat n =>
proxy n -> Natural
natVal

-- | 'getFinite' and 'fromIntegral'.
finiteInt :: KnownNat n => Finite n -> Int
{-# INLINE finiteInt #-}
finiteInt :: Finite n -> Int
finiteInt = Integer -> Int
forall a b. (Integral a, Num b) => a -> b
fromIntegral (Integer -> Int) -> (Finite n -> Integer) -> Finite n -> Int
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Finite n -> Integer
forall (n :: Nat). Finite n -> Integer
getFinite

ratVal0 :: (KnownNat n, KnownNat d) => Proxy n -> Proxy d -> Proxy (n / d) -> Rational
{-# INLINE ratVal0 #-}
ratVal0 :: Proxy n -> Proxy d -> Proxy (n / d) -> Rational
ratVal0 Proxy n
prxyn Proxy d
prxyd Proxy (n / d)
_ = Natural -> Integer
forall a b. (Integral a, Num b) => a -> b
fromIntegral (Proxy n -> Natural
forall (n :: Nat) (proxy :: Nat -> Type).
KnownNat n =>
proxy n -> Natural
natVal Proxy n
prxyn) Integer -> Integer -> Rational
forall a. Integral a => a -> a -> Ratio a
% Natural -> Integer
forall a b. (Integral a, Num b) => a -> b
fromIntegral (Proxy d -> Natural
forall (n :: Nat) (proxy :: Nat -> Type).
KnownNat n =>
proxy n -> Natural
natVal Proxy d
prxyd)