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