[Solved] Strange behaviour of Random() wrt Double

[MathMan]

Please consider the following small procedure under FPC 3.2.2

Code: Pascal  [Select][+][-]

1. procedure Test;
2.  
3. var
4.   Count1, Count2: NativeInt;
5.  
6.   Temp: UInt64;
7.  
8.   Mult: Double;
9.   Fp64: Double;
10.  
11. begin
12.   // I take random Double values across the complete binary exponent range
13.   Temp := UInt64( 1 ) shl 52;
14.   Mult := pDouble( @Temp )^;
15.  
16.   for Count1 := 1 to 2046 do begin
17.     for Count2 := 4095 downto 0 do begin
18.       Fp64 := Mult * ( 1.0+Random() ); // <= this is the point I stumble
19.     end;
20.  
21.     if( Count1<2046 ) then begin
22.       Mult := Double( 2.0 ) * Mult;
23.     end;
24.   end;
25. end;
26.  

When I look at Fp64 on the bit level I find that the low 16 bit of the mantissa are always 0 and not random as expected. Random() is declared Extended in the RTL, Fp64 and Mult are declared Double. I tried all variants of the marked line I could think of, explicitely converting to Double or shuffle the order, but the result stays the same.

Question: Did I miss something or is Random really not delivering back a 'fully' random Double?

Cheers,
MathMan

[Thaddy]

Question: Did I miss something or is Random really not delivering back a 'fully' random Double?

Depends a bit on the exact fpc version, because some higher "gods" decided to change the prng.
- code for 80 bit random never existed, double reduction and a shr to 32 bit
- mersenne twister (upto 3.2.4): 53 bit precision from 64 bit? or is it also a shr?
- xoshiro128**(from 3.3.1):64 bit precision from 128 bit

There are some higher precision randoms available for 3.2.2 using 128 bit to return full 64 bit double random values.
I have to check system.pas on what is different. (not algorithm but resolution)

I don't think it is the random itself, but how it (a) random is handled.

[MathMan]

Question: Did I miss something or is Random really not delivering back a 'fully' random Double?

Depends a bit on the exact fpc version, because some higher "gods" decided to change the prng.
- code for 80 bit random never existed, double reduction and a shr to 32 bit
- mersenne twister (upto 3.2.4): 53 bit precision from 64 bit? or is it also a shr?
- xoshiro128**(from 3.3.1):64 bit precision from 128 bit

There are some higher precision randoms available for 3.2.2 using 128 bit to return full 64 bit double random values.
I have to check system.pas on what is different. (not algorithm but resolution)

Thanks Thaddy - that already explains everything. On a closer look I detected, that actually the low 20 bit are 0 - which clarifies point 2 from your list. No need to further investigate - I'll implement my own Random for Double now. Will of course be based on RTL random but use bit-wizardry to have the full mantissa randomised.