RND breadth

General FreeBASIC programming questions.
albert
Posts: 6000
Joined: Sep 28, 2006 2:41
Location: California, USA

Re: RND breadth

Post by albert »

The problem is in the RND Function... it only generates upper values... > .0010
Here's a test to see: it only generates out of 1,000,000 loops , it only generates 1,000 below .0010

Code: Select all


screen 19

dim as double l0 = .00000 , h0 = .0001
dim as double l1 = .00011 , h1 = .0002
dim as double l2 = .00021 , h2 = .0003
dim as double l3 = .00031 , h3 = .0004
dim as double l4 = .00041 , h4 = .0005
dim as double l5 = .00051 , h5 = .0006
dim as double l6 = .00061 , h6 = .0007
dim as double l7 = .00071 , h7 = .0008
dim as double l8 = .00081 , h8 = .0009
dim as double l9 = .00091 , h9 = .0010

dim as integer array(0 to 9)

dim as double n = 0
for a as longint = 1 to 1000000
    
    n = rnd
    
    select case n
    case l0 to h0 : array(0)+=1
    case l1 to h1 : array(1)+=1
    case l2 to h2 : array(2)+=1
    case l3 to h3 : array(3)+=1
    case l4 to h4 : array(4)+=1
    case l5 to h5 : array(5)+=1
    case l6 to h6 : array(6)+=1
    case l7 to h7 : array(7)+=1
    case l8 to h8 : array(8)+=1
    case l9 to h9 : array(9)+=1
    end select
    
next

dim as longint total = 0
for a as longint = lbound(array) to ubound(array) step 1
    print array(a)
    total+=array(a)
next
print "total = " ; total

sleep
end

dodicat
Posts: 7983
Joined: Jan 10, 2006 20:30
Location: Scotland

Re: RND breadth

Post by dodicat »

.001 is one thousandth of the way from 0 to 1.
So each leg you expect 100, and you have 10 legs up to .001.(as you see printed)
You expect about 1000 below .001 because you have 1000 of these intervals between 0 and one, 1000*1000 = 1,000,000
Here is the independent random generator and I have closed up your gaps so the intervals butt on to each other.

Code: Select all

 
 
'screen 19
' ulongint generator

type Rand
   a as ulongint
   b as ulongint
   c as ulongint
   d as ulongint
end type

function ranulong(x as rand) as ulongint
   dim e as ulongint = x.a - ((x.b shl 7) or (x.b shr (64 - 7)))
   x.a = x.b xor ((x.c shl 13) or (x.c shr (64 - 13)))
   x.b = x.c + ((x.d shl 37) or (x.d shr (64 - 37)))
   x.c = x.d + e
   x.d = e + x.a
   return x.d
end function

function randouble(x as rand) as double 
    return ranulong(x)/18446744073709551615ull
    end function

 sub init(x as rand, byval seed as ulongint=100)
   dim i as ulongint
    x=type(4058668781,seed,seed,seed)
    for i as ulongint=1 to 20
        ranulong(x)
        next
end sub

'=========
dim shared as rand z
init(z)
'========
randomize ,0

screen 19

dim as double l0 = .00000 , h0 = .0001
dim as double l1 = .0001 , h1 = .0002
dim as double l2 = .0002 , h2 = .0003
dim as double l3 = .0003 , h3 = .0004
dim as double l4 = .0004 , h4 = .0005
dim as double l5 = .0005 , h5 = .0006
dim as double l6 = .0006 , h6 = .0007
dim as double l7 = .0007 , h7 = .0008
dim as double l8 = .0008 , h8 = .0009
dim as double l9 = .0009 , h9 = .0010

dim as integer array(0 to 9)

dim as double n = 0
for a as longint = 1 to 1000000
    
    n = randouble(z)
    
    select case n
    case l0 to h0 : array(0)+=1
    case l1 to h1 : array(1)+=1
    case l2 to h2 : array(2)+=1
    case l3 to h3 : array(3)+=1
    case l4 to h4 : array(4)+=1
    case l5 to h5 : array(5)+=1
    case l6 to h6 : array(6)+=1
    case l7 to h7 : array(7)+=1
    case l8 to h8 : array(8)+=1
    case l9 to h9 : array(9)+=1
    end select
    
next

dim as longint total = 0
for a as longint = lbound(array) to ubound(array) step 1
    print array(a)
    total+=array(a)
next
print "total = " ; total

sleep
end
 
albert
Posts: 6000
Joined: Sep 28, 2006 2:41
Location: California, USA

Re: RND breadth

Post by albert »

The RND stuff is for my program Vari_Cph...

I needed to generate a 64 bit random binary string...
But it only generated upper values ,mostly "1's"

Right( string( 64 , "0" ) + bin( int( rnd * (2^64) ) ) , 64 ) , was generating mostly 1's..

I got it returning a full plethora of bits...With the following code..

Code: Select all


screen 19

do
    
    dim as string*64 str1 = right(string(64,"1"),64) ' + bin(int(rnd*2^64)),64)
    
    for a as longint = 1 to 48
        mid(str1,int(rnd*len(str1))+1,1) = "0"
    next
    
    print str1 , valulng("&B"+str1)

    sleep

loop until inkey = chr(27)

end
   
sancho3
Posts: 358
Joined: Sep 30, 2017 3:22

Re: RND breadth

Post by sancho3 »

Like Dodicat said, you have isolated a small range of numbers. If you isolate a same size range anywhere else you get similar amounts of hits.
Here I add the next range .00100 to .00200 to your code and we see similar amounts of hits:

Code: Select all

screen 19

dim as double l0 = .00000 , h0 = .0001
dim as double l1 = .00011 , h1 = .0002
dim as double l2 = .00021 , h2 = .0003
dim as double l3 = .00031 , h3 = .0004
dim as double l4 = .00041 , h4 = .0005
dim as double l5 = .00051 , h5 = .0006
dim as double l6 = .00061 , h6 = .0007
dim as double l7 = .00071 , h7 = .0008
dim as double l8 = .00081 , h8 = .0009
dim as double l9 = .00091 , h9 = .0010

dim As Double ll = .00100, hh = .0020

dim as integer array(0 to 10)

dim as double n = 0
for a as longint = 1 to 1000000
   
    n = rnd
   
    select case n
    case l0 to h0 : array(0)+=1
    case l1 to h1 : array(1)+=1
    case l2 to h2 : array(2)+=1
    case l3 to h3 : array(3)+=1
    case l4 to h4 : array(4)+=1
    case l5 to h5 : array(5)+=1
    case l6 to h6 : array(6)+=1
    case l7 to h7 : array(7)+=1
    case l8 to h8 : array(8)+=1
    case l9 to h9 : array(9)+=1
	 Case ll To hh : array(10)+= 1
    end select
   
next

dim as longint total = 0
for a as longint = lbound(array) to ubound(array) step 1
    print array(a)
    total+=array(a)
next
print "total = " ; total
Print "next range: "; array(10)
sleep
end

deltarho[1859]
Posts: 4308
Joined: Jan 02, 2017 0:34
Location: UK
Contact:

Re: RND breadth

Post by deltarho[1859] »

Oh dear!

dodicat and sancho3 beat me to it.

Lets look at the first two lines of the Select construct.

Code: Select all

Case l0 To h0 => 0 To .0001
Case l1 To h1 => .00011 To .0002
What happens to the values > .0001 and < .00011 ie a span of .00001. Simple - they don't get counted.

The Mersenne Twister passes a PractRand test to 256GB. If your 'problem', albert, really did exist PractRand would fail long before 256GB.

I am definitely out of this thread now. <smile>
jj2007
Posts: 2326
Joined: Oct 23, 2016 15:28
Location: Roma, Italia
Contact:

Re: RND breadth

Post by jj2007 »

albert wrote:But it only generated upper values ,mostly "1's"
If you use a 2^64 range, i.e. 0 ... 18446744073709551615, then the highest number in the first 0.1% of that range is 18446744073709551.615 (but there are probably a dozen different ways to explain the phenomenon - unfortunately it is not very intuitive).
dodicat
Posts: 7983
Joined: Jan 10, 2006 20:30
Location: Scotland

Re: RND breadth

Post by dodicat »

Albert.
I would create the zeros and ones this way:

Code: Select all

randomize
screen 19

#define range(f,l) Int(Rnd*(((l)+1)-(f)))+(f)

do
    dim as integer countzero,countone

    dim as string str1 = string(64,0)
    
    for a as longint = 0 to len(str1)-1
        var z=range(48,49)
        if z=48 then countzero+=1
        if z=49 then countone+=1
       
        str1[a]=z
    next
    
    print str1 , valulng("&B"+str1)


print "zeros   ";countzero
print "ones    ";countone
 if countone+countzero<>64 then print "Error!"
    sleep

loop until inkey = chr(27)

end
     
counting_pine
Site Admin
Posts: 6323
Joined: Jul 05, 2005 17:32
Location: Manchester, Lancs

Re: RND breadth

Post by counting_pine »

albert wrote:The RND stuff is for my program Vari_Cph...

I needed to generate a 64 bit random binary string...
But it only generated upper values ,mostly "1's"
I guess the problem you're having is that Rnd can't generate 64 bits of randomness in each number.
It should be able to manage up to 53 bits precision because it's a Double, but in practice, it seems to only return up to 32 bits of randomness.
Either way, a full 64 bits isn't possible.

The easiest way is probably just to generate two 32-bit numbers, and combine them:

Code: Select all

function rnd32() as ulong
    return int(rnd * 2^32)
end function

function rnd64() as ulongint
    return culngint(rnd32()) shl 32 + rnd32()
end function

randomize
print bin(rnd64, 64)
deltarho[1859]
Posts: 4308
Joined: Jan 02, 2017 0:34
Location: UK
Contact:

Re: RND breadth

Post by deltarho[1859] »

I have popped back to see what counting_pine posted. <smile>
it seems to only return up to 32 bits of randomness.
That is because the underlying FB generators are only 32 bit.

Since we are only requesting a couple of random numbers speed is not an issue so I would recommend 'Randomize , 5' before calling rnd32(). If our application uses random numbers for other purposes then we can always invoke another Randomize using one of the other PRNGs after calling rnd64().

I actually use something very similar to counting_pine's code.

Code: Select all

#define Get64Bit (Cast( Ulongint, Rnd*(2^32) ) Shl 32) Or Cast( Ulongint, Rnd*(2^32) )

Code: Select all

Randomize , 5
Dim As Ulongint x = Get64Bit
Print Bin(x,64)
Sleep
speedfixer
Posts: 606
Joined: Nov 28, 2012 1:27
Location: CA, USA moving to WA, USA
Contact:

Re: RND breadth

Post by speedfixer »

Just stumbled across this.
I know this is old and stale, but I hate to think that Albert never got comfortable with anyone's reason why he did not see what he expected.

Albert:

2^32 is not half of 2^64.

2^63 is half

2^62 is 1/4
2^61 is 1/8
...
...
...
2^32 is VERY small next to 2^64.

If we allow a perfect rnd generator, anything less than 2^32 simply will not be chosen frequently in the 2^64 range.

I hope this helps someone.

david
deltarho[1859]
Posts: 4308
Joined: Jan 02, 2017 0:34
Location: UK
Contact:

Re: RND breadth

Post by deltarho[1859] »

2^k * 2^n = 2^(k+n)

So, 2^64 * 1/2 = 2^64 * 2^(-1) = 2^63
MrSwiss
Posts: 3910
Joined: Jun 02, 2013 9:27
Location: Switzerland

Re: RND breadth

Post by MrSwiss »

Instead of a lot of words, some code to visualize things:

Code: Select all

' ULI_vs_UL_visualized.bas -- (c) 200-07-20, MrSwiss

Width 120, 25
Dim As ULongInt ULI = &hFFFFFFFFFFFFFFFF
Dim As ulong    UL  = &hFFFFFFFF

Print "64 bits vs. 32 bits visualized"
Print "~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~"
Print

Print "64 bits max:  "; Tab(15); Bin(ULI, 64), ULI
ULI = BitReset(ULI, 63)
Print "halve value:  "; Tab(15); Bin(ULI, 64), ULI
ULI = BitReset(ULI, 62)
Print " 4th value:   "; Tab(15); Bin(ULI, 64), ULI
ULI = BitReset(ULI, 61)
Print " 8th value:   "; Tab(15); Bin(ULI, 64), ULI
ULI = BitReset(ULI, 60)
Print "16th value:   "; Tab(15); Bin(ULI, 64), ULI
ULI = BitReset(ULI, 59)
Print "32th value:   "; Tab(15); Bin(ULI, 64), ULI
ULI = BitReset(ULI, 58)
Print "64th value:   "; Tab(15); Bin(ULI, 64), ULI
Print
Print "32 bits max:  "; Tab(47); Bin(UL, 32), UL

Sleep
Result:

Code: Select all

64 bits vs. 32 bits visualized
~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

64 bits max:  1111111111111111111111111111111111111111111111111111111111111111      18446744073709551615
halve value:  0111111111111111111111111111111111111111111111111111111111111111      9223372036854775807
 4th value:   0011111111111111111111111111111111111111111111111111111111111111      4611686018427387903
 8th value:   0001111111111111111111111111111111111111111111111111111111111111      2305843009213693951
16th value:   0000111111111111111111111111111111111111111111111111111111111111      1152921504606846975
32th value:   0000011111111111111111111111111111111111111111111111111111111111      576460752303423487
64th value:   0000001111111111111111111111111111111111111111111111111111111111      288230376151711743

32 bits max:                                  11111111111111111111111111111111      4294967295
deltarho[1859]
Posts: 4308
Joined: Jan 02, 2017 0:34
Location: UK
Contact:

Re: RND breadth

Post by deltarho[1859] »

For the life of me I cannot see how that will help anyone.

2^k * 2^n = 2^(k+n), says it all.

To keep in the spirit of the last post, whatever that is, just replace all the BitReset() with 'ULI Shr 1'.
Post Reply