system-random-effect 0.3.1 → 0.4.0
raw patch · 5 files changed
+429/−123 lines, 5 filesdep +bytestringdep +crypto-apidep +primitive
Dependencies added: bytestring, crypto-api, primitive
Files
- bench/Bench.hs +23/−2
- src/System/Random/Effect.hs +140/−99
- src/System/Random/Effect/Raw.hs +146/−0
- system-random-effect.cabal +9/−2
- test/Test.hs +111/−20
bench/Bench.hs view
@@ -3,6 +3,8 @@ {-# LANGUAGE FlexibleInstances #-} module Main ( main ) where +import Data.Word+ import Control.Eff import Control.Eff.State.Strict import Control.Applicative@@ -41,9 +43,28 @@ main = defaultMainWith config (return ()) [ bcompare [- bench "fastIntDist" (fastIntDist 1 1000 :: Eff (State Random :> ()) Int)- , bench "uniformIntDist" (uniformIntDist 1 1000 :: Eff (State Random :> ()) Integer)+ bench "fastIntDist" (fastIntDist 1 1000 :: Eff (State Random :> ()) Int)+ , bench "uniformIntDist" (uniformIntDist 1 1000 :: Eff (State Random :> ()) Integer)+ , bench "uniformIntegralDist" (uniformIntegralDist 1 1000 :: Eff (State Random :> ()) Word64) ]+ , bcompare [+ bench "normalDist" (normalDist 0 10 :: Eff (State Random :> ()) Double)+ , bench "lognormalDist" (lognormalDist 0 10 :: Eff (State Random :> ()) Double)+ ]+ , bcompare [+ bench "randomWord32" (randomWord :: Eff (State Random :> ()) Word)+ , bench "randomWord64" (randomWord64 :: Eff (State Random :> ()) Word64)+ , bench "randomBits32" (randomBits :: Eff (State Random :> ()) Word32)+ , bench "randomBits64" (randomBits :: Eff (State Random :> ()) Word64)+ , bench "randomBitList64" (randomBitList 64 :: Eff (State Random :> ()) [Bool])+ , bench "randomBitList128" (randomBitList 128 :: Eff (State Random :> ()) [Bool])+ ]+ , bcompare [+ bench "randomDouble" (randomDouble :: Eff (State Random :> ()) Double)+ , bench "uniformRealDist" (uniformRealDist 0 10 :: Eff (State Random :> ()) Double)+ , bench "linearRealDist" (linearRealDist 0 10 :: Eff (State Random :> ()) Double)+ ]+ , bench "binomialDist" (binomialDist 1000 0.5 :: Eff (State Random :> ()) Int) ] where config = defaultConfig {
src/System/Random/Effect.hs view
@@ -3,6 +3,8 @@ {-# LANGUAGE FlexibleContexts #-} {-# LANGUAGE ScopedTypeVariables #-} {-# LANGUAGE BangPatterns #-}+{-# LANGUAGE FlexibleInstances #-}+{-# LANGUAGE UndecidableInstances #-} -- | A random number effect, using a pure mersenne twister under -- the hood. This should be plug-and-play with any application -- making use of extensible effects.@@ -13,12 +15,16 @@ -- * Seeding , mkRandom , mkRandomIO+ , mkSecureRandomIO -- * Running+ , forRandEff , runRandomState -- * Uniform Distributions , uniformIntDist , uniformIntegralDist , uniformRealDist+ -- * Linear Distributions+ , linearRealDist -- * Bernoulli Distributions , bernoulliDist , binomialDist@@ -42,6 +48,10 @@ , buildDDH , discreteDist , piecewiseConstantDist+ , piecewiseLinearDist+ -- * Shuffling+ , knuthShuffle+ , knuthShuffleM -- * Raw Generators , randomInt , randomInt64@@ -54,15 +64,17 @@ import Control.Applicative import Control.Monad+import Control.Monad.Primitive import Control.Monad.ST import Data.Bits-import Data.Int import Data.List import Data.Ratio import Data.Typeable import Data.Vector.Algorithms.Search-import Data.Vector ( Vector )+import Data.Vector ( Vector, (//), (!) ) import qualified Data.Vector as V+import Data.Vector.Mutable (MVector)+import qualified Data.Vector.Mutable as M import Data.Word import Statistics.Distribution @@ -75,79 +87,23 @@ --import qualified Statistics.Distribution.Poisson as DP import qualified Statistics.Distribution.StudentT as DS -import qualified System.Random.Mersenne.Pure64 as SR- import Control.Eff import Control.Eff.Lift import Control.Eff.State.Strict --- | A pure mersenne twister pseudo-random number generator.-data Random = Random {-# UNPACK #-} !SR.PureMT- deriving Typeable---- | Create a random number generator from a 'Word64' seed.-mkRandom :: Word64 -> Random-mkRandom = Random . SR.pureMT-{-# INLINE mkRandom #-}---- | Create a new random number generator, using the clocktime as the base for--- the seed. This must be called from a computation with a lifted base effect--- of 'IO'.-mkRandomIO :: SetMember Lift (Lift IO) r- => Eff r Random-mkRandomIO = lift (fmap Random SR.newPureMT)-{-# INLINE mkRandomIO #-}---- | Runs an effectful random computation, returning the computation's result.-runRandomState :: Random- -> Eff (State Random :> r) w- -> Eff r w-runRandomState seed computation =- snd <$> runState seed computation-{-# INLINE runRandomState #-}---- | A generalized form of generating a random number of the correct type--- from System.Random.Mersenne.Pure64.-randomF :: Member (State Random) r- => (SR.PureMT -> (a, SR.PureMT))- -> Eff r a-randomF f = do- (Random old) <- get- let (val, new) = f old- put (Random new)- return val-{-# INLINE randomF #-}---- | Yield a new 'Int' value from the generator. The full 64 bits will be used--- on a 64 bit machine.-randomInt :: Member (State Random) r => Eff r Int-randomInt = randomF SR.randomInt-{-# INLINE randomInt #-}---- | Yield a new 'Word' value from the generator.-randomWord :: Member (State Random) r => Eff r Word-randomWord = randomF SR.randomWord-{-# INLINE randomWord #-}---- | Yield a new 'Int64' value from the generator.-randomInt64 :: Member (State Random) r => Eff r Int64-randomInt64 = randomF SR.randomInt64-{-# INLINE randomInt64 #-}---- | Yield a new 'Word64' value from the generator.-randomWord64 :: Member (State Random) r => Eff r Word64-randomWord64 = randomF SR.randomWord64-{-# INLINE randomWord64 #-}+import System.Random.Effect.Raw --- | Yield a new 53-bit precise 'Double' value from the generator.--- The returned number will be in the range [0, 1).-randomDouble :: Member (State Random) r => Eff r Double-randomDouble = randomF SR.randomDouble+randomDouble :: Member (State Random) r+ => Eff r Double+randomDouble = do+ r <- randomWord64+ return (fromIntegral (r `div` 2048) / 9007199254740992) {-# INLINE randomDouble #-} -- | Yields a set of random from the internal generator, -- using 'randomWord64' internally.-randomBits :: (Bits x, Member (State Random) r)+randomBits :: ( Member (State Random) r+ , Bits x) => Eff r x randomBits = do let z = clearBit (bit 0) 0 -- zero, so we can get the number of bits@@ -172,14 +128,10 @@ return $ take k (concatMap breakBits word64s) -- | Returns the maximum set bit in an integer.-maxBit :: Integer -> Int-maxBit x- | x == 0 = 0- | otherwise =- let loop y !k- | y == 1 = k- | otherwise = loop (y `unsafeShiftR` 1) (k+1)- in loop x 0+maxBit :: Num a => Integer -> a+maxBit = let loop !k 0 = k+ loop !k x = loop (k+1) (x `unsafeShiftR` 1)+ in loop 0 -- | Repeat a computation until it succeeds a test. loopUntil :: Monad m => (a -> Bool) -> m a -> m a@@ -238,8 +190,8 @@ -- since it relaxes type constraints, but passing in -- constant bounds such as @uniformIntegralDist 0 10@ -- will warn with -Wall.-uniformIntegralDist :: (Member (State Random) r- , Integral a)+uniformIntegralDist :: ( Member (State Random) r+ , Integral a ) => a -- ^ a -> a -- ^ b -> Eff r a@@ -260,7 +212,7 @@ return (d * range + a) -- | Generates a uniformly distributed random number in--- the inclusive range [a, b].+-- the range [a, b). -- -- NOTE: This code might not be correct, in that the -- returned value may not be perfectly uniformly@@ -280,6 +232,36 @@ in uniformRealDist' a range {-# INLINE uniformRealDist #-} +-- | Generates a linearly-distributed random number in the range @[a, b)@;+-- @a@ with a probability of 0.+-- This code is not guaranteed to be correct.+linearRealDist :: Member (State Random) r+ => Double+ -> Double+ -> Eff r Double+linearRealDist a' b' = do+ let (a, b) = (min a' b', max a' b')+ range = b - a+ -- P(x) =+ -- { 0 x < a+ -- { k*(x-a) a <= x < b+ -- { 0 b <= x+ -- CDF(x) =+ -- { 0 x < a+ -- { k*(1/2)*(x-a)^2 a <= x < b+ -- { 1 b <= x+ -- Choose k such that CDF is continuous:+ -- lim+ -- x->b- CDF(x) = k*(1/2)*(b-a)^2 = 1+ -- k = 2/(b-a)^2 = 2/range^2+ -- Thus CDF(x, a<=x<b) = ((x-a)/range)^2+ -- CDF^-1 :: (A, 0<=A<1) -> (x, a<=x<b)+ -- CDF^1(A) = sqrt(A) * range + a+ -- Given a uniformly-distributed number in [0, 1),+ -- we can make it linearly-distributed by applying CDF^-1.+ r <- randomDouble+ return $ sqrt(r) * range + a+ -- | Samples a continuous distribution. -- -- Generates random numbers as if they were sampled by the given@@ -287,8 +269,8 @@ -- -- This is implemented with the inverse transform rule: -- <http://en.wikipedia.org/wiki/Inverse_transform_sampling>.-sampleContDist :: (Member (State Random) r- , ContDistr d)+sampleContDist :: ( Member (State Random) r+ , ContDistr d ) => d -- ^ The distribution to sample. -> Eff r Double sampleContDist d =@@ -335,7 +317,7 @@ -- | The value represents the number of failures in a series of -- independent yes/no trials (each succeeds with probability p),--- before exactly k successes occur. +-- before exactly k successes occur. -- -- p must be in the range (0, 1] -- k must be >= 0@@ -406,15 +388,16 @@ gammaDist alpha beta = sampleContDist (DG.gammaDistr alpha beta) --- | ???------ Warning: NOT IMPLEMENTED!+-- | Generates random numbers as sampled from a Weibull+-- distribution. It was originally identified to describe+-- particle size distribution. weibullDist :: Member (State Random) r => Double -- ^ α. The shape parameter. -> Double -- ^ β. The scale parameter. -> Eff r Double-weibullDist =- error "system-random-effect: weibullDist: TODO: Patches welcome!"+weibullDist alpha beta = do -- lambda = beta, k = alpha+ r <- randomDouble+ return $ beta * (-log (1-r)) ** (1/alpha) -- | ??? --@@ -426,6 +409,7 @@ extremeValueDist = error "system-random-effect: extremeValueDist: TODO: Patches welcome!" + -- | Generates random numbers as sampled from the -- normal distribution. normalDist :: Member (State Random) r@@ -485,7 +469,7 @@ -- measurements, each with additive errors of unknown standard -- deviation, as in physical measurements. Or, alternatively, -- when estimating the unknown mean of a normal distribution--- with unknown standard deviation, given n+1 samples. +-- with unknown standard deviation, given n+1 samples. studentTDist :: Member (State Random) r => Double -- ^ The number of degrees of freedom -> Eff r Double@@ -526,30 +510,87 @@ mv <- V.unsafeThaw xs binarySearch mv y +-- Select an interval at random from a collection of intervals and+-- their weights, then run a provided function between its endpoints.+piecewiseDist :: Member (State Random) r+ => (Double -> Double -> Eff r Double)+ -> [Double]+ -> DiscreteDistributionHelper+ -> Eff r Double+piecewiseDist f intervals weights@(DDH rs)+ | V.length rs + 1 /= length intervals =+ error $ "system-random-effect: piecewiseDist:"+ ++ " Incongruent parameter lengths."+ ++ " intervals=" ++ show intervals+ ++ " weights=" ++ show rs+ | otherwise = do+ idx <- discreteDist weights++ let vints = V.fromList intervals+ (l, r) = (vints V.! idx, vints V.! (idx+1))+ f l r+{-# INLINE piecewiseDist #-}+ -- | This function produces random floating-point numbers, which -- are uniformly distributed within each of the several subintervals -- [b_i, b_(i+1)), each with its own weight w_i. The set of interval -- boundaries and the set of weights are the parameters of this -- distribution. ----- For example, `piecewiseConstantDist [ 0, 1, 10, 15 ]--- (buildDDH [ 1, 0, 1 ])`+-- For example, @piecewiseConstantDist [ 0, 1, 10, 15 ]+-- (buildDDH [ 1, 0, 1 ])@ -- will produce values between 0 and 1 half the time, and values -- between 10 and 15 the other half of the time. piecewiseConstantDist :: Member (State Random) r => [Double] -- ^ Intervals -> DiscreteDistributionHelper -- ^ Weights -> Eff r Double-piecewiseConstantDist intervals weights@(DDH rs)- | V.length rs + 1 /= length intervals =- error $ "system-random-effect: piecewiseConstantDist:"- ++ " Incongruent parameter lengths."- ++ " intervals=" ++ show intervals- ++ " weights=" ++ show rs- | otherwise = do- idx <- discreteDist weights+piecewiseConstantDist = piecewiseDist uniformRealDist - let vints = V.fromList intervals- (l, r) = (vints V.! idx, vints V.! (idx+1))+-- | This function produces random floating-point numbers, which+-- are distributed with linearly-increasing probability within+-- each of the several subintervals [b_i, b_(i+1)), each with its own+-- weight w_i. The set of interval boundaries and the set of weights+-- are the parameters of this distribution.+--+-- For example, `piecewiseLinearDist [ 0, 1, 10, 15 ]+-- (buildDDH [ 1, 0, 1 ])`+-- will produce values between 0 and 1 half the time, and values+-- between 10 and 15 the other half of the time.+piecewiseLinearDist :: Member (State Random) r+ => [Double] -- ^ Intervals+ -> DiscreteDistributionHelper -- ^ Weights+ -> Eff r Double+piecewiseLinearDist = piecewiseDist linearRealDist - uniformRealDist l r+-- | Shuffle a mutable vector.+knuthShuffleM :: ( PrimMonad m+ , Applicative m+ , Typeable1 m+ , Member (State Random) r+ , SetMember Lift (Lift m) r+ )+ => MVector (PrimState m) a+ -> Eff r ()+knuthShuffleM v = forM_ [0..iLast]+ $ \i -> uniformIntegralDist 0 iLast+ >>= lift . swap i+ where+ iLast = M.length v - 1++ swap i j = do+ [vi, vj] <- mapM (M.read v) [i, j]+ mapM_ (uncurry $ M.write v) [(i, vj), (j, vi)]++-- | Shuffle an immutable vector.+knuthShuffle :: Member (State Random) r+ => Vector a+ -> Eff r (Vector a)+knuthShuffle v0 = foldM swap v0 [0..iLast]+ where+ iLast = V.length v0 - 1++ swap :: Member (State Random) r => Vector a -> Int -> Eff r (Vector a)+ swap v i = do+ j <- uniformIntegralDist 0 iLast+ return $ v // zip [j, i] ((v !) <$> [i, j])
+ src/System/Random/Effect/Raw.hs view
@@ -0,0 +1,146 @@+-- | Contains a random number generator for the deterministic,+-- fast, but NOT SECURE mt19937 random number generator. This+-- is a very versatile (and pure!) generator that works for+-- pretty much everything EXCEPT security.+--+-- Seriously. Don't use this for crypto.+--+-- This is the minimal definition for a random number generator.+-- new ones can be implemented by adding to the sum type 'Random'.+{-# LANGUAGE DeriveDataTypeable #-}+{-# LANGUAGE TypeOperators #-}+{-# LANGUAGE FlexibleContexts #-}+{-# LANGUAGE BangPatterns #-}+module System.Random.Effect.Raw ( Random+ -- * Seeding+ , mkRandom+ , mkRandomIO+ , mkSecureRandomIO+ -- * Running+ , forRandEff+ , runRandomState+ -- * Raw Generators+ , randomInt+ , randomWord+ , randomInt64+ , randomWord64+ ) where++import Control.Applicative+import Control.Arrow ( second )+import Data.Bits+import qualified Data.ByteString as B+import Data.Typeable+import Data.Int+import Data.Word++import qualified Crypto.Random as CR+import qualified System.Random.Mersenne.Pure64 as SR++import Control.Eff+import Control.Eff.Lift+import Control.Eff.State.Strict++-- | A random number generator. Either a fast, insecure mersenne+-- twister or a secure one, depending on which smart constructor+-- is used to construct this type.+data Random = FastRandom {-# UNPACK #-} !SR.PureMT+ | SecureRandom {-# UNPACK #-} !CR.SystemRandom+ deriving Typeable++randomInt :: Member (State Random) r+ => Eff r Int+randomInt = randomF SR.randomInt srandomBits+{-# INLINE randomInt #-}++randomWord :: Member (State Random) r+ => Eff r Word+randomWord = randomF SR.randomWord srandomBits+{-# INLINE randomWord #-}++randomInt64 :: Member (State Random) r+ => Eff r Int64+randomInt64 = randomF SR.randomInt64 srandomBits+{-# INLINE randomInt64 #-}++randomWord64 :: Member (State Random) r+ => Eff r Word64+randomWord64 = randomF SR.randomWord64 srandomBits+{-# INLINE randomWord64 #-}++-- | Create a random number generator from a 'Word64' seed.+-- This uses the insecure (but fast) mersenne twister.+mkRandom :: Word64 -> Random+mkRandom = FastRandom . SR.pureMT+{-# INLINE mkRandom #-}++-- | Create a new random number generator, using the clocktime as the base for+-- the seed. This must be called from a computation with a lifted base effect+-- of 'IO'.+--+-- This is just a conveniently seeded mersenne twister.+mkRandomIO :: SetMember Lift (Lift IO) r+ => Eff r Random+mkRandomIO =+ FastRandom <$> lift SR.newPureMT+{-# INLINE mkRandomIO #-}++-- | Creates a new random number generator, using the system entropy source+-- as a seed. The random number generator returned from this function is+-- cryptographically secure, but not nearly as fast as the one returned+-- by 'mkRandom' and 'mkRandomIO'.+mkSecureRandomIO :: SetMember Lift (Lift IO) r+ => Eff r Random+mkSecureRandomIO = do+ SecureRandom <$> lift CR.newGenIO+{-# INLINE mkSecureRandomIO #-}++-- | Use a non-random effect as the Random source for running a random effect.+forRandEff :: Eff r Random -> Eff (State Random :> r) w -> Eff r w+forRandEff rndgen e = rndgen >>= (`runRandomState` e)+{-# INLINE forRandEff #-}++-- | Runs an effectful random computation, returning the computation's result.+runRandomState :: Random+ -> Eff (State Random :> r) w+ -> Eff r w+runRandomState seed computation =+ snd <$> runState seed computation+{-# INLINE runRandomState #-}++foldBits :: (Bits a, Num a)+ => B.ByteString+ -> a+foldBits bs =+ let addByte byte bits =+ (bits `unsafeShiftL` 8) .|. fromIntegral byte+ in B.foldr' addByte 0 bs+{-# INLINE foldBits #-}++-- | Securely generate some random bits.+srandomBits :: ( Bits a+ , Num a )+ => CR.SystemRandom+ -> (a, CR.SystemRandom)+srandomBits sr =+ let z = clearBit (bit 0) 0+ nBytes = bitSize z `div` 8+ in case CR.genBytes nBytes sr of+ Left err -> error $ "system-random-effect: System.Random.Effect.Secure: genBytes: " ++ show err+ Right (bs, sr') -> (z .|. foldBits bs, sr')+{-# INLINE srandomBits #-}++-- | A generalized form of generating a random number of the correct type+-- from System.Random.Mersenne.Pure64.+randomF :: Member (State Random) r+ => (SR.PureMT -> (a, SR.PureMT))+ -> (CR.SystemRandom -> (a, CR.SystemRandom))+ -> Eff r a+randomF f s = do+ old <- get+ let (val, new) = case old of+ (FastRandom r) -> second FastRandom (f r)+ (SecureRandom r) -> second SecureRandom (s r)+ put new+ return val+{-# INLINE randomF #-}
system-random-effect.cabal view
@@ -1,5 +1,5 @@ name: system-random-effect-version: 0.3.1+version: 0.4.0 synopsis: Random number generation for extensible effects. homepage: https://github.com/wowus/system-random-effect license: BSD3@@ -14,6 +14,7 @@ hs-source-dirs: src/ ghc-options: -Wall exposed-modules: System.Random.Effect+ other-modules: System.Random.Effect.Raw build-depends: base == 4.6.* , extensible-effects == 1.2.*@@ -21,6 +22,9 @@ , statistics == 0.10.* , vector == 0.10.* , vector-algorithms == 0.5.*+ , primitive == 0.5.*+ , crypto-api == 0.12.*+ , bytestring == 0.10.* default-language: Haskell2010 @@ -29,9 +33,10 @@ main-is: Test.hs hs-source-dirs: test/ - ghc-options: -rtsopts=all -threaded+ ghc-options: -rtsopts=all -threaded -O2 build-depends: base == 4.6.*+ , vector == 0.10.* , QuickCheck == 2.* , HUnit == 1.2.* , test-framework == 0.8.*@@ -47,11 +52,13 @@ hs-source-dirs: bench/ ghc-options: -Wall -rtsopts=all -threaded -fno-warn-orphans -O2 main-is: Bench.hs+ x-uses-tf: true build-depends: base == 4.6.* , deepseq == 1.3.* , criterion == 0.8.* , extensible-effects == 1.2.*+ , vector == 0.10.* , system-random-effect default-language: Haskell2010
test/Test.hs view
@@ -1,22 +1,28 @@ {-# LANGUAGE TypeOperators #-} {-# LANGUAGE FlexibleContexts #-}+{-# LANGUAGE BangPatterns #-} module Main ( main ) where +import Prelude hiding (all)+ import Control.Eff+import Control.Eff.Lift import Control.Eff.State.Strict import System.Random.Effect -import Control.Monad (void)+import Control.Applicative+import Control.Monad.ST+import Data.Vector ( Vector )+import qualified Data.Vector as V+import qualified Data.Vector.Mutable as MV import Data.Word-import Data.Typeable -import Test.Framework (defaultMain, testGroup)-import Test.Framework.Providers.HUnit+import Test.Framework ( defaultMain, Test ) import Test.Framework.Providers.QuickCheck2 -import Test.HUnit hiding (State) import Test.QuickCheck+import Test.QuickCheck.Property ( morallyDubiousIOProperty ) main :: IO () main = defaultMain tests@@ -24,6 +30,9 @@ runWithSeed :: Word64 -> Eff (State Random :> ()) a -> a runWithSeed seed = run . runRandomState (mkRandom seed) +runIOWithSeed :: Word64 -> Eff (State Random :> Lift IO :> ()) a -> IO a+runIOWithSeed seed = runLift . runRandomState (mkRandom seed)+ checkRange :: (Integer, Integer) -> Integer -> Bool checkRange (low, high) x = x >= low && x <= high@@ -36,33 +45,37 @@ in checkRange (low, high) . runWithSeed seed $ do uniformIntDist a b -testDiscreteDistributionInRange :: [Word64] -> Word64 -> Bool-testDiscreteDistributionInRange xs seed =+newtype DiscreteWeights = DW [Word64]+ deriving Show++instance Arbitrary DiscreteWeights where+ arbitrary = DW <$> suchThat (listOf arbitrary) ((> 0) . sum)+ shrink (DW xs) = map DW (shrink xs)++testDiscreteDistributionInRange :: DiscreteWeights -> Word64 -> Bool+testDiscreteDistributionInRange (DW xs) seed = let ddh = buildDDH xs minVal = 0 maxVal = length xs - 1- in sum xs == 0 ||- ((\x -> x >= minVal && x <= maxVal) . runWithSeed seed $- discreteDist ddh)+ in (\x -> x >= minVal && x <= maxVal) . runWithSeed seed $+ discreteDist ddh -testNoZeroDiscreteDistributionPick :: [Word64] -> Word64 -> Bool-testNoZeroDiscreteDistributionPick xs seed =+testNoZeroDiscreteDistributionPick :: DiscreteWeights -> Word64 -> Bool+testNoZeroDiscreteDistributionPick (DW xs) seed = let ddh = buildDDH xs- in sum xs == 0 ||- ((\x -> (xs !! x) /= 0) . runWithSeed seed $- discreteDist ddh)+ in (\x -> (xs !! x) /= 0) . runWithSeed seed $+ discreteDist ddh -testUnsafeThaw :: [Word64] -> Word64 -> Bool-testUnsafeThaw xs seed =+testUnsafeThaw :: DiscreteWeights -> Word64 -> Bool+testUnsafeThaw (DW xs) seed = let ddh = buildDDH xs- in sum xs == 0 ||- (runWithSeed seed $ do+ in runWithSeed seed $ do _ <- discreteDist ddh _ <- discreteDist ddh _ <- discreteDist ddh _ <- discreteDist ddh _ <- discreteDist ddh- return True)+ return True testUniformIntegralDist :: Integer -> Integer -> Word64 -> Bool testUniformIntegralDist a b seed =@@ -70,10 +83,88 @@ r2 = runWithSeed seed $ uniformIntegralDist a b in r1 == r2 +testKnuthShuffle :: [Int] -> Word64 -> Bool+testKnuthShuffle xs' seed =+ let xs = V.fromList xs'+ countIf f = V.length . V.filter f+ shuffled = runWithSeed seed (knuthShuffle xs)+ sameCount v1 v2 = V.all id+ $ V.map (\x -> countIf (== x) v1+ == countIf (== x) v2) v1+ in sameCount xs shuffled++testKnuthShuffleM :: [Int] -> Word64 -> IO Bool+testKnuthShuffleM xs' seed = do+ let xs = V.fromList xs'+ countIf f = V.length . V.filter f+ shuffled = do+ vs <- V.thaw xs+ runIOWithSeed seed (knuthShuffleM vs)+ V.freeze vs+ sameCount v1 v2 = V.all id+ $ V.map (\x -> countIf (== x) v1+ == countIf (== x) v2) v1+ shuf <- shuffled+ return (sameCount xs shuf)++testKnuthShuffleEquivalence :: [Int] -> Word64 -> IO Bool+testKnuthShuffleEquivalence xs seed = do+ let vs = V.fromList xs+ ks1 = runWithSeed seed (knuthShuffle vs)++ xs' <- V.thaw vs+ runIOWithSeed seed (knuthShuffleM xs')+ ks2 <- V.freeze xs'++ return (ks1 == ks2)++testSecureRandom :: Integer -> Integer -> IO Bool+testSecureRandom a b = do+ let low = min a b+ high = max a b++ runLift $ do+ rng <- mkSecureRandomIO+ return $ run $ runRandomState rng $+ checkRange (low, high) <$> uniformIntDist a b++(|>) :: b -> (b -> c) -> c+(|>) = flip ($)++histogram :: Vector Integer -> Vector Int+histogram v = runST $ do+ mv <- MV.replicate (fromIntegral (V.maximum v + 1)) 0+ V.forM_ v $ \i' -> do+ let i = (fromIntegral i') :: Int+ !x <- MV.read mv i+ MV.write mv i (x+1)+ V.unsafeFreeze mv++-- check if all uniformly distributed numbers are within 10% of the mean.+-- 10% was a number chosen arbitrarily.+simpleUniformIntDistTest :: Word64 -> Bool+simpleUniformIntDistTest seed =+ let nBuckets = 5 :: Int+ samplesPerBucket = 4000 :: Int+ nSamples = nBuckets * samplesPerBucket+ maxDelta = (fromIntegral samplesPerBucket) `div` 10+ nums = uniformIntDist 0 (fromIntegral (nBuckets - 1))+ |> V.replicateM nSamples+ |> runWithSeed seed+ hist = histogram nums+ in V.all (\x -> x >= samplesPerBucket - maxDelta+ && x <= samplesPerBucket + maxDelta) hist++tests :: [Test] tests = [ testProperty "random range" testUniformRandom , testProperty "discrete dist range" testDiscreteDistributionInRange , testProperty "no non-zero discrete dist pick" testNoZeroDiscreteDistributionPick , testProperty "unsafeThaw is okay to use" testUnsafeThaw , testProperty "testUniformIntegralDist == testUniformIntDist" testUniformIntegralDist+ , testProperty "knuth shuffle" testKnuthShuffle+ , testProperty "knuth shuffle (monadic)" (\xs seed -> morallyDubiousIOProperty $ testKnuthShuffleM xs seed)+ , testProperty "knuth shuffle equivalence" (\xs seed -> morallyDubiousIOProperty $ testKnuthShuffleEquivalence xs seed)+ , testProperty "secure random" (\a b -> morallyDubiousIOProperty $ testSecureRandom a b)+ , testProperty "uniformIntDist is uniform-ish" simpleUniformIntDistTest ]