OpenEuphoria

Re: Algorithm needed for random number filtering

• Posted by petelomax Dec 04, 2016
• 1235 views

Here is a very simple demonstration of Von Neuman's technique.
To keep things simple I am using 4 bit numbers, it should not be difficult to extend to any (even) power of 2.
First, I generate 2,000,000 random numbers and record the initial distribution (idist).
I average each consecutive pair to create a bell curve distribution (bdist).
To remove bias, sample each (non overlapping) pair of bits: make 01 -> 1, 10 -> 0, and ignore 00 and 11. (rdist)

sequence idist = repeat(0,16), 
         bdist = repeat(0,16), 
         rdist = repeat(0,16) 
 
integer bits = 0 
integer r1, r2, r3 = 0 
for i=1 to 1000000 do 
    -- initial distribution 
    r1 = rand(16) 
    idist[r1] += 1 
    r2 = rand(16) 
    idist[r2] += 1 
    -- make a bell curve distribution 
    r1 = round((r1+r2)/2) 
    bdist[r1] += 1 
    -- remove bias 
    for j=1 to 1 do 
        r2 = and_bits(r1,3) 
        if r2=1 
        or r2=2 then 
            r3 = r3*2+and_bits(r2,1) 
            bits += 1 
            if bits=4 then 
                rdist[r3+1] += 1 
                r3 = 0 
                bits = 0 
            end if 
        end if 
        -- (must use 0..15 not 1..16 here, hence the -1/+1) 
        r1 = floor((r1-1)/4)+1 
    end for 
end for 
  
procedure showdist(sequence d) 
integer di, maxdi = 0, tdi = 0 
    for i=1 to length(d) do 
        di = d[i] 
        if di>maxdi then 
            maxdi = di 
        end if 
        tdi += di 
    end for 
    printf(1,"(%d samples)\n",tdi) 
    for i=1 to length(d) do 
        printf(1,"%s\n",{repeat('*',round((d[i]/maxdi)*80))}) 
    end for 
end procedure 
 
puts(1,"initial\n") 
showdist(idist) 
puts(1,"biased\n") 
showdist(bdist) 
puts(1,"recovered\n") 
showdist(rdist) 

output:

initial 
(2000000 samples) 
******************************************************************************** 
******************************************************************************* 
******************************************************************************** 
******************************************************************************** 
******************************************************************************** 
******************************************************************************** 
******************************************************************************** 
******************************************************************************* 
******************************************************************************** 
******************************************************************************** 
******************************************************************************* 
******************************************************************************* 
******************************************************************************** 
******************************************************************************** 
******************************************************************************* 
******************************************************************************** 
biased 
(1000000 samples) 
*** 
************* 
*********************** 
********************************** 
******************************************** 
******************************************************* 
***************************************************************** 
*************************************************************************** 
******************************************************************************** 
********************************************************************** 
************************************************************ 
************************************************* 
*************************************** 
**************************** 
****************** 
******** 
recovered 
(125088 samples) 
******************************************************************************** 
****************************************************************************** 
******************************************************************************** 
****************************************************************************** 
******************************************************************************* 
***************************************************************************** 
******************************************************************************* 
******************************************************************************** 
******************************************************************************* 
******************************************************************************* 
****************************************************************************** 
****************************************************************************** 
******************************************************************************* 
******************************************************************************* 
****************************************************************************** 
****************************************************************************