Range only PRNG

General FreeBASIC programming questions.
jj2007
Posts: 2326
Joined: Oct 23, 2016 15:28
Location: Roma, Italia
Contact:

Re: Range only PRNG

Post by jj2007 »

deltarho[1859] wrote:I get something like this '184467440737.0945'. This only happens when First < 0. I am certain that I have seen this before
Looks like the unsigned representation of a -1 QWORD, divided by 100,000,000. The exact value is 184467440737.09551615
deltarho[1859]
Posts: 4308
Joined: Jan 02, 2017 0:34
Location: UK
Contact:

Re: Range only PRNG

Post by deltarho[1859] »

Oh, boy Jochen what a catch. Image

I was overflowing Ulongint, not something we ever expect to happen. Sometimes the random numbers generated were safe, and sometimes they were not. This is the result of using ' k +=' as opposed to 'k ='.

It was easily fixed by using k as a Double and not a Ulongint.

We can now change SFC32.RangeN back to:

Code: Select all

Function SFC32.RangeN ( ByVal First As Long, ByVal Last As long ) As Long
Dim As Ulong dw_tmp
Dim As Long Adj
This.s(3) += 1
dw_tmp = This.s(0) + This.s(1) + This.s(3)
This.s(0) = This.s(1) Xor ( This.s(1) shr 9 )
This.s(1) = This.s(2) + ( This.s(2) shl 3 )
This.s(2) = rotl(This.s(2), 21) + dw_tmp
Return Clng( This.s(2) Mod ( Last - First + 1) ) + First ' By dodicat   
End Function
I think that SFC32.bi is now finished - it is fairly comprehensive.

It was good fun writing it but if anyone is using PCG32II then, for most, I don't think switching to SFC32 will make much difference to their applications. If anyone finds the decided throughput increase beneficial, then let me know - I'd be interested in what they are doing.

Thanks also for srvaldez's help during this project.
deltarho[1859]
Posts: 4308
Joined: Jan 02, 2017 0:34
Location: UK
Contact:

Re: Range only PRNG

Post by deltarho[1859] »

I have decided to change how the generator is initialized. As is the state vector is populated via Windows BCryptGenRandom and then we employ a warm-up.

I don't know if SFC332 needs a warm-up and, if so, for how long.

Another way, avoiding a warm-up, is to seed using a balanced Hamming weight mentioned in the Seeding FB generators thread.

Now, whilst I was at it, I did a bit of tidying up. There are 88 instances of 'This'. I have a habit of using that whether needed or not. Anyway we don't need them here, so I ripped them out. The generator engine could be put into a function, but that would slow things down, so I used a macro.

The new way of seeding uses the asm instruction popcnt which was introduced with Intel's Nehalem in 2008 and AMD's Barcelona in 2007. If you don't have that D.J.Peters wrote a nice piece of code in the above link.

The Sub SFC32.Initialize() is small now but, more to the point, SFC32.bi is now Windows independent - no more BCryptGenRandom. The code, in general, seems to me to be a little more readable as well.

SFC32.bi

Code: Select all

#Include Once "windows.bi"
#Inclib "bcrypt"
#Include Once "win/wincrypt.bi"
 
#define rotl(x,k) ( (x Shl k) Or (x Shr(32-k)) ) ' Note the extra parentheses
 
#macro Engine
  s(3) += 1
  dw_tmp = s(0) + s(1) + s(3)
  s(0) = s(1) Xor ( s(1) shr 9 )
  s(1) = s(2) + ( s(2) shl 3 )
  s(2) = rotl(s(2), 21) + dw_tmp
#endmacro
 
TYPE SFC32
  Declare Constructor
  Declare Sub Initialize()
  Declare Function Rand() As Ulong
  Declare Function Rnd() As Double ' 32-bit granularity
  Declare Function RndD() As Double ' 53-bit granularity
  Declare Function RangeN( ByVal One As Long, ByVal Two As Long ) As Long
  Declare Function RangeF( ByVal One as Long, ByVal Two as Long ) As Double
  Declare Sub GetSnapshot()
  Declare Sub SetSnapShot()
  As Ulong s(0 To 3)
  As Ulong snaps(0 to 3)
END TYPE
 
Function SFC32.Rand( ) As Ulong
Dim As Ulong dw_tmp
Engine
Return s(2)
End Function
 
Function SFC32.Rnd( ) As Double
Dim As Ulong dw_tmp
Engine
Return s(2)/2^32
End Function
 
Function SFC32.RndD( ) As Double
Dim As Ulong dw_tmp, TempVar1, TempVar2
Engine
TempVar1 = s(2)
Engine
TempVar2 = s(2)
Asm
  mov eax, dword ptr [TempVar1]
  movd xmm0, eax
  mov eax, dword ptr [TempVar2]
  movd xmm1, eax
  punpckldq xmm0, xmm1
  psrlq xmm0, 12
  mov eax, 1
  cvtsi2sd xmm1, eax
  por xmm0, xmm1
  subsd xmm0, xmm1
  movq [Function], xmm0
End Asm
End Function
 
Function SFC32.RangeN ( ByVal First As Long, ByVal Last As long ) As Long
Dim As Ulong dw_tmp
Engine
Return Clng( s(2) Mod ( Last - First + 1) ) + First ' By dodicat
End Function
 
Function SFC32.RangeF ( ByVal First As Long, ByVal Last As long ) As Double
Dim As Ulong dw_tmp
Engine
Return s(2)/2^32*( Last - First ) + First
End Function
 
Sub SFC32.GetSnapshot ( )
  For i As Ulong = 0 to 3
    snaps(i) = s(i)
  Next
End Sub
 
Sub SFC32.SetSnapshot ( )
  For i As Ulong = 0 To 3
    s(i)= snaps(i)
  Next
End Sub
 
Function HammingSeed() As Ulong
Dim As Ulong Seed, CopySeed, numBits = 0
While numBits <> 16
  Asm
    rdtsc    ' Get a fast changing number from the processor
    mov Dword Ptr [Seed], eax
  End Asm
  CopySeed = Seed : numBits = 0
  While CopySeed <> 0
    CopySeed = CopySeed And CopySeed - 1
    numBits = numBits + 1
  Wend
Wend
Asm
  mov eax, Dword Ptr [Seed]
  bswap eax    ' Reverse byte order to break sequence
  mov [Function], eax
End Asm
End Function
 
Sub SFC32.Initialize( )
s(0) = HammingSeed
s(1) = HammingSeed
s(2) = HammingSeed
s(3) = 1
GetSnapshot( )
End Sub
 
Constructor SFC32
  Initialize
  GetSnapshot
End Constructor
deltarho[1859]
Posts: 4308
Joined: Jan 02, 2017 0:34
Location: UK
Contact:

Re: Range only PRNG

Post by deltarho[1859] »

I have just replaced the Function HammingSeed() by another which does not rely on your CPU having the popcnt asm instruction by incorporating code by adeyblue. It takes a little longer to come up with a 50% Hamming weight seed, but I very much doubt that you will perceive a wait when initializing the generator. Actually, you won't - I will put money on it.
Post Reply