<!DOCTYPE HTML PUBLIC "-//W3C//DTD HTML 4.01 Transitional//EN">
<html>
<head>
<meta content="text/html; charset=ISO-8859-1"
http-equiv="Content-Type">
</head>
<body bgcolor="#ffffff" text="#000000">
On 03/02/2012 09:02 PM, Richard O'Keefe wrote:
<blockquote
cite="mid:209543D1-42F8-44CB-B1EB-6DA3F71F9209@cs.otago.ac.nz"
type="cite">
<pre wrap="">
On 3/03/2012, at 1:32 PM, Michael Truog wrote:
</pre>
<blockquote type="cite">
<pre wrap="">If you use this implementation, then it would provide 408 bits of randomness, uniformly distributed.
</pre>
</blockquote>
<pre wrap="">
I don't see how it can. Just because all the numbers are in a given range does not mean that
all the numbers in that range are possible. The Wichmann-Hill generator has a state with
4 31-bit numbers, meaning that you can't possibly get more than 124 bits out of it. If you
really get more, you have a different algorithm.
Using 64-bit integer arithmetic, the algorithm looks like
a = (a * 11600LL) % 2147483579;
b = (b * 47003LL) % 2147483543;
c = (c * 23000LL) % 2147483423;
d = (d * 33000LL) % 2147483123;
w = a/2147483579.0 + b/2147483543.0
+ c/2147483423.0 + d/2147483123.0;
if (w >= 2.0) w -= 2.0;
if (w >= 1.0) w -= 1.0;
return w;
</pre>
</blockquote>
I agree. I think the random_wh06 module had precision problems with
the floating point values it was using for R. If I modify the
module to fit the 64bit implementation, it does match the 124 bit
maximum with the range [0...21267638781707063560975648195455661513)
for the uniform distribution (coming from "w" in the 64bit version,
without division, or "I" in the random_wh06 module), which is then
input for the modulus N to have the algorithm use only a subset of
the total range for a uniform distribution. The Erlang code
equivalent is:<br>
<pre><div style="background-color: transparent;" class="line" id="LC99"> <span class="nv">B1</span> <span class="o">=</span> <span class="p">(</span><span class="mi">11600</span> <span class="o">*</span> <span class="nv">A1</span><span class="p">)</span> <span class="ow">rem</span> <span class="mi">2147483579</span><span class="p">,</span></div><div style="background-color: transparent;" class="line" id="LC100"> <span class="nv">B2</span> <span class="o">=</span> <span class="p">(</span><span class="mi">47003</span> <span class="o">*</span> <span class="nv">A2</span><span class="p">)</span> <span class="ow">rem</span> <span class="mi">2147483543</span><span class="p">,</span></div><div style="background-color: transparent;" class="line" id="LC101"> <span class="nv">B3</span> <span class="o">=</span> <span class="p">(</span><span class="mi">23000</span> <span class="o">*</span> <span class="nv">A3</span><
span class="p">)</span> <span class="ow">rem</span> <span class="mi">2147483423</span><span class="p">,</span></div><div style="background-color: transparent;" class="line" id="LC102"> <span class="nv">B4</span> <span class="o">=</span> <span class="p">(</span><span class="mi">33000</span> <span class="o">*</span> <span class="nv">A4</span><span class="p">)</span> <span class="ow">rem</span> <span class="mi">2147483123</span><span class="p">,</span></div><div style="background-color: transparent;" class="line" id="LC103">
</div><div style="background-color: transparent;" class="line" id="LC104"> <span class="nb">put</span><span class="p">(</span><span class="n">random_wh06_seed</span><span class="p">,</span> <span class="p">{</span><span class="nv">B1</span><span class="p">,</span> <span class="nv">B2</span><span class="p">,</span> <span class="nv">B3</span><span class="p">,</span> <span class="nv">B4</span><span class="p">}),</span></div><div style="background-color: transparent;" class="line" id="LC105">
</div><div style="background-color: transparent;" class="line" id="LC106"> <span class="nv">I</span> <span class="o">=</span> <span class="p">((</span><span class="nv">B1</span> <span class="o">*</span> <span class="mi">9903516371291919229607132747</span><span class="p">)</span> <span class="o">+</span></div><div style="background-color: transparent;" class="line" id="LC107"> <span class="p">(</span><span class="nv">B2</span> <span class="o">*</span> <span class="mi">9903516537312557910938853791</span><span class="p">)</span> <span class="o">+</span></div><div style="background-color: transparent;" class="line" id="LC108"> <span class="p">(</span><span class="nv">B3</span> <span class="o">*</span> <span class="mi">9903517090714727049595319831</span><span class="p">)</span> <span class="o">+</span></div><div style="background-color: transparent;" cla
ss="line" id="LC109"> <span class="p">(</span><span class="nv">B4</span> <span class="o">*</span> <span class="mi">9903518474220420479167438931</span><span class="p">))</span></div><div style="background-color: transparent;" class="line" id="LC110"> <span class="ow">rem</span> <span class="mi">21267638781707063560975648195455661513</span><span class="p">,</span></div></pre>
I updated the implementation with this, to correct the problem:<br>
<a class="moz-txt-link-freetext" href="https://github.com/okeuday/erlbench/blob/master/src/random_wh06.erl">https://github.com/okeuday/erlbench/blob/master/src/random_wh06.erl</a><br>
<br>
Thank you for finding the 64bit C implementation. Please tell me if
you think my approach with the recent changes is erroneous.<br>
<br>
Thanks,<br>
Michael<br>
</body>
</html>