Here is xoroshiro128.bas using Sebastian Vigna's xoroshiro128**. It should be noted that is was written with David Blackman. As a 'by the way' O'Neill and Blackman get on quite well.
Benefits:
Passes PractRand to 1TB - will test beyond that as time permits.
Outputs Doubles with 53-bit granularity.
It is marginally faster in 32-bit mode than PCG32II which outputs with 32-bit granularity and much faster in 64-bit mode.
Has a sequence repeat, period, after 2^128 randm number requests.
Is thread safe.
If we did nothing but request random numbers then with my machine, running at 3.9GHz, the sequence would repeat after 2.833 trillion x The estimated current age of the universe (13.5 billion years). Even with a 64-bit period we get 659 x The estimated current age of the universe. PowerBASIC's Rnd has a period of 32-bits and will repeat after 54 minutes.
The state vector is two Ulongints. MyRandomize("") gives a random state vector and MyRandomize( <Val(String)> ) is used for fixed seeding; MyRandomize("123456789") for example. Don't worry about the quality of fixed seeding as the generator is warmed up using 8K of 64-bit values, in any event, which should prove more than adequate. On my machine warming up takes 27 microseconds in 32-bit mode and 7 microseconds in 64-bt mode.
The function rand returns a Ulongint. There is no randse, which my other generators use, because the output is 64-bit giving 53-bit granularity and called randD.
There are two range functions: One is discrete outputting Longs and the other is continuous outputting Doubles. The discrete function uses dodicat's world famous MOD method.
We have GetSnapshot and SetSnapshot with regard the state vector. The Constructor employs MyRandomize( "" ) followed GetSnapshot. We can then repeat a sequence from the outset, after the warm up, without calling GetSnapshot. These are very useful procedures and PCG32II has had them for a while. I never got around to including them with my other generators.
At the bottom of xoroshiro128.bas is: ... [1]
Code: Select all
Dim Shared As xoroshiro128 x128
#undef Rnd
#define Rnd x128.randD
so we will use x128 in our code.
On a Copy&Paste from the forum any use of Rnd will get replaced with x128.randD when we include xoroshiro128.bas. For example, UEZ's Fireworks has had ' #include "E:\FreeBASIC\xoroshiro128\xoroshiro128.bas" ' put at the head - job done. The Fireworks experience has not been enhanced, I don't think, but I wanted to test the above idea and, of course, I did it because I could. Having said that because Fireworks uses dodicat's Regulate procedure then Fireworks may have a little more idle time. We can, of course, develop using Rnd and that will also get replaced with x128.randD when we include xoroshiro128.bas. A random range is often used in applications and Rnd's use in them will also get replaced but it is better to use one of the range functions in this package; especially the discrete version.
[1] Note that the above tip assumes a primary thread only - you will need to use another 'Dim Shared As xoroshiro128 ...' otherwise you will not have thread safety; the x128's will collide.
Following xoroshiro128.bas is a usage example and typical output.
xoroshiro128.bas
Code: Select all
' http://prng.di.unimi.it/xoroshiro128starstar.c
/'
static inline uint64_t rotl(const uint64_t x, int k) {
return (x << k) | (x >> (64 - k));
}
static uint64_t s[2];
uint64_t next(void) {
const uint64_t s0 = s[0];
uint64_t s1 = s[1];
const uint64_t result = rotl(s0 * 5, 7) * 9;
s1 ^= s0;
s[0] = rotl(s0, 24) ^ s1 ^ (s1 << 16); // a, b
s[1] = rotl(s1, 37); // c
return result;
}
'/
#define Make64Bit ( Cast( Ulongint, Rnd*(2^32) ) Shl 32) Or Cast( Ulongint, Rnd*(2^32) )
#define rotl(x,k) ( (x Shl k) Or (x Shr(64-k)) ) ' Note the extra parentheses
#macro Engine
s0 = This.s(0)
s1 = This.s(1)
result = rotl(s0*5, 7)*9
s1 Xor= s0
This.s(0) = rotl(s0,24) Xor s1 Xor (s1 shl 16)
This.s(1) = rotl(s1,37)
#endmacro
' Generator is Sebastian Vigna's xoroshiro128**
Type xoroshiro128
Public:
Declare Constructor
Declare Sub MyRandomize( seed as string )
Declare Function rand() As Ulongint
Declare Function randD() As Double
Declare Function range overload( Byval One As Long, Byval Two As Long ) As Long
Declare Function range overload ( byval One as double, Byval Two as double ) as double
Declare Sub GetSnapshot()
Declare Sub SetSnapshot()
Private:
As Ulongint s(0 To 1)
As Ulongint snapshot(0 to 1)
End Type
Private Function Get64Bit() As Ulongint
Dim As Ulongint dummy = Make64Bit
Return rotl(dummy, 16)
End Function
Private Sub xoroshiro128.MyRandomize( seed as string )
Dim As Ulong i
' If seed = "" then use cryptographic random numbers else use Mersenne Twister
If seed = "" Then Randomize , 5 Else Randomize Val(seed), 3
This.s(0) = Get64Bit
This.s(1) = Get64Bit
' Warm up generator - essential
For i = 1 To 8192
this.rand()
Next i
End Sub
Private Function xoroshiro128.rand() As ulongint
Dim As Ulongint s0, s1, result
Engine
Return result
End Function
Private Function xoroshiro128.randD() As Double
Dim As ulongint s0, s1, result
Engine
Return result/2^64
End Function
Private Function xoroshiro128.range( Byval One As Long, Byval Two As Long ) As Long
Dim As ulongint s0, s1, result
Engine
' Edit: Changed 'result' to Culng(result) - MOD is slow when working with 64-bit.
Return Clng(Culng(result) Mod ( Two-One+1 )) + One ' By dodicat
End Function
Private Function xoroshiro128.range( Byval One As Double, Byval Two As Double ) As Double
Dim As Ulongint s0, s1, result
Engine
Return result/2^64*( Two - One ) + One
End Function
Private Sub xoroshiro128.GetSnapshot()
This.snapshot(0) = This.s(0)
This.snapshot(1) = This.s(1)
End Sub
Private Sub xoroshiro128.SetSnapshot()
This.s(0) = This.snapshot(0)
This.s(1) = This.snapshot(1)
End Sub
Constructor xoroshiro128
This.MyRandomize( "" )
This.GetSnapshot
End Constructor
Dim Shared As xoroshiro128 x128
#undef Rnd
#define Rnd x128.randD
Code: Select all
'#Console On
#Include "xoroshiro128.bas"
Dim As ulong i
Dim As Double tot
Print "Random floats [0,1) (53-bit granularity)"
Print
For i = 1 To 6
Print Rnd
Next
Print
x128.SetSnapshot
Print "SetSnapshot executed - no GetSnapshot executed by user"
Print
For i = 1 To 6
Print Rnd
Next
Print
For i = 1 To 10^8
tot += Rnd
Next
Print "Average [0,1)";
Print tot/10^8
Print
Print "Random Ulongints"
Print
For i = 1 To 6
Print x128.rand
Next
Print
Print "Discrete range, 5 to 20"
Print
For i = 1 To 6
Print x128.range(5,20)
Next
Print
Print "Continuous range, 1. to 10."
Print
For i = 1 To 6
Print x128.range(1.,10.)
Next
Print
Print "Randomize - random seed - GetSnapshot executed"
Print
x128.MyRandomize( "" )
x128.GetSnapshot
For i = 1 To 6
Print Rnd
Next
Print
x128.SetSnapshot
Print "SetSnapshot executed"
Print
For i = 1 To 6
Print Rnd
Next
Print
Print "Randomize - fixed seed - GetSnapshot executed"
Print
x128.MyRandomize( "123456" )
x128.GetSnapshot
For i = 1 To 6
Print Rnd
Next
Print
x128.SetSnapshot
Print "SetSnapshot executed"
Print
For i = 1 To 6
Print Rnd
Next
Print
Print " Done"
Sleep
Code: Select all
Random floats [0,1) (53-bit granularity)
0.4049058391992926
0.002008626534225989
0.456443287069138
0.04185532376693672
0.2414107387053619
0.7075873397474699
SetSnapshot executed - no GetSnapshot executed by user
0.4049058391992926
0.002008626534225989
0.456443287069138
0.04185532376693672
0.2414107387053619
0.7075873397474699
Average [0,1) 0.5000014459389598
Random Ulongints
1601524761088817538
12336677007064841681
14695406431682103560
12754510237835497810
4165814302991093027
17288837923636644936
Discrete range, 5 to 20
9
14
19
10
11
20
Continuous range, 1. to 10.
8.902743816548306
5.189807367551428
4.003102644336549
9.649630648651838
7.006182164074057
6.712164673793058
Randomize - random seed - GetSnapshot executed
0.1005232132552173
0.8435787285616101
0.214314369648843
0.2620323823819769
0.7587725589312045
0.389518875631241
SetSnapshot executed
0.1005232132552173
0.8435787285616101
0.214314369648843
0.2620323823819769
0.7587725589312045
0.389518875631241
Randomize - fixed seed - GetSnapshot executed
0.2030120460316923
0.7836776707789495
0.06840766063610897
0.1860217678276363
0.8058282182269503
0.1017877019893167
SetSnapshot executed
0.2030120460316923
0.7836776707789495
0.06840766063610897
0.1860217678276363
0.8058282182269503
0.1017877019893167
Done
