From 2879ad46210ce13282ca0de7bc5ed4bafcd2a95b Mon Sep 17 00:00:00 2001 From: Alexey Kuleshevich Date: Wed, 12 Feb 2020 04:48:15 +0300 Subject: [PATCH] parent cfdfe6f09ad414fde5b855cc5f90207533413241 author Alexey Kuleshevich 1581472095 +0300 committer Leonhard Markert 1590493894 +0200 This patch is mostly backwards compatible. See "Breaking Changes" below for the full list of backwards incompatible changes. This patch fixes quality and performance issues, addresses additional miscellaneous issues, and introduces a monadic API. Issues addressed ================ Priority issues fixed in this patch: - Title: "The seeds generated by split are not independent" Link: https://github.com/haskell/random/issues/25 Fixed: changed algorithm to SplitMix, which provides a robust 'split' operation - Title: "Very low throughput" Link: https://github.com/haskell/random/issues/51 Fixed: see "Performance" below Additional issues addressed in this patch: - Title: "Add Random instances for tuples" Link: https://github.com/haskell/random/issues/26 Addressed: added 'Uniform' instances for up to 6-tuples - Title: "Add Random instance for Natural" Link: https://github.com/haskell/random/issues/44 Addressed: added 'UniformRange' instance for 'Natural' - Title: "incorrect distribution of randomR for floating-point numbers" Link: https://github.com/haskell/random/issues/53 Addressed: see "Regarding floating-point numbers" below - Title: "System/Random.hs:43:1: warning: [-Wtabs]" Link: https://github.com/haskell/random/issues/55 Fixed: no more tabs - Title: "Why does random for Float and Double produce exactly 24 or 53 bits?" Link: https://github.com/haskell/random/issues/58 Fixed: see "Regarding floating-point numbers" below - Title: "read :: StdGen fails for strings longer than 6" Link: https://github.com/haskell/random/issues/59 Addressed: 'StdGen' is no longer an instance of 'Read' Regarding floating-point numbers: with this patch, the relevant instances for 'Float' and 'Double' sample more bits than before but do not sample every possible representable value. The documentation now clearly spells out what this means for users. Quality (issue 25) ================== The algorithm [1] in version 1.1 of this library fails empirical PRNG tests when used to generate "split sequences" as proposed in [3]. SplitMix [2] passes the same tests. This patch changes 'StdGen' to use the SplitMix implementation provided by the splitmix package. Test batteries used: dieharder, TestU1, PractRand. [1]: P. L'Ecuyer, "Efficient and portable combined random number generators". https://doi.org/10.1145/62959.62969 [2]: G. L. Steele, D. Lea, C. H. Flood, "Fast splittable pseudorandom number generators". https://doi.org/10.1145/2714064.2660195 [3]: H. G. Schaathun, "Evaluation of splittable pseudo-random generators". https://doi.org/10.1017/S095679681500012X Performance (issue 51) ====================== The "improvement" column in the following table is a multiplier: the improvement for 'random' for type 'Float' is 1038, so this operation is 1038 times faster with this patch. | Name | Mean (1.1) | Mean (patch) | Improvement| | ----------------------- | ---------- | ------------ | ---------- | | pure/random/Float | 30 | 0.03 | 1038| | pure/random/Double | 52 | 0.03 | 1672| | pure/random/Integer | 43 | 0.33 | 131| | pure/uniform/Word8 | 14 | 0.03 | 422| | pure/uniform/Word16 | 13 | 0.03 | 375| | pure/uniform/Word32 | 21 | 0.03 | 594| | pure/uniform/Word64 | 42 | 0.03 | 1283| | pure/uniform/Word | 44 | 0.03 | 1491| | pure/uniform/Int8 | 15 | 0.03 | 511| | pure/uniform/Int16 | 15 | 0.03 | 507| | pure/uniform/Int32 | 22 | 0.03 | 749| | pure/uniform/Int64 | 44 | 0.03 | 1405| | pure/uniform/Int | 43 | 0.03 | 1512| | pure/uniform/Char | 17 | 0.49 | 35| | pure/uniform/Bool | 18 | 0.03 | 618| | pure/uniform/CChar | 14 | 0.03 | 485| | pure/uniform/CSChar | 14 | 0.03 | 455| | pure/uniform/CUChar | 13 | 0.03 | 448| | pure/uniform/CShort | 14 | 0.03 | 473| | pure/uniform/CUShort | 13 | 0.03 | 457| | pure/uniform/CInt | 21 | 0.03 | 737| | pure/uniform/CUInt | 21 | 0.03 | 742| | pure/uniform/CLong | 43 | 0.03 | 1544| | pure/uniform/CULong | 42 | 0.03 | 1460| | pure/uniform/CPtrdiff | 43 | 0.03 | 1494| | pure/uniform/CSize | 43 | 0.03 | 1475| | pure/uniform/CWchar | 22 | 0.03 | 785| | pure/uniform/CSigAtomic | 21 | 0.03 | 749| | pure/uniform/CLLong | 43 | 0.03 | 1554| | pure/uniform/CULLong | 42 | 0.03 | 1505| | pure/uniform/CIntPtr | 43 | 0.03 | 1476| | pure/uniform/CUIntPtr | 42 | 0.03 | 1463| | pure/uniform/CIntMax | 43 | 0.03 | 1535| | pure/uniform/CUIntMax | 42 | 0.03 | 1493| API changes =========== MonadRandom ----------- This patch adds a class 'MonadRandom': -- | 'MonadRandom' is an interface to monadic pseudo-random number generators. class Monad m => MonadRandom g s m | g m -> s where {-# MINIMAL freezeGen,thawGen,(uniformWord32|uniformWord64) #-} type Frozen g = (f :: Type) | f -> g freezeGen :: g s -> m (Frozen g) thawGen :: Frozen g -> m (g s) uniformWord32 :: g s -> m Word32 -- default implementation in terms of uniformWord64 uniformWord64 :: g s -> m Word64 -- default implementation in terms of uniformWord32 -- plus methods for other word sizes and for byte strings -- all have default implementations so the MINIMAL pragma holds Conceptually, in 'MonadRandom g s m', 'g s' is the type of the generator, 's' is the state type, and 'm' the underlying monad. Via the functional dependency 'g m -> s', the state type is determined by the generator and monad. 'Frozen' is the type of the generator's state "at rest". It is defined as an injective type family via 'f -> g', so there is no ambiguity as to which 'g' any 'Frozen g' belongs to. This definition is generic enough to accommodate, for example, the 'Gen' type from the 'mwc-random' package, which itself abstracts over the underlying primitive monad and state token. The documentation shows the full 'MonadRandom Gen' instance. Four 'MonadRandom' instances ("monadic adapters") are provided for pure generators to enable their use in monadic code. The documentation describes them in detail. 'Uniform' and 'UniformRange' ---------------------------- The 'Random' typeclass has conceptually been split into 'Uniform' and 'UniformRange'. The 'Random' typeclass is still included for backwards compatibility. 'Uniform' is for types where it is possible to sample from the type's entire domain; 'UniformRange' is for types where one can sample from a specified range. Breaking Changes ================ This patch introduces these breaking changes: * requires 'base >= 4.10' (GHC-8.2) * 'StdGen' is no longer an instance of 'Read' * 'randomIO' and 'randomRIO' where extracted from the 'Random' class into separate functions In addition, there may be import clashes with new functions, e.g. 'uniform' and 'uniformR'. Deprecations ============ This patch introduces 'genWord64', 'genWord32' and similar methods to the 'RandomGen' class. The significantly slower method 'next' and its companion 'genRange' are now deprecated. Co-authored-by: Alexey Kuleshevich Co-authored-by: idontgetoutmuch Co-authored-by: Leonhard Markert --- .darcs-boring | 5 - .gitignore | 15 +- .travis.yml | 173 ++- Benchmark/BinSearch.hs | 142 --- Benchmark/Makefile | 24 - Benchmark/SimpleRNGBench.hs | 322 ------ CHANGELOG.md | 110 +- DEVLOG.md | 196 ---- README.md | 19 +- Setup.hs | 27 + System/Random.hs | 609 ----------- bench-legacy/BinSearch.hs | 149 +++ bench-legacy/SimpleRNGBench.hs | 269 +++++ bench/Main.hs | 235 ++++ cabal.project | 2 + prologue.txt | 1 - random.cabal | 230 ++-- src/System/Random.hs | 467 ++++++++ src/System/Random/Internal.hs | 999 ++++++++++++++++++ src/System/Random/Monad.hs | 555 ++++++++++ stack-coveralls.yaml | 8 + stack-old.yaml | 25 + stack.yaml | 7 + test-legacy/Legacy.hs | 15 + .../Random1283.hs | 20 +- .../rangeTest.hs => test-legacy/RangeTest.hs | 73 +- {tests => test-legacy}/T7936.hs | 3 +- {tests => test-legacy}/TestRandomIOs.hs | 3 +- {tests => test-legacy}/TestRandomRs.hs | 5 +- test/Spec.hs | 195 ++++ test/Spec/Range.hs | 34 + test/Spec/Run.hs | 12 + test/doctests.hs | 12 + tests/Makefile | 14 - tests/all.T | 10 - tests/rangeTest.stdout | 67 -- tests/slowness.hs | 55 - 37 files changed, 3533 insertions(+), 1574 deletions(-) delete mode 100644 .darcs-boring delete mode 100644 Benchmark/BinSearch.hs delete mode 100644 Benchmark/Makefile delete mode 100644 Benchmark/SimpleRNGBench.hs delete mode 100644 DEVLOG.md delete mode 100644 System/Random.hs create mode 100644 bench-legacy/BinSearch.hs create mode 100644 bench-legacy/SimpleRNGBench.hs create mode 100644 bench/Main.hs create mode 100644 cabal.project delete mode 100644 prologue.txt create mode 100644 src/System/Random.hs create mode 100644 src/System/Random/Internal.hs create mode 100644 src/System/Random/Monad.hs create mode 100644 stack-coveralls.yaml create mode 100644 stack-old.yaml create mode 100644 stack.yaml create mode 100644 test-legacy/Legacy.hs rename tests/random1283.hs => test-legacy/Random1283.hs (59%) rename tests/rangeTest.hs => test-legacy/RangeTest.hs (82%) rename {tests => test-legacy}/T7936.hs (89%) rename {tests => test-legacy}/TestRandomIOs.hs (93%) rename {tests => test-legacy}/TestRandomRs.hs (89%) create mode 100644 test/Spec.hs create mode 100644 test/Spec/Range.hs create mode 100644 test/Spec/Run.hs create mode 100644 test/doctests.hs delete mode 100644 tests/Makefile delete mode 100644 tests/all.T delete mode 100644 tests/rangeTest.stdout delete mode 100644 tests/slowness.hs diff --git a/.darcs-boring b/.darcs-boring deleted file mode 100644 index 6682606a0..000000000 --- a/.darcs-boring +++ /dev/null @@ -1,5 +0,0 @@ -^dist(/|$) -^setup(/|$) -^GNUmakefile$ -^Makefile.local$ -^.depend(.bak)?$ diff --git a/.gitignore b/.gitignore index 41c6d8c61..89f91002f 100644 --- a/.gitignore +++ b/.gitignore @@ -1,12 +1,3 @@ -*~ - -Thumbs.db -.DS_Store - -GNUmakefile -dist-install/ -ghc.mk - -dist -.cabal-sandbox -cabal.sandbox.config +*.lock +.stack-work/ +cabal.project.local diff --git a/.travis.yml b/.travis.yml index 6a03bcb12..057d0c7e6 100644 --- a/.travis.yml +++ b/.travis.yml @@ -1,5 +1,168 @@ -language: haskell -ghc: - - 7.4 - - 7.6 - - 7.8 +############################################################################### +# To see the local changes made to this file, run: +# +# $ curl -sS https://raw.githubusercontent.com/commercialhaskell/stack/0cbcce0b1b89e5431121b745afdd1dfbd3ba8a43/doc/travis-complex.yml | diff -u - .travis.yml +############################################################################### + +# This is the complex Travis configuration, which is intended for use +# on open source libraries which need compatibility across multiple GHC +# versions, must work with cabal-install, and should be +# cross-platform. For more information and other options, see: +# +# https://docs.haskellstack.org/en/stable/travis_ci/ +# +# Copy these contents into the root directory of your Github project in a file +# named .travis.yml + +# Run jobs on Linux unless "os" is specified explicitly. +os: linux + +# Do not choose a language; we provide our own build tools. +language: generic + +# Caching so the next build will be fast too. +cache: + directories: + - $HOME/.ghc + - $HOME/.cabal + - $HOME/.stack + - $TRAVIS_BUILD_DIR/.stack-work + +# The different configurations we want to test. We have BUILD=cabal which uses +# cabal-install, and BUILD=stack which uses Stack. More documentation on each +# of those below. +# +# We set the compiler values here to tell Travis to use a different +# cache file per set of arguments. +# +# If you need to have different apt packages for each combination in the +# job matrix, you can use a line such as: +# addons: {apt: {packages: [libfcgi-dev,libgmp-dev]}} +jobs: + include: + # We grab the appropriate GHC and cabal-install versions from hvr's PPA. See: + # https://github.com/hvr/multi-ghc-travis + - env: BUILD=stack ARGS="--resolver lts-11.22 --stack-yaml stack-old.yaml" + compiler: ": #stack 8.2.2" + addons: {apt: {packages: [libgmp-dev]}} + + - env: BUILD=stack ARGS="--resolver lts-12.26 --stack-yaml stack-old.yaml" + compiler: ": #stack 8.4.4" + addons: {apt: {packages: [libgmp-dev]}} + + - env: BUILD=stack ARGS="--resolver lts-14.27" COVERALLS_STACK_YAML="stack-coveralls.yaml" + compiler: ": #stack 8.6.5" + addons: {apt: {packages: [libgmp-dev]}} + + - env: BUILD=stack ARGS="--resolver lts-15" + compiler: ": #stack 8.8.3" + addons: {apt: {packages: [libgmp-dev]}} + + - env: BUILD=stack ARGS="--resolver nightly" + compiler: ": #stack" + addons: {apt: {packages: [libgmp-dev]}} + + allow_failures: + - env: BUILD=cabal GHCVER=head CABALVER=head HAPPYVER=1.19.5 ALEXVER=3.1.7 + - env: BUILD=stack ARGS="--resolver nightly" + +before_install: +# Using compiler above sets CC to an invalid value, so unset it +- unset CC + +# We want to always allow newer versions of packages when building on GHC HEAD +- CABALARGS="" +- if [ "x$GHCVER" = "xhead" ]; then CABALARGS=--allow-newer; fi + +# Download and unpack the stack executable +- export PATH=/opt/ghc/$GHCVER/bin:/opt/cabal/$CABALVER/bin:$HOME/.local/bin:/opt/alex/$ALEXVER/bin:/opt/happy/$HAPPYVER/bin:$HOME/.cabal/bin:$PATH +- mkdir -p ~/.local/bin +- | + if [ `uname` = "Darwin" ] + then + travis_retry curl --insecure -L https://get.haskellstack.org/stable/osx-x86_64.tar.gz | tar xz --strip-components=1 --include '*/stack' -C ~/.local/bin + else + travis_retry curl -L https://get.haskellstack.org/stable/linux-x86_64.tar.gz | tar xz --wildcards --strip-components=1 -C ~/.local/bin '*/stack' + fi + + # Use the more reliable S3 mirror of Hackage + mkdir -p $HOME/.cabal + echo 'remote-repo: hackage.haskell.org:http://hackage.fpcomplete.com/' > $HOME/.cabal/config + echo 'remote-repo-cache: $HOME/.cabal/packages' >> $HOME/.cabal/config + + +install: +- echo "$(ghc --version) [$(ghc --print-project-git-commit-id 2> /dev/null || echo '?')]" +- if [ -f configure.ac ]; then autoreconf -i; fi +- | + set -ex + case "$BUILD" in + stack) + # Add in extra-deps for older snapshots, as necessary + # + # This is disabled by default, as relying on the solver like this can + # make builds unreliable. Instead, if you have this situation, it's + # recommended that you maintain multiple stack-lts-X.yaml files. + + #stack --no-terminal --install-ghc $ARGS test --bench --dry-run || ( \ + # stack --no-terminal $ARGS build cabal-install && \ + # stack --no-terminal $ARGS solver --update-config) + + # Build the dependencies + # stack --no-terminal --install-ghc $ARGS test --bench --only-dependencies + ;; + cabal) + cabal --version + travis_retry cabal update + + # Get the list of packages from the stack.yaml file. Note that + # this will also implicitly run hpack as necessary to generate + # the .cabal files needed by cabal-install. + PACKAGES=$(stack --install-ghc query locals | grep '^ *path' | sed 's@^ *path:@@') + + cabal install --only-dependencies --enable-tests --enable-benchmarks --force-reinstalls --ghc-options=-O0 --reorder-goals --max-backjumps=-1 $CABALARGS $PACKAGES + ;; + esac + set +ex + +script: +- | + set -ex + case "$BUILD" in + stack) + BUILD_ARGS="--bench --no-run-benchmarks --haddock --no-haddock-deps" + if [ -n "${COVERALLS_STACK_YAML}" ] && [ -n "${COVERALLS_REPO_TOKEN}" ]; then + stack $ARGS --stack-yaml="$COVERALLS_STACK_YAML" test --coverage $BUILD_ARGS + stack $ARGS --stack-yaml="$COVERALLS_STACK_YAML" hpc report --all + travis_retry curl -L https://github.com/lehins/stack-hpc-coveralls/releases/download/0.0.5.0/shc.tar.gz | tar xz shc + STACK_YAML="$COVERALLS_STACK_YAML" ./shc --repo-token=$COVERALLS_REPO_TOKEN combined custom + else + stack --no-terminal $ARGS test $BUILD_ARGS + fi + ;; + cabal) + cabal install --enable-tests --enable-benchmarks --force-reinstalls --ghc-options=-O0 --reorder-goals --max-backjumps=-1 $CABALARGS $PACKAGES + + ORIGDIR=$(pwd) + for dir in $PACKAGES + do + cd $dir + cabal check || [ "$CABALVER" == "1.16" ] + cabal sdist + PKGVER=$(cabal info . | awk '{print $2;exit}') + SRC_TGZ=$PKGVER.tar.gz + cd dist + tar zxfv "$SRC_TGZ" + cd "$PKGVER" + cabal configure --enable-tests --ghc-options -O0 + cabal build + if [ "$CABALVER" = "1.16" ] || [ "$CABALVER" = "1.18" ]; then + cabal test + else + cabal test --show-details=streaming --log=/dev/stdout + fi + cd $ORIGDIR + done + ;; + esac + set +ex diff --git a/Benchmark/BinSearch.hs b/Benchmark/BinSearch.hs deleted file mode 100644 index f61164855..000000000 --- a/Benchmark/BinSearch.hs +++ /dev/null @@ -1,142 +0,0 @@ - -{- - Binary search over benchmark input sizes. - - There are many good ways to measure the time it takes to perform a - certain computation on a certain input. However, frequently, it's - challenging to pick the right input size for all platforms and all - compilataion modes. - - Sometimes for linear-complexity benchmarks it is better to measure - /throughput/, i.e. elements processed per second. That is, fixing - the time of execution and measuring the amount of work done (rather - than the reverse). This library provides a simple way to search for - an appropriate input size that results in the desired execution time. - - An alternative approach is to kill the computation after a certain - amount of time and observe how much work it has completed. - -} -module BinSearch - ( - binSearch - ) -where - -import Control.Monad -import Data.Time.Clock -- Not in 6.10 -import Data.List -import System.IO -import Prelude hiding (min,max,log) - - - --- | Binary search for the number of inputs to a computation that --- results in a specified amount of execution time in seconds. For example: --- --- > binSearch verbose N (min,max) kernel --- --- ... will find the right input size that results in a time --- between min and max, then it will then run for N trials and --- return the median (input,time-in-seconds) pair. -binSearch :: Bool -> Integer -> (Double,Double) -> (Integer -> IO ()) -> IO (Integer, Double) -binSearch verbose trials (min,max) kernel = - do - when(verbose)$ putStrLn$ "[binsearch] Binary search for input size resulting in time in range "++ show (min,max) - - let desired_exec_length = 1.0 - good_trial t = (toRational t <= toRational max) && (toRational t >= toRational min) - - -- At some point we must give up... - loop n | n > ((2::Integer) ^ (100::Integer)) = error "ERROR binSearch: This function doesn't seem to scale in proportion to its last argument." - - -- Not allowed to have "0" size input, bump it back to one: - loop 0 = loop 1 - - loop n = - do - when(verbose)$ putStr$ "[binsearch:"++ show n ++ "] " - time <- timeit$ kernel n - when(verbose)$ putStrLn$ "Time consumed: "++ show time - let rate = fromIntegral n / time - - -- [2010.06.09] Introducing a small fudge factor to help our guess get over the line: - let initial_fudge_factor = 1.10 - fudge_factor = 1.01 -- Even in the steady state we fudge a little - guess = desired_exec_length * rate - - -- TODO: We should keep more history here so that we don't re-explore input space we have already explored. - -- This is a balancing act because of randomness in execution time. - - if good_trial time - then do - when(verbose)$ putStrLn$ "[binsearch] Time in range. LOCKING input size and performing remaining trials." - print_trial 1 n time - lockin (trials-1) n [time] - - -- Here we're still in the doubling phase: - else if time < 0.100 - then loop (2*n) - - else do when(verbose)$ - putStrLn$ "[binsearch] Estimated rate to be " - ++show (round rate::Integer)++" per second. Trying to scale up..." - - -- Here we've exited the doubling phase, but we're making our first guess as to how big a real execution should be: - if time > 0.100 && time < 0.33 * desired_exec_length - then do when(verbose)$ putStrLn$ "[binsearch] (Fudging first guess a little bit extra)" - loop (round$ guess * initial_fudge_factor) - else loop (round$ guess * fudge_factor) - - -- Termination condition: Done with all trials. - lockin 0 n log = do when(verbose)$ putStrLn$ "[binsearch] Time-per-unit for all trials: "++ - (concat $ intersperse " " (map (show . (/ toDouble n) . toDouble) $ sort log)) - return (n, log !! ((length log) `quot` 2)) -- Take the median - - lockin trials_left n log = - do when(verbose)$ putStrLn$ "[binsearch]------------------------------------------------------------" - time <- timeit$ kernel n - -- hFlush stdout - print_trial (trials - trials_left +1 ) n time - -- when(verbose)$ hFlush stdout - lockin (trials_left - 1) n (time : log) - - print_trial :: Integer -> Integer -> NominalDiffTime -> IO () - print_trial trialnum n time = - let rate = fromIntegral n / time - timeperunit = time / fromIntegral n - in - when(verbose)$ putStrLn$ "[binsearch] TRIAL: "++show trialnum ++ - " secPerUnit: "++ showTime timeperunit ++ - " ratePerSec: "++ show (rate) ++ - " seconds: "++showTime time - - - - (n,t) <- loop 1 - return (n, fromRational$ toRational t) - -showTime :: NominalDiffTime -> String -showTime t = show ((fromRational $ toRational t) :: Double) - -toDouble :: Real a => a -> Double -toDouble = fromRational . toRational - - --- Could use cycle counters here.... but the point of this is to time --- things on the order of a second. -timeit :: IO () -> IO NominalDiffTime -timeit io = - do strt <- getCurrentTime - io - end <- getCurrentTime - return (diffUTCTime end strt) -{- -test :: IO (Integer,Double) -test = - binSearch True 3 (1.0, 1.05) - (\n -> - do v <- newIORef 0 - forM_ [1..n] $ \i -> do - old <- readIORef v - writeIORef v (old+i)) --} diff --git a/Benchmark/Makefile b/Benchmark/Makefile deleted file mode 100644 index 8a84e6479..000000000 --- a/Benchmark/Makefile +++ /dev/null @@ -1,24 +0,0 @@ - - -#OPTS= -O2 -Wall -XCPP -OPTS= -O2 -Wall -XCPP -Werror - -all: lib bench - -lib: - (cd .. && ghc $(OPTS) --make System/Random.hs) - -bench: - ghc $(OPTS) -i.. --make SimpleRNGBench.hs - -# When benchmarking against other packages installed via cabal it is -# necessary to IGNORE the System/Random.hs files in the current directory. -# (Otherwise instances of RandomGen are confused.) -bench2: - ghc $(OPTS) -DTEST_COMPETITORS --make SimpleRNGBench.hs - -clean: - rm -f *.o *.hi SimpleRNGBench -# cabal clean -# (cd Benchmark/ && rm -f *.o *.hi SimpleRNGBench) -# (cd System/ && rm -f *.o *.hi SimpleRNGBench) diff --git a/Benchmark/SimpleRNGBench.hs b/Benchmark/SimpleRNGBench.hs deleted file mode 100644 index c25b75d47..000000000 --- a/Benchmark/SimpleRNGBench.hs +++ /dev/null @@ -1,322 +0,0 @@ -{-# LANGUAGE BangPatterns, ScopedTypeVariables, ForeignFunctionInterface #-} -{-# OPTIONS_GHC -fwarn-unused-imports #-} - --- | A simple script to do some very basic timing of the RNGs. - -module Main where - -import System.Exit (exitSuccess, exitFailure) -import System.Environment -import System.Random -import System.CPUTime (getCPUTime) -import System.CPUTime.Rdtsc -import System.Console.GetOpt - -import GHC.Conc -import Control.Concurrent -import Control.Monad -import Control.Exception - -import Data.IORef -import Data.Word -import Data.List hiding (last,sum) -import Data.Int -import Data.List.Split hiding (split) -import Text.Printf - -import Foreign.Ptr -import Foreign.C.Types -import Foreign.ForeignPtr -import Foreign.Storable (peek,poke) - -import Prelude hiding (last,sum) -import BinSearch - -#ifdef TEST_COMPETITORS -import System.Random.Mersenne.Pure64 -import System.Random.MWC -import Control.Monad.Primitive --- import System.IO.Unsafe -import GHC.IO -#endif - ----------------------------------------------------------------------------------------------------- --- Miscellaneous helpers: - --- Readable large integer printing: -commaint :: Integral a => a -> String -commaint n = - reverse $ concat $ - intersperse "," $ - chunk 3 $ - reverse (show n) - -padleft :: Int -> String -> String -padleft n str | length str >= n = str -padleft n str | otherwise = take (n - length str) (repeat ' ') ++ str - -padright :: Int -> String -> String -padright n str | length str >= n = str -padright n str | otherwise = str ++ take (n - length str) (repeat ' ') - -fmt_num :: (RealFrac a, PrintfArg a) => a -> String -fmt_num n = if n < 100 - then printf "%.2f" n - else commaint (round n :: Integer) - - --- Measure clock frequency, spinning rather than sleeping to try to --- stay on the same core. -measureFreq :: IO Int64 -measureFreq = do - let second = 1000 * 1000 * 1000 * 1000 -- picoseconds are annoying - t1 <- rdtsc - start <- getCPUTime - let loop !n !last = - do t2 <- rdtsc - when (t2 < last) $ - putStrLn$ "COUNTERS WRAPPED "++ show (last,t2) - cput <- getCPUTime - if (cput - start < second) - then loop (n+1) t2 - else return (n,t2) - (n,t2) <- loop 0 t1 - putStrLn$ " Approx getCPUTime calls per second: "++ commaint (n::Int64) - when (t2 < t1) $ - putStrLn$ "WARNING: rdtsc not monotonically increasing, first "++show t1++" then "++show t2++" on the same OS thread" - - return$ fromIntegral (t2 - t1) - ----------------------------------------------------------------------------------------------------- - --- Test overheads without actually generating any random numbers: -data NoopRNG = NoopRNG -instance RandomGen NoopRNG where - next g = (0,g) -#ifdef ENABLE_SPLITTABLEGEN - genRange _ = (0,0) -instance SplittableGen NoopRNG where -#endif - split g = (g,g) - --- An RNG generating only 0 or 1: -data BinRNG = BinRNG StdGen -instance RandomGen BinRNG where - next (BinRNG g) = (x `mod` 2, BinRNG g') - where (x,g') = next g -#ifdef ENABLE_SPLITTABLEGEN - genRange _ = (0,1) -instance SplittableGen BinRNG where -#endif - split (BinRNG g) = (BinRNG g1, BinRNG g2) - where (g1,g2) = split g - - - -#ifdef TEST_COMPETITORS -data MWCRNG = MWCRNG (Gen (PrimState IO)) --- data MWCRNG = MWCRNG GenIO -instance RandomGen MWCRNG where - -- For testing purposes we hack this to be non-monadic: --- next g@(MWCRNG gen) = unsafePerformIO $ - next g@(MWCRNG gen) = unsafeDupablePerformIO $ - do v <- uniform gen - return (v, g) -#endif - ----------------------------------------------------------------------------------------------------- --- Drivers to get random numbers repeatedly. - -type Kern = Int -> Ptr Int -> IO () - --- [2011.01.28] Changing this to take "count" and "accumulator ptr" as arguments: --- foreign import ccall "cbits/c_test.c" blast_rands :: Kern --- foreign import ccall "cbits/c_test.c" store_loop :: Kern --- foreign import ccall unsafe "stdlib.hs" rand :: IO Int - -{-# INLINE timeit #-} -timeit :: (Random a, RandomGen g) => - Int -> Int64 -> String -> g -> (g -> (a,g)) -> IO () -timeit numthreads freq msg gen nxt = - do - counters <- forM [1..numthreads] (const$ newIORef (1::Int64)) - tids <- forM counters $ \counter -> - forkIO $ infloop counter (nxt gen) - threadDelay (1000*1000) -- One second - mapM_ killThread tids - - finals <- mapM readIORef counters - let mean :: Double = fromIntegral (foldl1 (+) finals) / fromIntegral numthreads - cycles_per :: Double = fromIntegral freq / mean - printResult (round mean :: Int64) msg cycles_per - - where - infloop !counter !(!_,!g) = - do incr counter - infloop counter (nxt g) - - incr !counter = - do -- modifyIORef counter (+1) -- Not strict enough! - c <- readIORef counter - let c' = c+1 - _ <- evaluate c' - writeIORef counter c' - - --- This function times an IO function on one or more threads. Rather --- than running a fixed number of iterations, it uses a binary search --- to find out how many iterations can be completed in a second. -timeit_foreign :: Int -> Int64 -> String -> (Int -> Ptr Int -> IO ()) -> IO Int64 -timeit_foreign numthreads freq msg ffn = do - ptr :: ForeignPtr Int <- mallocForeignPtr - - let kern = if numthreads == 1 - then ffn - else replicate_kernel numthreads ffn - wrapped n = withForeignPtr ptr (kern$ fromIntegral n) - (n,t) <- binSearch False 1 (1.0, 1.05) wrapped - - let total_per_second = round $ fromIntegral n * (1 / t) - cycles_per = fromIntegral freq * t / fromIntegral n - printResult total_per_second msg cycles_per - return total_per_second - - where - -- This lifts a C kernel to operate simultaneously on N threads. - replicate_kernel :: Int -> Kern -> Kern - replicate_kernel nthreads kern n ptr = do - ptrs <- forM [1..nthreads] - (const mallocForeignPtr) - tmpchan <- newChan - -- let childwork = ceiling$ fromIntegral n / fromIntegral nthreads - let childwork = n -- Keep it the same.. interested in per-thread throughput. - -- Fork/join pattern: - _ <- forM ptrs $ \pt -> forkIO $ - withForeignPtr pt $ \p -> do - kern (fromIntegral childwork) p - result <- peek p - writeChan tmpchan result - - results <- forM [1..nthreads] $ \_ -> - readChan tmpchan - -- Meaningless semantics here... sum the child ptrs and write to the input one: - poke ptr (foldl1 (+) results) - return () - - -printResult :: Int64 -> String -> Double -> IO () -printResult total msg cycles_per = - putStrLn$ " "++ padleft 11 (commaint total) ++" randoms generated "++ padright 27 ("["++msg++"]") ++" ~ " - ++ fmt_num cycles_per ++" cycles/int" - ----------------------------------------------------------------------------------------------------- --- Main Script - -data Flag = NoC | Help - deriving (Show, Eq) - -options :: [OptDescr Flag] -options = - [ Option ['h'] ["help"] (NoArg Help) "print program help" - , Option [] ["noC"] (NoArg NoC) "omit C benchmarks, haskell only" - ] - -main :: IO () -main = do - argv <- getArgs - let (opts,_,other) = getOpt Permute options argv - - when (not$ null other) $ do - putStrLn$ "ERROR: Unrecognized options: " - mapM_ putStr other - exitFailure - - when (Help `elem` opts) $ do - putStr$ usageInfo "Benchmark random number generation" options - exitSuccess - - putStrLn$ "\nHow many random numbers can we generate in a second on one thread?" - - t1 <- rdtsc - t2 <- rdtsc - putStrLn (" Cost of rdtsc (ffi call): " ++ show (t2 - t1)) - - freq <- measureFreq - putStrLn$ " Approx clock frequency: " ++ commaint freq - - let - randInt = random :: RandomGen g => g -> (Int,g) - randWord16 = random :: RandomGen g => g -> (Word16,g) - randFloat = random :: RandomGen g => g -> (Float,g) - randCFloat = random :: RandomGen g => g -> (CFloat,g) - randDouble = random :: RandomGen g => g -> (Double,g) - randCDouble = random :: RandomGen g => g -> (CDouble,g) - randInteger = random :: RandomGen g => g -> (Integer,g) - randBool = random :: RandomGen g => g -> (Bool,g) - randChar = random :: RandomGen g => g -> (Char,g) - - gen = mkStdGen 23852358661234 - gamut th = do - putStrLn$ " First, timing System.Random.next:" - timeit th freq "constant zero gen" NoopRNG next - timeit th freq "System.Random stdGen/next" gen next - - putStrLn$ "\n Second, timing System.Random.random at different types:" - timeit th freq "System.Random Ints" gen randInt - timeit th freq "System.Random Word16" gen randWord16 - timeit th freq "System.Random Floats" gen randFloat - timeit th freq "System.Random CFloats" gen randCFloat - timeit th freq "System.Random Doubles" gen randDouble - timeit th freq "System.Random CDoubles" gen randCDouble - timeit th freq "System.Random Integers" gen randInteger - timeit th freq "System.Random Bools" gen randBool - timeit th freq "System.Random Chars" gen randChar - -#ifdef TEST_COMPETITORS - putStrLn$ "\n Next test other RNG packages on Hackage:" - let gen_mt = pureMT 39852 - randInt2 = random :: RandomGen g => g -> (Int,g) - randFloat2 = random :: RandomGen g => g -> (Float,g) - timeit th freq "System.Random.Mersenne.Pure64 next" gen_mt next - timeit th freq "System.Random.Mersenne.Pure64 Ints" gen_mt randInt2 - timeit th freq "System.Random.Mersenne.Pure64 Floats" gen_mt randFloat2 - --- gen_mwc <- create - withSystemRandom $ \ gen_mwc -> do - let randInt3 = random :: RandomGen g => g -> (Int,g) - randFloat3 = random :: RandomGen g => g -> (Float,g) - - timeit th freq "System.Random.MWC next" (MWCRNG gen_mwc) next - timeit th freq "System.Random.MWC Ints" (MWCRNG gen_mwc) randInt3 - timeit th freq "System.Random.MWC Floats" (MWCRNG gen_mwc) randFloat3 - -#endif - - putStrLn$ "\n Next timing range-restricted System.Random.randomR:" - timeit th freq "System.Random Ints" gen (randomR (-100, 100::Int)) - timeit th freq "System.Random Word16s" gen (randomR (-100, 100::Word16)) - timeit th freq "System.Random Floats" gen (randomR (-100, 100::Float)) - timeit th freq "System.Random CFloats" gen (randomR (-100, 100::CFloat)) - timeit th freq "System.Random Doubles" gen (randomR (-100, 100::Double)) - timeit th freq "System.Random CDoubles" gen (randomR (-100, 100::CDouble)) - timeit th freq "System.Random Integers" gen (randomR (-100, 100::Integer)) - timeit th freq "System.Random Bools" gen (randomR (False, True::Bool)) - timeit th freq "System.Random Chars" gen (randomR ('a', 'z')) - timeit th freq "System.Random BIG Integers" gen (randomR (0, (2::Integer) ^ (5000::Int))) - - -- when (not$ NoC `elem` opts) $ do - -- putStrLn$ " Comparison to C's rand():" - -- timeit_foreign th freq "ptr store in C loop" store_loop - -- timeit_foreign th freq "rand/store in C loop" blast_rands - -- timeit_foreign th freq "rand in Haskell loop" (\n ptr -> forM_ [1..n]$ \_ -> rand ) - -- timeit_foreign th freq "rand/store in Haskell loop" (\n ptr -> forM_ [1..n]$ \_ -> do n <- rand; poke ptr n ) - -- return () - - -- Test with 1 thread and numCapabilities threads: - gamut 1 - when (numCapabilities > 1) $ do - putStrLn$ "\nNow "++ show numCapabilities ++" threads, reporting mean randoms-per-second-per-thread:" - gamut numCapabilities - return () - - putStrLn$ "Finished." diff --git a/CHANGELOG.md b/CHANGELOG.md index 15c882af9..0e6fd8f09 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -1,3 +1,112 @@ +# 1.2 + +1. Breaking change which mostly maintains backwards compatibility, see + "Breaking Changes" below. +2. Support for monadic generators e.g. [mwc-random](https://hackage.haskell.org/package/mwc-random). +3. Monadic adapters for pure generators (providing a uniform monadic + interface to pure and monadic generators). +4. Faster in all cases except one by more than x18 (N.B. x18 not 18%) and + some cases (depending on the type) faster by more than x1000 - see + below for benchmarks. +5. Passes a large number of random number test suites: + * [dieharder](http://webhome.phy.duke.edu/~rgb/General/dieharder.php "venerable") + * [TestU01 (SmallCrush, Crush, BigCrush)](http://simul.iro.umontreal.ca/testu01/tu01.html "venerable") + * [PractRand](http://pracrand.sourceforge.net/ "active") + * [gjrand](http://gjrand.sourceforge.net/ "active") + * See [random-quality](https://github.com/tweag/random-quality) + for details on how to do this yourself. +6. Better quality split as judged by these + [tests](https://www.cambridge.org/core/journals/journal-of-functional-programming/article/evaluation-of-splittable-pseudorandom-generators/3EBAA9F14939C5BB5560E32D1A132637). Again + see [random-quality](https://github.com/tweag/random-quality) for + details on how to do this yourself. +7. Unbiased generation of ranges. +8. Updated tests and benchmarks. +9. [Continuous integration](https://travis-ci.org/github/idontgetoutmuch/random). +10. Fully documented - for more details see the [haddock](https://htmlpreview.github.io/?https://github.com/idontgetoutmuch/random/blob/release-notes/docs/System-Random.html). + +### Breaking Changes + +Version 1.2 introduces these breaking changes: + +* requires `base >= 4.10` (GHC-8.2) +* `StdGen` is no longer an instance of `Read` +* `randomIO` and `randomRIO` where extracted from the `Random` class into + separate functions + +In addition, there may be import clashes with new functions, e.g. `uniform` and +`uniformR`. + +### Deprecations + +Version 1.2 introduces `genWord64`, `genWord32` and similar methods to the +`RandomGen` class. The significantly slower method `next` and its companion +`genRange` are now deprecated. + +### Issues Addressed + + Issue Number | Description | Comment +--------------|-------------|-------- + [25](https://github.com/haskell/random/issues/25) | The seeds generated by split are not independent | Fixed: changed algorithm to SplitMix, which provides a robust split operation + [26](https://github.com/haskell/random/issues/26) | Add Random instances for tuples | Addressed: added `Uniform` instances for up to 6-tuples + [44](https://github.com/haskell/random/issues/44) | Add Random instance for Natural | Addressed: added UniformRange instance for Natural + [51](https://github.com/haskell/random/issues/51) | Very low throughput | Fixed: see benchmarks below + [53](https://github.com/haskell/random/issues/53) | incorrect distribution of randomR for floating-point numbers | (\*) + [55](https://github.com/haskell/random/issues/55) | System/Random.hs:43:1: warning: [-Wtabs] | Fixed: No more tabs + [58](https://github.com/haskell/random/issues/58) | Why does random for Float and Double produce exactly 24 or 53 bits? | (\*) + [59](https://github.com/haskell/random/issues/59) | read :: StdGen fails for strings longer than 6 | Addressed: StdGen is no longer an instance of Read + +#### Comments + +(\*) 1.2 samples more bits but does not sample every `Float` or +`Double`. There are methods to do this but they have some downsides; +see [here](https://github.com/idontgetoutmuch/random/issues/105) for a +fuller discussion. + +## Benchmarks + +Here are some benchmarks run on a 3.1 GHz Intel Core i7. The full +benchmarks can be run using e.g. `stack bench`. The benchmarks are +measured in milliseconds per 100,000 generations. In some cases, the +performance is over x1000 times better; the minimum performance +increase for the types listed below is more than x36. + + Name | 1.1 Mean | 1.2 Mean +------------|----------|---------- + Float | 27.819 | 0.305 + Double | 50.644 | 0.328 + Integer | 42.332 | 0.332 + Word8 | 12.591 | 0.028 + Word16 | 12.726 | 0.028 + Word32 | 20.429 | 0.027 + Word64 | 42.299 | 0.028 + Word | 40.739 | 0.027 + Int8 | 13.479 | 0.027 + Int16 | 13.218 | 0.027 + Int32 | 20.562 | 0.027 + Int64 | 46.513 | 0.029 + Int | 43.847 | 0.028 + Char | 17.009 | 0.462 + Bool | 17.542 | 0.027 + CChar | 13.276 | 0.027 + CSChar | 13.287 | 0.027 + CUChar | 13.409 | 0.027 + CShort | 13.158 | 0.027 + CUShort | 12.865 | 0.027 + CInt | 20.705 | 0.028 + CUInt | 19.895 | 0.027 + CLong | 41.679 | 0.027 + CULong | 40.806 | 0.027 + CPtrdiff | 41.878 | 0.027 + CSize | 40.739 | 0.027 + CWchar | 20.718 | 0.027 + CSigAtomic | 20.768 | 0.029 + CLLong | 42.011 | 0.028 + CULLong | 41.428 | 0.027 + CIntPtr | 45.385 | 0.027 + CUIntPtr | 40.797 | 0.027 + CIntMax | 41.778 | 0.027 + CUIntMax | 40.467 | 0.027 + # 1.1 * breaking change to `randomIValInteger` to improve RNG quality and performance see https://github.com/haskell/random/pull/4 and @@ -23,4 +132,3 @@ bump for bug fixes, # 1.0.0.4 bumped version for float/double range bugfix - diff --git a/DEVLOG.md b/DEVLOG.md deleted file mode 100644 index 6e0b28dc4..000000000 --- a/DEVLOG.md +++ /dev/null @@ -1,196 +0,0 @@ - - -DEVLOG: A collection of notes accumulated during development. -============================================================= - - -[2011.06.24] (transient) Regression in stdGen performance. ----------------------------------------------------------- - -I just added a simple benchmark to make sure that whatever fix I -introduce for trac ticket #5133 does not regress performance. Yet in -doing so I discovered that I'm getting much worse performance out of -rev 130e421e912d than I'm seeing in my installed random-1.0.0.3 package. - -Current version: - How many random numbers can we generate in a second on one thread? - Cost of rdtsc (ffi call): 100 - Approx getCPUTime calls per second: 234,553 - Approx clock frequency: 3,335,220,196 - First, timing with System.Random interface: - 68,550,189 random ints generated [constant zero gen] ~ 48.65 cycles/int - 900,889 random ints generated [System.Random stdGen] ~ 3,702 cycles/int - -random-1.0.0.3 version: - How many random numbers can we generate in a second on one thread? - Cost of rdtsc (ffi call): 75 - Approx getCPUTime calls per second: 215,332 - Approx clock frequency: 3,334,964,738 - First, timing with System.Random interface: - 71,683,748 random ints generated [constant zero gen] ~ 46.52 cycles/int - 13,609,559 random ints generated [System.Random stdGen] ~ 245 cycles/int - -A >13X difference!! -Both are compiled with the same options. The only difference is which -System.Random is used. - -When did the regression occur? - - * e059ed15172585310f9c -- 10/13/2010 perf still good - * 6c43f80f48178ac617 -- SplittableGen introduced, still good perf - * 130e421e912d394653a4 -- most recent, bad performance - -Ok... this is very odd. It was a heisenbug becuase it's disappeared -now! I'll leave this note here to help remember to look for it in the -future. - -Ryan - - -[2011.06.24] Timing non-int types ---------------------------------- - -The results are highly uneven: - - Cost of rdtsc (ffi call): 84 - Approx getCPUTime calls per second: 220,674 - Approx clock frequency: 3,336,127,273 - First, timing with System.Random interface: - 112,976,933 randoms generated [constant zero gen] ~ 29.53 cycles/int - 14,415,176 randoms generated [System.Random stdGen] ~ 231 cycles/int - 70,751 randoms generated [System.Random Floats] ~ 47,153 cycles/int - 70,685 randoms generated [System.Random CFloats] ~ 47,197 cycles/int - 2,511,635 randoms generated [System.Random Doubles] ~ 1,328 cycles/int - 70,494 randoms generated [System.Random CDoubles] ~ 47,325 cycles/int - 858,012 randoms generated [System.Random Integers] ~ 3,888 cycles/int - 4,756,213 randoms generated [System.Random Bools] ~ 701 cycles/int - -As you can see, all the types that use the generic randomIvalFrac / -randomFrac definitions perform badly. What's more, the above results -INCLUDE an attempt to inline: - - {-# INLINE randomIvalFrac #-} - {-# INLINE randomFrac #-} - {-# INLINE randomIvalDouble #-} - -After reimplementing random/Float these are the new results: - - Cost of rdtsc (ffi call): 100 - Approx getCPUTime calls per second: 200,582 - Approx clock frequency: 3,334,891,942 - First, timing with System.Random interface: - 105,266,949 randoms generated [constant zero gen] ~ 31.68 cycles/int - 13,593,392 randoms generated [System.Random stdGen] ~ 245 cycles/int - 10,962,597 randoms generated [System.Random Floats] ~ 304 cycles/int - 11,926,573 randoms generated [System.Random CFloats] ~ 280 cycles/int - 2,421,520 randoms generated [System.Random Doubles] ~ 1,377 cycles/int - 2,535,087 randoms generated [System.Random CDoubles] ~ 1,315 cycles/int - 856,276 randoms generated [System.Random Integers] ~ 3,895 cycles/int - 4,976,373 randoms generated [System.Random Bools] ~ 670 cycles/int - -(But I still need to propagate these changes throughout all types / API calls.) - - - -[2011.06.28] Integer Generation via random and randomR -------------------------------------------------------- - -Back on the master branch I notice that while randomIvalInteger does -well for small ranges, it's advantage doesn't scale to larger ranges: - - range (-100,100): - 5,105,290 randoms generated [System.Random Integers] ~ 653 cycles/int - - range (0,2^5000): - 8,969 randoms generated [System.Random BIG Integers] ~ 371,848 cycles/int - - - -[2011.08.25] Validating release version 1.0.1.0 rev 40bbfd2867 --------------------------------------------------------------- - -This is a bugfix release without SplittableGen. It passed (cd tests; -make test) on my Mac Os 10.6 machine. - -I ran GHC validate using the following fingerprint - - .|c5056b932a06b4adce5167a5cb69f1f0768d28ec - ghc-tarballs|e7b7b152083f7c3e3559e557a239757d41ac02a6 - libraries/Cabal|3dcc425495523ab6142027097cb598a4d2ad810a - libraries/Win32|085b11285b6adbc6484d9c21f5e0b8105556869c - libraries/array|fa295423e7404d3d1d3d82655b2b44d50f045a44 - libraries/base|a57369f54bd25a1de5d477f3c363b3bafd17d168 - libraries/binary|9065df2120254601c68c3a28264fd062abde9354 - libraries/bytestring|caad22630f73e0e9b1b61f4da34f8aefcc6001d8 - libraries/containers|667591b168c804d3eeae503dff1c848ed6852412 - libraries/directory|d44f52926278319286804d8333149dd13d4ecc82 - libraries/dph|b23b45a9e8fce985327b076997d61ab0ddc7b2f7 - libraries/extensible-exceptions|e77722871a5812d52c467e3a8fd9c7b97cdec521 - libraries/filepath|fd381017dca45de5c94dac85a6233516a6b6963d - libraries/ghc-prim|0a84a755e1248b4d50f6535a0ce75af0bb21b0ad - libraries/haskeline|8787a64417500efffc9c48032ee7c37315fb2547 - libraries/haskell2010|470b34b6c0336339eac9fbcfb6020e46b5154bfe - libraries/haskell98|5590c0f042d6d07352e0bf49cedeef5ba0821f23 - libraries/hoopl|b98db91cd0c53ddb2c275c04823f9c379774104b - libraries/hpc|7c726abec939b11af1ecf89740ca8d04e6a1360d - libraries/integer-gmp|65c73842bca2f075b65f418a5ff34753b185e0d7 - libraries/integer-simple|10202111c59f0695ef782d1ec9e6fc992933fc9a - libraries/mtl|a41470c1890073e678f0dca2a9ef4c02d55bf7c6 - libraries/old-locale|104e7e5a7b33424f34e98825a0d4ccb7614ca7c2 - libraries/old-time|81e0c8a4b98d4b084eefe75bedf91a44edd31888 - libraries/pretty|036fb8dfbb9d4a787fcd150c2756b4899be4e942 - libraries/primitive|74479e07b92b8859eae473e5cc86b40decae1d6e - libraries/process|68ba490d6691f55eab19a249379144831055e2ac - libraries/random|3fb0e9e42b54d7b01b794fc27d4d678d7d74ff0e - libraries/template-haskell|02362d12e5ae0af20d637eec97db51f6827a1625 - libraries/terminfo|baec6aff59d13ba294b370f9563e8068706392ce - libraries/unix|f55638fb5c6badd385c51a41de7ff96ef106de42 - libraries/utf8-string|ec2b85940a256dbc8771e5e2065ca8f76cc802d0 - libraries/vector|1e72d46bdc8488a84558b64ac63632cef1d8a695 - libraries/xhtml|cb2cf6c34d815fdf4ed74efeb65e1993e7bda514 - testsuite|26c608a0c31d56917099e4f48bf58c1d1e92e61c - utils/haddock|d54959189f33105ed09a59efee5ba34f53369282 - utils/hsc2hs|f8cbf37ab28ab4512d932678c08c263aa412e008 - - - -First validating in the context of a slightly stale GHC head -(7.3.20110727) on a mac. - - -[2011.09.30] Redoing timings after bugfix in version 1.0.1.1 ------------------------------------------------------------- - -It looks like there has been serious performance regression (3.33ghz -nehalem still). - - How many random numbers can we generate in a second on one thread? - Cost of rdtsc (ffi call): 38 - Approx getCPUTime calls per second: 7,121 - Approx clock frequency: 96,610,524 - First, timing System.Random.next: - 148,133,038 randoms generated [constant zero gen] ~ 0.65 cycles/int - 12,656,455 randoms generated [System.Random stdGen/next] ~ 7.63 cycles/int - - Second, timing System.Random.random at different types: - 676,066 randoms generated [System.Random Ints] ~ 143 cycles/int - 3,917,247 randoms generated [System.Random Word16] ~ 24.66 cycles/int - 2,231,460 randoms generated [System.Random Floats] ~ 43.29 cycles/int - 2,269,993 randoms generated [System.Random CFloats] ~ 42.56 cycles/int - 686,363 randoms generated [System.Random Doubles] ~ 141 cycles/int - 2,165,679 randoms generated [System.Random CDoubles] ~ 44.61 cycles/int - 713,702 randoms generated [System.Random Integers] ~ 135 cycles/int - 3,647,551 randoms generated [System.Random Bools] ~ 26.49 cycles/int - 4,296,919 randoms generated [System.Random Chars] ~ 22.48 cycles/int - - Next timing range-restricted System.Random.randomR: - 4,307,214 randoms generated [System.Random Ints] ~ 22.43 cycles/int - 4,068,982 randoms generated [System.Random Word16s] ~ 23.74 cycles/int - 2,059,264 randoms generated [System.Random Floats] ~ 46.92 cycles/int - 1,960,359 randoms generated [System.Random CFloats] ~ 49.28 cycles/int - 678,978 randoms generated [System.Random Doubles] ~ 142 cycles/int - 2,009,665 randoms generated [System.Random CDoubles] ~ 48.07 cycles/int - 4,296,452 randoms generated [System.Random Integers] ~ 22.49 cycles/int - 3,689,999 randoms generated [System.Random Bools] ~ 26.18 cycles/int - 4,367,577 randoms generated [System.Random Chars] ~ 22.12 cycles/int - 6,650 randoms generated [System.Random BIG Integers] ~ 14,528 cycles/int - diff --git a/README.md b/README.md index 9d5bb51b2..98005110b 100644 --- a/README.md +++ b/README.md @@ -1,6 +1,19 @@ -The Haskell Standard Library -- Random Number Generation -======================================================== -[![Build Status](https://secure.travis-ci.org/haskell/random.svg?branch=master)](http://travis-ci.org/haskell/random) +# The Haskell Standard Library + +## Random Number Generation + +### Status + +| Language | Travis | Coveralls | +|:--------:|:------:|:---------:| +| ![GitHub top language](https://img.shields.io/github/languages/top/idontgetoutmuch/random.svg) | [![Build Status](https://secure.travis-ci.org/idontgetoutmuch/random.svg?v1.2-proposal)](http://travis-ci.org/idontgetoutmuch/random) | [![Coverage Status](https://coveralls.io/repos/github/idontgetoutmuch/random/badge.svg?branch=v1.2-proposal)](https://coveralls.io/github/idontgetoutmuch/random?branch=v1.2-proposal) + +| Package | Hackage | Nightly | LTS | +|:-------------------|:-------:|:-------:|:---:| +| [`random`](https://github.com/idontgetoutmuch/random)| [![Hackage](https://img.shields.io/hackage/v/random.svg)](https://hackage.haskell.org/package/random)| [![Nightly](https://www.stackage.org/package/random/badge/nightly)](https://www.stackage.org/nightly/package/random)| [![Nightly](https://www.stackage.org/package/random/badge/lts)](https://www.stackage.org/lts/package/random) + + +### Description This library provides a basic interface for (splittable) random number generators. diff --git a/Setup.hs b/Setup.hs index 6fa548caf..8ec54a08d 100644 --- a/Setup.hs +++ b/Setup.hs @@ -1,6 +1,33 @@ +{-# LANGUAGE CPP #-} +{-# OPTIONS_GHC -Wall #-} module Main (main) where +#ifndef MIN_VERSION_cabal_doctest +#define MIN_VERSION_cabal_doctest(x,y,z) 0 +#endif + +#if MIN_VERSION_cabal_doctest(1,0,0) + +import Distribution.Extra.Doctest ( defaultMainWithDoctests ) +main :: IO () +main = defaultMainWithDoctests "doctests" + +#else + +#ifdef MIN_VERSION_Cabal +-- If the macro is defined, we have new cabal-install, +-- but for some reason we don't have cabal-doctest in package-db +-- +-- Probably we are running cabal sdist, when otherwise using new-build +-- workflow +#warning You are configuring this package without cabal-doctest installed. \ + The doctests test-suite will not work as a result. \ + To fix this, install cabal-doctest before configuring. +#endif + import Distribution.Simple main :: IO () main = defaultMain + +#endif diff --git a/System/Random.hs b/System/Random.hs deleted file mode 100644 index ab7727405..000000000 --- a/System/Random.hs +++ /dev/null @@ -1,609 +0,0 @@ -#if __GLASGOW_HASKELL__ >= 701 -{-# LANGUAGE Trustworthy #-} -#endif - ------------------------------------------------------------------------------ --- | --- Module : System.Random --- Copyright : (c) The University of Glasgow 2001 --- License : BSD-style (see the file LICENSE in the 'random' repository) --- --- Maintainer : libraries@haskell.org --- Stability : stable --- Portability : portable --- --- This library deals with the common task of pseudo-random number --- generation. The library makes it possible to generate repeatable --- results, by starting with a specified initial random number generator, --- or to get different results on each run by using the system-initialised --- generator or by supplying a seed from some other source. --- --- The library is split into two layers: --- --- * A core /random number generator/ provides a supply of bits. --- The class 'RandomGen' provides a common interface to such generators. --- The library provides one instance of 'RandomGen', the abstract --- data type 'StdGen'. Programmers may, of course, supply their own --- instances of 'RandomGen'. --- --- * The class 'Random' provides a way to extract values of a particular --- type from a random number generator. For example, the 'Float' --- instance of 'Random' allows one to generate random values of type --- 'Float'. --- --- This implementation uses the Portable Combined Generator of L'Ecuyer --- ["System.Random\#LEcuyer"] for 32-bit computers, transliterated by --- Lennart Augustsson. It has a period of roughly 2.30584e18. --- ------------------------------------------------------------------------------ - -#include "MachDeps.h" - -module System.Random - ( - - -- $intro - - -- * Random number generators - -#ifdef ENABLE_SPLITTABLEGEN - RandomGen(next, genRange) - , SplittableGen(split) -#else - RandomGen(next, genRange, split) -#endif - -- ** Standard random number generators - , StdGen - , mkStdGen - - -- ** The global random number generator - - -- $globalrng - - , getStdRandom - , getStdGen - , setStdGen - , newStdGen - - -- * Random values of various types - , Random ( random, randomR, - randoms, randomRs, - randomIO, randomRIO ) - - -- * References - -- $references - - ) where - -import Prelude - -import Data.Bits -import Data.Int -import Data.Word -import Foreign.C.Types - -#ifdef __NHC__ -import CPUTime ( getCPUTime ) -import Foreign.Ptr ( Ptr, nullPtr ) -import Foreign.C ( CTime, CUInt ) -#else -import System.CPUTime ( getCPUTime ) -import Data.Time ( getCurrentTime, UTCTime(..) ) -import Data.Ratio ( numerator, denominator ) -#endif -import Data.Char ( isSpace, chr, ord ) -import System.IO.Unsafe ( unsafePerformIO ) -import Data.IORef ( IORef, newIORef, readIORef, writeIORef ) -#if MIN_VERSION_base (4,6,0) -import Data.IORef ( atomicModifyIORef' ) -#else -import Data.IORef ( atomicModifyIORef ) -#endif -import Numeric ( readDec ) - -#ifdef __GLASGOW_HASKELL__ -import GHC.Exts ( build ) -#else --- | A dummy variant of build without fusion. -{-# INLINE build #-} -build :: ((a -> [a] -> [a]) -> [a] -> [a]) -> [a] -build g = g (:) [] -#endif - -#if !MIN_VERSION_base (4,6,0) -atomicModifyIORef' :: IORef a -> (a -> (a,b)) -> IO b -atomicModifyIORef' ref f = do - b <- atomicModifyIORef ref - (\x -> let (a, b) = f x - in (a, a `seq` b)) - b `seq` return b -#endif - --- The standard nhc98 implementation of Time.ClockTime does not match --- the extended one expected in this module, so we lash-up a quick --- replacement here. -#ifdef __NHC__ -foreign import ccall "time.h time" readtime :: Ptr CTime -> IO CTime -getTime :: IO (Integer, Integer) -getTime = do CTime t <- readtime nullPtr; return (toInteger t, 0) -#else -getTime :: IO (Integer, Integer) -getTime = do - utc <- getCurrentTime - let daytime = toRational $ utctDayTime utc - return $ quotRem (numerator daytime) (denominator daytime) -#endif - --- | The class 'RandomGen' provides a common interface to random number --- generators. --- -#ifdef ENABLE_SPLITTABLEGEN --- Minimal complete definition: 'next'. -#else --- Minimal complete definition: 'next' and 'split'. -#endif - -class RandomGen g where - - -- |The 'next' operation returns an 'Int' that is uniformly distributed - -- in the range returned by 'genRange' (including both end points), - -- and a new generator. - next :: g -> (Int, g) - - -- |The 'genRange' operation yields the range of values returned by - -- the generator. - -- - -- It is required that: - -- - -- * If @(a,b) = 'genRange' g@, then @a < b@. - -- - -- * 'genRange' always returns a pair of defined 'Int's. - -- - -- The second condition ensures that 'genRange' cannot examine its - -- argument, and hence the value it returns can be determined only by the - -- instance of 'RandomGen'. That in turn allows an implementation to make - -- a single call to 'genRange' to establish a generator's range, without - -- being concerned that the generator returned by (say) 'next' might have - -- a different range to the generator passed to 'next'. - -- - -- The default definition spans the full range of 'Int'. - genRange :: g -> (Int,Int) - - -- default method - genRange _ = (minBound, maxBound) - -#ifdef ENABLE_SPLITTABLEGEN --- | The class 'SplittableGen' proivides a way to specify a random number --- generator that can be split into two new generators. -class SplittableGen g where -#endif - -- |The 'split' operation allows one to obtain two distinct random number - -- generators. This is very useful in functional programs (for example, when - -- passing a random number generator down to recursive calls), but very - -- little work has been done on statistically robust implementations of - -- 'split' (["System.Random\#Burton", "System.Random\#Hellekalek"] - -- are the only examples we know of). - split :: g -> (g, g) - -{- | -The 'StdGen' instance of 'RandomGen' has a 'genRange' of at least 30 bits. - -The result of repeatedly using 'next' should be at least as statistically -robust as the /Minimal Standard Random Number Generator/ described by -["System.Random\#Park", "System.Random\#Carta"]. -Until more is known about implementations of 'split', all we require is -that 'split' deliver generators that are (a) not identical and -(b) independently robust in the sense just given. - -The 'Show' and 'Read' instances of 'StdGen' provide a primitive way to save the -state of a random number generator. -It is required that @'read' ('show' g) == g@. - -In addition, 'reads' may be used to map an arbitrary string (not necessarily one -produced by 'show') onto a value of type 'StdGen'. In general, the 'Read' -instance of 'StdGen' has the following properties: - -* It guarantees to succeed on any string. - -* It guarantees to consume only a finite portion of the string. - -* Different argument strings are likely to result in different results. - --} - -data StdGen - = StdGen !Int32 !Int32 - -instance RandomGen StdGen where - next = stdNext - genRange _ = stdRange - -#ifdef ENABLE_SPLITTABLEGEN -instance SplittableGen StdGen where -#endif - split = stdSplit - -instance Show StdGen where - showsPrec p (StdGen s1 s2) = - showsPrec p s1 . - showChar ' ' . - showsPrec p s2 - -instance Read StdGen where - readsPrec _p = \ r -> - case try_read r of - r'@[_] -> r' - _ -> [stdFromString r] -- because it shouldn't ever fail. - where - try_read r = do - (s1, r1) <- readDec (dropWhile isSpace r) - (s2, r2) <- readDec (dropWhile isSpace r1) - return (StdGen s1 s2, r2) - -{- - If we cannot unravel the StdGen from a string, create - one based on the string given. --} -stdFromString :: String -> (StdGen, String) -stdFromString s = (mkStdGen num, rest) - where (cs, rest) = splitAt 6 s - num = foldl (\a x -> x + 3 * a) 1 (map ord cs) - - -{- | -The function 'mkStdGen' provides an alternative way of producing an initial -generator, by mapping an 'Int' into a generator. Again, distinct arguments -should be likely to produce distinct generators. --} -mkStdGen :: Int -> StdGen -- why not Integer ? -mkStdGen s = mkStdGen32 $ fromIntegral s - -{- -From ["System.Random\#LEcuyer"]: "The integer variables s1 and s2 ... must be -initialized to values in the range [1, 2147483562] and [1, 2147483398] -respectively." --} -mkStdGen32 :: Int32 -> StdGen -mkStdGen32 sMaybeNegative = StdGen (s1+1) (s2+1) - where - -- We want a non-negative number, but we can't just take the abs - -- of sMaybeNegative as -minBound == minBound. - s = sMaybeNegative .&. maxBound - (q, s1) = s `divMod` 2147483562 - s2 = q `mod` 2147483398 - -createStdGen :: Integer -> StdGen -createStdGen s = mkStdGen32 $ fromIntegral s - -{- | -With a source of random number supply in hand, the 'Random' class allows the -programmer to extract random values of a variety of types. - -Minimal complete definition: 'randomR' and 'random'. - --} - -class Random a where - -- | Takes a range /(lo,hi)/ and a random number generator - -- /g/, and returns a random value uniformly distributed in the closed - -- interval /[lo,hi]/, together with a new generator. It is unspecified - -- what happens if /lo>hi/. For continuous types there is no requirement - -- that the values /lo/ and /hi/ are ever produced, but they may be, - -- depending on the implementation and the interval. - randomR :: RandomGen g => (a,a) -> g -> (a,g) - - -- | The same as 'randomR', but using a default range determined by the type: - -- - -- * For bounded types (instances of 'Bounded', such as 'Char'), - -- the range is normally the whole type. - -- - -- * For fractional types, the range is normally the semi-closed interval - -- @[0,1)@. - -- - -- * For 'Integer', the range is (arbitrarily) the range of 'Int'. - random :: RandomGen g => g -> (a, g) - - -- | Plural variant of 'randomR', producing an infinite list of - -- random values instead of returning a new generator. - {-# INLINE randomRs #-} - randomRs :: RandomGen g => (a,a) -> g -> [a] - randomRs ival g = build (\cons _nil -> buildRandoms cons (randomR ival) g) - - -- | Plural variant of 'random', producing an infinite list of - -- random values instead of returning a new generator. - {-# INLINE randoms #-} - randoms :: RandomGen g => g -> [a] - randoms g = build (\cons _nil -> buildRandoms cons random g) - - -- | A variant of 'randomR' that uses the global random number generator - -- (see "System.Random#globalrng"). - randomRIO :: (a,a) -> IO a - randomRIO range = getStdRandom (randomR range) - - -- | A variant of 'random' that uses the global random number generator - -- (see "System.Random#globalrng"). - randomIO :: IO a - randomIO = getStdRandom random - --- | Produce an infinite list-equivalent of random values. -{-# INLINE buildRandoms #-} -buildRandoms :: RandomGen g - => (a -> as -> as) -- ^ E.g. '(:)' but subject to fusion - -> (g -> (a,g)) -- ^ E.g. 'random' - -> g -- ^ A 'RandomGen' instance - -> as -buildRandoms cons rand = go - where - -- The seq fixes part of #4218 and also makes fused Core simpler. - go g = x `seq` (x `cons` go g') where (x,g') = rand g - - -instance Random Integer where - randomR ival g = randomIvalInteger ival g - random g = randomR (toInteger (minBound::Int), toInteger (maxBound::Int)) g - -instance Random Int where randomR = randomIvalIntegral; random = randomBounded -instance Random Int8 where randomR = randomIvalIntegral; random = randomBounded -instance Random Int16 where randomR = randomIvalIntegral; random = randomBounded -instance Random Int32 where randomR = randomIvalIntegral; random = randomBounded -instance Random Int64 where randomR = randomIvalIntegral; random = randomBounded - -#ifndef __NHC__ --- Word is a type synonym in nhc98. -instance Random Word where randomR = randomIvalIntegral; random = randomBounded -#endif -instance Random Word8 where randomR = randomIvalIntegral; random = randomBounded -instance Random Word16 where randomR = randomIvalIntegral; random = randomBounded -instance Random Word32 where randomR = randomIvalIntegral; random = randomBounded -instance Random Word64 where randomR = randomIvalIntegral; random = randomBounded - -instance Random CChar where randomR = randomIvalIntegral; random = randomBounded -instance Random CSChar where randomR = randomIvalIntegral; random = randomBounded -instance Random CUChar where randomR = randomIvalIntegral; random = randomBounded -instance Random CShort where randomR = randomIvalIntegral; random = randomBounded -instance Random CUShort where randomR = randomIvalIntegral; random = randomBounded -instance Random CInt where randomR = randomIvalIntegral; random = randomBounded -instance Random CUInt where randomR = randomIvalIntegral; random = randomBounded -instance Random CLong where randomR = randomIvalIntegral; random = randomBounded -instance Random CULong where randomR = randomIvalIntegral; random = randomBounded -instance Random CPtrdiff where randomR = randomIvalIntegral; random = randomBounded -instance Random CSize where randomR = randomIvalIntegral; random = randomBounded -instance Random CWchar where randomR = randomIvalIntegral; random = randomBounded -instance Random CSigAtomic where randomR = randomIvalIntegral; random = randomBounded -instance Random CLLong where randomR = randomIvalIntegral; random = randomBounded -instance Random CULLong where randomR = randomIvalIntegral; random = randomBounded -instance Random CIntPtr where randomR = randomIvalIntegral; random = randomBounded -instance Random CUIntPtr where randomR = randomIvalIntegral; random = randomBounded -instance Random CIntMax where randomR = randomIvalIntegral; random = randomBounded -instance Random CUIntMax where randomR = randomIvalIntegral; random = randomBounded - -instance Random Char where - randomR (a,b) g = - case (randomIvalInteger (toInteger (ord a), toInteger (ord b)) g) of - (x,g') -> (chr x, g') - random g = randomR (minBound,maxBound) g - -instance Random Bool where - randomR (a,b) g = - case (randomIvalInteger (bool2Int a, bool2Int b) g) of - (x, g') -> (int2Bool x, g') - where - bool2Int :: Bool -> Integer - bool2Int False = 0 - bool2Int True = 1 - - int2Bool :: Int -> Bool - int2Bool 0 = False - int2Bool _ = True - - random g = randomR (minBound,maxBound) g - -{-# INLINE randomRFloating #-} -randomRFloating :: (Fractional a, Num a, Ord a, Random a, RandomGen g) => (a, a) -> g -> (a, g) -randomRFloating (l,h) g - | l>h = randomRFloating (h,l) g - | otherwise = let (coef,g') = random g in - (2.0 * (0.5*l + coef * (0.5*h - 0.5*l)), g') -- avoid overflow - -instance Random Double where - randomR = randomRFloating - random rng = - case random rng of - (x,rng') -> - -- We use 53 bits of randomness corresponding to the 53 bit significand: - ((fromIntegral (mask53 .&. (x::Int64)) :: Double) - / fromIntegral twoto53, rng') - where - twoto53 = (2::Int64) ^ (53::Int64) - mask53 = twoto53 - 1 - -instance Random Float where - randomR = randomRFloating - random rng = - -- TODO: Faster to just use 'next' IF it generates enough bits of randomness. - case random rng of - (x,rng') -> - -- We use 24 bits of randomness corresponding to the 24 bit significand: - ((fromIntegral (mask24 .&. (x::Int32)) :: Float) - / fromIntegral twoto24, rng') - -- Note, encodeFloat is another option, but I'm not seeing slightly - -- worse performance with the following [2011.06.25]: --- (encodeFloat rand (-24), rng') - where - mask24 = twoto24 - 1 - twoto24 = (2::Int32) ^ (24::Int32) - --- CFloat/CDouble are basically the same as a Float/Double: -instance Random CFloat where - randomR = randomRFloating - random rng = case random rng of - (x,rng') -> (realToFrac (x::Float), rng') - -instance Random CDouble where - randomR = randomRFloating - -- A MYSTERY: - -- Presently, this is showing better performance than the Double instance: - -- (And yet, if the Double instance uses randomFrac then its performance is much worse!) - random = randomFrac - -- random rng = case random rng of - -- (x,rng') -> (realToFrac (x::Double), rng') - -mkStdRNG :: Integer -> IO StdGen -mkStdRNG o = do - ct <- getCPUTime - (sec, psec) <- getTime - return (createStdGen (sec * 12345 + psec + ct + o)) - -randomBounded :: (RandomGen g, Random a, Bounded a) => g -> (a, g) -randomBounded = randomR (minBound, maxBound) - --- The two integer functions below take an [inclusive,inclusive] range. -randomIvalIntegral :: (RandomGen g, Integral a) => (a, a) -> g -> (a, g) -randomIvalIntegral (l,h) = randomIvalInteger (toInteger l, toInteger h) - -{-# SPECIALIZE randomIvalInteger :: (Num a) => - (Integer, Integer) -> StdGen -> (a, StdGen) #-} - -randomIvalInteger :: (RandomGen g, Num a) => (Integer, Integer) -> g -> (a, g) -randomIvalInteger (l,h) rng - | l > h = randomIvalInteger (h,l) rng - | otherwise = case (f 1 0 rng) of (v, rng') -> (fromInteger (l + v `mod` k), rng') - where - (genlo, genhi) = genRange rng - b = fromIntegral genhi - fromIntegral genlo + 1 - - -- Probabilities of the most likely and least likely result - -- will differ at most by a factor of (1 +- 1/q). Assuming the RandomGen - -- is uniform, of course - - -- On average, log q / log b more random values will be generated - -- than the minimum - q = 1000 - k = h - l + 1 - magtgt = k * q - - -- generate random values until we exceed the target magnitude - f mag v g | mag >= magtgt = (v, g) - | otherwise = v' `seq`f (mag*b) v' g' where - (x,g') = next g - v' = (v * b + (fromIntegral x - fromIntegral genlo)) - - --- The continuous functions on the other hand take an [inclusive,exclusive) range. -randomFrac :: (RandomGen g, Fractional a) => g -> (a, g) -randomFrac = randomIvalDouble (0::Double,1) realToFrac - -randomIvalDouble :: (RandomGen g, Fractional a) => (Double, Double) -> (Double -> a) -> g -> (a, g) -randomIvalDouble (l,h) fromDouble rng - | l > h = randomIvalDouble (h,l) fromDouble rng - | otherwise = - case (randomIvalInteger (toInteger (minBound::Int32), toInteger (maxBound::Int32)) rng) of - (x, rng') -> - let - scaled_x = - fromDouble (0.5*l + 0.5*h) + -- previously (l+h)/2, overflowed - fromDouble ((0.5*h - 0.5*l) / (0.5 * realToFrac int32Count)) * -- avoid overflow - fromIntegral (x::Int32) - in - (scaled_x, rng') - -int32Count :: Integer -int32Count = toInteger (maxBound::Int32) - toInteger (minBound::Int32) + 1 -- GHC ticket #3982 - -stdRange :: (Int,Int) -stdRange = (1, 2147483562) - -stdNext :: StdGen -> (Int, StdGen) --- Returns values in the range stdRange -stdNext (StdGen s1 s2) = (fromIntegral z', StdGen s1'' s2'') - where z' = if z < 1 then z + 2147483562 else z - z = s1'' - s2'' - - k = s1 `quot` 53668 - s1' = 40014 * (s1 - k * 53668) - k * 12211 - s1'' = if s1' < 0 then s1' + 2147483563 else s1' - - k' = s2 `quot` 52774 - s2' = 40692 * (s2 - k' * 52774) - k' * 3791 - s2'' = if s2' < 0 then s2' + 2147483399 else s2' - -stdSplit :: StdGen -> (StdGen, StdGen) -stdSplit std@(StdGen s1 s2) - = (left, right) - where - -- no statistical foundation for this! - left = StdGen new_s1 t2 - right = StdGen t1 new_s2 - - new_s1 | s1 == 2147483562 = 1 - | otherwise = s1 + 1 - - new_s2 | s2 == 1 = 2147483398 - | otherwise = s2 - 1 - - StdGen t1 t2 = snd (next std) - --- The global random number generator - -{- $globalrng #globalrng# - -There is a single, implicit, global random number generator of type -'StdGen', held in some global variable maintained by the 'IO' monad. It is -initialised automatically in some system-dependent fashion, for example, by -using the time of day, or Linux's kernel random number generator. To get -deterministic behaviour, use 'setStdGen'. --} - --- |Sets the global random number generator. -setStdGen :: StdGen -> IO () -setStdGen sgen = writeIORef theStdGen sgen - --- |Gets the global random number generator. -getStdGen :: IO StdGen -getStdGen = readIORef theStdGen - -theStdGen :: IORef StdGen -theStdGen = unsafePerformIO $ do - rng <- mkStdRNG 0 - newIORef rng - --- |Applies 'split' to the current global random generator, --- updates it with one of the results, and returns the other. -newStdGen :: IO StdGen -newStdGen = atomicModifyIORef' theStdGen split - -{- |Uses the supplied function to get a value from the current global -random generator, and updates the global generator with the new generator -returned by the function. For example, @rollDice@ gets a random integer -between 1 and 6: - -> rollDice :: IO Int -> rollDice = getStdRandom (randomR (1,6)) - --} - -getStdRandom :: (StdGen -> (a,StdGen)) -> IO a -getStdRandom f = atomicModifyIORef' theStdGen (swap . f) - where swap (v,g) = (g,v) - -{- $references - -1. FW #Burton# Burton and RL Page, /Distributed random number generation/, -Journal of Functional Programming, 2(2):203-212, April 1992. - -2. SK #Park# Park, and KW Miller, /Random number generators - -good ones are hard to find/, Comm ACM 31(10), Oct 1988, pp1192-1201. - -3. DG #Carta# Carta, /Two fast implementations of the minimal standard -random number generator/, Comm ACM, 33(1), Jan 1990, pp87-88. - -4. P #Hellekalek# Hellekalek, /Don\'t trust parallel Monte Carlo/, -Department of Mathematics, University of Salzburg, -, 1998. - -5. Pierre #LEcuyer# L'Ecuyer, /Efficient and portable combined random -number generators/, Comm ACM, 31(6), Jun 1988, pp742-749. - -The Web site is a great source of information. - --} diff --git a/bench-legacy/BinSearch.hs b/bench-legacy/BinSearch.hs new file mode 100644 index 000000000..81a57930e --- /dev/null +++ b/bench-legacy/BinSearch.hs @@ -0,0 +1,149 @@ + +{- + Binary search over benchmark input sizes. + + There are many good ways to measure the time it takes to perform a + certain computation on a certain input. However, frequently, it's + challenging to pick the right input size for all platforms and all + compilataion modes. + + Sometimes for linear-complexity benchmarks it is better to measure + /throughput/, i.e. elements processed per second. That is, fixing + the time of execution and measuring the amount of work done (rather + than the reverse). This library provides a simple way to search for + an appropriate input size that results in the desired execution time. + + An alternative approach is to kill the computation after a certain + amount of time and observe how much work it has completed. + -} +module BinSearch + ( + binSearch + ) +where + +import Control.Monad +import Data.Time.Clock -- Not in 6.10 +import Data.List +import System.IO +import Prelude hiding (min,max,log) + + + +-- | Binary search for the number of inputs to a computation that +-- results in a specified amount of execution time in seconds. For example: +-- +-- > binSearch verbose N (min,max) kernel +-- +-- ... will find the right input size that results in a time +-- between min and max, then it will then run for N trials and +-- return the median (input,time-in-seconds) pair. +binSearch :: Bool -> Integer -> (Double,Double) -> (Integer -> IO ()) -> IO (Integer, Double) +binSearch verbose trials (min, max) kernel = do + when verbose $ + putStrLn $ + "[binsearch] Binary search for input size resulting in time in range " ++ + show (min, max) + let desired_exec_length = 1.0 + good_trial t = + (toRational t <= toRational max) && (toRational t >= toRational min) + -- At some point we must give up... + loop n + | n > ((2 :: Integer) ^ (100 :: Integer)) = + error + "ERROR binSearch: This function doesn't seem to scale in proportion to its last argument." + -- Not allowed to have "0" size input, bump it back to one: + loop 0 = loop 1 + loop n = do + when verbose $ putStr $ "[binsearch:" ++ show n ++ "] " + time <- timeit $ kernel n + when verbose $ putStrLn $ "Time consumed: " ++ show time + let rate = fromIntegral n / time + -- [2010.06.09] Introducing a small fudge factor to help our guess get over the line: + let initial_fudge_factor = 1.10 + fudge_factor = 1.01 -- Even in the steady state we fudge a little + guess = desired_exec_length * rate + -- TODO: We should keep more history here so that we don't re-explore input space we + -- have already explored. This is a balancing act because of randomness in + -- execution time. + if good_trial time + then do + when verbose $ + putStrLn + "[binsearch] Time in range. LOCKING input size and performing remaining trials." + print_trial 1 n time + lockin (trials - 1) n [time] + else if time < 0.100 + then loop (2 * n) + else do + when verbose $ + putStrLn $ + "[binsearch] Estimated rate to be " ++ + show (round rate :: Integer) ++ + " per second. Trying to scale up..." + -- Here we've exited the doubling phase, but we're making our + -- first guess as to how big a real execution should be: + if time > 0.100 && time < 0.33 * desired_exec_length + then do + when verbose $ + putStrLn + "[binsearch] (Fudging first guess a little bit extra)" + loop (round $ guess * initial_fudge_factor) + else loop (round $ guess * fudge_factor) + -- Termination condition: Done with all trials. + lockin 0 n log = do + when verbose $ + putStrLn $ + "[binsearch] Time-per-unit for all trials: " ++ + concat + (intersperse " " (map (show . (/ toDouble n) . toDouble) $ sort log)) + return (n, log !! (length log `quot` 2)) -- Take the median + lockin trials_left n log = do + when verbose $ + putStrLn + "[binsearch]------------------------------------------------------------" + time <- timeit $ kernel n + -- hFlush stdout + print_trial (trials - trials_left + 1) n time + -- whenverbose$ hFlush stdout + lockin (trials_left - 1) n (time : log) + print_trial :: Integer -> Integer -> NominalDiffTime -> IO () + print_trial trialnum n time = + let rate = fromIntegral n / time + timeperunit = time / fromIntegral n + in when verbose $ + putStrLn $ + "[binsearch] TRIAL: " ++ + show trialnum ++ + " secPerUnit: " ++ + showTime timeperunit ++ + " ratePerSec: " ++ show rate ++ " seconds: " ++ showTime time + (n, t) <- loop 1 + return (n, fromRational $ toRational t) + + +showTime :: NominalDiffTime -> String +showTime t = show ((fromRational $ toRational t) :: Double) + +toDouble :: Real a => a -> Double +toDouble = fromRational . toRational + + +-- Could use cycle counters here.... but the point of this is to time +-- things on the order of a second. +timeit :: IO () -> IO NominalDiffTime +timeit io = do + strt <- getCurrentTime + io + end <- getCurrentTime + return (diffUTCTime end strt) +{- +test :: IO (Integer,Double) +test = + binSearch True 3 (1.0, 1.05) + (\n -> + do v <- newIORef 0 + forM_ [1..n] $ \i -> do + old <- readIORef v + writeIORef v (old+i)) +-} diff --git a/bench-legacy/SimpleRNGBench.hs b/bench-legacy/SimpleRNGBench.hs new file mode 100644 index 000000000..61df451f5 --- /dev/null +++ b/bench-legacy/SimpleRNGBench.hs @@ -0,0 +1,269 @@ +{-# LANGUAGE BangPatterns, ScopedTypeVariables, ForeignFunctionInterface #-} +{-# OPTIONS_GHC -fwarn-unused-imports #-} + +-- | A simple script to do some very basic timing of the RNGs. + +module Main where + +import System.Exit (exitSuccess, exitFailure) +import System.Environment +import System.Random +import System.CPUTime (getCPUTime) +import System.CPUTime.Rdtsc +import System.Console.GetOpt + +import GHC.Conc +import Control.Concurrent +import Control.Monad +import Control.Exception + +import Data.IORef +import Data.Word +import Data.List hiding (last,sum) +import Data.Int +import Data.List.Split hiding (split) +import Text.Printf + +import Foreign.Ptr +import Foreign.C.Types +import Foreign.ForeignPtr +import Foreign.Storable (peek,poke) + +import Prelude hiding (last,sum) +import BinSearch + +---------------------------------------------------------------------------------------------------- +-- Miscellaneous helpers: + +-- Readable large integer printing: +commaint :: Show a => a -> String +commaint n = reverse $ concat $ intersperse "," $ chunk 3 $ reverse (show n) + +padleft :: Int -> String -> String +padleft n str | length str >= n = str +padleft n str | otherwise = take (n - length str) (repeat ' ') ++ str + +padright :: Int -> String -> String +padright n str | length str >= n = str +padright n str | otherwise = str ++ take (n - length str) (repeat ' ') + +fmt_num :: (RealFrac a, PrintfArg a) => a -> String +fmt_num n = + if n < 100 + then printf "%.2f" n + else commaint (round n :: Integer) + + +-- Measure clock frequency, spinning rather than sleeping to try to +-- stay on the same core. +measureFreq :: IO Int64 +measureFreq = do + let second = 1000 * 1000 * 1000 * 1000 -- picoseconds are annoying + t1 <- rdtsc + start <- getCPUTime + let loop !n !last = do + t2 <- rdtsc + when (t2 < last) $ putStrLn $ "COUNTERS WRAPPED " ++ show (last, t2) + cput <- getCPUTime + if cput - start < second + then loop (n + 1) t2 + else return (n, t2) + (n, t2) <- loop 0 t1 + putStrLn $ " Approx getCPUTime calls per second: " ++ commaint (n :: Int64) + when (t2 < t1) $ + putStrLn $ + "WARNING: rdtsc not monotonically increasing, first " ++ + show t1 ++ " then " ++ show t2 ++ " on the same OS thread" + return $ fromIntegral (t2 - t1) + +---------------------------------------------------------------------------------------------------- + +-- Test overheads without actually generating any random numbers: +data NoopRNG = NoopRNG +instance RandomGen NoopRNG where + next g = (0, g) + genRange _ = (0, 0) + split g = (g, g) + +-- An RNG generating only 0 or 1: +data BinRNG = BinRNG StdGen +instance RandomGen BinRNG where + next (BinRNG g) = (x `mod` 2, BinRNG g') + where + (x, g') = next g + genRange _ = (0, 1) + split (BinRNG g) = (BinRNG g1, BinRNG g2) + where + (g1, g2) = split g + + +---------------------------------------------------------------------------------------------------- +-- Drivers to get random numbers repeatedly. + +type Kern = Int -> Ptr Int -> IO () + +-- [2011.01.28] Changing this to take "count" and "accumulator ptr" as arguments: +-- foreign import ccall "cbits/c_test.c" blast_rands :: Kern +-- foreign import ccall "cbits/c_test.c" store_loop :: Kern +-- foreign import ccall unsafe "stdlib.hs" rand :: IO Int + +{-# INLINE timeit #-} +timeit :: (Random a, RandomGen g) => Int -> Int64 -> String -> g -> (g -> (a,g)) -> IO () +timeit numthreads freq msg gen nxt = do + counters <- forM [1 .. numthreads] (const $ newIORef (1 :: Int64)) + tids <- forM counters $ \counter -> forkIO $ infloop counter (nxt gen) + threadDelay (1000 * 1000) -- One second + mapM_ killThread tids + finals <- mapM readIORef counters + let mean :: Double = + fromIntegral (foldl1 (+) finals) / fromIntegral numthreads + cycles_per :: Double = fromIntegral freq / mean + printResult (round mean :: Int64) msg cycles_per + where + infloop !counter (!_, !g) = do + incr counter + infloop counter (nxt g) + incr !counter + -- modifyIORef counter (+1) -- Not strict enough! + = do + c <- readIORef counter + let c' = c + 1 + _ <- evaluate c' + writeIORef counter c' + + +-- This function times an IO function on one or more threads. Rather +-- than running a fixed number of iterations, it uses a binary search +-- to find out how many iterations can be completed in a second. +timeit_foreign :: Int -> Int64 -> String -> (Int -> Ptr Int -> IO ()) -> IO Int64 +timeit_foreign numthreads freq msg ffn = do + ptr :: ForeignPtr Int <- mallocForeignPtr + let kern = + if numthreads == 1 + then ffn + else replicate_kernel numthreads ffn + wrapped n = withForeignPtr ptr (kern $ fromIntegral n) + (n, t) <- binSearch False 1 (1.0, 1.05) wrapped + let total_per_second = round $ fromIntegral n * (1 / t) + cycles_per = fromIntegral freq * t / fromIntegral n + printResult total_per_second msg cycles_per + return total_per_second + -- This lifts a C kernel to operate simultaneously on N threads. + where + replicate_kernel :: Int -> Kern -> Kern + replicate_kernel nthreads kern n ptr = do + ptrs <- forM [1 .. nthreads] (const mallocForeignPtr) + tmpchan <- newChan + -- let childwork = ceiling$ fromIntegral n / fromIntegral nthreads + let childwork = n -- Keep it the same.. interested in per-thread throughput. + -- Fork/join pattern: + forM_ ptrs $ \pt -> + forkIO $ + withForeignPtr pt $ \p -> do + kern (fromIntegral childwork) p + result <- peek p + writeChan tmpchan result + results <- forM [1 .. nthreads] $ \_ -> readChan tmpchan + -- Meaningless semantics here... sum the child ptrs and write to the input one: + poke ptr (foldl1 (+) results) + + +printResult :: Int64 -> String -> Double -> IO () +printResult total msg cycles_per = + putStrLn $ + " " ++ + padleft 11 (commaint total) ++ + " randoms generated " ++ + padright 27 ("[" ++ msg ++ "]") ++ + " ~ " ++ fmt_num cycles_per ++ " cycles/int" + +---------------------------------------------------------------------------------------------------- +-- Main Script + +data Flag = NoC | Help + deriving (Show, Eq) + +options :: [OptDescr Flag] +options = + [ Option ['h'] ["help"] (NoArg Help) "print program help" + , Option [] ["noC"] (NoArg NoC) "omit C benchmarks, haskell only" + ] + +main :: IO () +main = do + argv <- getArgs + let (opts,_,other) = getOpt Permute options argv + + unless (null other) $ do + putStrLn "ERROR: Unrecognized options: " + mapM_ putStr other + exitFailure + + when (Help `elem` opts) $ do + putStr $ usageInfo "Benchmark random number generation" options + exitSuccess + + putStrLn "\nHow many random numbers can we generate in a second on one thread?" + + t1 <- rdtsc + t2 <- rdtsc + putStrLn (" Cost of rdtsc (ffi call): " ++ show (t2 - t1)) + + freq <- measureFreq + putStrLn $ " Approx clock frequency: " ++ commaint freq + + let randInt = random :: RandomGen g => g -> (Int,g) + randWord16 = random :: RandomGen g => g -> (Word16,g) + randFloat = random :: RandomGen g => g -> (Float,g) + randCFloat = random :: RandomGen g => g -> (CFloat,g) + randDouble = random :: RandomGen g => g -> (Double,g) + randCDouble = random :: RandomGen g => g -> (CDouble,g) + randInteger = random :: RandomGen g => g -> (Integer,g) + randBool = random :: RandomGen g => g -> (Bool,g) + randChar = random :: RandomGen g => g -> (Char,g) + + gen = mkStdGen 23852358661234 + gamut th = do + putStrLn " First, timing System.Random.next:" + timeit th freq "constant zero gen" NoopRNG next + timeit th freq "System.Random stdGen/next" gen next + + putStrLn "\n Second, timing System.Random.random at different types:" + timeit th freq "System.Random Ints" gen randInt + timeit th freq "System.Random Word16" gen randWord16 + timeit th freq "System.Random Floats" gen randFloat + timeit th freq "System.Random CFloats" gen randCFloat + timeit th freq "System.Random Doubles" gen randDouble + timeit th freq "System.Random CDoubles" gen randCDouble + timeit th freq "System.Random Integers" gen randInteger + timeit th freq "System.Random Bools" gen randBool + timeit th freq "System.Random Chars" gen randChar + + putStrLn "\n Next timing range-restricted System.Random.randomR:" + timeit th freq "System.Random Ints" gen (randomR (-100, 100::Int)) + timeit th freq "System.Random Word16s" gen (randomR ( 100, 300::Word16)) + timeit th freq "System.Random Floats" gen (randomR (-100, 100::Float)) + timeit th freq "System.Random CFloats" gen (randomR (-100, 100::CFloat)) + timeit th freq "System.Random Doubles" gen (randomR (-100, 100::Double)) + timeit th freq "System.Random CDoubles" gen (randomR (-100, 100::CDouble)) + timeit th freq "System.Random Integers" gen (randomR (-100, 100::Integer)) + timeit th freq "System.Random Bools" gen (randomR (False, True::Bool)) + timeit th freq "System.Random Chars" gen (randomR ('a', 'z')) + timeit th freq "System.Random BIG Integers" gen (randomR (0, (2::Integer) ^ (5000::Int))) + + -- when (not$ NoC `elem` opts) $ do + -- putStrLn$ " Comparison to C's rand():" + -- timeit_foreign th freq "ptr store in C loop" store_loop + -- timeit_foreign th freq "rand/store in C loop" blast_rands + -- timeit_foreign th freq "rand in Haskell loop" (\n ptr -> forM_ [1..n]$ \_ -> rand ) + -- timeit_foreign th freq "rand/store in Haskell loop" (\n ptr -> forM_ [1..n]$ \_ -> do n <- rand; poke ptr n ) + -- return () + + -- Test with 1 thread and numCapabilities threads: + gamut 1 + when (numCapabilities > 1) $ do + putStrLn $ "\nNow "++ show numCapabilities ++" threads, reporting mean randoms-per-second-per-thread:" + void $ gamut numCapabilities + + putStrLn "Finished." + diff --git a/bench/Main.hs b/bench/Main.hs new file mode 100644 index 000000000..38904c5c7 --- /dev/null +++ b/bench/Main.hs @@ -0,0 +1,235 @@ +{-# LANGUAGE AllowAmbiguousTypes #-} +{-# LANGUAGE BangPatterns #-} +{-# LANGUAGE ExplicitForAll #-} +{-# LANGUAGE FlexibleContexts #-} +{-# LANGUAGE ScopedTypeVariables #-} +{-# LANGUAGE TypeApplications #-} + +module Main (main) where + +import Control.Monad +import Data.Int +import Data.Proxy +import Data.Typeable +import Data.Word +import Foreign.C.Types +import Gauge.Main +import Numeric.Natural (Natural) +import System.Random.SplitMix as SM + +import System.Random.Monad + +main :: IO () +main = do + let !sz = 100000 + genLengths = + -- create 5000 small lengths that are needed for ShortByteString generation + runStateGen (mkStdGen 2020) $ \g -> replicateM 5000 (uniformRM (16 + 1, 16 + 7) g) + defaultMain + [ bgroup "baseline" + [ let !smGen = SM.mkSMGen 1337 in bench "nextWord32" $ nf (genMany SM.nextWord32 smGen) sz + , let !smGen = SM.mkSMGen 1337 in bench "nextWord64" $ nf (genMany SM.nextWord64 smGen) sz + , let !smGen = SM.mkSMGen 1337 in bench "nextInt" $ nf (genMany SM.nextInt smGen) sz + , let !smGen = SM.mkSMGen 1337 in bench "split" $ nf (genMany SM.splitSMGen smGen) sz + ] + , bgroup "pure" + [ bgroup "random" + [ pureRandomBench @Float sz + , pureRandomBench @Double sz + , pureRandomBench @Integer sz + ] + , bgroup "uniform" + [ pureUniformBench @Word8 sz + , pureUniformBench @Word16 sz + , pureUniformBench @Word32 sz + , pureUniformBench @Word64 sz + , pureUniformBench @Word sz + , pureUniformBench @Int8 sz + , pureUniformBench @Int16 sz + , pureUniformBench @Int32 sz + , pureUniformBench @Int64 sz + , pureUniformBench @Int sz + , pureUniformBench @Char sz + , pureUniformBench @Bool sz + , pureUniformBench @CBool sz + , pureUniformBench @CChar sz + , pureUniformBench @CSChar sz + , pureUniformBench @CUChar sz + , pureUniformBench @CShort sz + , pureUniformBench @CUShort sz + , pureUniformBench @CInt sz + , pureUniformBench @CUInt sz + , pureUniformBench @CLong sz + , pureUniformBench @CULong sz + , pureUniformBench @CPtrdiff sz + , pureUniformBench @CSize sz + , pureUniformBench @CWchar sz + , pureUniformBench @CSigAtomic sz + , pureUniformBench @CLLong sz + , pureUniformBench @CULLong sz + , pureUniformBench @CIntPtr sz + , pureUniformBench @CUIntPtr sz + , pureUniformBench @CIntMax sz + , pureUniformBench @CUIntMax sz + ] + , bgroup "uniformR" + [ bgroup "full" + [ pureUniformRFullBench @Word8 sz + , pureUniformRFullBench @Word16 sz + , pureUniformRFullBench @Word32 sz + , pureUniformRFullBench @Word64 sz + , pureUniformRFullBench @Word sz + , pureUniformRFullBench @Int8 sz + , pureUniformRFullBench @Int16 sz + , pureUniformRFullBench @Int32 sz + , pureUniformRFullBench @Int64 sz + , pureUniformRFullBench @Int sz + , pureUniformRFullBench @Char sz + , pureUniformRFullBench @Bool sz + , pureUniformRFullBench @CBool sz + , pureUniformRFullBench @CChar sz + , pureUniformRFullBench @CSChar sz + , pureUniformRFullBench @CUChar sz + , pureUniformRFullBench @CShort sz + , pureUniformRFullBench @CUShort sz + , pureUniformRFullBench @CInt sz + , pureUniformRFullBench @CUInt sz + , pureUniformRFullBench @CLong sz + , pureUniformRFullBench @CULong sz + , pureUniformRFullBench @CPtrdiff sz + , pureUniformRFullBench @CSize sz + , pureUniformRFullBench @CWchar sz + , pureUniformRFullBench @CSigAtomic sz + , pureUniformRFullBench @CLLong sz + , pureUniformRFullBench @CULLong sz + , pureUniformRFullBench @CIntPtr sz + , pureUniformRFullBench @CUIntPtr sz + , pureUniformRFullBench @CIntMax sz + , pureUniformRFullBench @CUIntMax sz + ] + , bgroup "excludeMax" + [ pureUniformRExcludeMaxBench @Word8 sz + , pureUniformRExcludeMaxBench @Word16 sz + , pureUniformRExcludeMaxBench @Word32 sz + , pureUniformRExcludeMaxBench @Word64 sz + , pureUniformRExcludeMaxBench @Word sz + , pureUniformRExcludeMaxBench @Int8 sz + , pureUniformRExcludeMaxBench @Int16 sz + , pureUniformRExcludeMaxBench @Int32 sz + , pureUniformRExcludeMaxBench @Int64 sz + , pureUniformRExcludeMaxBench @Int sz + , pureUniformRExcludeMaxEnumBench @Char sz + , pureUniformRExcludeMaxEnumBench @Bool sz + , pureUniformRExcludeMaxBench @CBool sz + , pureUniformRExcludeMaxBench @CChar sz + , pureUniformRExcludeMaxBench @CSChar sz + , pureUniformRExcludeMaxBench @CUChar sz + , pureUniformRExcludeMaxBench @CShort sz + , pureUniformRExcludeMaxBench @CUShort sz + , pureUniformRExcludeMaxBench @CInt sz + , pureUniformRExcludeMaxBench @CUInt sz + , pureUniformRExcludeMaxBench @CLong sz + , pureUniformRExcludeMaxBench @CULong sz + , pureUniformRExcludeMaxBench @CPtrdiff sz + , pureUniformRExcludeMaxBench @CSize sz + , pureUniformRExcludeMaxBench @CWchar sz + , pureUniformRExcludeMaxBench @CSigAtomic sz + , pureUniformRExcludeMaxBench @CLLong sz + , pureUniformRExcludeMaxBench @CULLong sz + , pureUniformRExcludeMaxBench @CIntPtr sz + , pureUniformRExcludeMaxBench @CUIntPtr sz + , pureUniformRExcludeMaxBench @CIntMax sz + , pureUniformRExcludeMaxBench @CUIntMax sz + ] + , bgroup "includeHalf" + [ pureUniformRIncludeHalfBench @Word8 sz + , pureUniformRIncludeHalfBench @Word16 sz + , pureUniformRIncludeHalfBench @Word32 sz + , pureUniformRIncludeHalfBench @Word64 sz + , pureUniformRIncludeHalfBench @Word sz + , pureUniformRIncludeHalfBench @Int8 sz + , pureUniformRIncludeHalfBench @Int16 sz + , pureUniformRIncludeHalfBench @Int32 sz + , pureUniformRIncludeHalfBench @Int64 sz + , pureUniformRIncludeHalfBench @Int sz + , pureUniformRIncludeHalfEnumBench @Char sz + , pureUniformRIncludeHalfEnumBench @Bool sz + , pureUniformRIncludeHalfBench @CBool sz + , pureUniformRIncludeHalfBench @CChar sz + , pureUniformRIncludeHalfBench @CSChar sz + , pureUniformRIncludeHalfBench @CUChar sz + , pureUniformRIncludeHalfBench @CShort sz + , pureUniformRIncludeHalfBench @CUShort sz + , pureUniformRIncludeHalfBench @CInt sz + , pureUniformRIncludeHalfBench @CUInt sz + , pureUniformRIncludeHalfBench @CLong sz + , pureUniformRIncludeHalfBench @CULong sz + , pureUniformRIncludeHalfBench @CPtrdiff sz + , pureUniformRIncludeHalfBench @CSize sz + , pureUniformRIncludeHalfBench @CWchar sz + , pureUniformRIncludeHalfBench @CSigAtomic sz + , pureUniformRIncludeHalfBench @CLLong sz + , pureUniformRIncludeHalfBench @CULLong sz + , pureUniformRIncludeHalfBench @CIntPtr sz + , pureUniformRIncludeHalfBench @CUIntPtr sz + , pureUniformRIncludeHalfBench @CIntMax sz + , pureUniformRIncludeHalfBench @CUIntMax sz + ] + , bgroup "unbounded" + [ pureUniformRBench @Float (1.23e-4, 5.67e8) sz + , pureUniformRBench @Double (1.23e-4, 5.67e8) sz + , let !i = (10 :: Integer) ^ (100 :: Integer) + !range = (-i - 1, i + 1) + in pureUniformRBench @Integer range sz + , let !n = (10 :: Natural) ^ (100 :: Natural) + !range = (1, n - 1) + in pureUniformRBench @Natural range sz + ] + , bgroup "ShortByteString" + [ env (pure genLengths) $ \ ~(ns, gen) -> + bench "genShortByteString" $ + nfIO $ runStateGenT_ gen $ \g -> mapM (`uniformShortByteString` g) ns + ] + ] + ] + ] + +pureRandomBench :: forall a. (Typeable a, Random a) => Int -> Benchmark +pureRandomBench = let !stdGen = mkStdGen 1337 in pureBench @a (genMany (random @a) stdGen) + +pureUniformBench :: forall a. (Typeable a, Uniform a) => Int -> Benchmark +pureUniformBench = let !stdGen = mkStdGen 1337 in pureBench @a (genMany (uniform @_ @a) stdGen) + +pureUniformRFullBench :: forall a. (Typeable a, UniformRange a, Bounded a) => Int -> Benchmark +pureUniformRFullBench = let !range = (minBound @a, maxBound @a) in pureUniformRBench range + +pureUniformRExcludeMaxBench :: forall a. (Typeable a, UniformRange a, Bounded a, Num a) => Int -> Benchmark +pureUniformRExcludeMaxBench = let !range = (minBound @a, maxBound @a - 1) in pureUniformRBench range + +pureUniformRExcludeMaxEnumBench :: forall a. (Typeable a, UniformRange a, Bounded a, Enum a) => Int -> Benchmark +pureUniformRExcludeMaxEnumBench = let !range = (minBound @a, pred (maxBound @a)) in pureUniformRBench range + +pureUniformRIncludeHalfBench :: forall a. (Typeable a, UniformRange a, Bounded a, Integral a) => Int -> Benchmark +pureUniformRIncludeHalfBench = let !range = (minBound @a, (maxBound @a `div` 2) + 1) in pureUniformRBench range + +pureUniformRIncludeHalfEnumBench :: forall a. (Typeable a, UniformRange a, Bounded a, Enum a) => Int -> Benchmark +pureUniformRIncludeHalfEnumBench = + let !range = (succ (minBound @a), toEnum ((fromEnum (maxBound @a) `div` 2) + 1)) + in pureUniformRBench range + +pureUniformRBench :: forall a. (Typeable a, UniformRange a) => (a, a) -> Int -> Benchmark +pureUniformRBench range = + let !stdGen = mkStdGen 1337 + in pureBench @a (genMany (uniformR range) stdGen) + +pureBench :: forall a. (Typeable a) => (Int -> ()) -> Int -> Benchmark +pureBench f sz = bench (showsTypeRep (typeRep (Proxy :: Proxy a)) "") $ nf f sz + +genMany :: (g -> (a, g)) -> g -> Int -> () +genMany f g0 n = go g0 0 + where + go g i + | i < n = + case f g of + (x, g') -> x `seq` go g' (i + 1) + | otherwise = g `seq` () diff --git a/cabal.project b/cabal.project new file mode 100644 index 000000000..3f330dd76 --- /dev/null +++ b/cabal.project @@ -0,0 +1,2 @@ +packages: . +constraints: splitmix -random diff --git a/prologue.txt b/prologue.txt deleted file mode 100644 index ea0344b9c..000000000 --- a/prologue.txt +++ /dev/null @@ -1 +0,0 @@ -Random number library. diff --git a/random.cabal b/random.cabal index fd29840fb..6d6f69e9e 100644 --- a/random.cabal +++ b/random.cabal @@ -1,70 +1,174 @@ -name: random -version: 1.1 - - - - -license: BSD3 -license-file: LICENSE -maintainer: core-libraries-committee@haskell.org -bug-reports: https://github.com/haskell/random/issues -synopsis: random number library -category: System +cabal-version: >=1.10 +name: random +version: 1.2 +license: BSD3 +license-file: LICENSE +maintainer: core-libraries-committee@haskell.org +bug-reports: https://github.com/haskell/random/issues +synopsis: Pseudo-random number generation description: - This package provides a basic random number generation - library, including the ability to split random number - generators. - + This package provides basic pseudo-random number generation, including the + ability to split random number generators. + . + == "System.Random": pure pseudo-random number interface + . + In pure code, use 'System.Random.uniform' and 'System.Random.uniformR' from + "System.Random" to generate pseudo-random numbers with a pure pseudo-random + number generator like 'System.Random.StdGen'. + . + As an example, here is how you can simulate rolls of a six-sided die using + 'System.Random.uniformR': + . + >>> let roll = uniformR (1, 6) :: RandomGen g => g -> (Word8, g) + >>> let rolls = unfoldr (Just . roll) :: RandomGen g => g -> [Word8] + >>> let pureGen = mkStdGen 42 + >>> take 10 (rolls pureGen) :: [Word8] + [1,1,3,2,4,5,3,4,6,2] + . + See "System.Random" for more details. + . + == "System.Random.Monad": monadic pseudo-random number interface + . + In monadic code, use 'System.Random.Monad.uniformM' and + 'System.Random.Monad.uniformRM' from "System.Random.Monad" to generate + pseudo-random numbers with a monadic pseudo-random number generator, or + using a monadic adapter. + . + As an example, here is how you can simulate rolls of a six-sided die using + 'System.Random.Monad.uniformRM': + . + >>> let rollM = uniformRM (1, 6) :: MonadRandom g s m => g s -> m Word8 + >>> let pureGen = mkStdGen 42 + >>> runGenState_ pureGen (replicateM 10 . rollM) :: m [Word8] + [1,1,3,2,4,5,3,4,6,2] + . + The monadic adapter 'System.Random.Monad.runGenState_' is used here to lift + the pure pseudo-random number generator @pureGen@ into the + 'System.Random.Monad.MonadRandom' context. + . + The monadic interface can also be used with existing monadic pseudo-random + number generators. In this example, we use the one provided in the + package: + . + >>> import System.Random.MWC as MWC + >>> let rollM = uniformRM (1, 6) :: MonadRandom g s m => g s -> m Word8 + >>> monadicGen <- MWC.create + >>> (replicateM 10 . rollM) monadicGen :: m [Word8] + [2,3,6,6,4,4,3,1,5,4] + . + See "System.Random.Monad" for more details. + +category: System +build-type: Custom extra-source-files: - .travis.yml - README.md - CHANGELOG.md - .gitignore - .darcs-boring - - - -build-type: Simple --- cabal-version 1.8 needed because "the field 'build-depends: random' refers --- to a library which is defined within the same package" -cabal-version: >= 1.8 - - - -Library - exposed-modules: - System.Random - extensions: CPP - GHC-Options: -O2 - build-depends: base >= 3 && < 5, time + README.md + CHANGELOG.md source-repository head type: git - location: http://git.haskell.org/packages/random.git + location: https://github.com/haskell/random.git --- To run the Test-Suite: --- $ cabal configure --enable-tests --- $ cabal test --show-details=always --test-options="+RTS -M1M -RTS" +custom-setup + setup-depends: + base >=4.10 && <5, + Cabal >=1.10 && <3.3, + cabal-doctest >=1.0.6 && <1.1 -Test-Suite T7936 - type: exitcode-stdio-1.0 - main-is: T7936.hs - hs-source-dirs: tests - build-depends: base >= 3 && < 5, random - ghc-options: -rtsopts -O2 - -Test-Suite TestRandomRs - type: exitcode-stdio-1.0 - main-is: TestRandomRs.hs - hs-source-dirs: tests - build-depends: base >= 3 && < 5, random - ghc-options: -rtsopts -O2 - -- TODO. Why does the following not work? - --test-options: +RTS -M1M -RTS - -Test-Suite TestRandomIOs - type: exitcode-stdio-1.0 - main-is: TestRandomIOs.hs - hs-source-dirs: tests - build-depends: base >= 3 && < 5, random - ghc-options: -rtsopts -O2 +library + exposed-modules: + System.Random + System.Random.Internal + System.Random.Monad + + hs-source-dirs: src + default-language: Haskell2010 + ghc-options: + -Weverything -Wno-implicit-prelude -Wno-missing-import-lists + -Wno-missing-local-signatures -Wno-redundant-constraints + -Wno-unsafe + + build-depends: + base >=4.10 && <5, + bytestring >=0.10 && <0.11, + deepseq >=1.1 && <2, + mtl >=2.2 && <2.3, + splitmix >=0.0.3 && <0.1 + +test-suite legacy-test + type: exitcode-stdio-1.0 + main-is: Legacy.hs + hs-source-dirs: test-legacy + other-modules: + T7936 + TestRandomIOs + TestRandomRs + Random1283 + RangeTest + + default-language: Haskell2010 + ghc-options: -with-rtsopts=-M4M -Wno-deprecations + build-depends: + base >=4.10 && <5, + containers >=0.5 && <0.7, + random -any + +test-suite doctests + type: exitcode-stdio-1.0 + main-is: doctests.hs + hs-source-dirs: test + default-language: Haskell2010 + build-depends: + base >=4.10 && <5, + doctest >=0.15 && <0.17, + mwc-random >=0.13 && <0.15, + primitive >=0.6 && <0.8, + random -any, + unliftio >=0.2 && <0.3 + +test-suite spec + type: exitcode-stdio-1.0 + main-is: Spec.hs + hs-source-dirs: test + other-modules: + Spec.Range + Spec.Run + + default-language: Haskell2010 + ghc-options: -Wall + build-depends: + base >=4.10 && <5, + bytestring >=0.10 && <0.11, + random -any, + smallcheck >=1.1 && <1.2, + tasty >=1.0 && <1.3, + tasty-smallcheck >=0.8 && <0.9, + tasty-expected-failure >=0.11 && <0.12, + tasty-hunit >=0.10 && <0.11 + +benchmark legacy-bench + type: exitcode-stdio-1.0 + main-is: SimpleRNGBench.hs + hs-source-dirs: bench-legacy + other-modules: BinSearch + default-language: Haskell2010 + ghc-options: + -Wall -O2 -threaded -rtsopts -with-rtsopts=-N -Wno-deprecations + + build-depends: + base >=4.10 && <5, + random -any, + rdtsc -any, + split >=0.2 && <0.3, + time >=1.8 && <1.11 + +benchmark bench + type: exitcode-stdio-1.0 + main-is: Main.hs + hs-source-dirs: bench + default-language: Haskell2010 + ghc-options: -Wall -O2 + build-depends: + base >=4.10 && <5, + gauge >=0.2.3 && <0.3, + random -any, + splitmix >=0.0.3 && <0.1 diff --git a/src/System/Random.hs b/src/System/Random.hs new file mode 100644 index 000000000..120dd37da --- /dev/null +++ b/src/System/Random.hs @@ -0,0 +1,467 @@ +{-# LANGUAGE DefaultSignatures #-} +{-# LANGUAGE Trustworthy #-} + +-- | +-- Module : System.Random +-- Copyright : (c) The University of Glasgow 2001 +-- License : BSD-style (see the file LICENSE in the 'random' repository) +-- Maintainer : libraries@haskell.org +-- Stability : stable +-- +-- This library deals with the common task of pseudo-random number generation. +module System.Random + ( + -- * Introduction + -- $introduction + + -- * Usage + -- $usagepure + + -- * Pure number generator interface + -- $interfaces + RandomGen(..) + , uniform + , uniformR + , genByteString + , Random(..) + , Uniform + , UniformRange + -- ** Standard pseudo-random number generator + , StdGen + , mkStdGen + + -- ** Global standard pseudo-random number generator + -- $globalstdgen + , getStdRandom + , getStdGen + , setStdGen + , newStdGen + , randomIO + , randomRIO + + -- * Compatibility and reproducibility + -- ** Backwards compatibility and deprecations + -- $deprecations + + -- ** Reproducibility + -- $reproducibility + + -- * Notes for pseudo-random number generator implementors + -- ** How to implement 'RandomGen' + -- $implementrandomgen + + -- * References + -- $references + ) where + +import Control.Arrow +import Control.Monad.IO.Class +import Data.ByteString (ByteString) +import Data.Int +import Data.IORef +import Data.Word +import Foreign.C.Types +import GHC.Exts +import System.IO.Unsafe (unsafePerformIO) +import System.Random.Internal +import qualified System.Random.SplitMix as SM + +-- $introduction +-- +-- This module provides type classes and instances for the following concepts: +-- +-- [Pure pseudo-random number generators] 'RandomGen' is an interface to pure +-- pseudo-random number generators. +-- +-- 'StdGen', the standard pseudo-random number generator provided in this +-- library, is an instance of 'RandomGen'. It uses the SplitMix +-- implementation provided by the +-- package. +-- Programmers may, of course, supply their own instances of 'RandomGen'. +-- +-- $usagepure +-- +-- In pure code, use 'uniform' and 'uniformR' to generate pseudo-random values +-- with a pure pseudo-random number generator like 'StdGen'. +-- +-- >>> :{ +-- let rolls :: RandomGen g => Int -> g -> [Word8] +-- rolls n = take n . unfoldr (Just . uniformR (1, 6)) +-- pureGen = mkStdGen 137 +-- in +-- rolls 10 pureGen :: [Word8] +-- :} +-- [1,2,6,6,5,1,4,6,5,4] +-- +-- To run use a /monadic/ pseudo-random computation in pure code with a pure +-- pseudo-random number generator, use 'runGenState' and its variants. +-- +-- >>> :{ +-- let rollsM :: MonadRandom g s m => Int -> g s -> m [Word8] +-- rollsM n = replicateM n . uniformRM (1, 6) +-- pureGen = mkStdGen 137 +-- in +-- runStateGen_ pureGen (rollsM 10) :: [Word8] +-- :} +-- [1,2,6,6,5,1,4,6,5,4] + +------------------------------------------------------------------------------- +-- Pseudo-random number generator interfaces +------------------------------------------------------------------------------- + +-- $interfaces +-- +-- Pseudo-random number generators come in two flavours: /pure/ and /monadic/. +-- +-- ['RandomGen': pure pseudo-random number generators] These generators produce +-- a new pseudo-random value together with a new instance of the +-- pseudo-random number generator. +-- +-- Pure pseudo-random number generators should implement 'split' if they +-- are /splittable/, that is, if there is an efficient method to turn one +-- instance of the generator into two such that the pseudo-random numbers +-- produced by the two resulting generators are not correlated. See [1] for +-- some background on splittable pseudo-random generators. +-- +-- ['System.Random.Monad.MonadRandom': monadic pseudo-random number generators] +-- See "System.Random.Monad" module +-- + +-- | Generates a value uniformly distributed over all possible values of that +-- type. +-- +-- This is a pure version of 'System.Random.Monad.uniformM'. +-- +-- @since 1.2 +uniform :: (RandomGen g, Uniform a) => g -> (a, g) +uniform g = runStateGen g uniformM + +-- | Generates a value uniformly distributed over the provided range, which +-- is interpreted as inclusive in the lower and upper bound. +-- +-- * @uniformR (1 :: Int, 4 :: Int)@ generates values uniformly from the set +-- \(\{1,2,3,4\}\) +-- +-- * @uniformR (1 :: Float, 4 :: Float)@ generates values uniformly from the +-- set \(\{x\;|\;1 \le x \le 4\}\) +-- +-- The following law should hold to make the function always defined: +-- +-- > uniformR (a, b) = uniformR (b, a) +-- +-- This is a pure version of 'System.Random.Monad.uniformRM'. +-- +-- @since 1.2 +uniformR :: (RandomGen g, UniformRange a) => (a, a) -> g -> (a, g) +uniformR r g = runStateGen g (uniformRM r) + +-- | Generates a 'ByteString' of the specified size using a pure pseudo-random +-- number generator. See 'uniformByteString' for the monadic version. +-- +-- @since 1.2 +genByteString :: RandomGen g => Int -> g -> (ByteString, g) +genByteString n g = runStateGenST g (uniformByteString n) +{-# INLINE genByteString #-} + +-- | The class of types for which uniformly distributed values can be +-- generated. +-- +-- 'Random' exists primarily for backwards compatibility with version 1.1 of +-- this library. In new code, use the better specified 'Uniform' and +-- 'UniformRange' instead. +class Random a where + + -- | Takes a range /(lo,hi)/ and a pseudo-random number generator + -- /g/, and returns a pseudo-random value uniformly distributed over the + -- closed interval /[lo,hi]/, together with a new generator. It is unspecified + -- what happens if /lo>hi/. For continuous types there is no requirement + -- that the values /lo/ and /hi/ are ever produced, but they may be, + -- depending on the implementation and the interval. + {-# INLINE randomR #-} + randomR :: RandomGen g => (a, a) -> g -> (a, g) + default randomR :: (RandomGen g, UniformRange a) => (a, a) -> g -> (a, g) + randomR r g = runStateGen g (uniformRM r) + + -- | The same as 'randomR', but using a default range determined by the type: + -- + -- * For bounded types (instances of 'Bounded', such as 'Char'), + -- the range is normally the whole type. + -- + -- * For fractional types, the range is normally the semi-closed interval + -- @[0,1)@. + -- + -- * For 'Integer', the range is (arbitrarily) the range of 'Int'. + {-# INLINE random #-} + random :: RandomGen g => g -> (a, g) + default random :: (RandomGen g, Uniform a) => g -> (a, g) + random g = runStateGen g uniformM + + -- | Plural variant of 'randomR', producing an infinite list of + -- pseudo-random values instead of returning a new generator. + {-# INLINE randomRs #-} + randomRs :: RandomGen g => (a,a) -> g -> [a] + randomRs ival g = build (\cons _nil -> buildRandoms cons (randomR ival) g) + + -- | Plural variant of 'random', producing an infinite list of + -- pseudo-random values instead of returning a new generator. + {-# INLINE randoms #-} + randoms :: RandomGen g => g -> [a] + randoms g = build (\cons _nil -> buildRandoms cons random g) + + +-- | Produce an infinite list-equivalent of pseudo-random values. +{-# INLINE buildRandoms #-} +buildRandoms :: RandomGen g + => (a -> as -> as) -- ^ E.g. '(:)' but subject to fusion + -> (g -> (a,g)) -- ^ E.g. 'random' + -> g -- ^ A 'RandomGen' instance + -> as +buildRandoms cons rand = go + where + -- The seq fixes part of #4218 and also makes fused Core simpler. + go g = x `seq` (x `cons` go g') where (x,g') = rand g + +-- Generate values in the Int range +instance Random Integer where + random = first (toInteger :: Int -> Integer) . random +instance Random Int8 +instance Random Int16 +instance Random Int32 +instance Random Int64 +instance Random Int +instance Random Word +instance Random Word8 +instance Random Word16 +instance Random Word32 +instance Random Word64 +instance Random CBool +instance Random CChar +instance Random CSChar +instance Random CUChar +instance Random CShort +instance Random CUShort +instance Random CInt +instance Random CUInt +instance Random CLong +instance Random CULong +instance Random CPtrdiff +instance Random CSize +instance Random CWchar +instance Random CSigAtomic +instance Random CLLong +instance Random CULLong +instance Random CIntPtr +instance Random CUIntPtr +instance Random CIntMax +instance Random CUIntMax +instance Random CFloat where + randomR (CFloat l, CFloat h) = first CFloat . randomR (l, h) + random = first CFloat . random +instance Random CDouble where + randomR (CDouble l, CDouble h) = first CDouble . randomR (l, h) + random = first CDouble . random + +instance Random Char +instance Random Bool +instance Random Double where + randomR r g = runStateGen g (uniformRM r) + random g = runStateGen g (uniformRM (0, 1)) +instance Random Float where + randomR r g = runStateGen g (uniformRM r) + random g = runStateGen g (uniformRM (0, 1)) + +------------------------------------------------------------------------------- +-- Global pseudo-random number generator +------------------------------------------------------------------------------- + +-- $globalstdgen +-- +-- There is a single, implicit, global pseudo-random number generator of type +-- 'StdGen', held in a global variable maintained by the 'IO' monad. It is +-- initialised automatically in some system-dependent fashion. To get +-- deterministic behaviour, use 'setStdGen'. +-- +-- Note that 'mkStdGen' also gives deterministic behaviour without requiring an +-- 'IO' context. + +-- |Sets the global pseudo-random number generator. +setStdGen :: MonadIO m => StdGen -> m () +setStdGen = liftIO . writeIORef theStdGen + +-- |Gets the global pseudo-random number generator. +getStdGen :: MonadIO m => m StdGen +getStdGen = liftIO $ readIORef theStdGen + +theStdGen :: IORef StdGen +theStdGen = unsafePerformIO $ SM.initSMGen >>= newIORef . StdGen +{-# NOINLINE theStdGen #-} + +-- |Applies 'split' to the current global pseudo-random generator, +-- updates it with one of the results, and returns the other. +newStdGen :: MonadIO m => m StdGen +newStdGen = liftIO $ atomicModifyIORef' theStdGen split + +{- |Uses the supplied function to get a value from the current global +random generator, and updates the global generator with the new generator +returned by the function. For example, @rollDice@ gets a pseudo-random integer +between 1 and 6: + +> rollDice :: IO Int +> rollDice = getStdRandom (randomR (1,6)) + +-} +getStdRandom :: MonadIO m => (StdGen -> (a, StdGen)) -> m a +getStdRandom f = liftIO $ atomicModifyIORef' theStdGen (swap . f) + where swap (v, g) = (g, v) + + +-- | A variant of 'randomR' that uses the global pseudo-random number +-- generator. +randomRIO :: (Random a, MonadIO m) => (a, a) -> m a +randomRIO range = liftIO $ getStdRandom (randomR range) + +-- | A variant of 'random' that uses the global pseudo-random number +-- generator. +randomIO :: (Random a, MonadIO m) => m a +randomIO = liftIO $ getStdRandom random + +------------------------------------------------------------------------------- +-- Notes +------------------------------------------------------------------------------- + +-- $implementrandomgen +-- +-- Consider these points when writing a 'RandomGen' instance for a given pure +-- pseudo-random number generator: +-- +-- * If the pseudo-random number generator has a power-of-2 modulus, that is, +-- it natively outputs @2^n@ bits of randomness for some @n@, implement +-- 'genWord8', 'genWord16', 'genWord32' and 'genWord64'. See below for more +-- details. +-- +-- * If the pseudo-random number generator does not have a power-of-2 +-- modulus, implement 'next' and 'genRange'. See below for more details. +-- +-- * If the pseudo-random number generator is splittable, implement 'split'. +-- If there is no suitable implementation, 'split' should fail with a +-- helpful error message. +-- +-- === How to implement 'RandomGen' for a pseudo-random number generator with power-of-2 modulus +-- +-- Suppose you want to implement a [permuted congruential +-- generator](https://en.wikipedia.org/wiki/Permuted_congruential_generator). +-- +-- >>> data PCGen = PCGen !Word64 !Word64 +-- +-- It produces a full 'Word32' of randomness per iteration. +-- +-- >>> import Data.Bits +-- >>> :{ +-- let stepGen :: PCGen -> (Word32, PCGen) +-- stepGen (PCGen state inc) = let +-- newState = state * 6364136223846793005 + (inc .|. 1) +-- xorShifted = fromIntegral (((state `shiftR` 18) `xor` state) `shiftR` 27) :: Word32 +-- rot = fromIntegral (state `shiftR` 59) :: Word32 +-- out = (xorShifted `shiftR` (fromIntegral rot)) .|. (xorShifted `shiftL` fromIntegral ((-rot) .&. 31)) +-- in (out, PCGen newState inc) +-- :} +-- +-- >>> fst $ stepGen $ snd $ stepGen (PCGen 17 29) +-- 3288430965 +-- +-- You can make it an instance of 'RandomGen' as follows: +-- +-- >>> :{ +-- instance RandomGen PCGen where +-- genWord32 = stepGen +-- split _ = error "PCG is not splittable" +-- :} +-- +-- +-- === How to implement 'RandomGen' for a pseudo-random number generator without a power-of-2 modulus +-- +-- __We do not recommend you implement any new pseudo-random number generators without a power-of-2 modulus.__ +-- +-- Pseudo-random number generators without a power-of-2 modulus perform +-- /significantly worse/ than pseudo-random number generators with a power-of-2 +-- modulus with this library. This is because most functionality in this +-- library is based on generating and transforming uniformly pseudo-random +-- machine words, and generating uniformly pseudo-random machine words using a +-- pseudo-random number generator without a power-of-2 modulus is expensive. +-- +-- The pseudo-random number generator from +-- natively +-- generates an integer value in the range @[1, 2147483562]@. This is the +-- generator used by this library before it was replaced by SplitMix in version +-- 1.2. +-- +-- >>> data LegacyGen = LegacyGen !Int32 !Int32 +-- >>> :{ +-- let legacyNext :: LegacyGen -> (Int, LegacyGen) +-- legacyNext (LegacyGen s1 s2) = (fromIntegral z', LegacyGen s1'' s2'') where +-- z' = if z < 1 then z + 2147483562 else z +-- z = s1'' - s2'' +-- k = s1 `quot` 53668 +-- s1' = 40014 * (s1 - k * 53668) - k * 12211 +-- s1'' = if s1' < 0 then s1' + 2147483563 else s1' +-- k' = s2 `quot` 52774 +-- s2' = 40692 * (s2 - k' * 52774) - k' * 3791 +-- s2'' = if s2' < 0 then s2' + 2147483399 else s2' +-- :} +-- +-- You can make it an instance of 'RandomGen' as follows: +-- +-- >>> :{ +-- instance RandomGen LegacyGen where +-- next = legacyNext +-- genRange _ = (1, 2147483562) +-- split _ = error "Not implemented" +-- :} +-- +-- $deprecations +-- +-- Version 1.2 mostly maintains backwards compatibility with version 1.1. This +-- has a few consequences users should be aware of: +-- +-- * The type class 'Random' is only provided for backwards compatibility. +-- New code should use 'Uniform' and 'UniformRange' instead. +-- +-- * The methods 'next' and 'genRange' in 'RandomGen' are deprecated and only +-- provided for backwards compatibility. New instances of 'RandomGen' should +-- implement word-based methods instead. See below for more information +-- about how to write a 'RandomGen' instance. +-- +-- * This library provides instances for 'Random' for some unbounded types +-- for backwards compatibility. For an unbounded type, there is no way +-- to generate a value with uniform probability out of its entire domain, so +-- the 'random' implementation for unbounded types actually generates a +-- value based on some fixed range. +-- +-- For 'Integer', 'random' generates a value in the 'Int' range. For 'Float' +-- and 'Double', 'random' generates a floating point value in the range @[0, +-- 1)@. +-- +-- This library does not provide 'Uniform' instances for any unbounded +-- types. +-- +-- $reproducibility +-- +-- If you have two builds of a particular piece of code against this library, +-- any deterministic function call should give the same result in the two +-- builds if the builds are +-- +-- * compiled against the same major version of this library +-- * on the same architecture (32-bit or 64-bit) +-- +-- $references +-- +-- 1. Guy L. Steele, Jr., Doug Lea, and Christine H. Flood. 2014. Fast +-- splittable pseudorandom number generators. In Proceedings of the 2014 ACM +-- International Conference on Object Oriented Programming Systems Languages & +-- Applications (OOPSLA '14). ACM, New York, NY, USA, 453-472. DOI: +-- + +-- $setup +-- +-- >>> import Control.Monad (replicateM) +-- >>> import Data.List (unfoldr) diff --git a/src/System/Random/Internal.hs b/src/System/Random/Internal.hs new file mode 100644 index 000000000..f91120f79 --- /dev/null +++ b/src/System/Random/Internal.hs @@ -0,0 +1,999 @@ +{-# LANGUAGE BangPatterns #-} +{-# LANGUAGE CPP #-} +{-# LANGUAGE DefaultSignatures #-} +{-# LANGUAGE DerivingStrategies #-} +{-# LANGUAGE FlexibleContexts #-} +{-# LANGUAGE FlexibleInstances #-} +{-# LANGUAGE FunctionalDependencies #-} +{-# LANGUAGE GHCForeignImportPrim #-} +{-# LANGUAGE GeneralizedNewtypeDeriving #-} +{-# LANGUAGE MagicHash #-} +{-# LANGUAGE RankNTypes #-} +{-# LANGUAGE ScopedTypeVariables #-} +{-# LANGUAGE Trustworthy #-} +{-# LANGUAGE TypeFamilyDependencies #-} +{-# LANGUAGE UnboxedTuples #-} +{-# LANGUAGE UndecidableInstances #-} +{-# LANGUAGE UnliftedFFITypes #-} +{-# OPTIONS_HADDOCK hide, not-home #-} +#include "MachDeps.h" + +-- | +-- Module : System.Random.Internal +-- Copyright : (c) The University of Glasgow 2001 +-- License : BSD-style (see the file LICENSE in the 'random' repository) +-- Maintainer : libraries@haskell.org +-- Stability : stable +-- +-- This library deals with the common task of pseudo-random number generation. +module System.Random.Internal + (-- * Pure and monadic pseudo-random number generator interfaces + RandomGen(..) + , MonadRandom(..) + + -- ** Standard pseudo-random number generator + , StdGen(..) + , mkStdGen + + -- * Monadic adapters for pure pseudo-random number generators + -- ** Pure adapter + , StateGen(..) + , StateGenM(..) + , splitGen + , runStateGen + , runStateGen_ + , runStateGenT + , runStateGenT_ + , runStateGenST + + -- * Pseudo-random values of various types + , Uniform(..) + , UniformRange(..) + , uniformByteString + + -- * Generators for sequences of pseudo-random bytes + , genShortByteStringIO + , genShortByteStringST + ) where + +import Control.Arrow +import Control.DeepSeq (NFData) +import Control.Monad.IO.Class +import Control.Monad.ST +import Control.Monad.ST.Unsafe +import Control.Monad.State.Strict +import Data.Bits +import Data.ByteString.Builder.Prim (word64LE) +import Data.ByteString.Builder.Prim.Internal (runF) +import Data.ByteString.Internal (ByteString(PS)) +import Data.ByteString.Short.Internal (ShortByteString(SBS), fromShort) +import Data.Int +import Data.Kind (Type) +import Data.Word +import Foreign.C.Types +import Foreign.Ptr (plusPtr) +import Foreign.Storable (pokeByteOff) +import GHC.Exts +import GHC.ForeignPtr +import GHC.IO (IO(..)) +import GHC.Word +import Numeric.Natural (Natural) +import System.IO.Unsafe (unsafePerformIO) +import qualified System.Random.SplitMix as SM +import qualified System.Random.SplitMix32 as SM32 + +-- | 'RandomGen' is an interface to pure pseudo-random number generators. +-- +-- 'StdGen' is the standard 'RandomGen' instance provided by this library. +{-# DEPRECATED next "No longer used" #-} +{-# DEPRECATED genRange "No longer used" #-} +class RandomGen g where + {-# MINIMAL split,(genWord32|genWord64|(next,genRange)) #-} + -- | Returns an 'Int' that is uniformly distributed over the range returned by + -- 'genRange' (including both end points), and a new generator. Using 'next' + -- is inefficient as all operations go via 'Integer'. See + -- [here](https://alexey.kuleshevi.ch/blog/2019/12/21/random-benchmarks) for + -- more details. It is thus deprecated. + next :: g -> (Int, g) + next g = runStateGen g (uniformRM (genRange g)) + + -- | Returns a 'Word8' that is uniformly distributed over the entire 'Word8' + -- range. + -- + -- @since 1.2 + genWord8 :: g -> (Word8, g) + genWord8 = first fromIntegral . genWord32 + + -- | Returns a 'Word16' that is uniformly distributed over the entire 'Word16' + -- range. + -- + -- @since 1.2 + genWord16 :: g -> (Word16, g) + genWord16 = first fromIntegral . genWord32 + + -- | Returns a 'Word32' that is uniformly distributed over the entire 'Word32' + -- range. + -- + -- @since 1.2 + genWord32 :: g -> (Word32, g) + genWord32 = randomIvalIntegral (minBound, maxBound) + -- Once `next` is removed, this implementation should be used instead: + -- first fromIntegral . genWord64 + + -- | Returns a 'Word64' that is uniformly distributed over the entire 'Word64' + -- range. + -- + -- @since 1.2 + genWord64 :: g -> (Word64, g) + genWord64 g = + case genWord32 g of + (l32, g') -> + case genWord32 g' of + (h32, g'') -> + ((fromIntegral h32 `unsafeShiftL` 32) .|. fromIntegral l32, g'') + + -- | @genWord32R upperBound g@ returns a 'Word32' that is uniformly + -- distributed over the range @[0, upperBound]@. + -- + -- @since 1.2 + genWord32R :: Word32 -> g -> (Word32, g) + genWord32R m g = runStateGen g (unbiasedWordMult32 m) + + -- | @genWord64R upperBound g@ returns a 'Word64' that is uniformly + -- distributed over the range @[0, upperBound]@. + -- + -- @since 1.2 + genWord64R :: Word64 -> g -> (Word64, g) + genWord64R m g = runStateGen g (unsignedBitmaskWithRejectionM uniformWord64 m) + + -- | @genShortByteString n g@ returns a 'ShortByteString' of length @n@ + -- filled with pseudo-random bytes. + -- + -- @since 1.2 + genShortByteString :: Int -> g -> (ShortByteString, g) + genShortByteString n g = + unsafePerformIO $ runStateGenT g (genShortByteStringIO n . uniformWord64) + {-# INLINE genShortByteString #-} + + -- | Yields the range of values returned by 'next'. + -- + -- It is required that: + -- + -- * If @(a, b) = 'genRange' g@, then @a < b@. + -- * 'genRange' must not examine its argument so the value it returns is + -- determined only by the instance of 'RandomGen'. + -- + -- The default definition spans the full range of 'Int'. + genRange :: g -> (Int, Int) + genRange _ = (minBound, maxBound) + + -- | Returns two distinct pseudo-random number generators. + -- + -- Implementations should take care to ensure that the resulting generators + -- are not correlated. Some pseudo-random number generators are not + -- splittable. In that case, the 'split' implementation should fail with a + -- descriptive 'error' message. + split :: g -> (g, g) + + +-- | 'MonadRandom' is an interface to monadic pseudo-random number generators. +class Monad m => MonadRandom g s m | g m -> s where + -- | Represents the state of the pseudo-random number generator for use with + -- 'thawGen' and 'freezeGen'. + -- + -- @since 1.2 + type Frozen g = (f :: Type) | f -> g + {-# MINIMAL freezeGen,thawGen,(uniformWord32|uniformWord64) #-} + + -- | Restores the pseudo-random number generator from its 'Frozen' + -- representation. + -- + -- @since 1.2 + thawGen :: Frozen g -> m (g s) + + -- | Saves the state of the pseudo-random number generator to its 'Frozen' + -- representation. + -- + -- @since 1.2 + freezeGen :: g s -> m (Frozen g) + + -- | @uniformWord32R upperBound g@ generates a 'Word32' that is uniformly + -- distributed over the range @[0, upperBound]@. + -- + -- @since 1.2 + uniformWord32R :: Word32 -> g s -> m Word32 + uniformWord32R = unsignedBitmaskWithRejectionM uniformWord32 + + -- | @uniformWord64R upperBound g@ generates a 'Word64' that is uniformly + -- distributed over the range @[0, upperBound]@. + -- + -- @since 1.2 + uniformWord64R :: Word64 -> g s -> m Word64 + uniformWord64R = unsignedBitmaskWithRejectionM uniformWord64 + + -- | Generates a 'Word8' that is uniformly distributed over the entire 'Word8' + -- range. + -- + -- The default implementation extracts a 'Word8' from 'uniformWord32'. + -- + -- @since 1.2 + uniformWord8 :: g s -> m Word8 + uniformWord8 = fmap fromIntegral . uniformWord32 + + -- | Generates a 'Word16' that is uniformly distributed over the entire + -- 'Word16' range. + -- + -- The default implementation extracts a 'Word16' from 'uniformWord32'. + -- + -- @since 1.2 + uniformWord16 :: g s -> m Word16 + uniformWord16 = fmap fromIntegral . uniformWord32 + + -- | Generates a 'Word32' that is uniformly distributed over the entire + -- 'Word32' range. + -- + -- The default implementation extracts a 'Word32' from 'uniformWord64'. + -- + -- @since 1.2 + uniformWord32 :: g s -> m Word32 + uniformWord32 = fmap fromIntegral . uniformWord64 + + -- | Generates a 'Word64' that is uniformly distributed over the entire + -- 'Word64' range. + -- + -- The default implementation combines two 'Word32' from 'uniformWord32' into + -- one 'Word64'. + -- + -- @since 1.2 + uniformWord64 :: g s -> m Word64 + uniformWord64 g = do + l32 <- uniformWord32 g + h32 <- uniformWord32 g + pure (unsafeShiftL (fromIntegral h32) 32 .|. fromIntegral l32) + + -- | @uniformShortByteString n g@ generates a 'ShortByteString' of length @n@ + -- filled with pseudo-random bytes. + -- + -- @since 1.2 + uniformShortByteString :: Int -> g s -> m ShortByteString + default uniformShortByteString :: MonadIO m => Int -> g s -> m ShortByteString + uniformShortByteString n = genShortByteStringIO n . uniformWord64 + {-# INLINE uniformShortByteString #-} + + + +data MBA s = MBA (MutableByteArray# s) + + +-- | Efficiently generates a sequence of pseudo-random bytes in a platform +-- independent manner. The allocated memory is be pinned, so it is safe to use +-- with FFI calls. +-- +-- @since 1.2 +genShortByteStringIO :: MonadIO m => Int -> m Word64 -> m ShortByteString +genShortByteStringIO n0 gen64 = do + let !n@(I# n#) = max 0 n0 + (n64, nrem64) = n `quotRem` 8 + MBA mba# <- + liftIO $ + IO $ \s# -> + case newPinnedByteArray# n# s# of + (# s'#, mba# #) -> (# s'#, MBA mba# #) + let go i ptr + | i < n64 = do + w64 <- gen64 + -- Writing 8 bytes at a time in a Little-endian order gives us + -- platform portability + liftIO $ runF word64LE w64 ptr + go (i + 1) (ptr `plusPtr` 8) + | otherwise = return ptr + ptr <- go 0 (Ptr (byteArrayContents# (unsafeCoerce# mba#))) + when (nrem64 > 0) $ do + w64 <- gen64 + -- In order to not mess up the byte order we write generated Word64 into a + -- temporary pointer and then copy only the missing bytes over to the array. + -- It is tempting to simply generate as many bytes as we still need using + -- smaller generators (eg. uniformWord8), but that would result in + -- inconsistent tail when total length is slightly varied. + liftIO $ do + let goRem64 z i = + when (i < nrem64) $ do + pokeByteOff ptr i (fromIntegral z :: Word8) + goRem64 (z `unsafeShiftR` 8) (i + 1) + goRem64 w64 0 + liftIO $ + IO $ \s# -> + case unsafeFreezeByteArray# mba# s# of + (# s'#, ba# #) -> (# s'#, SBS ba# #) +{-# INLINE genShortByteStringIO #-} + +-- | Same as 'genShortByteStringIO', but runs in 'ST'. +-- +-- @since 1.2 +genShortByteStringST :: Int -> ST s Word64 -> ST s ShortByteString +genShortByteStringST n action = + unsafeIOToST (genShortByteStringIO n (unsafeSTToIO action)) + +pinnedByteArrayToByteString :: ByteArray# -> ByteString +pinnedByteArrayToByteString ba# = + PS (pinnedByteArrayToForeignPtr ba#) 0 (I# (sizeofByteArray# ba#)) +{-# INLINE pinnedByteArrayToByteString #-} + +pinnedByteArrayToForeignPtr :: ByteArray# -> ForeignPtr a +pinnedByteArrayToForeignPtr ba# = + ForeignPtr (byteArrayContents# ba#) (PlainPtr (unsafeCoerce# ba#)) +{-# INLINE pinnedByteArrayToForeignPtr #-} + + +-- | Generates a pseudo-random 'ByteString' of the specified size. +-- +-- @since 1.2 +uniformByteString :: MonadRandom g s m => Int -> g s -> m ByteString +uniformByteString n g = do + ba@(SBS ba#) <- uniformShortByteString n g + pure $ + if isTrue# (isByteArrayPinned# ba#) + then pinnedByteArrayToByteString ba# + else fromShort ba +{-# INLINE uniformByteString #-} + + +-- | Opaque data type that carries the type of a pure pseudo-random number +-- generator. +-- +-- @since 1.2 +data StateGenM g s = StateGenM +newtype StateGen g = StateGen g + +instance (RandomGen g, MonadState g m) => MonadRandom (StateGenM g) g m where + type Frozen (StateGenM g) = StateGen g + thawGen (StateGen g) = StateGenM <$ put g + freezeGen _ = fmap StateGen get + uniformWord32R r _ = state (genWord32R r) + uniformWord64R r _ = state (genWord64R r) + uniformWord8 _ = state genWord8 + uniformWord16 _ = state genWord16 + uniformWord32 _ = state genWord32 + uniformWord64 _ = state genWord64 + uniformShortByteString n _ = state (genShortByteString n) + + + +-- | Splits a pseudo-random number generator into two. Updates the state with +-- one of the resulting generators and returns the other. +-- +-- @since 1.2 +splitGen :: (MonadState g m, RandomGen g) => m g +splitGen = state split + +-- | Runs a monadic generating action in the `State` monad using a pure +-- pseudo-random number generator. +-- +-- @since 1.2 +runStateGen :: RandomGen g => g -> (StateGenM g g -> State g a) -> (a, g) +runStateGen g f = runState (f StateGenM) g + +-- | Runs a monadic generating action in the `State` monad using a pure +-- pseudo-random number generator. Returns only the resulting pseudo-random +-- value. +-- +-- @since 1.2 +runStateGen_ :: RandomGen g => g -> (StateGenM g g -> State g a) -> a +runStateGen_ g = fst . runStateGen g + +-- | Runs a monadic generating action in the `StateT` monad using a pure +-- pseudo-random number generator. +-- +-- @since 1.2 +runStateGenT :: RandomGen g => g -> (StateGenM g g -> StateT g m a) -> m (a, g) +runStateGenT g f = runStateT (f StateGenM) g + +-- | Runs a monadic generating action in the `StateT` monad using a pure +-- pseudo-random number generator. Returns only the resulting pseudo-random +-- value. +-- +-- @since 1.2 +runStateGenT_ :: (RandomGen g, Functor f) => g -> (StateGenM g g -> StateT g f a) -> f a +runStateGenT_ g = fmap fst . runStateGenT g + +-- | Runs a monadic generating action in the `ST` monad using a pure +-- pseudo-random number generator. +-- +-- @since 1.2 +runStateGenST :: RandomGen g => g -> (forall s . StateGenM g g -> StateT g (ST s) a) -> (a, g) +runStateGenST g action = runST $ runStateGenT g action +{-# INLINE runStateGenST #-} + + +-- | The standard pseudo-random number generator. +newtype StdGen = StdGen { unStdGen :: SM.SMGen } + deriving newtype (RandomGen, NFData) + deriving stock (Show) + +instance Eq StdGen where + StdGen x1 == StdGen x2 = SM.unseedSMGen x1 == SM.unseedSMGen x2 + +instance RandomGen SM.SMGen where + next = SM.nextInt + genWord32 = SM.nextWord32 + genWord64 = SM.nextWord64 + split = SM.splitSMGen + +instance RandomGen SM32.SMGen where + next = SM32.nextInt + genWord32 = SM32.nextWord32 + genWord64 = SM32.nextWord64 + split = SM32.splitSMGen + +-- | Constructs a 'StdGen' deterministically. +mkStdGen :: Int -> StdGen +mkStdGen = StdGen . SM.mkSMGen . fromIntegral + +-- | The class of types for which a uniformly distributed value can be drawn +-- from all possible values of the type. +-- +-- @since 1.2 +class Uniform a where + -- | Generates a value uniformly distributed over all possible values of that + -- type. + -- + -- @since 1.2 + uniformM :: MonadRandom g s m => g s -> m a + +-- | The class of types for which a uniformly distributed value can be drawn +-- from a range. +-- +-- @since 1.2 +class UniformRange a where + -- | Generates a value uniformly distributed over the provided range, which + -- is interpreted as inclusive in the lower and upper bound. + -- + -- * @uniformRM (1 :: Int, 4 :: Int)@ generates values uniformly from the + -- set \(\{1,2,3,4\}\) + -- + -- * @uniformRM (1 :: Float, 4 :: Float)@ generates values uniformly from + -- the set \(\{x\;|\;1 \le x \le 4\}\) + -- + -- The following law should hold to make the function always defined: + -- + -- > uniformRM (a, b) = uniformRM (b, a) + -- + -- @since 1.2 + uniformRM :: MonadRandom g s m => (a, a) -> g s -> m a + +instance UniformRange Integer where + uniformRM = uniformIntegralM + +instance UniformRange Natural where + uniformRM = uniformIntegralM + +instance Uniform Int8 where + uniformM = fmap (fromIntegral :: Word8 -> Int8) . uniformWord8 +instance UniformRange Int8 where + uniformRM = signedBitmaskWithRejectionRM (fromIntegral :: Int8 -> Word8) fromIntegral + +instance Uniform Int16 where + uniformM = fmap (fromIntegral :: Word16 -> Int16) . uniformWord16 +instance UniformRange Int16 where + uniformRM = signedBitmaskWithRejectionRM (fromIntegral :: Int16 -> Word16) fromIntegral + {-# INLINE uniformRM #-} + +instance Uniform Int32 where + uniformM = fmap (fromIntegral :: Word32 -> Int32) . uniformWord32 +instance UniformRange Int32 where + uniformRM = signedBitmaskWithRejectionRM (fromIntegral :: Int32 -> Word32) fromIntegral + {-# INLINE uniformRM #-} + +instance Uniform Int64 where + uniformM = fmap (fromIntegral :: Word64 -> Int64) . uniformWord64 +instance UniformRange Int64 where + uniformRM = signedBitmaskWithRejectionRM (fromIntegral :: Int64 -> Word64) fromIntegral + {-# INLINE uniformRM #-} + +instance Uniform Int where +#if WORD_SIZE_IN_BITS < 64 + uniformM = fmap (fromIntegral :: Word32 -> Int) . uniformWord32 +#else + uniformM = fmap (fromIntegral :: Word64 -> Int) . uniformWord64 +#endif + {-# INLINE uniformM #-} +instance UniformRange Int where + uniformRM = signedBitmaskWithRejectionRM (fromIntegral :: Int -> Word) fromIntegral + {-# INLINE uniformRM #-} + +instance Uniform Word where +#if WORD_SIZE_IN_BITS < 64 + uniformM = fmap (fromIntegral :: Word32 -> Word) . uniformWord32 +#else + uniformM = fmap (fromIntegral :: Word64 -> Word) . uniformWord64 +#endif +instance UniformRange Word where + {-# INLINE uniformRM #-} + uniformRM = unsignedBitmaskWithRejectionRM + +instance Uniform Word8 where + {-# INLINE uniformM #-} + uniformM = uniformWord8 +instance UniformRange Word8 where + {-# INLINE uniformRM #-} + uniformRM = unbiasedWordMult32RM + +instance Uniform Word16 where + {-# INLINE uniformM #-} + uniformM = uniformWord16 +instance UniformRange Word16 where + {-# INLINE uniformRM #-} + uniformRM = unbiasedWordMult32RM + +instance Uniform Word32 where + {-# INLINE uniformM #-} + uniformM = uniformWord32 +instance UniformRange Word32 where + {-# INLINE uniformRM #-} + uniformRM = unbiasedWordMult32RM + +instance Uniform Word64 where + {-# INLINE uniformM #-} + uniformM = uniformWord64 +instance UniformRange Word64 where + {-# INLINE uniformRM #-} + uniformRM = unsignedBitmaskWithRejectionRM + +instance Uniform CBool where + uniformM = fmap CBool . uniformM +instance UniformRange CBool where + uniformRM (CBool b, CBool t) = fmap CBool . uniformRM (b, t) + {-# INLINE uniformRM #-} + +instance Uniform CChar where + uniformM = fmap CChar . uniformM +instance UniformRange CChar where + uniformRM (CChar b, CChar t) = fmap CChar . uniformRM (b, t) + {-# INLINE uniformRM #-} + +instance Uniform CSChar where + uniformM = fmap CSChar . uniformM +instance UniformRange CSChar where + uniformRM (CSChar b, CSChar t) = fmap CSChar . uniformRM (b, t) + {-# INLINE uniformRM #-} + +instance Uniform CUChar where + uniformM = fmap CUChar . uniformM +instance UniformRange CUChar where + uniformRM (CUChar b, CUChar t) = fmap CUChar . uniformRM (b, t) + {-# INLINE uniformRM #-} + +instance Uniform CShort where + uniformM = fmap CShort . uniformM +instance UniformRange CShort where + uniformRM (CShort b, CShort t) = fmap CShort . uniformRM (b, t) + {-# INLINE uniformRM #-} + +instance Uniform CUShort where + uniformM = fmap CUShort . uniformM +instance UniformRange CUShort where + uniformRM (CUShort b, CUShort t) = fmap CUShort . uniformRM (b, t) + {-# INLINE uniformRM #-} + +instance Uniform CInt where + uniformM = fmap CInt . uniformM +instance UniformRange CInt where + uniformRM (CInt b, CInt t) = fmap CInt . uniformRM (b, t) + {-# INLINE uniformRM #-} + +instance Uniform CUInt where + uniformM = fmap CUInt . uniformM +instance UniformRange CUInt where + uniformRM (CUInt b, CUInt t) = fmap CUInt . uniformRM (b, t) + {-# INLINE uniformRM #-} + +instance Uniform CLong where + uniformM = fmap CLong . uniformM +instance UniformRange CLong where + uniformRM (CLong b, CLong t) = fmap CLong . uniformRM (b, t) + {-# INLINE uniformRM #-} + +instance Uniform CULong where + uniformM = fmap CULong . uniformM +instance UniformRange CULong where + uniformRM (CULong b, CULong t) = fmap CULong . uniformRM (b, t) + {-# INLINE uniformRM #-} + +instance Uniform CPtrdiff where + uniformM = fmap CPtrdiff . uniformM +instance UniformRange CPtrdiff where + uniformRM (CPtrdiff b, CPtrdiff t) = fmap CPtrdiff . uniformRM (b, t) + {-# INLINE uniformRM #-} + +instance Uniform CSize where + uniformM = fmap CSize . uniformM +instance UniformRange CSize where + uniformRM (CSize b, CSize t) = fmap CSize . uniformRM (b, t) + {-# INLINE uniformRM #-} + +instance Uniform CWchar where + uniformM = fmap CWchar . uniformM +instance UniformRange CWchar where + uniformRM (CWchar b, CWchar t) = fmap CWchar . uniformRM (b, t) + {-# INLINE uniformRM #-} + +instance Uniform CSigAtomic where + uniformM = fmap CSigAtomic . uniformM +instance UniformRange CSigAtomic where + uniformRM (CSigAtomic b, CSigAtomic t) = fmap CSigAtomic . uniformRM (b, t) + {-# INLINE uniformRM #-} + +instance Uniform CLLong where + uniformM = fmap CLLong . uniformM +instance UniformRange CLLong where + uniformRM (CLLong b, CLLong t) = fmap CLLong . uniformRM (b, t) + {-# INLINE uniformRM #-} + +instance Uniform CULLong where + uniformM = fmap CULLong . uniformM +instance UniformRange CULLong where + uniformRM (CULLong b, CULLong t) = fmap CULLong . uniformRM (b, t) + {-# INLINE uniformRM #-} + +instance Uniform CIntPtr where + uniformM = fmap CIntPtr . uniformM +instance UniformRange CIntPtr where + uniformRM (CIntPtr b, CIntPtr t) = fmap CIntPtr . uniformRM (b, t) + {-# INLINE uniformRM #-} + +instance Uniform CUIntPtr where + uniformM = fmap CUIntPtr . uniformM +instance UniformRange CUIntPtr where + uniformRM (CUIntPtr b, CUIntPtr t) = fmap CUIntPtr . uniformRM (b, t) + {-# INLINE uniformRM #-} + +instance Uniform CIntMax where + uniformM = fmap CIntMax . uniformM +instance UniformRange CIntMax where + uniformRM (CIntMax b, CIntMax t) = fmap CIntMax . uniformRM (b, t) + {-# INLINE uniformRM #-} + +instance Uniform CUIntMax where + uniformM = fmap CUIntMax . uniformM +instance UniformRange CUIntMax where + uniformRM (CUIntMax b, CUIntMax t) = fmap CUIntMax . uniformRM (b, t) + {-# INLINE uniformRM #-} + +instance UniformRange CFloat where + uniformRM (CFloat l, CFloat h) = fmap CFloat . uniformRM (l, h) + {-# INLINE uniformRM #-} + +instance UniformRange CDouble where + uniformRM (CDouble l, CDouble h) = fmap CDouble . uniformRM (l, h) + {-# INLINE uniformRM #-} + + +-- The `chr#` and `ord#` are the prim functions that will be called, regardless of which +-- way you gonna do the `Char` conversion, so it is better to call them directly and +-- bypass all the hoops. Also because `intToChar` and `charToInt` are internal functions +-- and are called on valid character ranges it is impossible to generate an invalid +-- `Char`, therefore it is totally fine to omit all the unnecessary checks involved in +-- other paths of conversion. +word32ToChar :: Word32 -> Char +word32ToChar (W32# w#) = C# (chr# (word2Int# w#)) +{-# INLINE word32ToChar #-} + +charToWord32 :: Char -> Word32 +charToWord32 (C# c#) = W32# (int2Word# (ord# c#)) +{-# INLINE charToWord32 #-} + +instance Uniform Char where + uniformM g = word32ToChar <$> unbiasedWordMult32 (charToWord32 maxBound) g + {-# INLINE uniformM #-} +instance UniformRange Char where + uniformRM (l, h) g = + word32ToChar <$> unbiasedWordMult32RM (charToWord32 l, charToWord32 h) g + {-# INLINE uniformRM #-} + +instance Uniform Bool where + uniformM = fmap wordToBool . uniformWord8 + where wordToBool w = (w .&. 1) /= 0 +instance UniformRange Bool where + uniformRM (False, False) _g = return False + uniformRM (True, True) _g = return True + uniformRM _ g = uniformM g + +-- | See /Floating point number caveats/ in "System.Random.Monad". +instance UniformRange Double where + uniformRM (l, h) g = do + w64 <- uniformWord64 g + let x = word64ToDoubleInUnitInterval w64 + return $ (h - l) * x + l + +-- | Turns a given uniformly distributed 'Word64' value into a uniformly +-- distributed 'Double' value in the range [0, 1]. +word64ToDoubleInUnitInterval :: Word64 -> Double +word64ToDoubleInUnitInterval w64 = d / m + where + d = fromIntegral w64 :: Double + m = fromIntegral (maxBound :: Word64) :: Double +{-# INLINE word64ToDoubleInUnitInterval #-} + +-- | See /Floating point number caveats/ in "System.Random.Monad". +instance UniformRange Float where + uniformRM (l, h) g = do + w32 <- uniformWord32 g + let x = word32ToFloatInUnitInterval w32 + return $ (h - l) * x + l + +-- | Turns a given uniformly distributed 'Word32' value into a uniformly +-- distributed 'Float' value in the range [0,1]. +word32ToFloatInUnitInterval :: Word32 -> Float +word32ToFloatInUnitInterval w32 = f / m + where + f = fromIntegral w32 :: Float + m = fromIntegral (maxBound :: Word32) :: Float +{-# INLINE word32ToFloatInUnitInterval #-} + +-- The two integer functions below take an [inclusive,inclusive] range. +randomIvalIntegral :: (RandomGen g, Integral a) => (a, a) -> g -> (a, g) +randomIvalIntegral (l,h) = randomIvalInteger (toInteger l, toInteger h) + +{-# SPECIALIZE randomIvalInteger :: (Num a) => + (Integer, Integer) -> StdGen -> (a, StdGen) #-} + +randomIvalInteger :: (RandomGen g, Num a) => (Integer, Integer) -> g -> (a, g) +randomIvalInteger (l,h) rng + | l > h = randomIvalInteger (h,l) rng + | otherwise = case f 1 0 rng of (v, rng') -> (fromInteger (l + v `mod` k), rng') + where + (genlo, genhi) = genRange rng + b = fromIntegral genhi - fromIntegral genlo + 1 + + -- Probabilities of the most likely and least likely result + -- will differ at most by a factor of (1 +- 1/q). Assuming the RandomGen + -- is uniform, of course + + -- On average, log q / log b more pseudo-random values will be generated + -- than the minimum + q = 1000 :: Integer + k = h - l + 1 + magtgt = k * q + + -- generate pseudo-random values until we exceed the target magnitude + f mag v g | mag >= magtgt = (v, g) + | otherwise = v' `seq`f (mag*b) v' g' where + (x,g') = next g + v' = v * b + (fromIntegral x - fromIntegral genlo) + +-- | Generate an integral in the range @[l, h]@ if @l <= h@ and @[h, l]@ +-- otherwise. +uniformIntegralM :: (Bits a, Integral a, MonadRandom g s m) => (a, a) -> g s -> m a +uniformIntegralM (l, h) gen = case l `compare` h of + LT -> do + let limit = h - l + let limitAsWord64 :: Word64 = fromIntegral limit + bounded <- + if fromIntegral limitAsWord64 == limit + -- Optimisation: if 'limit' fits into 'Word64', generate a bounded + -- 'Word64' and then convert to 'Integer' + then fromIntegral <$> unsignedBitmaskWithRejectionM uniformWord64 limitAsWord64 gen + else boundedExclusiveIntegralM (limit + 1) gen + return $ l + bounded + GT -> uniformIntegralM (h, l) gen + EQ -> pure l +{-# INLINEABLE uniformIntegralM #-} + +-- | Generate an integral in the range @[0, s)@ using a variant of Lemire's +-- multiplication method. +-- +-- Daniel Lemire. 2019. Fast Random Integer Generation in an Interval. In ACM +-- Transactions on Modeling and Computer Simulation +-- https://doi.org/10.1145/3230636 +-- +-- PRECONDITION (unchecked): s > 0 +boundedExclusiveIntegralM :: (Bits a, Integral a, MonadRandom g s m) => a -> g s -> m a +boundedExclusiveIntegralM (s :: a) gen = go + where + n = integralWordSize s + -- We renamed 'L' from the paper to 'k' here because 'L' is not a valid + -- variable name in Haskell and 'l' is already used in the algorithm. + k = WORD_SIZE_IN_BITS * n + twoToK = (1::a) `shiftL` k + modTwoToKMask = twoToK - 1 + + t = (twoToK - s) `mod` s + go = do + x <- uniformIntegralWords n gen + let m = x * s + -- m .&. modTwoToKMask == m `mod` twoToK + let l = m .&. modTwoToKMask + if l < t + then go + -- m `shiftR` k == m `quot` twoToK + else return $ m `shiftR` k +{-# INLINE boundedExclusiveIntegralM #-} + +-- | @integralWordSize i@ returns that least @w@ such that +-- @i <= WORD_SIZE_IN_BITS^w@. +integralWordSize :: (Bits a, Num a) => a -> Int +integralWordSize = go 0 + where + go !acc i + | i == 0 = acc + | otherwise = go (acc + 1) (i `shiftR` WORD_SIZE_IN_BITS) +{-# INLINE integralWordSize #-} + +-- | @uniformIntegralWords n@ is a uniformly pseudo-random integral in the range +-- @[0, WORD_SIZE_IN_BITS^n)@. +uniformIntegralWords :: (Bits a, Integral a, MonadRandom g s m) => Int -> g s -> m a +uniformIntegralWords n gen = go 0 n + where + go !acc i + | i == 0 = return acc + | otherwise = do + (w :: Word) <- uniformM gen + go ((acc `shiftL` WORD_SIZE_IN_BITS) .|. fromIntegral w) (i - 1) +{-# INLINE uniformIntegralWords #-} + +-- | Uniformly generate an 'Integral' in an inclusive-inclusive range. +-- +-- Only use for integrals size less than or equal to that of 'Word32'. +unbiasedWordMult32RM :: (MonadRandom g s m, Integral a) => (a, a) -> g s -> m a +unbiasedWordMult32RM (b, t) g + | b <= t = (+b) . fromIntegral <$> unbiasedWordMult32 (fromIntegral (t - b)) g + | otherwise = (+t) . fromIntegral <$> unbiasedWordMult32 (fromIntegral (b - t)) g +{-# SPECIALIZE unbiasedWordMult32RM :: MonadRandom g s m => (Word8, Word8) -> g s -> m Word8 #-} + +-- | Uniformly generate Word32 in @[0, s]@. +unbiasedWordMult32 :: MonadRandom g s m => Word32 -> g s -> m Word32 +unbiasedWordMult32 s g + | s == maxBound = uniformWord32 g + | otherwise = unbiasedWordMult32Exclusive (s+1) g +{-# INLINE unbiasedWordMult32 #-} + +-- | See [Lemire's paper](https://arxiv.org/pdf/1805.10941.pdf), +-- [O\'Neill's +-- blogpost](https://www.pcg-random.org/posts/bounded-rands.html) and +-- more directly [O\'Neill's github +-- repo](https://github.com/imneme/bounded-rands/blob/3d71f53c975b1e5b29f2f3b05a74e26dab9c3d84/bounded32.cpp#L234). +-- N.B. The range is [0,t) **not** [0,t]. +unbiasedWordMult32Exclusive :: MonadRandom g s m => Word32 -> g s -> m Word32 +unbiasedWordMult32Exclusive r g = go + where + t :: Word32 + t = (-r) `mod` r -- Calculates 2^32 `mod` r!!! + go = do + x <- uniformWord32 g + let m :: Word64 + m = fromIntegral x * fromIntegral r + l :: Word32 + l = fromIntegral m + if l >= t then return (fromIntegral $ m `shiftR` 32) else go + +-- | This only works for unsigned integrals +unsignedBitmaskWithRejectionRM :: + (MonadRandom g s m, FiniteBits a, Num a, Ord a, Uniform a) + => (a, a) + -> g s + -> m a +unsignedBitmaskWithRejectionRM (bottom, top) gen + | bottom == top = pure top + | otherwise = (b +) <$> unsignedBitmaskWithRejectionM uniformM r gen + where + (b, r) = if bottom > top then (top, bottom - top) else (bottom, top - bottom) +{-# INLINE unsignedBitmaskWithRejectionRM #-} + +-- | This works for signed integrals by explicit conversion to unsigned and abusing overflow +signedBitmaskWithRejectionRM :: + (Num a, Num b, Ord b, Ord a, FiniteBits a, MonadRandom g s f, Uniform a) + => (b -> a) + -> (a -> b) + -> (b, b) + -> g s + -> f b +signedBitmaskWithRejectionRM toUnsigned fromUnsigned (bottom, top) gen + | bottom == top = pure top + | otherwise = + (b +) . fromUnsigned <$> unsignedBitmaskWithRejectionM uniformM r gen + -- This works in all cases, see Appendix 1 at the end of the file. + where + (b, r) = + if bottom > top + then (top, toUnsigned bottom - toUnsigned top) + else (bottom, toUnsigned top - toUnsigned bottom) +{-# INLINE signedBitmaskWithRejectionRM #-} + +unsignedBitmaskWithRejectionM :: (Ord a, FiniteBits a, Num a, MonadRandom g s m) => (g s -> m a) -> a -> g s -> m a +unsignedBitmaskWithRejectionM genUniformM range gen = go + where + mask = complement zeroBits `shiftR` countLeadingZeros (range .|. 1) + go = do + x <- genUniformM gen + let x' = x .&. mask + if x' > range + then go + else pure x' +{-# INLINE unsignedBitmaskWithRejectionM #-} + +------------------------------------------------------------------------------- +-- 'Uniform' instances for tuples +------------------------------------------------------------------------------- + +instance (Uniform a, Uniform b) => Uniform (a, b) where + uniformM g = (,) <$> uniformM g <*> uniformM g + +instance (Uniform a, Uniform b, Uniform c) => Uniform (a, b, c) where + uniformM g = (,,) <$> uniformM g <*> uniformM g <*> uniformM g + +instance (Uniform a, Uniform b, Uniform c, Uniform d) => Uniform (a, b, c, d) where + uniformM g = (,,,) <$> uniformM g <*> uniformM g <*> uniformM g <*> uniformM g + +instance (Uniform a, Uniform b, Uniform c, Uniform d, Uniform e) => Uniform (a, b, c, d, e) where + uniformM g = (,,,,) <$> uniformM g <*> uniformM g <*> uniformM g <*> uniformM g <*> uniformM g + +instance (Uniform a, Uniform b, Uniform c, Uniform d, Uniform e, Uniform f) => Uniform (a, b, c, d, e, f) where + uniformM g = (,,,,,) <$> uniformM g <*> uniformM g <*> uniformM g <*> uniformM g <*> uniformM g <*> uniformM g + +instance (Uniform a, Uniform b, Uniform c, Uniform d, Uniform e, Uniform f, Uniform g) => Uniform (a, b, c, d, e, f, g) where + uniformM g = (,,,,,,) <$> uniformM g <*> uniformM g <*> uniformM g <*> uniformM g <*> uniformM g <*> uniformM g <*> uniformM g + +-- Appendix 1. +-- +-- @top@ and @bottom@ are signed integers of bit width @n@. @toUnsigned@ +-- converts a signed integer to an unsigned number of the same bit width @n@. +-- +-- range = toUnsigned top - toUnsigned bottom +-- +-- This works out correctly thanks to modular arithmetic. Conceptually, +-- +-- toUnsigned x | x >= 0 = x +-- toUnsigned x | x < 0 = 2^n + x +-- +-- The following combinations are possible: +-- +-- 1. @bottom >= 0@ and @top >= 0@ +-- 2. @bottom < 0@ and @top >= 0@ +-- 3. @bottom < 0@ and @top < 0@ +-- +-- Note that @bottom >= 0@ and @top < 0@ is impossible because of the +-- invariant @bottom < top@. +-- +-- For any signed integer @i@ of width @n@, we have: +-- +-- -2^(n-1) <= i <= 2^(n-1) - 1 +-- +-- Considering each combination in turn, we have +-- +-- 1. @bottom >= 0@ and @top >= 0@ +-- +-- range = (toUnsigned top - toUnsigned bottom) `mod` 2^n +-- --^ top >= 0, so toUnsigned top == top +-- --^ bottom >= 0, so toUnsigned bottom == bottom +-- = (top - bottom) `mod` 2^n +-- --^ top <= 2^(n-1) - 1 and bottom >= 0 +-- --^ top - bottom <= 2^(n-1) - 1 +-- --^ 0 < top - bottom <= 2^(n-1) - 1 +-- = top - bottom +-- +-- 2. @bottom < 0@ and @top >= 0@ +-- +-- range = (toUnsigned top - toUnsigned bottom) `mod` 2^n +-- --^ top >= 0, so toUnsigned top == top +-- --^ bottom < 0, so toUnsigned bottom == 2^n + bottom +-- = (top - (2^n + bottom)) `mod` 2^n +-- --^ summand -2^n cancels out in calculation modulo 2^n +-- = (top - bottom) `mod` 2^n +-- --^ top <= 2^(n-1) - 1 and bottom >= -2^(n-1) +-- --^ top - bottom <= (2^(n-1) - 1) - (-2^(n-1)) = 2^n - 1 +-- --^ 0 < top - bottom <= 2^n - 1 +-- = top - bottom +-- +-- 3. @bottom < 0@ and @top < 0@ +-- +-- range = (toUnsigned top - toUnsigned bottom) `mod` 2^n +-- --^ top < 0, so toUnsigned top == 2^n + top +-- --^ bottom < 0, so toUnsigned bottom == 2^n + bottom +-- = ((2^n + top) - (2^n + bottom)) `mod` 2^n +-- --^ summand 2^n cancels out in calculation modulo 2^n +-- = (top - bottom) `mod` 2^n +-- --^ top <= -1 +-- --^ bottom >= -2^(n-1) +-- --^ top - bottom <= -1 - (-2^(n-1)) = 2^(n-1) - 1 +-- --^ 0 < top - bottom <= 2^(n-1) - 1 +-- = top - bottom diff --git a/src/System/Random/Monad.hs b/src/System/Random/Monad.hs new file mode 100644 index 000000000..a5f737905 --- /dev/null +++ b/src/System/Random/Monad.hs @@ -0,0 +1,555 @@ +{-# LANGUAGE BangPatterns #-} +{-# LANGUAGE DerivingStrategies #-} +{-# LANGUAGE FlexibleInstances #-} +{-# LANGUAGE MultiParamTypeClasses #-} +{-# LANGUAGE RankNTypes #-} +{-# LANGUAGE Safe #-} +{-# LANGUAGE TypeFamilies #-} +{-# LANGUAGE UndecidableInstances #-} + +-- | +-- Module : System.Random.Monad +-- Copyright : (c) The University of Glasgow 2001 +-- License : BSD-style (see the file LICENSE in the 'random' repository) +-- Maintainer : libraries@haskell.org +-- Stability : stable +-- +-- This library deals with the common task of pseudo-random number generation. +module System.Random.Monad + ( + -- * Pure Random Generator + module System.Random + -- * Monadic Random Generator + -- $introduction + + -- * Usage + -- $usagemonadic + + -- * Pure and monadic pseudo-random number generator interfaces + -- $interfaces + , MonadRandom(..) + , runGenM + , runGenM_ + , RandomGenM(..) + , randomM + , randomRM + , splitGenM + + -- * Monadic adapters for pure pseudo-random number generators + -- $monadicadapters + + -- ** Pure adapter + , StateGen(..) + , StateGenM(..) + , runStateGen + , runStateGen_ + , runStateGenT + , runStateGenT_ + , runStateGenST + -- ** Mutable adapter with atomic operations + , AtomicGen(..) + , AtomicGenM(..) + , applyAtomicGen + -- ** Mutable adapter in 'IO' + , IOGen(..) + , IOGenM(..) + , applyIOGen + -- ** Mutable adapter in 'ST' + , STGen(..) + , STGenM(..) + , applySTGen + , runSTGen + , runSTGen_ + + -- * Pseudo-random values of various types + -- $uniform + , Uniform(..) + , uniformListM + , UniformRange(..) + + -- * Generators for sequences of pseudo-random bytes + , genShortByteStringIO + , genShortByteStringST + , uniformByteString + + -- * Appendix + + -- ** How to implement 'MonadRandom' + -- $implementmonadrandom + + -- ** Floating point number caveats + -- $floating + + -- * References + -- $references + ) where + +import Control.Monad.IO.Class +import Control.Monad.ST +import Control.Monad.State.Strict +import Data.IORef +import Data.STRef +import System.Random +import System.Random.Internal + +-- $introduction +-- +-- This module provides type classes and instances for the following concepts: +-- +-- [Monadic pseudo-random number generators] 'MonadRandom' is an interface to +-- monadic pseudo-random number generators. +-- +-- [Monadic adapters] 'StateGenM', 'AtomicGenM', 'IOGenM' and 'STGenM' turn a +-- 'RandomGen' instance into a 'MonadRandom' instance. +-- +-- [Drawing from a range] 'UniformRange' is used to generate a value of a +-- type uniformly within a range. +-- +-- This library provides instances of 'UniformRange' for many common +-- numeric types. +-- +-- [Drawing from the entire domain of a type] 'Uniform' is used to generate a +-- value of a type uniformly over all possible values of that type. +-- +-- This library provides instances of 'Uniform' for many common bounded +-- numeric types. +-- +-- $usagemonadic +-- +-- In monadic code, use the relevant 'Uniform' and 'UniformRange' instances to +-- generate pseudo-random values via 'uniformM' and 'uniformRM', respectively. +-- +-- As an example, @rollsM@ generates @n@ pseudo-random values of @Word8@ in the +-- range @[1, 6]@ in a 'MonadRandom' context; given a /monadic/ pseudo-random +-- number generator, you can run this probabilistic computation as follows: +-- +-- >>> :{ +-- let rollsM :: MonadRandom g s m => Int -> g s -> m [Word8] +-- rollsM n = replicateM n . uniformRM (1, 6) +-- in do +-- monadicGen <- MWC.create +-- rollsM 10 monadicGen :: IO [Word8] +-- :} +-- [4,1,2,4,4,5,2,1,5,4] +-- +-- Given a /pure/ pseudo-random number generator, you can run the monadic +-- pseudo-random number computation @rollsM@ in an 'IO' or 'ST' context by +-- first applying a monadic adapter like 'AtomicGen', 'IOGen' or 'STGen' to the +-- pure pseudo-random number generator and then running it with 'runGenM'. +-- +-- >>> :{ +-- let rollsM :: MonadRandom g s m => Int -> g s -> m [Word8] +-- rollsM n = replicateM n . uniformRM (1, 6) +-- pureGen = mkStdGen 42 +-- in +-- runGenM_ (IOGen pureGen) (rollsM 10) :: IO [Word8] +-- :} +-- [5,1,4,3,3,2,5,2,2,4] + +------------------------------------------------------------------------------- +-- Pseudo-random number generator interfaces +------------------------------------------------------------------------------- + +-- $interfaces +-- +-- Pseudo-random number generators come in two flavours: /pure/ and /monadic/. +-- +-- ['System.Random.RandomGen': pure pseudo-random number generators] +-- See "System.Random" module. +-- +-- ['MonadRandom': monadic pseudo-random number generators] These generators +-- mutate their own state as they produce pseudo-random values. They +-- generally live in 'ST' or 'IO' or some transformer that implements +-- @PrimMonad@. +-- + +------------------------------------------------------------------------------- +-- Monadic adapters +------------------------------------------------------------------------------- + +-- $monadicadapters +-- +-- Pure pseudo-random number generators can be used in monadic code via the +-- adapters 'StateGenM', 'AtomicGenM', 'IOGenM' and 'STGenM'. +-- +-- * 'StateGenM' can be used in any state monad. With strict 'StateT' there is +-- no performance overhead compared to using the 'RandomGen' instance +-- directly. 'StateGenM' is /not/ safe to use in the presence of exceptions +-- and concurrency. +-- +-- * 'AtomicGenM' is safe in the presence of exceptions and concurrency since +-- it performs all actions atomically. +-- +-- * 'IOGenM' is a wrapper around an 'IORef' that holds a pure generator. +-- 'IOGenM' is safe in the presence of exceptions, but not concurrency. +-- +-- * 'STGenM' is a wrapper around an 'STRef' that holds a pure generator. +-- 'STGenM' is safe in the presence of exceptions, but not concurrency. + +-- | Interface to operations on 'RandomGen' wrappers like 'IOGenM' and 'StateGenM'. +-- +-- @since 1.2 +class (RandomGen r, MonadRandom (g r) s m) => RandomGenM g r s m where + applyRandomGenM :: (r -> (a, r)) -> g r s -> m a + +-- | Splits a pseudo-random number generator into two. Overwrites the mutable +-- wrapper with one of the resulting generators and returns the other. +-- +-- @since 1.2 +splitGenM :: RandomGenM g r s m => g r s -> m r +splitGenM = applyRandomGenM split + +instance (RandomGen r, MonadIO m) => RandomGenM IOGenM r RealWorld m where + applyRandomGenM = applyIOGen + +instance (RandomGen r, MonadIO m) => RandomGenM AtomicGenM r RealWorld m where + applyRandomGenM = applyAtomicGen + +instance (RandomGen r, MonadState r m) => RandomGenM StateGenM r r m where + applyRandomGenM f _ = state f + +instance RandomGen r => RandomGenM STGenM r s (ST s) where + applyRandomGenM = applySTGen + +-- | Runs a mutable pseudo-random number generator from its 'Frozen' state. +-- +-- >>> import Data.Int (Int8) +-- >>> runGenM (IOGen (mkStdGen 217)) (`uniformListM` 5) :: IO ([Int8], IOGen StdGen) +-- ([-74,37,-50,-2,3],IOGen {unIOGen = StdGen {unStdGen = SMGen 4273268533320920145 15251669095119325999}}) +-- +-- @since 1.2 +runGenM :: MonadRandom g s m => Frozen g -> (g s -> m a) -> m (a, Frozen g) +runGenM fg action = do + g <- thawGen fg + res <- action g + fg' <- freezeGen g + pure (res, fg') + +-- | Same as 'runGenM', but only returns the generated value. +-- +-- @since 1.2 +runGenM_ :: MonadRandom g s m => Frozen g -> (g s -> m a) -> m a +runGenM_ fg action = fst <$> runGenM fg action + +-- | Generates a list of pseudo-random values. +-- +-- @since 1.2 +uniformListM :: (MonadRandom g s m, Uniform a) => g s -> Int -> m [a] +uniformListM gen n = replicateM n (uniformM gen) + +-- | Generates a pseudo-random value using monadic interface and `Random` instance. +-- +-- @since 1.2 +randomM :: (RandomGenM g r s m, Random a) => g r s -> m a +randomM = applyRandomGenM random + +-- | Generates a pseudo-random value using monadic interface and `Random` instance. +-- +-- @since 1.2 +randomRM :: (RandomGenM g r s m, Random a) => (a, a) -> g r s -> m a +randomRM r = applyRandomGenM (randomR r) + +-- | Wraps an 'IORef' that holds a pure pseudo-random number generator. All +-- operations are performed atomically. +-- +-- * 'AtomicGenM' is safe in the presence of exceptions and concurrency. +-- * 'AtomicGenM' is the slowest of the monadic adapters due to the overhead +-- of its atomic operations. +-- +-- @since 1.2 +newtype AtomicGenM g s = AtomicGenM { unAtomicGenM :: IORef g} + +newtype AtomicGen g = AtomicGen { unAtomicGen :: g } + deriving stock (Eq, Show, Read) + +instance (RandomGen g, MonadIO m) => MonadRandom (AtomicGenM g) RealWorld m where + type Frozen (AtomicGenM g) = AtomicGen g + thawGen (AtomicGen g) = fmap AtomicGenM (liftIO $ newIORef g) + freezeGen (AtomicGenM gVar) = fmap AtomicGen (liftIO $ readIORef gVar) + uniformWord32R r = applyAtomicGen (genWord32R r) + {-# INLINE uniformWord32R #-} + uniformWord64R r = applyAtomicGen (genWord64R r) + {-# INLINE uniformWord64R #-} + uniformWord8 = applyAtomicGen genWord8 + {-# INLINE uniformWord8 #-} + uniformWord16 = applyAtomicGen genWord16 + {-# INLINE uniformWord16 #-} + uniformWord32 = applyAtomicGen genWord32 + {-# INLINE uniformWord32 #-} + uniformWord64 = applyAtomicGen genWord64 + {-# INLINE uniformWord64 #-} + uniformShortByteString n = applyAtomicGen (genShortByteString n) + +-- | Atomically applies a pure operation to the wrapped pseudo-random number +-- generator. +-- +-- @since 1.2 +applyAtomicGen :: MonadIO m => (g -> (a, g)) -> AtomicGenM g RealWorld -> m a +applyAtomicGen op (AtomicGenM gVar) = + liftIO $ atomicModifyIORef' gVar $ \g -> + case op g of + (a, g') -> (g', a) +{-# INLINE applyAtomicGen #-} + +-- | Wraps an 'IORef' that holds a pure pseudo-random number generator. +-- +-- * 'IOGenM' is safe in the presence of exceptions, but not concurrency. +-- * 'IOGenM' is slower than 'StateGenM' due to the extra pointer indirection. +-- * 'IOGenM' is faster than 'AtomicGenM' since the 'IORef' operations used by +-- 'IOGenM' are not atomic. +-- +-- An example use case is writing pseudo-random bytes into a file: +-- +-- >>> import UnliftIO.Temporary (withSystemTempFile) +-- >>> import Data.ByteString (hPutStr) +-- >>> let ioGen g = withSystemTempFile "foo.bin" $ \_ h -> uniformRM (0, 100) g >>= flip uniformByteString g >>= hPutStr h +-- +-- and then run it: +-- +-- >>> runGenM_ (IOGen (mkStdGen 1729)) ioGen +-- +-- @since 1.2 +newtype IOGenM g s = IOGenM { unIOGenM :: IORef g } + +newtype IOGen g = IOGen { unIOGen :: g } + deriving stock (Eq, Show, Read) + +instance (RandomGen g, MonadIO m) => MonadRandom (IOGenM g) RealWorld m where + type Frozen (IOGenM g) = IOGen g + thawGen (IOGen g) = fmap IOGenM (liftIO $ newIORef g) + freezeGen (IOGenM gVar) = fmap IOGen (liftIO $ readIORef gVar) + uniformWord32R r = applyIOGen (genWord32R r) + {-# INLINE uniformWord32R #-} + uniformWord64R r = applyIOGen (genWord64R r) + {-# INLINE uniformWord64R #-} + uniformWord8 = applyIOGen genWord8 + {-# INLINE uniformWord8 #-} + uniformWord16 = applyIOGen genWord16 + {-# INLINE uniformWord16 #-} + uniformWord32 = applyIOGen genWord32 + {-# INLINE uniformWord32 #-} + uniformWord64 = applyIOGen genWord64 + {-# INLINE uniformWord64 #-} + uniformShortByteString n = applyIOGen (genShortByteString n) + +-- | Applies a pure operation to the wrapped pseudo-random number generator. +-- +-- @since 1.2 +applyIOGen :: MonadIO m => (g -> (a, g)) -> IOGenM g RealWorld -> m a +applyIOGen f (IOGenM ref) = liftIO $ do + g <- readIORef ref + case f g of + (!a, !g') -> a <$ writeIORef ref g' +{-# INLINE applyIOGen #-} + +-- | Wraps an 'STRef' that holds a pure pseudo-random number generator. +-- +-- * 'STGenM' is safe in the presence of exceptions, but not concurrency. +-- * 'STGenM' is slower than 'StateGenM' due to the extra pointer indirection. +-- +-- @since 1.2 +newtype STGenM g s = STGenM { unSTGenM :: STRef s g } + +newtype STGen g = STGen { unSTGen :: g } + deriving stock (Eq, Show, Read) + +instance RandomGen g => MonadRandom (STGenM g) s (ST s) where + type Frozen (STGenM g) = STGen g + thawGen (STGen g) = fmap STGenM (newSTRef g) + freezeGen (STGenM gVar) = fmap STGen (readSTRef gVar) + uniformWord32R r = applySTGen (genWord32R r) + {-# INLINE uniformWord32R #-} + uniformWord64R r = applySTGen (genWord64R r) + {-# INLINE uniformWord64R #-} + uniformWord8 = applySTGen genWord8 + {-# INLINE uniformWord8 #-} + uniformWord16 = applySTGen genWord16 + {-# INLINE uniformWord16 #-} + uniformWord32 = applySTGen genWord32 + {-# INLINE uniformWord32 #-} + uniformWord64 = applySTGen genWord64 + {-# INLINE uniformWord64 #-} + uniformShortByteString n = applySTGen (genShortByteString n) + +-- | Applies a pure operation to the wrapped pseudo-random number generator. +-- +-- @since 1.2 +applySTGen :: (g -> (a, g)) -> STGenM g s -> ST s a +applySTGen f (STGenM ref) = do + g <- readSTRef ref + case f g of + (!a, !g') -> a <$ writeSTRef ref g' +{-# INLINE applySTGen #-} + +-- | Runs a monadic generating action in the `ST` monad using a pure +-- pseudo-random number generator. +-- +-- @since 1.2 +runSTGen :: RandomGen g => g -> (forall s . STGenM g s -> ST s a) -> (a, g) +runSTGen g action = unSTGen <$> runST (runGenM (STGen g) action) + +-- | Runs a monadic generating action in the `ST` monad using a pure +-- pseudo-random number generator. Returns only the resulting pseudo-random +-- value. +-- +-- @since 1.2 +runSTGen_ :: RandomGen g => g -> (forall s . STGenM g s -> ST s a) -> a +runSTGen_ g action = fst $ runSTGen g action + + +-- $uniform +-- +-- This library provides two type classes to generate pseudo-random values: +-- +-- * 'UniformRange' is used to generate a value of a type uniformly within a +-- range. +-- * 'Uniform' is used to generate a value of a type uniformly over all +-- possible values of that type. +-- +-- Types may have instances for both or just one of 'UniformRange' and +-- 'Uniform'. A few examples illustrate this: +-- +-- * 'Int', 'Word16' and 'Bool' are instances of both 'UniformRange' and +-- 'Uniform'. +-- * 'Integer', 'Float' and 'Double' each have an instance for 'UniformRange' +-- but no 'Uniform' instance. +-- * A hypothetical type @Radian@ representing angles by taking values in the +-- range @[0, 2π)@ has a trivial 'Uniform' instance, but no 'UniformRange' +-- instance: the problem is that two given @Radian@ values always span /two/ +-- ranges, one clockwise and one anti-clockwise. +-- * It is trivial to construct a @Uniform (a, b)@ instance given +-- @Uniform a@ and @Uniform b@ (and this library provides this tuple +-- instance). +-- * On the other hand, there is no correct way to construct a +-- @UniformRange (a, b)@ instance based on just @UniformRange a@ and +-- @UniformRange b@. + +------------------------------------------------------------------------------- +-- Notes +------------------------------------------------------------------------------- + +-- $floating +-- +-- The 'UniformRange' instances for 'Float' and 'Double' use the following +-- procedure to generate a random value in a range for @uniformRM (a, b) g@: +-- +-- 1. Generate \(x\) uniformly such that \(0 \leq x \leq 1\). +-- +-- The method by which \(x\) is sampled does not cover all representable +-- floating point numbers in the unit interval. The method never generates +-- denormal floating point numbers, for example. +-- +-- 2. Return \((b - a) * x + a\). +-- +-- Due to rounding errors, floating point operations are neither +-- commutative, nor associative, nor distributive the way the corresponding +-- operations on real numbers are. Additionally, floating point numbers +-- admit special values @NaN@ as well as negative and positive infinity. +-- +-- For pathological values, step 2 can yield surprising results. +-- +-- * The result may be greater than @max a b@. +-- +-- >>> :{ +-- let (a, b, x) = (-4.7021254e-38, -1.481e-42, 1.0) +-- result = (b - a) * x + a :: Float +-- in (result, result > max a b) +-- :} +-- (-1.48e-42,True) +-- +-- * The result may be @NaN@ even if \(a\) and \(b\) are not. +-- +-- >>> :{ +-- let (a, b, x) = (-1.7159568e38, 1.7159568e38, 0.0) +-- in (b - a) * x + a :: Float +-- :} +-- NaN +-- +-- What happens when @NaN@ or @Infinity@ are given to 'uniformRM'? We first +-- define them as constants: +-- +-- >>> nan = read "NaN" :: Float +-- >>> inf = read "Infinity" :: Float +-- +-- * If at least one of \(a\) or \(b\) is @NaN@, the result is @NaN@. +-- +-- >>> let (a, b, x) = (nan, 1, 0.5) in (b - a) * x + a +-- NaN +-- >>> let (a, b, x) = (-1, nan, 0.5) in (b - a) * x + a +-- NaN +-- +-- * Otherwise, if \(a\) is @Infinity@ or @-Infinity@, the result is @NaN@. +-- +-- >>> let (a, b, x) = (inf, 1, 0.5) in (b - a) * x + a +-- NaN +-- >>> let (a, b, x) = (-inf, 1, 0.5) in (b - a) * x + a +-- NaN +-- +-- * Otherwise, if \(b\) is @Infinity@ or @-Infinity@, the result is \(b\). +-- +-- >>> let (a, b, x) = (1, inf, 0.5) in (b - a) * x + a +-- Infinity +-- >>> let (a, b, x) = (1, -inf, 0.5) in (b - a) * x + a +-- -Infinity +-- +-- Note that the [GCC 10.1.0 C++ standard library](https://gcc.gnu.org/git/?p=gcc.git;a=blob;f=libstdc%2B%2B-v3/include/bits/random.h;h=19307fbc3ca401976ef6823e8fda893e4a263751;hb=63fa67847628e5f358e7e2e7edb8314f0ee31f30#l1859), +-- the [Java 10 standard library](https://docs.oracle.com/javase/10/docs/api/java/util/Random.html#doubles%28double,double%29) +-- and [CPython 3.8](https://github.com/python/cpython/blob/3.8/Lib/random.py#L417) +-- use the same procedure to generate floating point values in a range. +-- +-- $implementmonadrandom +-- +-- Typically, a monadic pseudo-random number generator has facilities to save +-- and restore its internal state in addition to generating pseudo-random +-- pseudo-random numbers. +-- +-- Here is an example instance for the monadic pseudo-random number generator +-- from the @mwc-random@ package: +-- +-- > instance (s ~ PrimState m, PrimMonad m) => MonadRandom MWC.Gen s m where +-- > type Frozen MWC.Gen = MWC.Seed +-- > thawGen = MWC.restore +-- > freezeGen = MWC.save +-- > uniformWord8 = MWC.uniform +-- > uniformWord16 = MWC.uniform +-- > uniformWord32 = MWC.uniform +-- > uniformWord64 = MWC.uniform +-- > uniformShortByteString n g = unsafeSTToPrim (genShortByteStringST n (MWC.uniform g)) +-- +-- $references +-- +-- 1. Guy L. Steele, Jr., Doug Lea, and Christine H. Flood. 2014. Fast +-- splittable pseudorandom number generators. In Proceedings of the 2014 ACM +-- International Conference on Object Oriented Programming Systems Languages & +-- Applications (OOPSLA '14). ACM, New York, NY, USA, 453-472. DOI: +-- + +-- $setup +-- >>> import Control.Arrow (first, second) +-- >>> import Control.Monad (replicateM) +-- >>> import Control.Monad.Primitive +-- >>> import Data.Bits +-- >>> import Data.Int (Int32) +-- >>> import Data.Word (Word8, Word16, Word32, Word64) +-- >>> import System.IO (IOMode(WriteMode), withBinaryFile) +-- >>> import System.Random.Monad +-- >>> import qualified System.Random.MWC as MWC +-- +-- >>> :set -XFlexibleContexts +-- >>> :set -XFlexibleInstances +-- >>> :set -XMultiParamTypeClasses +-- >>> :set -XTypeFamilies +-- >>> :set -XUndecidableInstances +-- +-- >>> :{ +-- instance (s ~ PrimState m, PrimMonad m) => MonadRandom MWC.Gen s m where +-- type Frozen MWC.Gen = MWC.Seed +-- thawGen = MWC.restore +-- freezeGen = MWC.save +-- uniformWord8 = MWC.uniform +-- uniformWord16 = MWC.uniform +-- uniformWord32 = MWC.uniform +-- uniformWord64 = MWC.uniform +-- uniformShortByteString n g = unsafeSTToPrim (genShortByteStringST n (MWC.uniform g)) +-- :} +-- diff --git a/stack-coveralls.yaml b/stack-coveralls.yaml new file mode 100644 index 000000000..5575c2ffb --- /dev/null +++ b/stack-coveralls.yaml @@ -0,0 +1,8 @@ +resolver: lts-14.27 +packages: +- . +extra-deps: +- rdtsc-1.3.0.1@sha256:0a6e8dc715ba82ad72c7e2b1c2f468999559bec059d50540719a80b00dcc4e66,1557 +flags: + splitmix: + random: false diff --git a/stack-old.yaml b/stack-old.yaml new file mode 100644 index 000000000..6f02f49c8 --- /dev/null +++ b/stack-old.yaml @@ -0,0 +1,25 @@ +resolver: lts-11.22 +packages: +- . +extra-deps: + # We require a splitmix version that supports the "random" flag. When this + # flag is set to false, splitmix does not depend on the random package. We + # require the flag to be set to false to avoid a circular dependency: + # random -> splitmix -> random. + # The oldest splitmix version to support the "random" flag is 0.0.3. LTS + # snapshots prior to lts-14 contain older versions of splitmix. To make sure + # the build goes through, we pin the splitmix version here. + # + # Note that although the resolver defined in this file comes with a splitmix + # >= 0.0.3 by default, this stack.yaml is also used in CI, where the resolver + # is overridden, making this explicit pinning of splitmix necessary. + # + # In addition, the doctests require 'SMGen' from the splitmix package to + # implement 'Prim'. So far, this is only the case on lehin's fork, so we use + # that. +- splitmix-0.0.4@sha256:fb9bb8b54a2e76c8a021fe5c4c3798047e1f60e168379a1f80693047fe00ad0e,4813 +- doctest-0.16.2@sha256:2f96e9bbe9aee11b47453c82c24b3dc76cdbb8a2a7c984dfd60b4906d08adf68,6942 +- cabal-doctest-1.0.8@sha256:34dff6369d417df2699af4e15f06bc181d495eca9c51efde173deae2053c197c,1491 +flags: + splitmix: + random: false diff --git a/stack.yaml b/stack.yaml new file mode 100644 index 000000000..396552681 --- /dev/null +++ b/stack.yaml @@ -0,0 +1,7 @@ +resolver: lts-15.1 +packages: +- . +extra-deps: [] +flags: + splitmix: + random: false diff --git a/test-legacy/Legacy.hs b/test-legacy/Legacy.hs new file mode 100644 index 000000000..f4660fdcd --- /dev/null +++ b/test-legacy/Legacy.hs @@ -0,0 +1,15 @@ +module Main (main) where + +import qualified Random1283 as Random1283 +import qualified RangeTest as RangeTest +import qualified T7936 as T7936 +import qualified TestRandomIOs as TestRandomIOs +import qualified TestRandomRs as TestRandomRs + +main :: IO () +main = do + Random1283.main + RangeTest.main + T7936.main + TestRandomIOs.main + TestRandomRs.main diff --git a/tests/random1283.hs b/test-legacy/Random1283.hs similarity index 59% rename from tests/random1283.hs rename to test-legacy/Random1283.hs index 33fc48862..239e29d11 100644 --- a/tests/random1283.hs +++ b/test-legacy/Random1283.hs @@ -1,19 +1,25 @@ +module Random1283 (main) where + import Control.Concurrent -import Control.Monad hiding (empty) -import Data.Sequence (ViewL(..), empty, fromList, viewl, (<|), (|>), (><)) +import Control.Monad +import Data.Sequence (Seq, ViewL(..), empty, fromList, viewl, (<|), (|>), (><)) import System.Random -- This test +threads, samples :: Int threads = 4 samples = 5000 +main :: IO () main = loopTest threads samples +loopTest :: Int -> Int -> IO () loopTest t s = do isClean <- testRace t s - when (not isClean) $ putStrLn "race condition!" + unless isClean $ putStrLn "race condition!" +testRace :: Int -> Int -> IO Bool testRace t s = do ref <- liftM (take (t*s) . randoms) getStdGen iss <- threadRandoms t s @@ -23,11 +29,12 @@ threadRandoms :: Random a => Int -> Int -> IO [[a]] threadRandoms t s = do vs <- sequence $ replicate t $ do v <- newEmptyMVar - forkIO (sequence (replicate s randomIO) >>= putMVar v) + _ <- forkIO (sequence (replicate s randomIO) >>= putMVar v) return v mapM takeMVar vs -isInterleavingOf xs yss = iio xs (viewl $ fromList yss) EmptyL where +isInterleavingOf :: Eq a => [a] -> [[a]] -> Bool +isInterleavingOf xs' yss' = iio xs' (viewl $ fromList yss') EmptyL where iio (x:xs) ((y:ys) :< yss) zss | x /= y = iio (x:xs) (viewl yss) (viewl (fromViewL zss |> (y:ys))) | x == y = iio xs (viewl ((ys <| yss) >< fromViewL zss)) EmptyL @@ -35,6 +42,7 @@ isInterleavingOf xs yss = iio xs (viewl $ fromList yss) EmptyL where iio [] EmptyL EmptyL = True iio _ _ _ = False -fromViewL (EmptyL) = empty +fromViewL :: ViewL a -> Seq a +fromViewL EmptyL = empty fromViewL (x :< xs) = x <| xs diff --git a/tests/rangeTest.hs b/test-legacy/RangeTest.hs similarity index 82% rename from tests/rangeTest.hs rename to test-legacy/RangeTest.hs index ac62c71dd..78cc17547 100644 --- a/tests/rangeTest.hs +++ b/test-legacy/RangeTest.hs @@ -1,25 +1,26 @@ -{-# LANGUAGE CPP #-} +module RangeTest (main) where + import Control.Monad import System.Random import Data.Int import Data.Word -import Data.Bits import Foreign.C.Types -- Take many measurements and record the max/min/average random values. -approxBounds :: (RandomGen g, Random a, Ord a, Num a) => - (g -> (a,g)) -> Int -> a -> (a,a) -> g -> ((a,a,a),g) +approxBounds :: + (RandomGen g, Random a, Ord a, Num a) => + (g -> (a,g)) -> Int -> a -> (a,a) -> g -> ((a,a,a),g) -- Here we do a little hack to essentiall pass in the type in the last argument: -approxBounds nxt iters unused (explo,exphi) initrng = - if False +approxBounds nxt iters unused (explo,exphi) initrng = + if False then ((unused,unused,unused),undefined) -- else loop initrng iters 100 (-100) 0 -- Oops, can't use minBound/maxBound here. - else loop initrng iters exphi explo 0 -- Oops, can't use minBound/maxBound here. - where - loop rng 0 mn mx sum = ((mn,mx,sum),rng) - loop rng n mn mx sum = - case nxt rng of - (x, rng') -> loop rng' (n-1) (min x mn) (max x mx) (x+sum) + else loop initrng iters exphi explo 0 + where + loop rng 0 mn mx sum' = ((mn,mx,sum'),rng) + loop rng n mn mx sum' = + case nxt rng of + (x, rng') -> loop rng' (n-1) (min x mn) (max x mx) (x+sum') -- We check that: @@ -27,34 +28,34 @@ approxBounds nxt iters unused (explo,exphi) initrng = -- (2) we get "close" to the bounds -- The with (2) is that we do enough trials to ensure that we can at -- least hit the 90% mark. -checkBounds:: (Real a, Show a, Ord a) => - String -> (Bool, a, a) -> ((a,a) -> StdGen -> ((a, a, t), StdGen)) -> IO () -checkBounds msg (exclusive,lo,hi) fun = - -- (lo,hi) is [inclusive,exclusive) - do putStr$ msg --- ++ ", expected range " ++ show (lo,hi) - ++ ": " - (mn,mx,sum) <- getStdRandom (fun (lo,hi)) - when (mn < lo)$ error$ "broke lower bound: " ++ show mn - when (mx > hi) $ error$ "broke upper bound: " ++ show mx - when (exclusive && mx >= hi)$ error$ "hit upper bound: " ++ show mx +checkBounds :: + (Real a, Show a, Ord a) => + String -> (Bool, a, a) -> ((a,a) -> StdGen -> ((a, a, t), StdGen)) -> IO () +checkBounds msg (exclusive,lo,hi) fun = do + -- (lo,hi) is [inclusive,exclusive) + putStr $ msg ++ ": " + (mn,mx,_) <- getStdRandom (fun (lo,hi)) + when (mn < lo) $ error $ "broke lower bound: " ++ show mn + when (mx > hi) $ error $ "broke upper bound: " ++ show mx + when (exclusive && mx >= hi)$ error$ "hit upper bound: " ++ show mx - let epsilon = 0.1 * (toRational hi - toRational lo) + let epsilon = 0.1 * (toRational hi - toRational lo) - when (toRational (hi - mx) > epsilon)$ error$ "didn't get close enough to upper bound: "++ show mx - when (toRational (mn - lo) > epsilon)$ error$ "didn't get close enough to lower bound: "++ show mn - putStrLn "Passed" + when (toRational (hi - mx) > epsilon) $ error $ "didn't get close enough to upper bound: "++ show mx + when (toRational (mn - lo) > epsilon) $ error $ "didn't get close enough to lower bound: "++ show mn + putStrLn "Passed" boundedRange :: (Num a, Bounded a) => (Bool, a, a) boundedRange = ( False, minBound, maxBound ) +trials :: Int trials = 5000 -- Keep in mind here that on some architectures (e.g. ARM) CChar, CWchar, and CSigAtomic -- are unsigned - -main = - do +main :: IO () +main = + do checkBounds "Int" boundedRange (approxBounds random trials (undefined::Int)) checkBounds "Integer" (False, fromIntegral (minBound::Int), fromIntegral (maxBound::Int)) (approxBounds random trials (undefined::Integer)) @@ -67,8 +68,8 @@ main = checkBounds "Word16" boundedRange (approxBounds random trials (undefined::Word16)) checkBounds "Word32" boundedRange (approxBounds random trials (undefined::Word32)) checkBounds "Word64" boundedRange (approxBounds random trials (undefined::Word64)) - checkBounds "Double" (True,0.0,1.0) (approxBounds random trials (undefined::Double)) - checkBounds "Float" (True,0.0,1.0) (approxBounds random trials (undefined::Float)) + checkBounds "Double" (False,0.0,1.0) (approxBounds random trials (undefined::Double)) + checkBounds "Float" (False,0.0,1.0) (approxBounds random trials (undefined::Float)) checkBounds "CChar" boundedRange (approxBounds random trials (undefined:: CChar)) checkBounds "CSChar" boundedRange (approxBounds random trials (undefined:: CSChar)) @@ -92,7 +93,9 @@ main = -- Then check all the range-restricted versions: checkBounds "Int R" (False,-100,100) (approxBounds (randomR (-100,100)) trials (undefined::Int)) - checkBounds "Integer R" (False,-100,100) (approxBounds (randomR (-100,100)) trials (undefined::Integer)) + checkBounds "Integer R" + (False,-100000000000000000000,100000000000000000000) + (approxBounds (randomR (-100000000000000000000,100000000000000000000)) trials (undefined::Integer)) checkBounds "Int8 R" (False,-100,100) (approxBounds (randomR (-100,100)) trials (undefined::Int8)) checkBounds "Int8 Rsmall" (False,-50,50) (approxBounds (randomR (-50,50)) trials (undefined::Int8)) checkBounds "Int8 Rmini" (False,3,4) (approxBounds (randomR (3,4)) trials (undefined::Int8)) @@ -106,8 +109,8 @@ main = checkBounds "Word16 R" (False,0,200) (approxBounds (randomR (0,200)) trials (undefined::Word16)) checkBounds "Word32 R" (False,0,200) (approxBounds (randomR (0,200)) trials (undefined::Word32)) checkBounds "Word64 R" (False,0,200) (approxBounds (randomR (0,200)) trials (undefined::Word64)) - checkBounds "Double R" (True,10.0,77.0) (approxBounds (randomR (10,77)) trials (undefined::Double)) - checkBounds "Float R" (True,10.0,77.0) (approxBounds (randomR (10,77)) trials (undefined::Float)) + checkBounds "Double R" (False,10.0,77.0) (approxBounds (randomR (10,77)) trials (undefined::Double)) + checkBounds "Float R" (False,10.0,77.0) (approxBounds (randomR (10,77)) trials (undefined::Float)) checkBounds "CChar R" (False,0,100) (approxBounds (randomR (0,100)) trials (undefined:: CChar)) checkBounds "CSChar R" (False,-100,100) (approxBounds (randomR (-100,100)) trials (undefined:: CSChar)) diff --git a/tests/T7936.hs b/test-legacy/T7936.hs similarity index 89% rename from tests/T7936.hs rename to test-legacy/T7936.hs index cfea911bc..47b30d1c5 100644 --- a/tests/T7936.hs +++ b/test-legacy/T7936.hs @@ -6,9 +6,10 @@ -- $ cabal test T7936 --test-options="+RTS -M1M -RTS" -- T7936: Heap exhausted; -module Main where +module T7936 where import System.Random (newStdGen) import Control.Monad (replicateM_) +main :: IO () main = replicateM_ 100000 newStdGen diff --git a/tests/TestRandomIOs.hs b/test-legacy/TestRandomIOs.hs similarity index 93% rename from tests/TestRandomIOs.hs rename to test-legacy/TestRandomIOs.hs index d8a00cc7c..4af2ddcd5 100644 --- a/tests/TestRandomIOs.hs +++ b/test-legacy/TestRandomIOs.hs @@ -6,7 +6,7 @@ -- $ cabal test TestRandomIOs --test-options="+RTS -M1M -RTS" -- TestRandomIOs: Heap exhausted; -module Main where +module TestRandomIOs where import Control.Monad (replicateM) import System.Random (randomIO) @@ -15,6 +15,7 @@ import System.Random (randomIO) -- the last one. -- Should use less than 1Mb of heap space, or we are generating a list of -- unevaluated thunks. +main :: IO () main = do rs <- replicateM 5000 randomIO :: IO [Int] print $ last rs diff --git a/tests/TestRandomRs.hs b/test-legacy/TestRandomRs.hs similarity index 89% rename from tests/TestRandomRs.hs rename to test-legacy/TestRandomRs.hs index cdae106a6..90145f219 100644 --- a/tests/TestRandomRs.hs +++ b/test-legacy/TestRandomRs.hs @@ -10,13 +10,14 @@ -- $ cabal test TestRandomRs --test-options="+RTS -M1M -RTS" -- TestRandomRs: Heap exhausted; -module Main where +module TestRandomRs where -import Control.Monad (liftM, replicateM) +import Control.Monad (liftM) import System.Random (randomRs, getStdGen) -- Return the five-thousandth random number: -- Should run in constant space (< 1Mb heap). +main :: IO () main = do n <- (last . take 5000 . randomRs (0, 1000000)) `liftM` getStdGen print (n::Integer) diff --git a/test/Spec.hs b/test/Spec.hs new file mode 100644 index 000000000..8fada9d6c --- /dev/null +++ b/test/Spec.hs @@ -0,0 +1,195 @@ +{-# LANGUAGE AllowAmbiguousTypes #-} +{-# LANGUAGE CPP #-} +{-# LANGUAGE FlexibleContexts #-} +{-# LANGUAGE FlexibleInstances #-} +{-# LANGUAGE MultiParamTypeClasses #-} +{-# LANGUAGE ScopedTypeVariables #-} +{-# LANGUAGE TypeApplications #-} +{-# OPTIONS_GHC -fno-warn-orphans #-} +module Main (main) where + +import Data.ByteString.Short as SBS +import Data.Coerce +import Data.Int +import Data.Typeable +import Data.Word +import Foreign.C.Types +import Numeric.Natural (Natural) +import System.Random +import Test.SmallCheck.Series as SC +import Test.Tasty +import Test.Tasty.HUnit +import Test.Tasty.SmallCheck as SC + +#include "HsBaseConfig.h" + +import qualified Spec.Range as Range +import qualified Spec.Run as Run + +main :: IO () +main = + defaultMain $ + testGroup + "Spec" + [ floatingSpec @Double + , floatingSpec @Float + , floatingSpec @CDouble + , floatingSpec @CFloat + , integralSpec @Word8 + , integralSpec @Word16 + , integralSpec @Word32 + , integralSpec @Word64 + , integralSpec @Word + , integralSpec @Int8 + , integralSpec @Int16 + , integralSpec @Int32 + , integralSpec @Int64 + , integralSpec @Int + , integralSpec @Char + , integralSpec @Bool + , integralSpec @CBool + , integralSpec @CChar + , integralSpec @CSChar + , integralSpec @CUChar + , integralSpec @CShort + , integralSpec @CUShort + , integralSpec @CInt + , integralSpec @CUInt + , integralSpec @CLong + , integralSpec @CULong + , integralSpec @CPtrdiff + , integralSpec @CSize + , integralSpec @CWchar + , integralSpec @CSigAtomic + , integralSpec @CLLong + , integralSpec @CULLong + , integralSpec @CIntPtr + , integralSpec @CUIntPtr + , integralSpec @CIntMax + , integralSpec @CUIntMax + , integralSpec @Integer + , integralSpec @Natural + -- , bitmaskSpec @Word8 + -- , bitmaskSpec @Word16 + -- , bitmaskSpec @Word32 + -- , bitmaskSpec @Word64 + -- , bitmaskSpec @Word + , runSpec + , floatTests + , byteStringSpec + ] + +floatTests :: TestTree +floatTests = testGroup "(Float)" + [ -- Check that https://github.com/haskell/random/issues/53 does not regress + + testCase "Subnormal generation not above upper bound" $ + [] @?= filter (>4.0e-45) (take 100000 $ randomRs (0, 4.0e-45::Float) $ mkStdGen 0) + + , testCase "Subnormal generation includes upper bound" $ + 1.0e-45 `elem` take 100 (randomRs (0, 1.0e-45::Float) $ mkStdGen 0) @? + "Does not contain 1.0e-45" + ] + +showsType :: forall t . Typeable t => ShowS +showsType = showsTypeRep (typeRep (Proxy :: Proxy t)) + +byteStringSpec :: TestTree +byteStringSpec = + testGroup + "ByteString" + [ SC.testProperty "genShortByteString" $ \(seed, n8) -> + let n = fromIntegral (n8 :: Word8) -- no need to generate huge collection of bytes + in SBS.length (fst (seeded (genShortByteString n) seed)) == n + , SC.testProperty "genByteString" $ \(seed, n8) -> + let n = fromIntegral (n8 :: Word8) + in SBS.toShort (fst (seeded (genByteString n) seed)) == + fst (seeded (genShortByteString n) seed) + ] + + +rangeSpec :: + forall a. + (SC.Serial IO a, Typeable a, Ord a, UniformRange a, Show a) + => TestTree +rangeSpec = + testGroup ("Range (" ++ showsType @a ")") + [ SC.testProperty "uniformR" $ seeded $ Range.uniformRangeWithin @_ @a + ] + +integralSpec :: + forall a. + (SC.Serial IO a, Typeable a, Ord a, UniformRange a, Show a) + => TestTree +integralSpec = + testGroup ("(" ++ showsType @a ")") + [ SC.testProperty "symmetric" $ seeded $ Range.symmetric @_ @a + , SC.testProperty "bounded" $ seeded $ Range.bounded @_ @a + , SC.testProperty "singleton" $ seeded $ Range.singleton @_ @a + , rangeSpec @a + -- TODO: Add more tests + ] + +floatingSpec :: + forall a. + (SC.Serial IO a, Typeable a, Num a, Ord a, Random a, UniformRange a, Show a) + => TestTree +floatingSpec = + testGroup ("(" ++ showsType @a ")") + [ SC.testProperty "uniformR" $ seeded $ Range.uniformRangeWithin @_ @a + -- TODO: Add more tests + ] + +runSpec :: TestTree +runSpec = testGroup "runGenState_ and runPrimGenIO_" + [ SC.testProperty "equal outputs" $ seeded $ \g -> monadic $ Run.runsEqual g ] + +-- | Create a StdGen instance from an Int and pass it to the given function. +seeded :: (StdGen -> a) -> Int -> a +seeded f = f . mkStdGen + + +instance Monad m => Serial m CFloat where + series = coerce <$> (series :: Series m Float) +instance Monad m => Serial m CDouble where + series = coerce <$> (series :: Series m Double) +instance Monad m => Serial m CBool where + series = coerce <$> (series :: Series m HTYPE_BOOL) +instance Monad m => Serial m CChar where + series = coerce <$> (series :: Series m HTYPE_CHAR) +instance Monad m => Serial m CSChar where + series = coerce <$> (series :: Series m HTYPE_SIGNED_CHAR) +instance Monad m => Serial m CUChar where + series = coerce <$> (series :: Series m HTYPE_UNSIGNED_CHAR) +instance Monad m => Serial m CShort where + series = coerce <$> (series :: Series m HTYPE_SHORT) +instance Monad m => Serial m CUShort where + series = coerce <$> (series :: Series m HTYPE_UNSIGNED_SHORT) +instance Monad m => Serial m CInt where + series = coerce <$> (series :: Series m HTYPE_INT) +instance Monad m => Serial m CUInt where + series = coerce <$> (series :: Series m HTYPE_UNSIGNED_INT) +instance Monad m => Serial m CLong where + series = coerce <$> (series :: Series m HTYPE_LONG) +instance Monad m => Serial m CULong where + series = coerce <$> (series :: Series m HTYPE_UNSIGNED_LONG) +instance Monad m => Serial m CPtrdiff where + series = coerce <$> (series :: Series m HTYPE_PTRDIFF_T) +instance Monad m => Serial m CSize where + series = coerce <$> (series :: Series m HTYPE_SIZE_T) +instance Monad m => Serial m CWchar where + series = coerce <$> (series :: Series m HTYPE_WCHAR_T) +instance Monad m => Serial m CSigAtomic where + series = coerce <$> (series :: Series m HTYPE_SIG_ATOMIC_T) +instance Monad m => Serial m CLLong where + series = coerce <$> (series :: Series m HTYPE_LONG_LONG) +instance Monad m => Serial m CULLong where + series = coerce <$> (series :: Series m HTYPE_UNSIGNED_LONG_LONG) +instance Monad m => Serial m CIntPtr where + series = coerce <$> (series :: Series m HTYPE_INTPTR_T) +instance Monad m => Serial m CUIntPtr where + series = coerce <$> (series :: Series m HTYPE_UINTPTR_T) +instance Monad m => Serial m CIntMax where + series = coerce <$> (series :: Series m HTYPE_INTMAX_T) +instance Monad m => Serial m CUIntMax where + series = coerce <$> (series :: Series m HTYPE_UINTMAX_T) diff --git a/test/Spec/Range.hs b/test/Spec/Range.hs new file mode 100644 index 000000000..ef0d346ff --- /dev/null +++ b/test/Spec/Range.hs @@ -0,0 +1,34 @@ +module Spec.Range + ( symmetric + , bounded + , singleton + , uniformRangeWithin + , uniformRangeWithinExcluded + ) where + +import System.Random.Monad + +symmetric :: (RandomGen g, UniformRange a, Eq a) => g -> (a, a) -> Bool +symmetric g (l, r) = fst (uniformR (l, r) g) == fst (uniformR (r, l) g) + +bounded :: (RandomGen g, UniformRange a, Ord a) => g -> (a, a) -> Bool +bounded g (l, r) = bottom <= result && result <= top + where + bottom = min l r + top = max l r + result = fst (uniformR (l, r) g) + +singleton :: (RandomGen g, UniformRange a, Eq a) => g -> a -> Bool +singleton g x = result == x + where + result = fst (uniformR (x, x) g) + +uniformRangeWithin :: (RandomGen g, UniformRange a, Ord a) => g -> (a, a) -> Bool +uniformRangeWithin gen (l, r) = + runStateGen_ gen $ \g -> + (\result -> min l r <= result && result <= max l r) <$> uniformRM (l, r) g + +uniformRangeWithinExcluded :: (RandomGen g, UniformRange a, Ord a) => g -> (a, a) -> Bool +uniformRangeWithinExcluded gen (l, r) = + runStateGen_ gen $ \g -> + (\result -> min l r <= result && (l == r || result < max l r)) <$> uniformRM (l, r) g diff --git a/test/Spec/Run.hs b/test/Spec/Run.hs new file mode 100644 index 000000000..c0e00edb8 --- /dev/null +++ b/test/Spec/Run.hs @@ -0,0 +1,12 @@ +module Spec.Run (runsEqual) where + +import Data.Word (Word64) +import System.Random.Monad + +runsEqual :: RandomGen g => g -> IO Bool +runsEqual g = do + let pureResult = runStateGen_ g uniformM :: Word64 + stResult = runSTGen_ g uniformM + ioResult <- runGenM_ (IOGen g) uniformM + atomicResult <- runGenM_ (AtomicGen g) uniformM + return $ all (pureResult ==) [stResult, ioResult, atomicResult] diff --git a/test/doctests.hs b/test/doctests.hs new file mode 100644 index 000000000..bbf7787f8 --- /dev/null +++ b/test/doctests.hs @@ -0,0 +1,12 @@ +module Main where + +import Build_doctests (flags, pkgs, module_sources) +import Data.Foldable (traverse_) +import Test.DocTest (doctest) + +main :: IO () +main = do + traverse_ putStrLn args + doctest args + where + args = flags ++ pkgs ++ module_sources diff --git a/tests/Makefile b/tests/Makefile deleted file mode 100644 index 39c71491b..000000000 --- a/tests/Makefile +++ /dev/null @@ -1,14 +0,0 @@ -# This Makefile runs the tests using GHC's testsuite framework. It -# assumes the package is part of a GHC build tree with the testsuite -# installed in ../../../testsuite. - -TOP=../../../testsuite -include $(TOP)/mk/boilerplate.mk -include $(TOP)/mk/test.mk - - -# Build tests locally without the central GHC testing infrastructure: -local: - ghc --make rangeTest.hs - ghc --make random1283.hs - diff --git a/tests/all.T b/tests/all.T deleted file mode 100644 index f1675ed5c..000000000 --- a/tests/all.T +++ /dev/null @@ -1,10 +0,0 @@ - -test('rangeTest', - normal, - compile_and_run, - ['']) - -test('random1283', - reqlib('containers'), - compile_and_run, - [' -package containers']) diff --git a/tests/rangeTest.stdout b/tests/rangeTest.stdout deleted file mode 100644 index 55ccaffb4..000000000 --- a/tests/rangeTest.stdout +++ /dev/null @@ -1,67 +0,0 @@ -Int: Passed -Integer: Passed -Int8: Passed -Int16: Passed -Int32: Passed -Int64: Passed -Word: Passed -Word8: Passed -Word16: Passed -Word32: Passed -Word64: Passed -Double: Passed -Float: Passed -CChar: Passed -CSChar: Passed -CUChar: Passed -CShort: Passed -CUShort: Passed -CInt: Passed -CUInt: Passed -CLong: Passed -CULong: Passed -CPtrdiff: Passed -CSize: Passed -CWchar: Passed -CSigAtomic: Passed -CLLong: Passed -CULLong: Passed -CIntPtr: Passed -CUIntPtr: Passed -CIntMax: Passed -CUIntMax: Passed -Int R: Passed -Integer R: Passed -Int8 R: Passed -Int8 Rsmall: Passed -Int8 Rmini: Passed -Int8 Rtrivial: Passed -Int16 R: Passed -Int32 R: Passed -Int64 R: Passed -Word R: Passed -Word8 R: Passed -Word16 R: Passed -Word32 R: Passed -Word64 R: Passed -Double R: Passed -Float R: Passed -CChar R: Passed -CSChar R: Passed -CUChar R: Passed -CShort R: Passed -CUShort R: Passed -CInt R: Passed -CUInt R: Passed -CLong R: Passed -CULong R: Passed -CPtrdiff R: Passed -CSize R: Passed -CWchar R: Passed -CSigAtomic R: Passed -CLLong R: Passed -CULLong R: Passed -CIntPtr R: Passed -CUIntPtr R: Passed -CIntMax R: Passed -CUIntMax R: Passed diff --git a/tests/slowness.hs b/tests/slowness.hs deleted file mode 100644 index 162a4cff0..000000000 --- a/tests/slowness.hs +++ /dev/null @@ -1,55 +0,0 @@ - --- http://hackage.haskell.org/trac/ghc/ticket/427 - --- Two (performance) problems in one: - -{-# OPTIONS -fffi #-} -module Main (main) where - -import Control.Monad -import System.Random - -foreign import ccall unsafe "random" _crandom :: IO Int -foreign import ccall unsafe "stdlib.hs" rand :: IO Int - -randomInt :: (Int, Int) -> IO Int -randomInt (min,max) = do --- n <- _crandom - n <- rand - return $ min + n `rem` range - where - range = max - min + 1 - -main = replicateM_ (5*10^6) $ do - x <- randomRIO (0::Int,1000) :: IO Int --- x <- randomInt (0::Int,1000) :: IO Int - x `seq` return () - return () - --- First, without the "seq" at the end, hardly anything is --- evaluated and we're building huge amounts of thunks. --- Three ideas about this one: --- - Blame the user :) --- - data StdGen = StdGen !Int !Int --- Use strict fields in StdGen. Doesn't actually help --- (at least in this example). --- - Force evaluation of the StdGen in getStdRandom. --- Does help in this example, but also changes behaviour --- of the library: --- x <- randomRIO undefined --- currently dies only when x (or the result of a later --- randomRIO) is evaluated. This change causes it to die --- immediately. - --- Second, even _with_ the "seq", replacing "randomRIO" by --- "randomInt" speeds the thing up with a factor of about --- 30. (2 to 3.6, in a "real world" university practicum --- exercise of 900 lines of code) --- Even given the fact that they're not really doing the --- same thing, this seems rather much :( - --------------------------------------------------------------------------------- - --- [2011.06.28] RRN: --- I'm currently seeing 1425 ms vs 43 ms for the above. 33X --- difference. If I use rand() instead it's about 52ms.