{-# LANGUAGE TypeFamilies #-}
module Simulation.Aivika.Distributed.Optimistic.Generator () where
import Control.Monad
import Control.Monad.Trans
import System.Random
import qualified System.Random.MWC as MWC
import Data.IORef
import Simulation.Aivika.Trans
import Simulation.Aivika.Trans.Generator.Primitive
import Simulation.Aivika.Distributed.Optimistic.Internal.DIO
import Simulation.Aivika.Distributed.Optimistic.Internal.IO
instance MonadGenerator DIO where
data Generator DIO =
Generator { Generator DIO -> DIO Double
generator01 :: DIO Double,
Generator DIO -> DIO Double
generatorNormal01 :: DIO Double,
Generator DIO -> DIO Int
generatorSequenceNo :: DIO Int
}
generateUniform :: Generator DIO -> Double -> Double -> DIO Double
generateUniform = DIO Double -> Double -> Double -> DIO Double
forall (m :: * -> *).
Monad m =>
m Double -> Double -> Double -> m Double
generateUniform01 (DIO Double -> Double -> Double -> DIO Double)
-> (Generator DIO -> DIO Double)
-> Generator DIO
-> Double
-> Double
-> DIO Double
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Generator DIO -> DIO Double
generator01
generateUniformInt :: Generator DIO -> Int -> Int -> DIO Int
generateUniformInt = DIO Double -> Int -> Int -> DIO Int
forall (m :: * -> *). Monad m => m Double -> Int -> Int -> m Int
generateUniformInt01 (DIO Double -> Int -> Int -> DIO Int)
-> (Generator DIO -> DIO Double)
-> Generator DIO
-> Int
-> Int
-> DIO Int
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Generator DIO -> DIO Double
generator01
generateTriangular :: Generator DIO -> Double -> Double -> Double -> DIO Double
generateTriangular = DIO Double -> Double -> Double -> Double -> DIO Double
forall (m :: * -> *).
Monad m =>
m Double -> Double -> Double -> Double -> m Double
generateTriangular01 (DIO Double -> Double -> Double -> Double -> DIO Double)
-> (Generator DIO -> DIO Double)
-> Generator DIO
-> Double
-> Double
-> Double
-> DIO Double
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Generator DIO -> DIO Double
generator01
generateNormal :: Generator DIO -> Double -> Double -> DIO Double
generateNormal = DIO Double -> Double -> Double -> DIO Double
forall (m :: * -> *).
Monad m =>
m Double -> Double -> Double -> m Double
generateNormal01 (DIO Double -> Double -> Double -> DIO Double)
-> (Generator DIO -> DIO Double)
-> Generator DIO
-> Double
-> Double
-> DIO Double
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Generator DIO -> DIO Double
generatorNormal01
generateLogNormal :: Generator DIO -> Double -> Double -> DIO Double
generateLogNormal = DIO Double -> Double -> Double -> DIO Double
forall (m :: * -> *).
Monad m =>
m Double -> Double -> Double -> m Double
generateLogNormal01 (DIO Double -> Double -> Double -> DIO Double)
-> (Generator DIO -> DIO Double)
-> Generator DIO
-> Double
-> Double
-> DIO Double
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Generator DIO -> DIO Double
generatorNormal01
generateExponential :: Generator DIO -> Double -> DIO Double
generateExponential = DIO Double -> Double -> DIO Double
forall (m :: * -> *). Monad m => m Double -> Double -> m Double
generateExponential01 (DIO Double -> Double -> DIO Double)
-> (Generator DIO -> DIO Double)
-> Generator DIO
-> Double
-> DIO Double
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Generator DIO -> DIO Double
generator01
generateErlang :: Generator DIO -> Double -> Int -> DIO Double
generateErlang = DIO Double -> Double -> Int -> DIO Double
forall (m :: * -> *).
Monad m =>
m Double -> Double -> Int -> m Double
generateErlang01 (DIO Double -> Double -> Int -> DIO Double)
-> (Generator DIO -> DIO Double)
-> Generator DIO
-> Double
-> Int
-> DIO Double
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Generator DIO -> DIO Double
generator01
generatePoisson :: Generator DIO -> Double -> DIO Int
generatePoisson = DIO Double -> Double -> DIO Int
forall (m :: * -> *). Monad m => m Double -> Double -> m Int
generatePoisson01 (DIO Double -> Double -> DIO Int)
-> (Generator DIO -> DIO Double)
-> Generator DIO
-> Double
-> DIO Int
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Generator DIO -> DIO Double
generator01
generateBinomial :: Generator DIO -> Double -> Int -> DIO Int
generateBinomial = DIO Double -> Double -> Int -> DIO Int
forall (m :: * -> *). Monad m => m Double -> Double -> Int -> m Int
generateBinomial01 (DIO Double -> Double -> Int -> DIO Int)
-> (Generator DIO -> DIO Double)
-> Generator DIO
-> Double
-> Int
-> DIO Int
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Generator DIO -> DIO Double
generator01
generateGamma :: Generator DIO -> Double -> Double -> DIO Double
generateGamma Generator DIO
g = DIO Double -> DIO Double -> Double -> Double -> DIO Double
forall (m :: * -> *).
Monad m =>
m Double -> m Double -> Double -> Double -> m Double
generateGamma01 (Generator DIO -> DIO Double
generatorNormal01 Generator DIO
g) (Generator DIO -> DIO Double
generator01 Generator DIO
g)
generateBeta :: Generator DIO -> Double -> Double -> DIO Double
generateBeta Generator DIO
g = DIO Double -> DIO Double -> Double -> Double -> DIO Double
forall (m :: * -> *).
Monad m =>
m Double -> m Double -> Double -> Double -> m Double
generateBeta01 (Generator DIO -> DIO Double
generatorNormal01 Generator DIO
g) (Generator DIO -> DIO Double
generator01 Generator DIO
g)
generateWeibull :: Generator DIO -> Double -> Double -> DIO Double
generateWeibull = DIO Double -> Double -> Double -> DIO Double
forall (m :: * -> *).
Monad m =>
m Double -> Double -> Double -> m Double
generateWeibull01 (DIO Double -> Double -> Double -> DIO Double)
-> (Generator DIO -> DIO Double)
-> Generator DIO
-> Double
-> Double
-> DIO Double
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Generator DIO -> DIO Double
generator01
generateDiscrete :: forall a. Generator DIO -> DiscretePDF a -> DIO a
generateDiscrete = DIO Double -> DiscretePDF a -> DIO a
forall (m :: * -> *) a. Monad m => m Double -> DiscretePDF a -> m a
generateDiscrete01 (DIO Double -> DiscretePDF a -> DIO a)
-> (Generator DIO -> DIO Double)
-> Generator DIO
-> DiscretePDF a
-> DIO a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Generator DIO -> DIO Double
generator01
generateSequenceNo :: Generator DIO -> DIO Int
generateSequenceNo = Generator DIO -> DIO Int
generatorSequenceNo
newGenerator :: GeneratorType DIO -> DIO (Generator DIO)
newGenerator GeneratorType DIO
tp =
case GeneratorType DIO
tp of
GeneratorType DIO
SimpleGenerator ->
do let g :: IO (IO Double)
g = Gen (PrimState IO) -> IO Double
forall a (m :: * -> *).
(Variate a, PrimMonad m) =>
Gen (PrimState m) -> m a
forall (m :: * -> *). PrimMonad m => Gen (PrimState m) -> m Double
MWC.uniform (Gen (PrimState IO) -> IO Double)
-> IO (Gen (PrimState IO)) -> IO (IO Double)
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
<$>
IO (Gen (PrimState IO))
MWC.createSystemRandom
IO Double
g' <- IO (IO Double) -> DIO (IO Double)
forall a. IO a -> DIO a
forall (m :: * -> *) a. MonadIOUnsafe m => IO a -> m a
liftIOUnsafe IO (IO Double)
g
DIO Double -> DIO (Generator DIO)
forall (m :: * -> *).
MonadGenerator m =>
m Double -> m (Generator m)
newRandomGenerator01 (IO Double -> DIO Double
forall a. IO a -> DIO a
forall (m :: * -> *) a. MonadIOUnsafe m => IO a -> m a
liftIOUnsafe IO Double
g')
SimpleGeneratorWithSeed Word32
x ->
[Char] -> DIO (Generator DIO)
forall a. HasCallStack => [Char] -> a
error [Char]
"Unsupported generator type SimpleGeneratorWithSeed: newGenerator"
CustomGenerator DIO (Generator DIO)
g ->
DIO (Generator DIO)
g
CustomGenerator01 DIO Double
g ->
DIO Double -> DIO (Generator DIO)
forall (m :: * -> *).
MonadGenerator m =>
m Double -> m (Generator m)
newRandomGenerator01 DIO Double
g
newRandomGenerator :: forall g. RandomGen g => g -> DIO (Generator DIO)
newRandomGenerator g
g =
do IORef g
r <- IO (IORef g) -> DIO (IORef g)
forall a. IO a -> DIO a
forall (m :: * -> *) a. MonadIOUnsafe m => IO a -> m a
liftIOUnsafe (IO (IORef g) -> DIO (IORef g)) -> IO (IORef g) -> DIO (IORef g)
forall a b. (a -> b) -> a -> b
$ g -> IO (IORef g)
forall a. a -> IO (IORef a)
newIORef g
g
let g01 :: DIO Double
g01 = do g
g <- IO g -> DIO g
forall a. IO a -> DIO a
forall (m :: * -> *) a. MonadIOUnsafe m => IO a -> m a
liftIOUnsafe (IO g -> DIO g) -> IO g -> DIO g
forall a b. (a -> b) -> a -> b
$ IORef g -> IO g
forall a. IORef a -> IO a
readIORef IORef g
r
let (Double
x, g
g') = g -> (Double, g)
forall g. RandomGen g => g -> (Double, g)
forall a g. (Random a, RandomGen g) => g -> (a, g)
random g
g
IO () -> DIO ()
forall a. IO a -> DIO a
forall (m :: * -> *) a. MonadIOUnsafe m => IO a -> m a
liftIOUnsafe (IO () -> DIO ()) -> IO () -> DIO ()
forall a b. (a -> b) -> a -> b
$ IORef g -> g -> IO ()
forall a. IORef a -> a -> IO ()
writeIORef IORef g
r g
g'
Double -> DIO Double
forall a. a -> DIO a
forall (m :: * -> *) a. Monad m => a -> m a
return Double
x
DIO Double -> DIO (Generator DIO)
forall (m :: * -> *).
MonadGenerator m =>
m Double -> m (Generator m)
newRandomGenerator01 DIO Double
g01
newRandomGenerator01 :: DIO Double -> DIO (Generator DIO)
newRandomGenerator01 DIO Double
g01 =
do DIO Double
gNormal01 <- DIO Double -> DIO (DIO Double)
newNormalGenerator01 DIO Double
g01
IORef Int
gSeqNoRef <- IO (IORef Int) -> DIO (IORef Int)
forall a. IO a -> DIO a
forall (m :: * -> *) a. MonadIOUnsafe m => IO a -> m a
liftIOUnsafe (IO (IORef Int) -> DIO (IORef Int))
-> IO (IORef Int) -> DIO (IORef Int)
forall a b. (a -> b) -> a -> b
$ Int -> IO (IORef Int)
forall a. a -> IO (IORef a)
newIORef Int
0
let gSeqNo :: DIO Int
gSeqNo =
do Int
x <- IO Int -> DIO Int
forall a. IO a -> DIO a
forall (m :: * -> *) a. MonadIOUnsafe m => IO a -> m a
liftIOUnsafe (IO Int -> DIO Int) -> IO Int -> DIO Int
forall a b. (a -> b) -> a -> b
$ IORef Int -> IO Int
forall a. IORef a -> IO a
readIORef IORef Int
gSeqNoRef
IO () -> DIO ()
forall a. IO a -> DIO a
forall (m :: * -> *) a. MonadIOUnsafe m => IO a -> m a
liftIOUnsafe (IO () -> DIO ()) -> IO () -> DIO ()
forall a b. (a -> b) -> a -> b
$ IORef Int -> (Int -> Int) -> IO ()
forall a. IORef a -> (a -> a) -> IO ()
modifyIORef' IORef Int
gSeqNoRef (Int -> Int -> Int
forall a. Num a => a -> a -> a
+Int
1)
Int -> DIO Int
forall a. a -> DIO a
forall (m :: * -> *) a. Monad m => a -> m a
return Int
x
Generator DIO -> DIO (Generator DIO)
forall a. a -> DIO a
forall (m :: * -> *) a. Monad m => a -> m a
return Generator { generator01 :: DIO Double
generator01 = DIO Double
g01,
generatorNormal01 :: DIO Double
generatorNormal01 = DIO Double
gNormal01,
generatorSequenceNo :: DIO Int
generatorSequenceNo = DIO Int
gSeqNo }
newNormalGenerator01 :: DIO Double
-> DIO (DIO Double)
newNormalGenerator01 :: DIO Double -> DIO (DIO Double)
newNormalGenerator01 DIO Double
g =
do IORef Double
nextRef <- IO (IORef Double) -> DIO (IORef Double)
forall a. IO a -> DIO a
forall (m :: * -> *) a. MonadIOUnsafe m => IO a -> m a
liftIOUnsafe (IO (IORef Double) -> DIO (IORef Double))
-> IO (IORef Double) -> DIO (IORef Double)
forall a b. (a -> b) -> a -> b
$ Double -> IO (IORef Double)
forall a. a -> IO (IORef a)
newIORef Double
0.0
IORef Bool
flagRef <- IO (IORef Bool) -> DIO (IORef Bool)
forall a. IO a -> DIO a
forall (m :: * -> *) a. MonadIOUnsafe m => IO a -> m a
liftIOUnsafe (IO (IORef Bool) -> DIO (IORef Bool))
-> IO (IORef Bool) -> DIO (IORef Bool)
forall a b. (a -> b) -> a -> b
$ Bool -> IO (IORef Bool)
forall a. a -> IO (IORef a)
newIORef Bool
False
IORef Double
xi1Ref <- IO (IORef Double) -> DIO (IORef Double)
forall a. IO a -> DIO a
forall (m :: * -> *) a. MonadIOUnsafe m => IO a -> m a
liftIOUnsafe (IO (IORef Double) -> DIO (IORef Double))
-> IO (IORef Double) -> DIO (IORef Double)
forall a b. (a -> b) -> a -> b
$ Double -> IO (IORef Double)
forall a. a -> IO (IORef a)
newIORef Double
0.0
IORef Double
xi2Ref <- IO (IORef Double) -> DIO (IORef Double)
forall a. IO a -> DIO a
forall (m :: * -> *) a. MonadIOUnsafe m => IO a -> m a
liftIOUnsafe (IO (IORef Double) -> DIO (IORef Double))
-> IO (IORef Double) -> DIO (IORef Double)
forall a b. (a -> b) -> a -> b
$ Double -> IO (IORef Double)
forall a. a -> IO (IORef a)
newIORef Double
0.0
IORef Double
psiRef <- IO (IORef Double) -> DIO (IORef Double)
forall a. IO a -> DIO a
forall (m :: * -> *) a. MonadIOUnsafe m => IO a -> m a
liftIOUnsafe (IO (IORef Double) -> DIO (IORef Double))
-> IO (IORef Double) -> DIO (IORef Double)
forall a b. (a -> b) -> a -> b
$ Double -> IO (IORef Double)
forall a. a -> IO (IORef a)
newIORef Double
0.0
let loop :: DIO ()
loop =
do Double
psi <- IO Double -> DIO Double
forall a. IO a -> DIO a
forall (m :: * -> *) a. MonadIOUnsafe m => IO a -> m a
liftIOUnsafe (IO Double -> DIO Double) -> IO Double -> DIO Double
forall a b. (a -> b) -> a -> b
$ IORef Double -> IO Double
forall a. IORef a -> IO a
readIORef IORef Double
psiRef
if (Double
psi Double -> Double -> Bool
forall a. Ord a => a -> a -> Bool
>= Double
1.0) Bool -> Bool -> Bool
|| (Double
psi Double -> Double -> Bool
forall a. Eq a => a -> a -> Bool
== Double
0.0)
then do Double
g1 <- DIO Double
g
Double
g2 <- DIO Double
g
let xi1 :: Double
xi1 = Double
2.0 Double -> Double -> Double
forall a. Num a => a -> a -> a
* Double
g1 Double -> Double -> Double
forall a. Num a => a -> a -> a
- Double
1.0
xi2 :: Double
xi2 = Double
2.0 Double -> Double -> Double
forall a. Num a => a -> a -> a
* Double
g2 Double -> Double -> Double
forall a. Num a => a -> a -> a
- Double
1.0
psi :: Double
psi = Double
xi1 Double -> Double -> Double
forall a. Num a => a -> a -> a
* Double
xi1 Double -> Double -> Double
forall a. Num a => a -> a -> a
+ Double
xi2 Double -> Double -> Double
forall a. Num a => a -> a -> a
* Double
xi2
IO () -> DIO ()
forall a. IO a -> DIO a
forall (m :: * -> *) a. MonadIOUnsafe m => IO a -> m a
liftIOUnsafe (IO () -> DIO ()) -> IO () -> DIO ()
forall a b. (a -> b) -> a -> b
$ IORef Double -> Double -> IO ()
forall a. IORef a -> a -> IO ()
writeIORef IORef Double
xi1Ref Double
xi1
IO () -> DIO ()
forall a. IO a -> DIO a
forall (m :: * -> *) a. MonadIOUnsafe m => IO a -> m a
liftIOUnsafe (IO () -> DIO ()) -> IO () -> DIO ()
forall a b. (a -> b) -> a -> b
$ IORef Double -> Double -> IO ()
forall a. IORef a -> a -> IO ()
writeIORef IORef Double
xi2Ref Double
xi2
IO () -> DIO ()
forall a. IO a -> DIO a
forall (m :: * -> *) a. MonadIOUnsafe m => IO a -> m a
liftIOUnsafe (IO () -> DIO ()) -> IO () -> DIO ()
forall a b. (a -> b) -> a -> b
$ IORef Double -> Double -> IO ()
forall a. IORef a -> a -> IO ()
writeIORef IORef Double
psiRef Double
psi
DIO ()
loop
else IO () -> DIO ()
forall a. IO a -> DIO a
forall (m :: * -> *) a. MonadIOUnsafe m => IO a -> m a
liftIOUnsafe (IO () -> DIO ()) -> IO () -> DIO ()
forall a b. (a -> b) -> a -> b
$ IORef Double -> Double -> IO ()
forall a. IORef a -> a -> IO ()
writeIORef IORef Double
psiRef (Double -> IO ()) -> Double -> IO ()
forall a b. (a -> b) -> a -> b
$ Double -> Double
forall a. Floating a => a -> a
sqrt (- Double
2.0 Double -> Double -> Double
forall a. Num a => a -> a -> a
* Double -> Double
forall a. Floating a => a -> a
log Double
psi Double -> Double -> Double
forall a. Fractional a => a -> a -> a
/ Double
psi)
DIO Double -> DIO (DIO Double)
forall a. a -> DIO a
forall (m :: * -> *) a. Monad m => a -> m a
return (DIO Double -> DIO (DIO Double)) -> DIO Double -> DIO (DIO Double)
forall a b. (a -> b) -> a -> b
$
do Bool
flag <- IO Bool -> DIO Bool
forall a. IO a -> DIO a
forall (m :: * -> *) a. MonadIOUnsafe m => IO a -> m a
liftIOUnsafe (IO Bool -> DIO Bool) -> IO Bool -> DIO Bool
forall a b. (a -> b) -> a -> b
$ IORef Bool -> IO Bool
forall a. IORef a -> IO a
readIORef IORef Bool
flagRef
if Bool
flag
then do IO () -> DIO ()
forall a. IO a -> DIO a
forall (m :: * -> *) a. MonadIOUnsafe m => IO a -> m a
liftIOUnsafe (IO () -> DIO ()) -> IO () -> DIO ()
forall a b. (a -> b) -> a -> b
$ IORef Bool -> Bool -> IO ()
forall a. IORef a -> a -> IO ()
writeIORef IORef Bool
flagRef Bool
False
IO Double -> DIO Double
forall a. IO a -> DIO a
forall (m :: * -> *) a. MonadIOUnsafe m => IO a -> m a
liftIOUnsafe (IO Double -> DIO Double) -> IO Double -> DIO Double
forall a b. (a -> b) -> a -> b
$ IORef Double -> IO Double
forall a. IORef a -> IO a
readIORef IORef Double
nextRef
else do IO () -> DIO ()
forall a. IO a -> DIO a
forall (m :: * -> *) a. MonadIOUnsafe m => IO a -> m a
liftIOUnsafe (IO () -> DIO ()) -> IO () -> DIO ()
forall a b. (a -> b) -> a -> b
$ IORef Double -> Double -> IO ()
forall a. IORef a -> a -> IO ()
writeIORef IORef Double
xi1Ref Double
0.0
IO () -> DIO ()
forall a. IO a -> DIO a
forall (m :: * -> *) a. MonadIOUnsafe m => IO a -> m a
liftIOUnsafe (IO () -> DIO ()) -> IO () -> DIO ()
forall a b. (a -> b) -> a -> b
$ IORef Double -> Double -> IO ()
forall a. IORef a -> a -> IO ()
writeIORef IORef Double
xi2Ref Double
0.0
IO () -> DIO ()
forall a. IO a -> DIO a
forall (m :: * -> *) a. MonadIOUnsafe m => IO a -> m a
liftIOUnsafe (IO () -> DIO ()) -> IO () -> DIO ()
forall a b. (a -> b) -> a -> b
$ IORef Double -> Double -> IO ()
forall a. IORef a -> a -> IO ()
writeIORef IORef Double
psiRef Double
0.0
DIO ()
loop
Double
xi1 <- IO Double -> DIO Double
forall a. IO a -> DIO a
forall (m :: * -> *) a. MonadIOUnsafe m => IO a -> m a
liftIOUnsafe (IO Double -> DIO Double) -> IO Double -> DIO Double
forall a b. (a -> b) -> a -> b
$ IORef Double -> IO Double
forall a. IORef a -> IO a
readIORef IORef Double
xi1Ref
Double
xi2 <- IO Double -> DIO Double
forall a. IO a -> DIO a
forall (m :: * -> *) a. MonadIOUnsafe m => IO a -> m a
liftIOUnsafe (IO Double -> DIO Double) -> IO Double -> DIO Double
forall a b. (a -> b) -> a -> b
$ IORef Double -> IO Double
forall a. IORef a -> IO a
readIORef IORef Double
xi2Ref
Double
psi <- IO Double -> DIO Double
forall a. IO a -> DIO a
forall (m :: * -> *) a. MonadIOUnsafe m => IO a -> m a
liftIOUnsafe (IO Double -> DIO Double) -> IO Double -> DIO Double
forall a b. (a -> b) -> a -> b
$ IORef Double -> IO Double
forall a. IORef a -> IO a
readIORef IORef Double
psiRef
IO () -> DIO ()
forall a. IO a -> DIO a
forall (m :: * -> *) a. MonadIOUnsafe m => IO a -> m a
liftIOUnsafe (IO () -> DIO ()) -> IO () -> DIO ()
forall a b. (a -> b) -> a -> b
$ IORef Bool -> Bool -> IO ()
forall a. IORef a -> a -> IO ()
writeIORef IORef Bool
flagRef Bool
True
IO () -> DIO ()
forall a. IO a -> DIO a
forall (m :: * -> *) a. MonadIOUnsafe m => IO a -> m a
liftIOUnsafe (IO () -> DIO ()) -> IO () -> DIO ()
forall a b. (a -> b) -> a -> b
$ IORef Double -> Double -> IO ()
forall a. IORef a -> a -> IO ()
writeIORef IORef Double
nextRef (Double -> IO ()) -> Double -> IO ()
forall a b. (a -> b) -> a -> b
$ Double
xi2 Double -> Double -> Double
forall a. Num a => a -> a -> a
* Double
psi
Double -> DIO Double
forall a. a -> DIO a
forall (m :: * -> *) a. Monad m => a -> m a
return (Double -> DIO Double) -> Double -> DIO Double
forall a b. (a -> b) -> a -> b
$ Double
xi1 Double -> Double -> Double
forall a. Num a => a -> a -> a
* Double
psi