System/Random.hs

   1 -----------------------------------------------------------------------------
   2 -- |
   3 -- Module      :  System.Random
   4 -- Copyright   :  (c) The University of Glasgow 2001
   5 -- License     :  BSD-style (see the file libraries/base/LICENSE)
   6 --
   7 -- Maintainer  :  libraries@haskell.org
   8 -- Stability   :  provisional
   9 -- Portability :  portable
  10 --
  11 -- Random numbers.
  12 --
  13 -----------------------------------------------------------------------------
  14
  15 module System.Random
  16         (
  17
  18         -- $intro
  19
  20         -- * The 'RandomGen' class, and the 'StdGen' generator
  21
  22           RandomGen(next, split, genRange)
  23         , StdGen
  24         , mkStdGen
  25
  26         -- * The 'Random' class
  27         , Random ( random,   randomR,
  28                    randoms,  randomRs,
  29                    randomIO, randomRIO )
  30
  31         -- * The global random number generator
  32
  33         -- $globalrng
  34
  35         , getStdRandom
  36         , getStdGen
  37         , setStdGen
  38         , newStdGen
  39
  40         -- * References
  41         -- $references
  42
  43         ) where
  44
  45 import Prelude
  46
  47 import System.CPUTime   ( getCPUTime )
  48 import System.Time      ( getClockTime, ClockTime(..) )
  49 import Data.Char        ( isSpace, chr, ord )
  50 import System.IO.Unsafe ( unsafePerformIO )
  51 import Data.IORef
  52 import Numeric          ( readDec )
  53
  54 #ifdef __GLASGOW_HASKELL__
  55 import GHC.IOBase       ( stToIO )
  56 #endif
  57
  58 {- $intro
  59
  60 The "Random" library deals with the common task of pseudo-random
  61 number generation. The library makes it possible to generate
  62 repeatable results, by starting with a specified initial random
  63 number generator; or to get different results on each run by using the
  64 system-initialised generator, or by supplying a seed from some other
  65 source.
  66
  67 The library is split into two layers:
  68
  69 * A core /random number generator/ provides a supply of bits. The class
  70 'RandomGen' provides a common interface to such generators.
  71
  72 * The class 'Random' provides a way to extract particular values from
  73 a random number generator. For example, the 'Float' instance of 'Random'
  74 allows one to generate random values of type 'Float'.
  75
  76 [Comment found in this file when merging with Library Report:]
  77
  78 The June 1988 (v31 \#6) issue of the Communications of the ACM has an
  79 article by Pierre L'Ecuyer called, /Efficient and Portable Combined
  80 Random Number Generators/.  Here is the Portable Combined Generator of
  81 L'Ecuyer for 32-bit computers.  It has a period of roughly 2.30584e18.
  82
  83 Transliterator: Lennart Augustsson
  84
  85 -}
  86
  87 -- |RandomGen
  88 -- The class 'RandomGen' provides a common interface to random number generators.
  89
  90 class RandomGen g where
  91
  92    -- |The 'next' operation allows one to extract at least 30 bits (one 'Int''s
  93    -- worth) from the generator, returning a new generator as well.  The
  94    -- integer returned may be positive or negative.
  95    next     :: g -> (Int, g)
  96
  97    -- |The 'split' operation allows one to obtain two distinct random number
  98    -- generators. This is very useful in functional programs (for example, when
  99    -- passing a random number generator down to recursive calls), but very
 100    -- little work has been done on statistically robust implementations of
 101    -- @split ([1,4]@ are the only examples we know of).
 102    split    :: g -> (g, g)
 103
 104    genRange :: g -> (Int,Int)
 105
 106    -- default mathod
 107    genRange g = (minBound,maxBound)
 108
 109 {- |The "Random" library provides one instance of 'RandomGen', the abstract data
 110 type 'StdGen'.
 111
 112 The result of repeatedly using next should be at least as statistically robust
 113 as the /Minimal Standard Random Number Generator/ described by
 114 ["System.Random\#Park", "System.Random\#Carta"]. Until more
 115 is known about implementations of 'split', all we require is that 'split' deliver
 116 generators that are (a) not identical and (b) independently robust in the sense
 117 just given.
 118
 119 The 'show'\/'Read' instances of 'StdGen' provide a primitive way to save the
 120 state of a random number generator. It is required that @read (show g) == g@.
 121
 122 In addition, 'read' may be used to map an arbitrary string (not necessarily one
 123 produced by 'show') onto a value of type 'StdGen'. In general, the 'read'
 124 instance of 'StdGen' has the following properties:
 125
 126 * It guarantees to succeed on any string.
 127
 128 *It guarantees to consume only a finite portion of the string.
 129
 130 * Different argument strings are likely to result in different results.
 131
 132 The function 'mkStdGen' provides an alternative way of producing an initial
 133 generator, by mapping an 'Int' into a generator. Again, distinct arguments
 134 should be likely to produce distinct generators.
 135
 136 Programmers may, of course, supply their own instances of 'RandomGen'.
 137
 138 -}
 139
 140 data StdGen
 141  = StdGen Int Int
 142
 143 instance RandomGen StdGen where
 144   next  = stdNext
 145   split = stdSplit
 146
 147 instance Show StdGen where
 148   showsPrec p (StdGen s1 s2) =
 149      showsPrec p s1 .
 150      showChar ' ' .
 151      showsPrec p s2
 152
 153 instance Read StdGen where
 154   readsPrec _p = \ r ->
 155      case try_read r of
 156        r@[_] -> r
 157        _   -> [stdFromString r] -- because it shouldn't ever fail.
 158     where
 159       try_read r = do
 160          (s1, r1) <- readDec (dropWhile isSpace r)
 161          (s2, r2) <- readDec (dropWhile isSpace r1)
 162          return (StdGen s1 s2, r2)
 163
 164 {-
 165  If we cannot unravel the StdGen from a string, create
 166  one based on the string given.
 167 -}
 168 stdFromString         :: String -> (StdGen, String)
 169 stdFromString s        = (mkStdGen num, rest)
 170         where (cs, rest) = splitAt 6 s
 171               num        = foldl (\a x -> x + 3 * a) 1 (map ord cs)
 172
 173
 174 mkStdGen :: Int -> StdGen -- why not Integer ?
 175 mkStdGen s
 176  | s < 0     = mkStdGen (-s)
 177  | otherwise = StdGen (s1+1) (s2+1)
 178       where
 179         (q, s1) = s `divMod` 2147483562
 180         s2      = q `mod` 2147483398
 181
 182 createStdGen :: Integer -> StdGen
 183 createStdGen s
 184  | s < 0     = createStdGen (-s)
 185  | otherwise = StdGen (fromInteger (s1+1)) (fromInteger (s2+1))
 186       where
 187         (q, s1) = s `divMod` 2147483562
 188         s2      = q `mod` 2147483398
 189
 190 -- FIXME: 1/2/3 below should be ** (vs@30082002) XXX
 191
 192 {- |The 'Random' class
 193 With a source of random number supply in hand, the 'Random' class allows the
 194 programmer to extract random values of a variety of types.
 195
 196 * 'randomR' takes a range /(lo,hi)/ and a random number generator /g/, and returns
 197 a random value uniformly distributed in the closed interval /[lo,hi]/, together
 198 with a new generator. It is unspecified what happens if /lo>hi/. For continuous
 199 types there is no requirement that the values /lo/ and /hi/ are ever produced,
 200 but they may be, depending on the implementation and the interval.
 201
 202 * 'random' does the same as 'randomR', but does not take a range.
 203
 204 (1) For bounded types (instances of 'Bounded', such as 'Char'), the range is
 205 normally the whole type.
 206
 207 (2) For fractional types, the range is normally the semi-closed interval @[0,1)@.
 208
 209 (3) For 'Integer', the range is (arbitrarily) the range of 'Int'.
 210
 211 * The plural versions, 'randomRs' and 'randoms', produce an infinite list of
 212 random values, and do not return a new generator.
 213
 214 * The 'IO' versions, 'randomRIO' and 'randomIO', use the global random number
 215 generator (see Section 17.3
 216 <http://www.haskell.org/onlinelibrary/random.html#global-rng>).
 217 -}
 218
 219 class Random a where
 220   -- |Minimal complete definition: 'random' and 'randomR'
 221   random  :: RandomGen g => g -> (a, g)
 222   randomR :: RandomGen g => (a,a) -> g -> (a,g)
 223
 224   -- |Default methods
 225   randoms  :: RandomGen g => g -> [a]
 226   randoms  g      = x : randoms g' where (x,g') = random g
 227
 228   randomRs :: RandomGen g => (a,a) -> g -> [a]
 229   randomRs ival g = x : randomRs ival g' where (x,g') = randomR ival g
 230
 231   randomIO  :: IO a
 232   randomIO         = getStdRandom random
 233
 234   randomRIO :: (a,a) -> IO a
 235   randomRIO range  = getStdRandom (randomR range)
 236
 237
 238 instance Random Int where
 239   randomR (a,b) g = randomIvalInteger (toInteger a, toInteger b) g
 240   random g        = randomR (minBound,maxBound) g
 241
 242 instance Random Char where
 243   randomR (a,b) g =
 244       case (randomIvalInteger (toInteger (ord a), toInteger (ord b)) g) of
 245         (x,g) -> (chr x, g)
 246   random g        = randomR (minBound,maxBound) g
 247
 248 instance Random Bool where
 249   randomR (a,b) g =
 250       case (randomIvalInteger (toInteger (bool2Int a), toInteger (bool2Int b)) g) of
 251         (x, g) -> (int2Bool x, g)
 252        where
 253          bool2Int False = 0
 254          bool2Int True  = 1
 255
 256          int2Bool 0     = False
 257          int2Bool _     = True
 258
 259   random g        = randomR (minBound,maxBound) g
 260
 261 instance Random Integer where
 262   randomR ival g = randomIvalInteger ival g
 263   random g       = randomR (toInteger (minBound::Int), toInteger (maxBound::Int)) g
 264
 265 instance Random Double where
 266   randomR ival g = randomIvalDouble ival id g
 267   random g       = randomR (0::Double,1) g
 268
 269 -- hah, so you thought you were saving cycles by using Float?
 270 instance Random Float where
 271   random g        = randomIvalDouble (0::Double,1) realToFrac g
 272   randomR (a,b) g = randomIvalDouble (realToFrac a, realToFrac b) realToFrac g
 273
 274 mkStdRNG :: Integer -> IO StdGen
 275 mkStdRNG o = do
 276     ct          <- getCPUTime
 277     (TOD sec _) <- getClockTime
 278     return (createStdGen (sec * 12345 + ct + o))
 279
 280 randomIvalInteger :: (RandomGen g, Num a) => (Integer, Integer) -> g -> (a, g)
 281 randomIvalInteger (l,h) rng
 282  | l > h     = randomIvalInteger (h,l) rng
 283  | otherwise = case (f n 1 rng) of (v, rng') -> (fromInteger (l + v `mod` k), rng')
 284      where
 285        k = h - l + 1
 286        b = 2147483561
 287        n = iLogBase b k
 288
 289        f 0 acc g = (acc, g)
 290        f n acc g =
 291           let
 292            (x,g')   = next g
 293           in
 294           f (n-1) (fromIntegral x + acc * b) g'
 295
 296 randomIvalDouble :: (RandomGen g, Fractional a) => (Double, Double) -> (Double -> a) -> g -> (a, g)
 297 randomIvalDouble (l,h) fromDouble rng
 298   | l > h     = randomIvalDouble (h,l) fromDouble rng
 299   | otherwise =
 300        case (randomIvalInteger (toInteger (minBound::Int), toInteger (maxBound::Int)) rng) of
 301          (x, rng') ->
 302             let
 303              scaled_x =
 304                 fromDouble ((l+h)/2) +
 305                 fromDouble ((h-l) / realToFrac intRange) *
 306                 fromIntegral (x::Int)
 307             in
 308             (scaled_x, rng')
 309
 310 intRange :: Integer
 311 intRange  = toInteger (maxBound::Int) - toInteger (minBound::Int)
 312
 313 iLogBase :: Integer -> Integer -> Integer
 314 iLogBase b i = if i < b then 1 else 1 + iLogBase b (i `div` b)
 315
 316 stdNext :: StdGen -> (Int, StdGen)
 317 stdNext (StdGen s1 s2) = (z', StdGen s1'' s2'')
 318         where   z'   = if z < 1 then z + 2147483562 else z
 319                 z    = s1'' - s2''
 320
 321                 k    = s1 `quot` 53668
 322                 s1'  = 40014 * (s1 - k * 53668) - k * 12211
 323                 s1'' = if s1' < 0 then s1' + 2147483563 else s1'
 324
 325                 k'   = s2 `quot` 52774
 326                 s2'  = 40692 * (s2 - k' * 52774) - k' * 3791
 327                 s2'' = if s2' < 0 then s2' + 2147483399 else s2'
 328
 329 stdSplit            :: StdGen -> (StdGen, StdGen)
 330 stdSplit std@(StdGen s1 s2)
 331                      = (left, right)
 332                        where
 333                         -- no statistical foundation for this!
 334                         left    = StdGen new_s1 t2
 335                         right   = StdGen t1 new_s2
 336
 337                         new_s1 | s1 == 2147483562 = 1
 338                                | otherwise        = s1 + 1
 339
 340                         new_s2 | s2 == 1          = 2147483398
 341                                | otherwise        = s2 - 1
 342
 343                         StdGen t1 t2 = snd (next std)
 344
 345 -- The global random number generator
 346
 347 {- $globalrng
 348
 349 There is a single, implicit, global random number generator of type
 350 'StdGen', held in some global variable maintained by the 'IO' monad. It is
 351 initialised automatically in some system-dependent fashion, for example, by
 352 using the time of day, or Linux's kernel random number generator. To get
 353 deterministic behaviour, use 'setStdGen'.
 354 -}
 355
 356 -- |'setStdGen' sets the global random number generator.
 357 setStdGen :: StdGen -> IO ()
 358 setStdGen sgen = writeIORef theStdGen sgen
 359
 360 -- |'getStdGen' gets the global random number generator.
 361 getStdGen :: IO StdGen
 362 getStdGen  = readIORef theStdGen
 363
 364 -- |'newStdGen' applies 'split' to the current global random generator, updates it
 365 -- with one of the results, and returns the other.
 366 theStdGen :: IORef StdGen
 367 theStdGen  = unsafePerformIO $ do
 368    rng <- mkStdRNG 0
 369    newIORef rng
 370
 371 newStdGen :: IO StdGen
 372 newStdGen = do
 373   rng <- getStdGen
 374   let (a,b) = split rng
 375   setStdGen a
 376   return b
 377
 378 {- |'getStdRandom' uses the supplied function to get a value from the current
 379 global random generator, and updates the global generator with the new generator
 380 returned by the function. For example, 'rollDice' gets a random integer between 1 and 6:
 381
 382 >  rollDice :: IO Int
 383 >  rollDice = getStdRandom (randomR (1,6))
 384
 385 -}
 386
 387 getStdRandom :: (StdGen -> (a,StdGen)) -> IO a
 388 getStdRandom f = do
 389    rng          <- getStdGen
 390    let (v, new_rng) = f rng
 391    setStdGen new_rng
 392    return v
 393
 394 {- $references
 395
 396 * [1] FW Burton and RL Page, /Distributed random number generation/,
 397 Journal of Functional Programming, 2(2):203-212, April 1992.
 398
 399 * [2] SK #Park# Park, and KW Miller, /Random number generators -
 400 good ones are hard to find/, Comm ACM 31(10), Oct 1988, pp1192-1201.
 401
 402 * [3] DG #Carta# Carta, /Two fast implementations of the minimal standard
 403 random number generator/, Comm ACM, 33(1), Jan 1990, pp87-88.
 404
 405 * [4] P Hellekalek, /Don\'t trust parallel Monte Carlo/,
 406 Department of Mathematics, University of Salzburg,
 407 <http://random.mat.sbg.ac.at/~peter/pads98.ps>, 1998.
 408
 409 The Web site <http://random.mat.sbg.ac.at/> is a great source of information.
 410
 411 -}