Control.Monad.Mersenne.Random

A fast random number generator monad.

Random monad

newtype Rand a

A basic random monad, for generating random numbers from pure mersenne twisters.

Constructors

Rand

Fields

runRand :: PureMT -> R a

Instances

Monad Rand

runRandom :: Rand a -> PureMT -> (a, PureMT)

Run a random computation using the generator g, returning the result and the updated generator.

evalRandom :: Rand a -> PureMT -> a

Evaluate a random computation using the mersenne generator g. Note that the generator g is not returned, so there's no way to recover the updated version of g.

Efficient generators (Word64 is the primitive getter).

getBool :: Rand Bool

getInt :: Rand Int

Yield a new Int value from the generator.

getWord :: Rand Word

Yield a new Word value from the generator.

getInt64 :: Rand Int64

Yield a new Int64 value from the generator.

getWord64 :: Rand Word64

Yield a new Word64 value from the generator.

getDouble :: Rand Double

Yield a new 53-bit precise Double value from the generator.

Internals

data R a

The state of a random monad, optimized for performance.

Constructors

R !a !PureMT

An example from a user on Stack Overflow -- taking a random walk, and printing a histogram.

 {-# LANGUAGE BangPatterns #-}

 import System.Environment
 import Text.Printf
 import Control.Monad.Mersenne.Random
 import System.Random.Mersenne.Pure64

 main :: IO ()
 main = do
   (size:iters:_) <- fmap (map read) getArgs
   let start = take size $ repeat 0
   rnd <- newPureMT
   let end = flip evalRandom rnd $ mapM (iterateM iters randStep) start
   putStr . unlines $ histogram "%.2g" end 13
 
 {-# INLINE iterateM #-}
 iterateM n f x = go n x
     where
         go 0 !x = return x
         go n !x = f x >>= go (n-1)
 
 randStep :: Double -> Rand Double
 randStep x = do
     v <- getBool
     return $! if v then x+1 else x-1
 

 histogram :: String -> [Double] -> Int -> [String]
 histogram _ _ 0 = []
 histogram fmt xs bins =
     let xmin = minimum xs
         xmax = maximum xs
         bsize = (xmax - xmin) / (fromIntegral bins)
         bs = take bins $ zip [xmin,xmin+bsize..] [xmin+bsize,xmin+2*bsize..]
         counts :: [Int]
         counts = let cs = map count bs
                  in  (init cs) ++ [last cs + (length $ filter (==xmax) xs)]
     in  map (format (maximum counts)) $ zip bs counts
   where
     toD :: (Real b) => b -> Double
     toD = fromRational . toRational
     count (xmin, xmax) = length $ filter (\x -> x >= xmin && x < xmax) xs
     format :: Int -> ((Double,Double), Int) -> String
     format maxc ((lo,hi), c) = 
         let cscale = 50.0 / toD maxc
             hashes = take (round $ (toD c)*cscale) $ repeat '#'
             label  = let los = printf fmt lo
                          his = printf fmt hi
                          l   = los ++ " .. " ++ his
                          pad = take (20 - (length l)) $ repeat ' '
                      in  pad ++ l
         in  label ++ ": " ++ hashes
 

Compiling this:

 $ ghc -O2 --make B.hs

And running it:

 $ time E 300 5000
   -194.00 .. -164.46: 
   -164.46 .. -134.92: #
   -134.92 .. -105.38: ####
    -105.38 .. -75.85: ###########
     -75.85 .. -46.31: #########################
     -46.31 .. -16.77: ##################################################
      -16.77 .. 12.77: #################################################
       12.77 .. 42.31: ###########################################
       42.31 .. 71.85: ###########################
      71.85 .. 101.38: ################
     101.38 .. 130.92: #######
     130.92 .. 160.46: #####
     160.46 .. 190.00: #
 ./E 500 3000  0.03s user 0.00s system 96% cpu 0.035 total