-- | Simple statistics code.
module ParkBench.Statistics
  ( Timed (..),
    Estimate (..),
    initialEstimate,
    updateEstimate,
    stdev,
    variance,
    goodness,
    Roll (..),
  )
where

import ParkBench.Prelude

-- | A value that took a certan time to compute.
data Timed a = Timed
  { Timed a -> Rational
nanoseconds :: {-# UNPACK #-} !Rational,
    Timed a -> a
value :: !a
  }
  deriving stock (a -> Timed b -> Timed a
(a -> b) -> Timed a -> Timed b
(forall a b. (a -> b) -> Timed a -> Timed b)
-> (forall a b. a -> Timed b -> Timed a) -> Functor Timed
forall a b. a -> Timed b -> Timed a
forall a b. (a -> b) -> Timed a -> Timed b
forall (f :: * -> *).
(forall a b. (a -> b) -> f a -> f b)
-> (forall a b. a -> f b -> f a) -> Functor f
<$ :: a -> Timed b -> Timed a
$c<$ :: forall a b. a -> Timed b -> Timed a
fmap :: (a -> b) -> Timed a -> Timed b
$cfmap :: forall a b. (a -> b) -> Timed a -> Timed b
Functor, Int -> Timed a -> ShowS
[Timed a] -> ShowS
Timed a -> String
(Int -> Timed a -> ShowS)
-> (Timed a -> String) -> ([Timed a] -> ShowS) -> Show (Timed a)
forall a. Show a => Int -> Timed a -> ShowS
forall a. Show a => [Timed a] -> ShowS
forall a. Show a => Timed a -> String
forall a.
(Int -> a -> ShowS) -> (a -> String) -> ([a] -> ShowS) -> Show a
showList :: [Timed a] -> ShowS
$cshowList :: forall a. Show a => [Timed a] -> ShowS
show :: Timed a -> String
$cshow :: forall a. Show a => Timed a -> String
showsPrec :: Int -> Timed a -> ShowS
$cshowsPrec :: forall a. Show a => Int -> Timed a -> ShowS
Show)

instance Monoid a => Monoid (Timed a) where
  mempty :: Timed a
mempty = Rational -> a -> Timed a
forall a. Rational -> a -> Timed a
Timed Rational
0 a
forall a. Monoid a => a
mempty
  mappend :: Timed a -> Timed a -> Timed a
mappend = Timed a -> Timed a -> Timed a
forall a. Semigroup a => a -> a -> a
(<>)

instance Semigroup a => Semigroup (Timed a) where
  Timed Rational
n0 a
x0 <> :: Timed a -> Timed a -> Timed a
<> Timed Rational
n1 a
x1 =
    Rational -> a -> Timed a
forall a. Rational -> a -> Timed a
Timed (Rational
n0 Rational -> Rational -> Rational
forall a. Num a => a -> a -> a
+ Rational
n1) (a
x0 a -> a -> a
forall a. Semigroup a => a -> a -> a
<> a
x1)

data Estimate a = Estimate
  { Estimate a -> Rational
kvariance :: {-# UNPACK #-} !Rational,
    Estimate a -> Timed a
mean :: {-# UNPACK #-} !(Timed a),
    Estimate a -> Word64
samples :: {-# UNPACK #-} !Word64
  }
  deriving stock (a -> Estimate b -> Estimate a
(a -> b) -> Estimate a -> Estimate b
(forall a b. (a -> b) -> Estimate a -> Estimate b)
-> (forall a b. a -> Estimate b -> Estimate a) -> Functor Estimate
forall a b. a -> Estimate b -> Estimate a
forall a b. (a -> b) -> Estimate a -> Estimate b
forall (f :: * -> *).
(forall a b. (a -> b) -> f a -> f b)
-> (forall a b. a -> f b -> f a) -> Functor f
<$ :: a -> Estimate b -> Estimate a
$c<$ :: forall a b. a -> Estimate b -> Estimate a
fmap :: (a -> b) -> Estimate a -> Estimate b
$cfmap :: forall a b. (a -> b) -> Estimate a -> Estimate b
Functor, Int -> Estimate a -> ShowS
[Estimate a] -> ShowS
Estimate a -> String
(Int -> Estimate a -> ShowS)
-> (Estimate a -> String)
-> ([Estimate a] -> ShowS)
-> Show (Estimate a)
forall a. Show a => Int -> Estimate a -> ShowS
forall a. Show a => [Estimate a] -> ShowS
forall a. Show a => Estimate a -> String
forall a.
(Int -> a -> ShowS) -> (a -> String) -> ([a] -> ShowS) -> Show a
showList :: [Estimate a] -> ShowS
$cshowList :: forall a. Show a => [Estimate a] -> ShowS
show :: Estimate a -> String
$cshow :: forall a. Show a => Estimate a -> String
showsPrec :: Int -> Estimate a -> ShowS
$cshowsPrec :: forall a. Show a => Int -> Estimate a -> ShowS
Show)

stdev :: Estimate a -> Double
stdev :: Estimate a -> Double
stdev =
  Double -> Double
forall a. Floating a => a -> a
sqrt (Double -> Double)
-> (Estimate a -> Double) -> Estimate a -> Double
forall b c a. (b -> c) -> (a -> b) -> a -> c
. (Fractional Double => Rational -> Double
forall a. Fractional a => Rational -> a
fromRational @Double) (Rational -> Double)
-> (Estimate a -> Rational) -> Estimate a -> Double
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Estimate a -> Rational
forall a. Estimate a -> Rational
variance

variance :: Estimate a -> Rational
variance :: Estimate a -> Rational
variance (Estimate Rational
kvariance Timed a
_ Word64
samples) =
  if Word64
samples Word64 -> Word64 -> Bool
forall a. Eq a => a -> a -> Bool
== Word64
1
    then Rational
0
    else Rational
kvariance Rational -> Rational -> Rational
forall a. Fractional a => a -> a -> a
/ Word64 -> Rational
w2r (Word64
samples Word64 -> Word64 -> Word64
forall a. Num a => a -> a -> a
- Word64
1)

-- | The "goodness" of an estimate, which is just how large its standard deviation is, relative to its mean.
--
-- Smaller is better, and the smallest possible value is 0.
goodness :: Estimate a -> Double
goodness :: Estimate a -> Double
goodness Estimate a
e =
  Estimate a -> Double
forall a. Estimate a -> Double
stdev Estimate a
e Double -> Double -> Double
forall a. Fractional a => a -> a -> a
/ Rational -> Double
r2d (Timed a -> Rational
forall a. Timed a -> Rational
nanoseconds (Estimate a -> Timed a
forall a. Estimate a -> Timed a
mean Estimate a
e))

-- | @initialEstimate v@ creates an estimate per thing-that-took-time @v@ that was a run of 1 iteration.
initialEstimate :: Timed a -> Estimate a
initialEstimate :: Timed a -> Estimate a
initialEstimate Timed a
mean =
  Estimate :: forall a. Rational -> Timed a -> Word64 -> Estimate a
Estimate
    { $sel:kvariance:Estimate :: Rational
kvariance = Rational
0,
      Timed a
mean :: Timed a
$sel:mean:Estimate :: Timed a
mean,
      $sel:samples:Estimate :: Word64
samples = Word64
1
    }

-- | @updateEstimate n v e@ updates estimate @e@ per thing-that-took-time @v@ that was a run of @n@ iterations.
updateEstimate :: Roll a => Word64 -> Timed a -> Estimate a -> Estimate a
updateEstimate :: Word64 -> Timed a -> Estimate a -> Estimate a
updateEstimate Word64
n (Timed Rational
tn a
value1) (Estimate Rational
kvariance (Timed Rational
mean a
value0) Word64
samples) =
  Rational -> Timed a -> Word64 -> Estimate a
forall a. Rational -> Timed a -> Word64 -> Estimate a
Estimate Rational
kvariance' (Rational -> a -> Timed a
forall a. Rational -> a -> Timed a
Timed Rational
mean' a
value') Word64
samples'
  where
    kvariance' :: Rational
kvariance' = Rational
kvariance Rational -> Rational -> Rational
forall a. Num a => a -> a -> a
+ Rational
nr Rational -> Rational -> Rational
forall a. Num a => a -> a -> a
* (Rational
t1 Rational -> Rational -> Rational
forall a. Num a => a -> a -> a
- Rational
mean) Rational -> Rational -> Rational
forall a. Num a => a -> a -> a
* (Rational
t1 Rational -> Rational -> Rational
forall a. Num a => a -> a -> a
- Rational
mean')
    mean' :: Rational
mean' = Rational -> Rational -> Rational
rollmean Rational
mean Rational
tn
    samples' :: Word64
samples' = Word64
samples Word64 -> Word64 -> Word64
forall a. Num a => a -> a -> a
+ Word64
n
    samplesr' :: Rational
samplesr' = Word64 -> Rational
w2r Word64
samples'
    t1 :: Rational
t1 = Rational
tn Rational -> Rational -> Rational
forall a. Fractional a => a -> a -> a
/ Rational
nr
    value' :: a
value' = (Rational -> Rational -> Rational) -> a -> a -> a
forall a.
Roll a =>
(Rational -> Rational -> Rational) -> a -> a -> a
roll Rational -> Rational -> Rational
rollmean a
value0 a
value1
    rollmean :: Rational -> Rational -> Rational
rollmean Rational
u0 Rational
u1 = Rational
u0 Rational -> Rational -> Rational
forall a. Num a => a -> a -> a
+ ((Rational
u1 Rational -> Rational -> Rational
forall a. Num a => a -> a -> a
- Rational
nr Rational -> Rational -> Rational
forall a. Num a => a -> a -> a
* Rational
u0) Rational -> Rational -> Rational
forall a. Fractional a => a -> a -> a
/ Rational
samplesr')
    nr :: Rational
nr = Word64 -> Rational
w2r Word64
n

class Roll a where
  roll :: (Rational -> Rational -> Rational) -> a -> a -> a