Firstly, credit to cbruce for bringing this PRNG to our attention and, secondly, credit to counting_pine for porting the minimal C implementation to FreeBASIC.

The wrapper overhead was expensive so I have used the age old method of code repetition by using the same pcg code in the head of each of the wrappers.

If we employ no seeding the first two random 32 bit values will be zero. All random number generators should be 'warmed up'. In An experimental exploration of Marsaglia’s xorshift generators,scrambled by Sebastiano Vigna he produces a graph showing the number of samples which are required to warm up some PRNGs, including his own. It is quite clear that the the larger the period of a generator the greater the number of samples should be executed before using a generator. Not surprisingly the Mersenne Twister needs a lot of warming up. Since PCG shares some similarities with Vigna's PRNGs and we have a period of 2^64 with this minimal PCG I have opted for a warm up of 200 samples. This should prove more than adequate and only takes one micro-second to complete, so you will not have time to put the kettle on. <smile>

The warm up is done in the function RandomizePCG32( Seed )

Seed is optional. If no seed is given a cryptographic random Ulongint is used from FB's PRNG option 5. To repeat a sequence a seed is required. Quite frankly, anything will do: 123456, say, is fine.

For single precision we simply execute PCG32S. For range, 0 to 100 for example, we simply execute PCG32R(0,100).

The code at the end of this post has the code to add to your apps and an example usage.

I now use these compiler option switches at all times: -gen gcc -Wc -O1

With no compiler optimization the throughput collapses.

This is a typical output using RandomizePCG32 - that is an unknown seed.

Code: Select all

` 0.3077029`

0.5957927

0.65415

0.7683138

0.6042935

0.1181773

0.07335913

0.6330329

0.7821538

0.4110711

116

204

86

7

211

37

128

147

240

157

Averge of 10^8 singles: 0.5000328068581712

Average of 10^8 range (0,100): 49.99446299

Estimate of 'raw' throughput (Millions per second)

Singles: 219.113

Range: 180.463

That 219.113 is twice as fast as FB's PRNG option 2.

In 64 bit mode I get 222.410 and 206.462 for singles and range respectively.

A period of 2^64 may be regarded as small but just to put that into perspective it would take over 3,000 years on my machine before a sequence repetition and that is just executing random numbers and doing nothing with them. Also, if we used 2^32 random numbers then we would be using 2^32/2^64 = 1/2^32 of those available, that is 23*10^-9% so the chance of 'landing' on numbers previously used is very small. Of course,the chance decreases as the number used decreases.

srvaldez has published the full implementation of PCG in [1] but I have not had a chance to look at it yet.

Some folk may find the following code adequate for their purposes.

Code: Select all

`'' *Really* minimal PCG32 code / (c) 2014 M.E. O'Neill / pcg-random.org`

'' Licensed under Apache License 2.0 (NO WARRANTY, etc. see website)

'' FreeBASIC translation by Matthew Fearnley 2017-06-04

' Include the code between ' ******/' ****** in your apps

' **********************************************

Type pcg32_random_t

As Ulongint state, inc

End Type

Dim Shared pcg32 As pcg32_random_t

Function pcg32_random_r(Byval rng As pcg32_random_t Ptr) As Ulong

Dim As Ulongint oldstate = rng->state

'' Advance internal state

rng->state = oldstate * 6364136223846793005ULL + (rng->inc Or 1)

'' Calculate output function (XSH RR), uses old state for max ILP

Dim As Ulong xorshifted = ((oldstate Shr 18u) xor oldstate) Shr 27u

Dim As Ulong rot = oldstate Shr 59u

Return (xorshifted Shr rot) Or (xorshifted Shl ((-rot) And 31))

End Function

Sub RandomizePCG32( Seed As Ulongint = 0 )

Dim i As Integer

If seed = 0 Then

Randomize , 5

pcg32.state = Rnd*(2^64-1)

Else

pcg32.state = seed

End If

For i As Integer = 1 To 200

pcg32_random_r(@pcg32)

Next i

End Sub

Function PCG32S() As Single

Dim tempvar As Ulong

Dim rng As pcg32_random_t Ptr = @pcg32

Dim As Ulongint oldstate = rng->state

'' Advance internal state

rng->state = oldstate * 6364136223846793005ULL + (rng->inc Or 1)

'' Calculate output function (XSH RR), uses old state for max ILP

Dim As Ulong xorshifted = ((oldstate Shr 18u) xor oldstate) Shr 27u

Dim As Ulong rot = oldstate Shr 59u

Tempvar = (xorshifted Shr rot) Or (xorshifted Shl ((-rot) And 31))

Asm

mov eax, dword Ptr [TempVar]

movd xmm0, eax

psrlq xmm0, 9

mov eax, 1

cvtsi2ss xmm1, eax

por xmm0, xmm1

subss xmm0, xmm1

movd [Function], xmm0

End Asm

End Function

Function PCG32R( Byval One As Long, Byval Two As Long ) As Long

Dim tempvar As Ulong

Dim rng As pcg32_random_t Ptr = @pcg32

Dim As Ulongint oldstate = rng->state

'' Advance internal state

rng->state = oldstate * 6364136223846793005ULL + (rng->inc Or 1)

'' Calculate output function (XSH RR), uses old state for max ILP

Dim As Ulong xorshifted = ((oldstate Shr 18u) xor oldstate) Shr 27u

Dim As Ulong rot = oldstate Shr 59u

Tempvar = (xorshifted Shr rot) Or (xorshifted Shl ((-rot) And 31))

Asm

mov edx, dword Ptr [TempVar]

mov ecx, [One]

mov eax, [Two]

cmp ecx, eax

jl Now1LT2

xchg eax, ecx

Now1LT2:

Sub eax, ecx

inc eax

jz doTheRnd

mul edx

Add edx, ecx

doTheRnd:

mov [Function], edx

End Asm

End Function

' *****************************************************

RandomizePCG32

Dim As Ulong I, N

' Knock out 10 singles

For i = 1 To 10

Print PCG32S

Next

Print

' Knock out 10 range

For i = 1 To 10

Print PCG32R(0,255)

Next

Print

Dim As Double tot

' Knock out average of 10^8 singles

For i = 1 To 10^8

tot += PCG32S

Next

tot = tot/10^8

Print "Averge of 10^8 singles: ";tot

Print

tot = 0

' Knock out average of 10^8 range (0,100)

For i = 1 To 10^8

tot += PCG32R(0,100)

Next

tot = tot/10^8

Print "Average of 10^8 range (0,100): ";tot

Print

Dim As Double t

Print "Estimate of 'raw' throughput (Millions per second)"

Print

t = Timer

For i = 1 To 10^8

PCG32S

Next

t = Timer - t

Print "Singles: ";Using "###.###";10^8/(t*1000000)

t = Timer

For i = 1 To 10^8

PCG32R(0,100)

Next

t = Timer - t

Print "Range: ";Using "###.###";10^8/(t*1000000)

Sleep