packages feed

mwc-random 0.13.1.0 → 0.13.1.1

raw patch · 4 files changed

+119/−26 lines, 4 files

Files

System/Random/MWC/Distributions.hs view
@@ -47,7 +47,20 @@ {-# INLINE normal #-} -- We express standard in terms of normal and not other way round -- because of bug in GHC. See bug #16 for more details.-normal m s gen = loop >>= (\x -> return $! m + s * x)+normal m s gen = do+  x <- standard gen+  return $! m + s * x++-- | Generate a normally distributed random variate with zero mean and+-- unit variance.+--+-- The implementation uses Doornik's modified ziggurat algorithm.+-- Compared to the ziggurat algorithm usually used, this is slower,+-- but generates more independent variates that pass stringent tests+-- of randomness.+standard :: PrimMonad m => Gen (PrimState m) -> m Double+{-# INLINE standard #-}+standard gen = loop   where     loop = do       u  <- (subtract 1 . (*2)) `liftM` uniform gen@@ -67,35 +80,35 @@              if e + c * (d - e) < 1                then return x                else loop-    blocks = (`I.snoc` 0) . I.cons (v/f) . I.cons r . I.unfoldrN 126 go $! T r f-      where-        go (T b g)   = let !u = T h (exp (-0.5 * h * h))-                           h  = sqrt (-2 * log (v / b + g))-                       in Just (h, u)-        v            = 9.91256303526217e-3-        f            = exp (-0.5 * r * r)-    {-# NOINLINE blocks #-}-    r                = 3.442619855899-    ratios           = I.zipWith (/) (I.tail blocks) blocks-    {-# NOINLINE ratios #-}     normalTail neg  = tailing       where tailing  = do-              x <- ((/r) . log) `liftM` uniform gen-              y <- log          `liftM` uniform gen+              x <- ((/rNorm) . log) `liftM` uniform gen+              y <- log              `liftM` uniform gen               if y * (-2) < x * x                 then tailing-                else return $! if neg then x - r else r - x+                else return $! if neg then x - rNorm else rNorm - x --- | Generate a normally distributed random variate with zero mean and--- unit variance.------ The implementation uses Doornik's modified ziggurat algorithm.--- Compared to the ziggurat algorithm usually used, this is slower,--- but generates more independent variates that pass stringent tests--- of randomness.-standard :: PrimMonad m => Gen (PrimState m) -> m Double-{-# INLINE standard #-}-standard = normal 0 1+-- Constants used by standard/normal. They are floated to the top+-- level to avoid performance regression (Bug #16) when blocks/ratios+-- are recalculated on each call to standard/normal. It's also+-- somewhat difficult to trigger reliably.+blocks :: I.Vector Double+blocks = (`I.snoc` 0) . I.cons (v/f) . I.cons rNorm . I.unfoldrN 126 go $! T rNorm f+  where+    go (T b g) = let !u = T h (exp (-0.5 * h * h))+                     h  = sqrt (-2 * log (v / b + g))+                 in Just (h, u)+    v = 9.91256303526217e-3+    f = exp (-0.5 * rNorm * rNorm)+{-# NOINLINE blocks #-}++rNorm :: Double+rNorm = 3.442619855899++ratios :: I.Vector Double+ratios = I.zipWith (/) (I.tail blocks) blocks+{-# NOINLINE ratios #-}+   -- | Generate an exponentially distributed random variate.
benchmarks/Benchmark.hs view
@@ -1,4 +1,6 @@+{-# LANGUAGE BangPatterns #-} import Control.Exception+import Control.Monad import Control.Monad.ST import Criterion.Main import Data.Int@@ -57,6 +59,19 @@       , bgroup "D"         [ bench "standard"    (standard      mwc :: IO Double)         , bench "normal"      (normal 1 3    mwc :: IO Double)+          -- Regression tests for #16. These functions should take 10x+          -- longer to execute.+          --+          -- N.B. Bang patterns are necessary to trigger the bug with+          --      GHC 7.6+        , bench "standard/N"  (replicateM_ 10 $ do+                                 !_ <- standard mwc :: IO Double+                                 return ()+                              )+        , bench "normal/N"    (replicateM_ 10 $ do+                                 !_ <- normal 1 3 mwc :: IO Double+                                 return ()+                              )         , bench "exponential" (exponential 3 mwc :: IO Double)         , bench "gamma,a<1"   (gamma 0.5 1   mwc :: IO Double)         , bench "gamma,a>1"   (gamma 2   1   mwc :: IO Double)
+ benchmarks/tsts.hs view
@@ -0,0 +1,65 @@+{-# LANGUAGE BangPatterns #-}+import Control.Monad+import System.Random.MWC+import System.Random.MWC.Distributions++main = do+  withSystemRandom $ \g -> replicateM_ (300*1000) $ do+    --+    !n <- normal 0 1  g+    !n <- normal 0 2  g+    !n <- normal 3 3  g+    !n <- normal 2 4  g+    !n <- normal 2 5  g+    !n <- normal 1 6  g+    !n <- normal 3 7  g+    !n <- normal 3 8  g+    !n <- normal 3 9  g+    !n <- normal 3 10 g+    --+    !n <- normal 0 1  g+    !n <- normal 0 2  g+    !n <- normal 3 3  g+    !n <- normal 2 4  g+    !n <- normal 2 5  g+    !n <- normal 1 6  g+    !n <- normal 3 7  g+    !n <- normal 3 8  g+    !n <- normal 3 9  g+    !n <- normal 3 10 g+    --+    !n <- normal 0 1  g+    !n <- normal 0 2  g+    !n <- normal 3 3  g+    !n <- normal 2 4  g+    !n <- normal 2 5  g+    !n <- normal 1 6  g+    !n <- normal 3 7  g+    !n <- normal 3 8  g+    !n <- normal 3 9  g+    !n <- normal 3 10 g+    --+    !n <- normal 0 1  g+    !n <- normal 0 2  g+    !n <- normal 3 3  g+    !n <- normal 2 4  g+    !n <- normal 2 5  g+    !n <- normal 1 6  g+    !n <- normal 3 7  g+    !n <- normal 3 8  g+    !n <- normal 3 9  g+    !n <- normal 3 10 g+    --+    !n <- normal 0 1  g+    !n <- normal 0 2  g+    !n <- normal 3 3  g+    !n <- normal 2 4  g+    !n <- normal 2 5  g+    !n <- normal 1 6  g+    !n <- normal 3 7  g+    !n <- normal 3 8  g+    !n <- normal 3 9  g+    !n <- normal 3 10 g+    --+    return () :: IO ()+    
mwc-random.cabal view
@@ -1,5 +1,5 @@ name:           mwc-random-version:        0.13.1.0+version:        0.13.1.1 synopsis:       Fast, high quality pseudo random number generation description:   This package contains code for generating high quality random