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