ghc/lib/std/Random.lhs

   1 %
   2 % (c) The GRASP/AQUA Project, Glasgow University, 1995-99
   3 %
   4
   5
   6 The June 1988 (v31 #6) issue of the Communications of the ACM has an
   7 article by Pierre L'Ecuyer called, "Efficient and Portable Combined
   8 Random Number Generators".  Here is the Portable Combined Generator of
   9 L'Ecuyer for 32-bit computers.  It has a period of roughly 2.30584e18.
  10
  11 Transliterator: Lennart Augustsson
  12
  13 sof 1/99 - code brought (kicking and screaming) into the new Random
  14 world..
  15
  16 \begin{code}
  17 module Random
  18         (
  19           RandomGen(next, split)
  20         , StdGen
  21         , mkStdGen
  22         , Random ( random,   randomR,
  23                    randoms,  randomRs,
  24                    randomIO, randomRIO )
  25         , getStdRandom
  26         , getStdGen
  27         , setStdGen
  28         , newStdGen
  29         ) where
  30
  31 #ifndef __HUGS__
  32 import PrelGHC          ( RealWorld )
  33 import PrelNum          ( fromInt )
  34 import PrelShow         ( showSignedInt, showSpace )
  35 import PrelRead         ( readDec )
  36 import PrelIOBase       ( unsafePerformIO, stToIO )
  37 import PrelArr          ( MutableVar, newVar, readVar, writeVar )
  38 import PrelReal         ( toInt )
  39 import PrelFloat        ( float2Double, double2Float )
  40 import Time             ( getClockTime, ClockTime(..) )
  41 #endif
  42 import CPUTime          ( getCPUTime )
  43 import Prelude
  44 import privileged Prelude
  45                         ( IORef
  46                         , newIORef
  47                         , readIORef
  48                         , writeIORef
  49                         , unsafePerformIO
  50                         )
  51
  52
  53 import Char ( isSpace, chr, ord )
  54 \end{code}
  55
  56 \begin{code}
  57 class RandomGen g where
  58    next  :: g -> (Int, g)
  59    split :: g -> (g, g)
  60
  61 \end{code}
  62
  63 \begin{code}
  64 data StdGen
  65  = StdGen Int Int
  66
  67 instance RandomGen StdGen where
  68   next  = stdNext
  69   split = stdSplit
  70
  71 #ifdef __HUGS__
  72 instance Show StdGen where
  73   showsPrec p (StdGen s1 s2) =
  74      showsPrec p s1 .
  75      showChar ' ' .
  76      showsPrec p s2
  77 #else
  78 instance Show StdGen where
  79   showsPrec p (StdGen s1 s2) =
  80      showSignedInt p s1 .
  81      showSpace          .
  82      showSignedInt p s2
  83 #endif
  84
  85 instance Read StdGen where
  86   readsPrec _p = \ r ->
  87      case try_read r of
  88        r@[_] -> r
  89        _   -> [stdFromString r] -- because it shouldn't ever fail.
  90     where
  91       try_read r = do
  92          (s1, r1) <- readDec (dropWhile isSpace r)
  93          (s2, r2) <- readDec (dropWhile isSpace r1)
  94          return (StdGen s1 s2, r2)
  95
  96 {-
  97  If we cannot unravel the StdGen from a string, create
  98  one based on the string given.
  99 -}
 100 stdFromString         :: String -> (StdGen, String)
 101 stdFromString s        = (mkStdGen num, rest)
 102         where (cs, rest) = splitAt 6 s
 103               num        = foldl (\a x -> x + 3 * a) 1 (map ord cs)
 104 \end{code}
 105
 106 \begin{code}
 107 mkStdGen :: Int -> StdGen -- why not Integer ?
 108 mkStdGen s
 109  | s < 0     = mkStdGen (-s)
 110  | otherwise = StdGen (s1+1) (s2+1)
 111       where
 112         (q, s1) = s `divMod` 2147483562
 113         s2      = q `mod` 2147483398
 114
 115 createStdGen :: Integer -> StdGen
 116 createStdGen s
 117  | s < 0     = createStdGen (-s)
 118  | otherwise = StdGen (toInt (s1+1)) (toInt (s2+1))
 119       where
 120         (q, s1) = s `divMod` 2147483562
 121         s2      = q `mod` 2147483398
 122
 123 \end{code}
 124
 125 The class definition - see library report for details.
 126
 127 \begin{code}
 128 class Random a where
 129   -- Minimal complete definition: random and randomR
 130   random  :: RandomGen g => g -> (a, g)
 131   randomR :: RandomGen g => (a,a) -> g -> (a,g)
 132
 133   randoms  :: RandomGen g => g -> [a]
 134   randoms  g      = x : randoms g' where (x,g') = random g
 135
 136   randomRs :: RandomGen g => (a,a) -> g -> [a]
 137   randomRs ival g = x : randomRs ival g' where (x,g') = randomR ival g
 138
 139   randomIO  :: IO a
 140   randomIO         = getStdRandom random
 141
 142   randomRIO :: (a,a) -> IO a
 143   randomRIO range  = getStdRandom (randomR range)
 144 \end{code}
 145
 146 \begin{code}
 147 instance Random Int where
 148   randomR (a,b) g = randomIvalInteger (toInteger a, toInteger b) g
 149   random g        = randomR (minBound,maxBound) g
 150
 151 instance Random Char where
 152   randomR (a,b) g =
 153       case (randomIvalInteger (toInteger (ord a), toInteger (ord b)) g) of
 154         (x,g) -> (chr x, g)
 155   random g        = randomR (minBound,maxBound) g
 156
 157 instance Random Bool where
 158   randomR (a,b) g =
 159       case (randomIvalInteger (toInteger (bool2Int a), toInteger (bool2Int b)) g) of
 160         (x, g) -> (int2Bool x, g)
 161        where
 162          bool2Int False = 0
 163          bool2Int True  = 1
 164
 165          int2Bool 0     = False
 166          int2Bool _     = True
 167
 168   random g        = randomR (minBound,maxBound) g
 169
 170 instance Random Integer where
 171   randomR ival g = randomIvalInteger ival g
 172   random g       = randomR (toInteger (minBound::Int), toInteger (maxBound::Int)) g
 173
 174 instance Random Double where
 175   randomR ival g = randomIvalDouble ival id g
 176   random g       = randomR (0::Double,1) g
 177
 178 -- hah, so you thought you were saving cycles by using Float?
 179
 180 #ifdef __HUGS__
 181 instance Random Float where
 182   random g        = randomIvalDouble (0::Double,1) realToFrac g
 183   randomR (a,b) g = randomIvalDouble (realToFrac a, realToFrac b) realToFrac g
 184 #else
 185 instance Random Float where
 186   randomR (a,b) g = randomIvalDouble (float2Double a, float2Double b) double2Float g
 187   random g        = randomIvalDouble (0::Double,1) double2Float g
 188 #endif
 189
 190 \end{code}
 191
 192
 193 \begin{code}
 194 #ifdef __HUGS__
 195 mkStdRNG :: Integer -> IO StdGen
 196 mkStdRNG o = do
 197     ct          <- getCPUTime
 198     return (createStdGen (ct + o))
 199 #else
 200 mkStdRNG :: Integer -> IO StdGen
 201 mkStdRNG o = do
 202     ct          <- getCPUTime
 203     (TOD sec _) <- getClockTime
 204     return (createStdGen (sec * 12345 + ct + o))
 205 #endif
 206
 207 randomIvalInteger :: (RandomGen g, Num a) => (Integer, Integer) -> g -> (a, g)
 208 randomIvalInteger (l,h) rng
 209  | l > h     = randomIvalInteger (h,l) rng
 210  | otherwise = case (f n 1 rng) of (v, rng') -> (fromInteger (l + v `mod` k), rng')
 211      where
 212        k = h - l + 1
 213        b = 2147483561
 214        n = iLogBase b k
 215
 216        f 0 acc g = (acc, g)
 217        f n acc g =
 218           let
 219            (x,g')   = next g
 220           in
 221           f (n-1) (fromInt x + acc * b) g'
 222
 223 randomIvalDouble :: (RandomGen g, Fractional a) => (Double, Double) -> (Double -> a) -> g -> (a, g)
 224 randomIvalDouble (l,h) fromDouble rng
 225   | l > h     = randomIvalDouble (h,l) fromDouble rng
 226   | otherwise =
 227        case (randomIvalInteger (toInteger (minBound::Int), toInteger (maxBound::Int)) rng) of
 228          (x, rng') ->
 229             let
 230              scaled_x =
 231                 fromDouble ((l+h)/2) +
 232                 fromDouble ((h-l) / realToFrac intRange) *
 233                 fromIntegral (x::Int)
 234             in
 235             (scaled_x, rng')
 236
 237 intRange :: Integer
 238 intRange  = toInteger (maxBound::Int) - toInteger (minBound::Int)
 239
 240 iLogBase :: Integer -> Integer -> Integer
 241 iLogBase b i = if i < b then 1 else 1 + iLogBase b (i `div` b)
 242
 243 stdNext :: StdGen -> (Int, StdGen)
 244 stdNext (StdGen s1 s2) = (z', StdGen s1'' s2'')
 245         where   z'   = if z < 1 then z + 2147483562 else z
 246                 z    = s1'' - s2''
 247
 248                 k    = s1 `quot` 53668
 249                 s1'  = 40014 * (s1 - k * 53668) - k * 12211
 250                 s1'' = if s1' < 0 then s1' + 2147483563 else s1'
 251
 252                 k'   = s2 `quot` 52774
 253                 s2'  = 40692 * (s2 - k' * 52774) - k' * 3791
 254                 s2'' = if s2' < 0 then s2' + 2147483399 else s2'
 255
 256 stdSplit            :: StdGen -> (StdGen, StdGen)
 257 stdSplit std@(StdGen s1 s2)
 258                      = (left, right)
 259                        where
 260                         -- no statistical foundation for this!
 261                         left    = StdGen new_s1 t2
 262                         right   = StdGen t1 new_s2
 263
 264                         new_s1 | s1 == 2147483562 = 1
 265                                | otherwise        = s1 + 1
 266
 267                         new_s2 | s2 == 1          = 2147483398
 268                                | otherwise        = s2 - 1
 269
 270                         StdGen t1 t2 = snd (next std)
 271 \end{code}
 272
 273
 274 \begin{code}
 275 #ifdef __HUGS__
 276
 277 setStdGen :: StdGen -> IO ()
 278 setStdGen sgen = writeIORef theStdGen sgen
 279
 280 getStdGen :: IO StdGen
 281 getStdGen  = readIORef theStdGen
 282
 283 theStdGen :: IORef StdGen
 284 theStdGen  = unsafePerformIO (newIORef (createStdGen 0))
 285
 286 #else
 287
 288 global_rng :: MutableVar RealWorld StdGen
 289 global_rng = unsafePerformIO $ do
 290    rng <- mkStdRNG 0
 291    stToIO (newVar rng)
 292
 293 setStdGen :: StdGen -> IO ()
 294 setStdGen sgen = stToIO (writeVar global_rng sgen)
 295
 296 getStdGen :: IO StdGen
 297 getStdGen = stToIO (readVar global_rng)
 298
 299 #endif
 300
 301
 302 newStdGen :: IO StdGen
 303 newStdGen = do
 304   rng <- getStdGen
 305   let (a,b) = split rng
 306   setStdGen a
 307   return b
 308
 309 getStdRandom :: (StdGen -> (a,StdGen)) -> IO a
 310 getStdRandom f = do
 311    rng          <- getStdGen
 312    let (v, new_rng) = f rng
 313    setStdGen new_rng
 314    return v
 315 \end{code}