RNGs
flip
there are many RNG's out there with huge periods.
For example:
1. Mersenne Twister (MT)
2. Add With Carry (AWC)
3. Subtract With Borrow (SWB)
4. KISS
There is also Marsaglia's DIEHARD battery of tests to test the quality of
these RNGs (in fact, 2 - 4 are his RNG's). Obviously, these RNGs are not
good for crypto applications, but are good for anything Monte Carlo.
Questions:
1. My question is this, which is the best overall RNG out there, where best
means that it passes all of the DIEHARD tests and has a huge period?
2. Is there a RNG which can be implemented in a FPGA to simulate "White
Noise"?
Thank you for any references (as I am looking through Google as we speak),
but cannot believe that many haven't been down this road.
Flip
Herman Rubin
>there are many RNG's out there with huge periods.
>For example:
>1. Mersenne Twister (MT)
>2. Add With Carry (AWC)
>3. Subtract With Borrow (SWB)
>4. KISS
>There is also Marsaglia's DIEHARD battery of tests to test the quality of
>these RNGs (in fact, 2 - 4 are his RNG's). Obviously, these RNGs are not
>good for crypto applications, but are good for anything Monte Carlo.
>Questions:
>1. My question is this, which is the best overall RNG out there, where best
>means that it passes all of the DIEHARD tests and has a huge period?
Is this reasonable? Every PRNG fails a test not too
difficult to state, namely, that it is produced the way
it is.
Periods much larger than the amount of computation which
will be carried out before the end of the functional
universe according to all theories are easy to get.
>2. Is there a RNG which can be implemented in a FPGA to simulate "White
>Noise"?
>Thank you for any references (as I am looking through Google as we speak),
>but cannot believe that many haven't been down this road.
>Flip
--
This address is for information only. I do not claim that these views
are those of the Statistics Department or of Purdue University.
Herman Rubin, Deptartment of Statistics, Purdue University
hru...@stat.purdue.edu Phone: (765)494-6054 FAX: (765)494-0558
S. Morgan
physical process".
For example, record N number of H/T of an un-biased coin flip. Store
the values as 0/1. When you want the RN sequence, just read out the
results from the stored table.
So, my answers to the question you asked (but not in terms of a computer
generated RNG that you might have been thinking about) are:
1. "best overall RNG" comes from a physical process 'deamed' to be
random
2. Sorry, this obviously does not support "field programmable gate
arrays". However, if you can control the hardware environment and don't
want reproduceable results (needed for comparison of different
techniques to the same PRN sequence), then why not add a tunnel diode
and generate your own RN sequence from the tunnel diode's output?
George Marsaglia
"flip" <flip_...@safebunch.com> wrote in message
news:10460956...@news-1.nethere.net...
It doesn't seem reasonable to consider that there is a best.
Just as the "Miss America Contest" judges have to view
talent, swim suit, evening gown and interview scores
before voting, for RNG's one should consider test results,
speed, simplicity, period length, size and possibly other
factors, depending on the application. (My own judging
tends to give extra weight to the equivalent of the
swim-suit competition---beauty of the underlying theory.)
Most of us are likely to give the most weight to
test results: the RNG should provide results
consistent with the underlying probability theory
for a wide variety of simulation problems, such
as those in the Diehard battery of tests.
The C procedure MWC256( ) listed below rates highly in all
categories. So does CMWC4096( ) below, with, so far,
record period.
The Mersenne Twister rates highly in test results and
period-length, but poorly in simplicity and swim-suit,
(beauty of underlying theory). C versions I have
seen take several dozens of statements.
Add-With-Carry (AWC) and Subtract-With-Borrow (SWB)
have been supplemented with Multiply-With-Carry (MWC)
and Complimentary-Multiply-With-Carry (CMWC).
I use SWB now only to provide IEEE standard 64-bit reals
in [0,1) directly, without the usual floating of integers. It was
the basis for the RNG in MATLAB.
The KISS generator rates highly in test results and
simplicity, as Keep It Simple Stupid was the theme for
its name. Its swim suit rating is high (in each of the
three pieces). The KISS RNG is used to provide randomness
by several gaming companies. The third component has had
several versions, depending on whether it is implemented in
assembler, C or Fortran. Unfortunately, Fortran does
not provide an easy way to get the top and bottom halves
of the 64-bit product of two 32-bit integers.
Here is the latest version of KISS, in C:
unsigned long KISS(){
static unsigned long x=123456789,y=362436000,z=521288629,c=7654321;
unsigned long long t, a=698769069LL;
x=69069*x+12345;
y^=(y<<13); y^=(y>>17); y^=(y<<5);
t=a*z+c; c=(t>>32);
return x+y+(z=t);
}
Choosing random seed values for x,y,z,and c provides
an excellent source of 32-bit random integers, with
period >2^125---adequate, but much shorter than the
remarkable 2^19937 of the Mersenne Twister or the
2^8222 of MWC256( ). For a really really long period,
consider CMWC4096( ) with period 2^131086,
listed after MWC256( ):
static unsigned long Q[256],c=362436;
unsigned long MWC256(void){
unsigned long long t,a=1540315826LL;
unsigned long x;
static unsigned char i=255;
t=a*Q[++i]+c; c=(t>>32);
x=t+c; if(x<c){x++;c++;}
return(Q[i]=x); }
Choosing random values for the static array Q[256] and c
puts you at a random starting place in the base b=2^32-1
expansion of 1/p, where p is the prime 1540315826*b^256-1.
The expansion of 1/p has cycles of more than 10^2475 base-b
'digits', and they are returned in reverse order, from a
random starting point in the cycle. One attraction of
MWC256 is that it allows us to index the Q-array with
a byte.
Now for the really really big:
static unsigned long Q[4096],c=362436;
unsigned long CMWC4096(void){
unsigned long long t, a=18782LL;
static unsigned long i=4095;
unsigned long x,r=0xfffffffe;
i=(i+1)&4095;
t=a*Q[i]+c;
c=(t>>32); x=t+c; if(x<c){x++;c++;}
return(Q[i]=r-x); }
This RNG has period > 2^131086, some 10^33459 times as
long as that of the Mersenne Twister and rates highly in
all categories. It provides the more than 10^39460 base-b digits
in the expansion of (p-1)/p, where p is the prime
p=18782*b^4096+1. , b=2^32-1. Those base-b 'digits' are
returned in reverse order from a random starting point
determined by the random choice of the initial values in Q[4096] and c.
The above are my own ratings. As a judge in the contest,
you may want to rate them yourself, and compare .
George Marsaglia
Bob Harris
shifts and xor's:
> y^=(y<<13); y^=(y>>17); y^=(y<<5);
This could also be considered as a PRNG. Where can I find the theory behind
that type of generator; i.e. in the first place how can I find shift counts
(in the example above, 13 left, 17 right, 5 left, like a combination lock)
that lead to long periods?
Bob H
George Marsaglia
"Bob Harris" <NspamI...@MINDnotSPRING.COM> wrote in message
news:BA810747.2E31D%NspamI...@MINDnotSPRING.COM...
Yes, it is somewhat like finding the numbers for
a combination lock.
Here is theory behind the method:
Rather than considering a computer word as a
32-bit integer, as we often do in order to use
number theory to develop RNGs, think of
y=(b1,b2,...,b32) as an element of the
vector space of 1x32 vectors with elements
in the field mod 2.
Then look for linear transformations, i.e, 32x32
matrices T, for which the sequence
yT,yT^2,yT^3,...
will be periodic with period 2^32-1.
It is not difficult to prove that this will happen
if and only if the matrix T is nonsingular with
order 2^32-1 in the group of nonsingular
32x32 binary matrices (elements in the field mod 2).
Now forming the matrix product yT might seem to be
too slow for practical use, and it is---for general T,
but not if we exploit certain operations common to
most CPU's: shift and exclusive-or (xor).
Let L be the 32x32 binary matrix that effects a
left shift of one position in the vector y, that is,
yL = (b1,b2,...,b32)L = (b2,b3,...,b32,0).
Then the matrix L^a will effect a left shift of 'a'
positions in y, and if R=L', (transpose), then
R^b will effect a right shift of b positions.
The matrices L and R are singular, since L^32=R^32=0,
but the matrices (I+L^a) and (I+R^b) are nonsingular
since their inverses can be expressed as finite sums
I+L^a+L^(2a)+... and I+R^b+R^(2b)+...
Now consider the nonsingular matrix T = (I+L^a)(I+R^b).
It certainly looks promising, as the matrix product
yT can be easily performed in a computer---for example,
in C, y^=(y<<a); y^=(y>>b);
But can we find 'a' and 'b' so that T=(I+l^a)(I+L^b)
has order 2^32-1 in the group of nonsingular binary
matrices?
The answer is no. There is no such pair a,b for
32x32 binary matrices.
But for T=(I+L^a)(I+R^b)(I+L^c) there are magic
numbers [a,b,c] such that the sequence
yT,yT^2,yT^3,...
has the full period 3^32-1 for any initial non-null
seed vector y.
For each choice of [a,b,c] with a<c there are 81 magic triples:
| 1, 3,10| 1, 5,16| 1, 5,19| 1, 9,29| 1,11, 6| 1,11,16| 1,19, 3| 1,21,20|
1,27,27|
| 2, 5,15| 2, 5,21| 2, 7, 7| 2, 7, 9| 2, 7,25| 2, 9,15| 2,15,17| 2,15,25|
2,21, 9|
| 3, 1,14| 3, 3,26| 3, 3,28| 3, 3,29| 3, 5,20| 3, 5,22| 3, 5,25| 3, 7,29|
3,13, 7|
| 3,23,25| 3,25,24| 3,27,11| 4, 3,17| 4, 3,27| 4, 5,15| 5, 3,21| 5, 7,22| 5,
9,7 |
| 5, 9,28| 5, 9,31| 5,13, 6| 5,15,17| 5,17,13| 5,21,12| 5,27, 8| 5,27,21|
5,27,25|
| 5,27,28| 6, 1,11| 6, 3,17| 6,17, 9| 6,21, 7| 6,21,13| 7, 1, 9| 7, 1,18| 7,
1,25|
| 7,13,25| 7,17,21| 7,25,12| 7,25,20| 8, 7,23| 8,9,23 | 9, 5,1 | 9, 5,25|
9,11,19|
| 9,21,16|10, 9,21|10, 9,25|11, 7,12|11, 7,16|11,17,13|11,21,13|12, 9,23|13,
3,17|
|13, 3,27|13, 5,19|13,17,15|14, 1,15|14,13,15|15,
1,29|17,15,20|17,15,23|17,15,26|
And for any choice of such a triple [a,b,c],
each one of these eight lines of C code:
y^=y<<a; y^=y>>b; y^=y<<c;
y^=y<<c; y^=y>>b; y^=y<<a;
y^=y>>a; y^=y<<b; y^=y>>c;
y^=y>>c; y^=y<<b; y^=y>>a;
y^=y<<a; y^=y<<c; y^=y>>b;
y^=y<<c; y^=y<<a; y^=y>>b;
y^=y>>a; y^=y>>c; y^=y<<b;
y^=y>>c; y^=y>>a; y^=y<<b;
may be used as the generating procedure for a RNG
with period 2^32-1, one that performs better in tests
of randomness than congruential RNG's.
Thus there are 8x81=648 magic combinations that
provide what I call full-period xorshift RNGs.
It is not difficult to write a C program that rapidly
forms T^2 or MT for two 32x32 binary matrices,
(multipling row x column is the 'and' operation)
and checking that T^(2^32)=T takes only a few seconds,
as does verifying that T^((2^32-1)/k) is not I for divisors
k of 2^32-1. So the entire process of finding all the magic triples
[a,b,c] takes only a few minutes. (But writing and verifying
the program took a few hours.)
Because of the impending dominance of 64-bit processors, I have
also found all the magic triples [a,b,c] for 64x64 binary
matrices. There are 275 basic triples, so that the corresponding
eight lines of C code provide 8x275=2200 full period xorshift RNGs
for 64-bit CPUs.
George Marsaglia
flip
news:lSidna-cF9d...@comcast.com...
Dr. Marsaglia,
two questions I received from these postings on KISS.
1. A question I have on KISS is wrt to the seed values. How critical is
the choice of the seed values to the properties of the numbers generated?
2. Finally, is there a method to generate 16 bit numbers with
period greater than 2^64 using something like KISS, it would really
help?
I suspect that it could utilize smaller constants to do this. I
am also looking at the Lagged Fibonacci approach to generating numbers
but have not spent too much time on it yet. While CMCW4096 looks awesome,
my implementation cannot afford an array of 4096 numbers.
Any thoughts?
Thanks, Flip
George Marsaglia
The KISS RNG uses four seeds, x, y, z, c,
with 0<=x<2^32, 0<y<2^32, 0<=z<2^32, 0<=c<698769069
except that the two pairs
z=0,c=0 and z=2^32-1,c=698769068
should be avoided.
That leaves
k=27681094672891588090390813844460011520,
some 10^37, possible choices for the seed set x,y,z,c.
A random choice of one of those k possible seed sets
will place you at a random position in a huge cycle
of over 10^37 32-bit integers, through which successive calls to
KISS() will move you one at time.
I think it would be hard to find a set of seeds that would not
produce a satisfactory starting place in that huge cycle.
There must be a few, but they should be extremely rare,
just as the monkey at a typewriter will type ilovesci.mathdon'tyou,
about once every 10^24 years. (60wpm,8hrday,noholidays).
Even the choice x=0,y=1,z=1,c=0 will produce a sequence of 32-bit
integers that pass all of the Diehard battery of tests of randomness,
with 269 satisfactory p-values.
(New edition at www.csis.hku.hk/~diehard/).
George Marsaglia
(The reason the seed pairs z=0,c=0 and z=2^32-1, c=698769068
should be avoided is that the third component of KISS is a
multiply-with-carry RNG, which produces the 'digits' in
the base-b expansion (b=2^32) of j/p, where p=a*b-1= 698769069*2^32-1
for some j. The choice z=0,c=0 produces 0/p=.00000000000...
and z=b-1,c=a-1 produces p/p=.ddddddddd... where d=b-1,
just as in base 10, p/p=.9999999999... .)