Add contrib/mersenne-random-pure64

contrib/mersenne-random-pure64-0.2.2.0/LICENSE

@@ -0,0 +1,30 @@
+Copyright (c) Don Stewart <dons@galois.com> 2008
+
+All rights reserved.
+
+Redistribution and use in source and binary forms, with or without
+modification, are permitted provided that the following conditions
+are met:
+
+1. Redistributions of source code must retain the above copyright
+   notice, this list of conditions and the following disclaimer.
+
+2. Redistributions in binary form must reproduce the above copyright
+   notice, this list of conditions and the following disclaimer in the
+   documentation and/or other materials provided with the distribution.
+
+3. Neither the name of the author nor the names of his contributors
+   may be used to endorse or promote products derived from this software
+   without specific prior written permission.
+
+THIS SOFTWARE IS PROVIDED BY THE AUTHORS ``AS IS'' AND
+ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
+IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
+ARE DISCLAIMED.  IN NO EVENT SHALL THE AUTHORS OR CONTRIBUTORS BE LIABLE
+FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
+DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS
+OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
+HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
+LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
+OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
+SUCH DAMAGE.

contrib/mersenne-random-pure64-0.2.2.0/Setup.lhs

@@ -0,0 +1,3 @@
+#!/usr/bin/env runhaskell
+> import Distribution.Simple
+> main = defaultMain

contrib/mersenne-random-pure64-0.2.2.0/System/Random/Mersenne/Pure64.hs

@@ -0,0 +1,126 @@
+{-# LANGUAGE CPP, ForeignFunctionInterface #-}
+{-# OPTIONS_GHC -fno-warn-orphans #-}
+--------------------------------------------------------------------
+-- |
+-- Module     : System.Random.Mersenne.Pure64
+-- Copyright  : Copyright (c) 2008, Don Stewart <dons@galois.com>
+-- License    : BSD3
+-- Maintainer : Don Stewart <dons@galois.com>
+-- Stability  : experimental
+-- Portability: CPP, FFI
+-- Tested with: GHC 6.8.3
+--
+-- A purely functional binding 64 bit binding to the classic mersenne
+-- twister random number generator. This is more flexible than the
+-- impure 'mersenne-random' library, at the cost of being a bit slower.
+-- This generator is however, many times faster than System.Random,
+-- and yields high quality randoms with a long period.
+--
+-- This generator may be used with System.Random, however, that is
+-- likely to be slower than using it directly.
+--
+module System.Random.Mersenne.Pure64 (
+
+    -- * The random number generator
+    PureMT          -- abstract: RandomGen
+
+    -- * Introduction
+    , pureMT        -- :: Word64 -> PureMT
+    , newPureMT     -- :: IO PureMT
+
+    -- $instance
+
+    -- * Low level access to the generator
+
+    -- $notes
+    , randomInt     -- :: PureMT -> (Int   ,PureMT)
+    , randomWord    -- :: PureMT -> (Word  ,PureMT)
+    , randomInt64   -- :: PureMT -> (Int64 ,PureMT)
+    , randomWord64  -- :: PureMT -> (Word64,PureMT)
+    , randomDouble  -- :: PureMT -> (Double,PureMT)
+
+    ) where
+
+------------------------------------------------------------------------
+
+import System.Random.Mersenne.Pure64.MTBlock
+import System.Random.Mersenne.Pure64.Internal
+import System.Random
+import Data.Word
+import Data.Int
+import Data.Time.Clock
+import Data.Time.Calendar
+import System.CPUTime
+
+-- | Create a PureMT generator from a 'Word64' seed.
+pureMT :: Word64 -> PureMT
+pureMT = mkPureMT . seedBlock . fromIntegral
+
+#if !MIN_VERSION_time(1,6,0)
+diffTimeToPicoseconds :: DiffTime -> Integer
+diffTimeToPicoseconds d =
+    round (1000 * 1000 * 1000 * 1000 * d)
+#endif
+
+-- | Create a new PureMT generator, using the clocktime as the base for the seed.
+newPureMT :: IO PureMT
+newPureMT = do
+    ct <- getCPUTime
+    t  <- getCurrentTime
+    let seed = toModifiedJulianDay (utctDay t) + diffTimeToPicoseconds (utctDayTime t) + ct
+    return $ pureMT $ fromIntegral seed
+
+------------------------------------------------------------------------
+-- System.Random interface.
+
+-- $instance
+--
+-- Being purely functional, the PureMT generator is an instance of
+-- RandomGen. However, it doesn't support 'split' yet.
+
+instance RandomGen PureMT where
+   next  = randomInt
+   split = error "System.Random.Mersenne.Pure: unable to split the mersenne twister"
+
+------------------------------------------------------------------------
+-- Direct access to Int, Word and Double types
+
+-- | Yield a new 'Int' value from the generator, returning a new
+-- generator and that 'Int'. The full 64 bits will be used on a 64 bit machine.
+randomInt :: PureMT -> (Int,PureMT)
+randomInt g = (fromIntegral i, g')
+        where (i, g') = randomWord64 g
+{-# INLINE randomInt #-}
+
+-- | Yield a new 'Word' value from the generator, returning a new
+-- generator and that 'Word'.
+randomWord :: PureMT -> (Word,PureMT)
+randomWord g = (fromIntegral i, g')
+        where (i, g') = randomWord64 g
+{-# INLINE randomWord #-}
+
+-- | Yield a new 'Int64' value from the generator, returning a new
+-- generator and that 'Int64'.
+randomInt64 :: PureMT -> (Int64,PureMT)
+randomInt64 g = (fromIntegral i, g')
+        where (i, g') = randomWord64 g
+{-# INLINE randomInt64 #-}
+
+-- | Efficiently yield a new 53-bit precise 'Double' value, and a new generator.
+randomDouble :: PureMT -> (Double,PureMT)
+randomDouble g = (fromIntegral (i `div` 2048) / 9007199254740992, g')
+        where (i, g') = randomWord64 g
+{-# INLINE randomDouble #-}
+
+-- | Yield a new 'Word64' value from the generator, returning a new
+-- generator and that 'Word64'.
+randomWord64 :: PureMT -> (Word64,PureMT)
+randomWord64 (PureMT block i nxt) = (mixWord64 (block `lookupBlock` i), mt)
+  where
+    mt | i < blockLen-1 = PureMT block (i+1) nxt
+       | otherwise      = mkPureMT nxt
+{-# INLINE randomWord64 #-}
+
+-- create a new PureMT from an MTBlock
+mkPureMT :: MTBlock -> PureMT
+mkPureMT block = PureMT block 0 (nextBlock block)

contrib/mersenne-random-pure64-0.2.2.0/System/Random/Mersenne/Pure64/Base.hsc

@@ -0,0 +1,78 @@
+{-# LANGUAGE EmptyDataDecls, CPP, ForeignFunctionInterface #-}
+--------------------------------------------------------------------
+-- |
+-- Module     : System.Random.Mersenne.Pure64.Base
+-- Copyright  : Copyright (c) 2008, Don Stewart <dons@galois.com>
+-- License    : BSD3
+-- Maintainer : Don Stewart <dons@galois.com>
+-- Stability  : experimental
+-- Portability: CPP, FFI, EmptyDataDecls
+-- Tested with: GHC 6.8.3
+--
+-- A purely functional binding 64 bit binding to the classic mersenne
+-- twister random number generator. This is more flexible than the
+-- impure 'mersenne-random' library, at the cost of being a bit slower.
+-- This generator is however, many times faster than System.Random,
+-- and yields high quality randoms with a long period.
+--
+module System.Random.Mersenne.Pure64.Base where
+
+#include "mt19937-64-block.h"
+#include "mt19937-64.h"
+#include "mt19937-64-unsafe.h"
+
+import Foreign.C.Types
+import Foreign
+
+data MTState
+
+type UInt64 = CULLong
+
+------------------------------------------------------------------------
+-- pure version:
+
+foreign import ccall unsafe "init_genrand64"
+    c_init_genrand64  :: Ptr MTState -> UInt64 -> IO ()
+
+foreign import ccall unsafe "genrand64_int64"
+    c_genrand64_int64 :: Ptr MTState -> IO UInt64
+
+foreign import ccall unsafe "genrand64_real2"
+    c_genrand64_real2 :: Ptr MTState -> IO CDouble
+
+sizeof_MTState :: Int
+sizeof_MTState = (# const sizeof (struct mt_state_t) ) -- 2504 bytes
+
+------------------------------------------------------------------------
+-- block based version:
+foreign import ccall unsafe "mix_bits"
+    c_mix_word64 :: Word64 -> Word64
+
+foreign import ccall unsafe "seed_genrand64_block"
+    c_seed_genrand64_block :: Ptr a -> Word64 -> IO ()
+
+foreign import ccall unsafe "next_genrand64_block"
+    c_next_genrand64_block :: Ptr a -> Ptr a -> IO ()
+
+-- | length of an MT block
+blockLen :: Int
+blockLen = (# const NN )
+
+-- | size of an MT block, in bytes
+blockSize :: Int
+blockSize = (# const sizeof (mt_block_struct) )
+
+------------------------------------------------------------------------
+-- model: (for testing purposes)
+
+foreign import ccall unsafe "init_genrand64_unsafe"
+    c_init_genrand64_unsafe :: UInt64 -> IO ()
+
+foreign import ccall unsafe "genrand64_int64_unsafe"
+    c_genrand64_int64_unsafe :: IO UInt64
+
+foreign import ccall unsafe "genrand64_real2_unsafe"
+    c_genrand64_real2_unsafe :: IO CDouble
+
+foreign import ccall unsafe "string.h memcpy"
+    c_memcpy :: Ptr Word8 -> Ptr Word8 -> CSize -> IO (Ptr Word8)

contrib/mersenne-random-pure64-0.2.2.0/System/Random/Mersenne/Pure64/Internal.hs

@@ -0,0 +1,20 @@
+{-# LANGUAGE MagicHash #-}
+module System.Random.Mersenne.Pure64.Internal
+  ( PureMT(..)
+
+  , blockLen
+  , blockSize
+  , MTBlock(..)
+  ) where
+
+import GHC.Exts
+import System.Random.Mersenne.Pure64.Base
+
+-- | 'PureMT', a pure mersenne twister pseudo-random number generator
+--
+data PureMT  = PureMT {-# UNPACK #-} !MTBlock {-# UNPACK #-} !Int MTBlock
+
+instance Show PureMT where
+    show _ = show "<PureMT>"
+
+data MTBlock = MTBlock ByteArray#

contrib/mersenne-random-pure64-0.2.2.0/System/Random/Mersenne/Pure64/MTBlock.hs

@@ -0,0 +1,97 @@
+{-# LANGUAGE MagicHash, UnboxedTuples, BangPatterns #-}
+--------------------------------------------------------------------
+-- |
+-- Module     : System.Random.Mersenne.Pure64
+-- Copyright  : Copyright (c) 2008, Bertram Felgenhauer <int-e@gmx.de>
+-- License    : BSD3
+-- Maintainer : Don Stewart <dons@galois.com>
+-- Stability  : experimental
+-- Portability:
+-- Tested with: GHC 6.8.3
+--
+-- A purely functional binding 64 bit binding to the classic mersenne
+-- twister random number generator. This is more flexible than the
+-- impure 'mersenne-random' library, at the cost of being a bit slower.
+-- This generator is however, many times faster than System.Random,
+-- and yields high quality randoms with a long period.
+--
+module System.Random.Mersenne.Pure64.MTBlock (
+    -- * Block type
+    MTBlock,
+
+    -- * Block functions
+    seedBlock,
+    nextBlock,
+    lookupBlock,
+
+    -- * Misc functions
+    blockLen,
+    mixWord64,
+) where
+
+import GHC.Exts
+#if __GLASGOW_HASKELL__ >= 706
+import GHC.IO
+#else
+import GHC.IOBase
+#endif
+import GHC.Word
+import System.Random.Mersenne.Pure64.Base
+import System.Random.Mersenne.Pure64.Internal
+
+allocateBlock :: IO MTBlock
+allocateBlock =
+    IO $ \s0 -> case newPinnedByteArray# blockSize# s0 of
+        (# s1, b0 #) -> case unsafeFreezeByteArray# b0 s1 of
+            (# s2, b1 #) -> (# s2, MTBlock b1 #)
+  where
+    !(I# blockSize#) = blockSize
+
+blockAsPtr :: MTBlock -> Ptr a
+blockAsPtr (MTBlock b) = Ptr (byteArrayContents# b)
+
+-- | create a new MT block, seeded with the given Word64 value
+seedBlock :: Word64 -> MTBlock
+seedBlock seed = unsafeDupablePerformIO $ do
+    b <- allocateBlock
+    c_seed_genrand64_block (blockAsPtr b) seed
+    c_next_genrand64_block (blockAsPtr b) (blockAsPtr b)
+    touch b
+    return b
+{-# NOINLINE seedBlock #-}
+
+-- | step: create a new MTBlock buffer from the previous one
+nextBlock :: MTBlock -> MTBlock
+nextBlock b = unsafeDupablePerformIO $ do
+    new <- allocateBlock
+    c_next_genrand64_block (blockAsPtr b) (blockAsPtr new)
+    touch b
+    touch new
+    return new
+{-# NOINLINE nextBlock #-}
+
+-- stolen from GHC.ForeignPtr - make sure the argument is still alive.
+touch :: a -> IO ()
+touch r = IO $ \s0 -> case touch# r s0 of s1 -> (# s1, () #)
+
+-- | look up an element of an MT block
+lookupBlock :: MTBlock -> Int -> Word64
+lookupBlock (MTBlock b) (I# i) = W64# (indexWord64Array# b i)
+
+-- | MT's word mix function.
+--
+-- (MT applies this function to each Word64 from the buffer before returning it)
+mixWord64 :: Word64 -> Word64
+mixWord64 = c_mix_word64
+
+-- Alternative implementation - it's probably faster on 64 bit machines, but
+-- on Athlon XP it loses.
+{-
+mixWord64 (W64# x0) = let
+    W64# x1 = W64# x0 `xor` (W64# (x0 `uncheckedShiftRL64#` 28#) .&. 0x5555555555555555)
+    W64# x2 = W64# x1 `xor` (W64# (x1 `uncheckedShiftL64#` 17#) .&. 0x71D67FFFEDA60000)
+    W64# x3 = W64# x2 `xor` (W64# (x2 `uncheckedShiftL64#` 37#) .&. 0xFFF7EEE000000000)
+    W64# x4 = W64# x3 `xor` (W64# (x3 `uncheckedShiftRL64#` 43#))
+  in
+    W64# x4
+-}