[Newbies] Re: What's wrong with this statement?
ncellier at ifrance.com
Fri Aug 1 01:00:20 UTC 2008
<johnps11 <at> bigpond.com> writes:
> The issue is almost certainly in atRandom:
> 30 timesRepeat: [ | a |
> a:= (2 raisedTo: 57) atRandom.
> Transcript show: a ; show: ' ' ; show: a even ; cr ].
> Gives me 30 falses.
> The odds of that are less than one in a billion (assuming uniformly
> distributed integers).
> This doesn't happen for 2^56 or lower, and does for all 2^n, n>= 57.
> My guess is that it's either a hardware issue or something about Random.
> I'm not clever enough to understand PRNGs, so I'll leave it others to work
> out what the answer is, although I suspect that for more than 56 bits you
> need a PRNG that uses LargeIntegers and not Floats.
Nothing surprising once again.
Float have a 53 bits integer mantissa m and and a biased exponent e:
f := m * (2 raisedTo: e).
m highBit = 53.
The PRNG has a uniform distribution on [0,1).
f < 1 means 100% of floats have e < -53.
Let us divide the uniform intervals in two recursively:
50% have e = -54 (0.5<=f<1)
50% have e < -54 (f<0.5)
25% have e = -55 (0.25<=f<0.5)
25% have e < -55 (f<0.25)
75% have e >=-55 (f>=0.25)
12.5% have e = -56 (0.125<=f<0.25)
12.5% have e < -56 (f<0.125)
87.5% have e >=-56 (f>=0.125)
6.25% have e = -57 (0.0625<=f<0.125)
6.25% have e < -57 (f<0.0625)
92.75% have e >=-57 (f>=0.0625)
3.125% have e = -58 (0.03125<=f<0.0625)
3.125% have e < -58 (f<0.03125)
96.875% have e >=-58 (f>=0.03125)
So when you multiply f by (2 raisedTo: 57),
that means (m bitShift: e+57),
87.5% have (e+57>0),
thus are shifted at least one bit to the left,
and you get 87.5% of even numbers.
since algorithm truncate then add 1, you then get 87.5% of odd results
I suspect that due to division by 16r7FFFFFFF, last bit of mantissa m is often
zero. In which case, you get 92.75% of odd results.
Or maybe last two bits, in which case you have 96.875% of odd results.
((1 to: 10000) count: [:i | (2 raisedTo: 57) atRandom odd]) / 10000.0 .
Yes, probably last two bits of these floats are set to 0...
The limit is by no way (2 raisedTo: 57).
The true limit at which you won't get any even number ?
well, smallest generated pseudo random float is 1/(2^31-1)
A little higher than 1.0 timesTwoPoser: -31.
smallest := ((1.0 timesTwoPower: 31) - 1) reciprocal.
smallest significandAsInteger printStringBase: 2.
The biased exponent is then (-31-52=-83).
(smallest significandAsInteger asFloat timesTwoPower: smallest exponent-52)
So, (2 raisedTo: 84) atRandom is guaranteed to be odd.
We conjectured that last two bits are zero,
maybe (2 raisedTo: 82) atRandom is already 100% odd?
(1 to: 10000000 by: 11) detect: [:e | ((e * smallest) significandAsInteger
bitAnd: 2r11) > 0] ifNone: [nil]
Hmm the conjecture is false...
Floats are tricky, but once you understand they are just rational with
denominator being a power of two and inexact operations truncating the mantissa
at 53 bits, the mystery vanishes a bit.
More information about the Beginners