On 02-Feb-2000, Simon Peyton-Jones <[EMAIL PROTECTED]> wrote:
> Guess what? I have been driven to the conclusion that
> the Haskell 98 Random library contains a fundamental, but
> easily fixed flaw. So I propose to declare a new 'typo'.
> (Yes, I am all too aware of the danger of typo-creep.)
>
> The library is currently based on the RandomGen class,
> whose signature is:
>
> class RandomGen g wher
> next :: g -> (Int, g)
> split :: g -> (g, g)
>
> The difficulty is that there is absolutely no robust way to
> use such a generator to generate random numbers uniformly
> distributed in a given range.
I think that is not true.
> Why not? Because there's no
> clue as to the precise range of Ints produced by next.
You can figure out the range as you go along, by seeing what `next'
returns. This makes things a bit more complex that it would otherwise
need to be, but it's still possible. See the code below.
> Suppose next produced output in the range 0..16. (It's specified
> to produce "at least 30 bits", but the argument still holds.)
> Suppose we want random numbers in the range 0..10. We can't just
> take "mod 11" because that would produce too many values in the
> range 0..5. Indeed, if we don't know the range of next's result,
> we simply cannot generate uniformly distributed random numbers.
I think the following function does the trick:
import Random
-- returns an infinite list of random integers in the range [0, n).
randsN :: RandomGen g => Integer -> g -> [Integer]
randsN n g = (x:xs) where
(x, g1) = randN n g
xs = randsN n g1
-- returns a random integer in the range [0, n),
-- together with the new state of the random number generator
randN :: RandomGen g => Integer -> g -> (Integer, g)
randN n _ | n <= 0 = error "randN: n <= 0"
randN n g = solve 0 1 n g (rands g)
--
-- `rands' is an infinite list containing random integers.
-- Each element of `rands' is a triple (r_i, max_i, g_i)
-- where max_i is an arbitrary but non-random positive integer,
-- r_i is a random integer evenly distributed in
-- the range [0, max_i), and g_i is the state of the
-- random number generator after having generated r_i.
--
--rands :: [(Integer, Integer, g)]
rands g = getRands leastMin leastMax g where
leastMin, leastMax :: Int
leastMin = -(2^29)
leastMax = 2^29 - 1
getRands min max g0 =
let (k, g1) = next g0 in
if k < min then getRands k max g1 -- use k as new min
else if k > max then getRands min k g1 -- use k as new max
else let min' = fromInt min
max' = fromInt max
k' = fromInt k
this = (k' - min', max' - min' + 1, g1)
rest = getRands min max g1
in (this : rest)
--
-- Compute a random integer in the range [0, n)
-- given a random integer `val' in the range [0, max0)
--
solve :: Integer -> Integer -> Integer -> g ->
[(Integer, Integer, g)] -> (Integer, g)
solve val max0 n g rands =
if max0 >= n then
let max0' = (max0 `div` n) * n in
if val < max0' then
(val `mod` n, g)
else
solve 0 1 n g rands
else
-- we need more random bits
let ((r, max1, g1):rest) = rands in
solve (val*max1 + r) (max0 * max1) n g1 rest
I don't object to improving the interface to Random...
but given the existence of solutions like the one above,
perhaps it should wait until the next revision of Haskell?
--
Fergus Henderson <[EMAIL PROTECTED]> | "I have always known that the pursuit
WWW: <http://www.cs.mu.oz.au/~fjh> | of excellence is a lethal habit"
PGP: finger [EMAIL PROTECTED] | -- the last words of T. S. Garp.
