Information on calculator Random Number Generators from PPC Journal articles

[Dieter]

(11-12-2016 07:00 PM)Namir Wrote:  Will do .. I am away from home, so give me a few days.

If I may chime in here: I did a few tests with 1 million random numbers and various seeds for the (x+pi)³ generator. The mean is around 0,499 and the variance near 0,0837. In this regard the R-D generator has very similar properties, but the numbers seem to be less evenly distributed.

Dieter

[ttw]

There are some more or less realiable methods of generating pseuo-random numbers on a calculator. Which particular calculator is used determines the suitable and programming of each method. The suggested ones in the calculator literature are not known to be very good. There are two easy-to implement classes that can be combined if desired, although detailed properties of the combination have not been investigated extensively.

The first is the LGC (Linear Congruential Generator), also called the Lehmer generator. It works on a simple recursive formula.

x(0): the seed
X(n): the nth number
A: the seed
M: the modulus
b: the addend

X(n)=A*X(n-1)+b Modulo M

This can be converted to a fraction between 0 and 1 by forming Y = X(n)/M

There are conditions on the parameters A, M, b, and X(0) to insure that the generator cycles over the numbers from 0 to M-1 or a large fraction thereof.

If M is a power of two, then the conditions are easy. The multiplier, A, must be congruent to either 3 or 5 modulo 8 and b must be odd; x(0) is irrelevent. X(n+1)=A*X(n)+b Modulo M.

Another possibility is to use the recursion X(n+1)=A*X(n) Modulo M with X(0) being odd. This is also sufficient the the period is 1/4 as long.

In general a long period is wanted. There are lots of choices (a few are given below). A really long (like 2^256 or so) generator is pretty good for numerical work but the generator is easily predictable. Note that the last bit cycles with period 2:...01010101...., the next with period 4, etc. The leading bit cycles with length M.

Several sets of parameters:

M = 256, A = 181, b odd (too small for real work, but one can step through it on a calculator to see what happens.

M = 2^32, A = 2891335453, b odd (A chosen by search)

This size has been popular on 32-bit computers for obvious reasons. One probably cannot cycle it on a calculator.

M = 2^64, A = 4976020386757901309, b odd (A chosen by search)

This generator is longer than most simulations are run for. I used it for caluclator games.

The second class of generators is the Shift Register Generator or Tausworthe Generator. It is based on shift a register over 1 bit and then XORing the result with a bit pattern (a primitive polynomial in the finite field the length of the generator to be technical) if the bit shifted out is a 1 bit. These generator generate good bit streams but not as good binary integers. There are lots of fixes. Similar generators are used for some encryption methods as (for example) using one generator to clock another or outputting one generator based on whether the other generator is a 1.

M is the size, P is the primitive polynomial (looked up in a table), T(0) must be non-zero. This genrator has a cycle 2^L-1 for length L; it generates all non-zero integers in range. The recurrence is given below (here, I'm shifting to the right to be easier to program). LSB stands for the Least Significant Bit (the units digit if read as a binary number), RS for right shift by one bit, end off, zero fill. The polynomials are given in hex.

IF LSB(T(n)) = 1
THEN T(n+1) = RS(T(n)) XOR P
ELSE T(n+1) = RS(T(n))
ENDIF

M = 2^8 (an 8 bit generator), P = B8
M = 2^32 (a 32 bit generator), P = 80200003
M = 2^64 (a 64 bit generator), P = D800000000000000

One way to create a longer generator is to combine two (or more) generator with differing lengths. With the proper combiner, the cycle length of the resulting generator is the Least Common Multiple of the cycles of the constituents. I've always used a shift register and a linear congruentialg with the same length (cycle 2^L and 2^L-1 so no thinking). The length L objects can be treated as integers or bit-strings. As bit strings the obvious operation is XOR; as integers, the obvious operations are Add and Subtract (modulo 2^L).

[Namir]

I tested the algorithms:

1) Frac(9821 * x + 0.211327)
2) Frac(R-D(x))
3) Frac((pi+x)^3

I used Excel VBA to automate the test. The initial seed that generated 5 sets of 10,000 numbers were calculated using:

initial_seed = 1234 * RND(1)

Regarding the maximum absolute auto-correlation values, Algorithm 2 was 0.0287, while Algorithms 1 and 3 where 0.036316 and 0.0375916, respectively. Thus Algorithm 1 seem to yield less correlated numbers.

Regarding the range of Chi-Square values, the results are:

Alg. 1 has the range (2.878, 16.79)
Alg. 2 has the range (7.48, 24.206)
Alg. 3 has the range (8.302, 14.976)

Thus Algorithm 1 seems to yield better distributed numbers. Algorithm 2 seems to yield the worst distributed numbers.

All three algorithms generated number with an average and standard deviation values that are close to the theoretical expected values.

[Paul Dale]

(11-13-2016 07:27 AM)Namir Wrote:  Regarding the maximum absolute auto-correlation values, Algorithm 2 was 0.0287, while Algorithms 1 and 3 where 0.036316 and 0.0375916, respectively. Thus Algorithm 1 seem to yield less correlated numbers.

Wouldn't this favour Algorithm 2? It's autocorrelation is the lowest of the three.

Earlier this year I created a "random" number generator that passes the NIST SP 800-90B tests for randomness:

\[ { x_0 = frac ( \pi ) } \] \[ { x_{n+1} = frac ( e^{2 x_n + 1} ) } \] \[ R_n = frac( 2^{20} x_n ) \]

These \( R_n \) are completely deterministic and are produced from a rather well behaved function. There are also no problems with floating point rounding errors where the exponential is evaluated.

The point is that it is impossible to demonstrate randomness. It is only possible to demonstrate non-randomness.


Pauli

[ As an aside, the TestU01 suite did pick this function as non-random.]

[Namir]

Pauli,

I ran the same test (just one run because you have a fixed way to initialize the random numbers). The maximum absolute auto-correlation value is 0.0281497. The Chi Square value is 21.612. The absolute auto-correlation value of your algorithm is better that the three algorithms I tested. The Chi Square value is higher than the Chi Square values or Algorithm 1 and Algorithm 2, but less than that of Algorithm 2. Your algorithm seems to generate numbers with better randomness, but a little bias in the distribution.

Overall, your algorithm is very nice and I consider it suitable for calculator games and application.

Have you tried to use a variable way to initialize x0?

Namir

[Namir]

I took Pauli's function and changed the way we calculated x0 as frac(1234*seed). I got five sets of 10,000 random numbers. Four sets generated auto-correlation values below 0.03 and one set generated above 0.03. The average absolute auto-correlation is very close to that of Pauli's original PRNG algorithm. The Chi Square values between 38 and 67, which is a range that is higher that the previous test I mentioned in this thread.

Overall, Pauli's original random number is a very nice and simple algorithm for calculators.

Namir

[Dieter]

(11-13-2016 07:27 AM)Namir Wrote:  Regarding the range of Chi-Square values, the results are:

Alg. 1 has the range (2.878, 16.79)
Alg. 2 has the range (7.48, 24.206)
Alg. 3 has the range (8.302, 14.976)

Thus Algorithm 1 seems to yield better distributed numbers. Algorithm 2 seems to yield the worst distributed numbers.

FTR, I have done a few hundred runs with 10.000 random numbers each, using the 9821 and (x+pi)³ RNG. The Chi² values came out very similar. For instance, a substantial part for the third RNG is below 5, values like 2 or 3 appear quite often, sometimes even less. That's why I do not think that we can determine this way if one of these two RNGs is "better" than the other.

For the R-D generator I got Chi² values that sometimes were significantly higher (like 30 or 40, in single cases even much more). So this seems to be not on par with the other two. But honestly... for calculator games: who cares?

Dieter

[Namir]

I finally tested frac(997*x) and got comparable results to the previous algorithms. Four sets of 10,000 random numbers had absolute auto-correlation below 0.03 and one above 0.03. The range of Chi Square was 5.41 to 10.944, which seems decent. The frac(997*x) sees to work well if random seeds (I used 1238*RND(1) in Excel VBA) and avoid "mine-field" numbers.

I think the general consensus is that the algorithms mentioned in this thread (with proper initialization when applicable) work fine with calculator applications, simulation, and games. In the case of PC use, it's a different ball game!!!

Namir

[Namir]

As my last very test, I tested frac((pi+x)^5) and was surprised that my (limited) test showed that the maximum absolute auto-correlations was 0.028144 (i.e. below 0.03) and the Chi Square is in the range of 3.47 and 11.222. This makes the frac((pi+x)^5) algorithm as slight favorite for calculators, over the other algorithms.

Namir

[Gene]

Ok, so for newseed = FRC (PI + Seed) ^ 5 how does this look to save as much of the stack as we can? Assumes seed is already in 00. FEEL free to condense. This preserves X, Y, Z in Y, Z and T. The original T and L are lost.

LBL R
PI
ST+ 00
X<> 00
STO 00
X^2
X^2
X<> 00
ST* 00
X<> 00
FRC
STO 00
END

[Dieter]

(11-13-2016 10:09 PM)Gene Wrote:  Ok, so for newseed = FRC (PI + Seed) ^ 5 how does this look to save as much of the stack as we can?

Do we have to?

Namir, could you do a few more runs (like 100 times a sequence of 10.000 random numbers)? According to my results the ^5 version is not better than the ^3 version, as far as Chi² and autocorrelation are concerned. I may be wrong, but just to be sure...

(11-13-2016 10:09 PM)Gene Wrote:  Assumes seed is already in 00. FEEL free to condense.

Simple:

Code:

LBL R
PI
ST+ 00
X<> 00
X^2
ST* L
ST* L
X<> L
FRC
STO 00
END

But now let's wait for an MCode version of the 9821 generator using R00 so that no buffer is required (cf. the other RNG thread). I think this will be the best solution as, unlike the user-code version, it does not sacrifice speed. RoberM already suggested a solution. Ángel, your turn now. ;-)

Dieter

[Gene]

We may not have / use the mcode generator.

This is the fallback option.

Code:

 01 LBL "S" 
 02 RCL 00 
 03 "SEED?" 
 04 PROMPT 
 05 ABS 
 06 FRC 
 07 STO 00 
 08 RTN 
 09 LBL "R" 
 10 PI 
 11 ST+ 00 
 12 X<> 00 
 13 X^2 
 14 ST* L 
 15 ST* L 
 16 X<> L 
 17 FRC 
 18 ABS 
 19 STO 00 
 20 END

[Namir]

(11-13-2016 10:46 PM)Dieter Wrote:  

(11-13-2016 10:09 PM)Gene Wrote:  Ok, so for newseed = FRC (PI + Seed) ^ 5 how does this look to save as much of the stack as we can?

Do we have to?

Namir, could you do a few more runs (like 100 times a sequence of 10.000 random numbers)? According to my results the ^5 version is not better than the ^3 version, as far as Chi² and autocorrelation are concerned. I may be wrong, but just to be sure...

(11-13-2016 10:09 PM)Gene Wrote:  Assumes seed is already in 00. FEEL free to condense.

Simple:

Code:

LBL R
PI
ST+ 00
X<> 00
X^2
ST* L
ST* L
X<> L
FRC
STO 00
END

But now let's wait for an MCode version of the 9821 generator using R00 so that no buffer is required (cf. the other RNG thread). I think this will be the best solution as, unlike the user-code version, it does not sacrifice speed. RoberM already suggested a solution. Ángel, your turn now. ;-)

Dieter

Sorry, I don't have time to do a 100 run. I don't doubt you results since 5 runs is below the minimum of 30 to begin to approach normal distribution.

Namir

[Dieter]

(11-14-2016 01:09 AM)Namir Wrote:  Sorry, I don't have time to do a 100 run. I don't doubt you results since 5 runs is below the minimum of 30 to begin to approach normal distribution.

In the meantime I've read a few papers on testing RNGs and found an interesting quote from Donald E. Knuth's well known book "The Art of Computer Programming" where a chapter deals with random number generation and testing. I do not own this book, but he seems to suggest the following way of evaluating the Chi² value obtained after a run with, say, 10.000 of random numbers that are grouped in, say, ten categories to check their distribution in the [0;1[ domain. According to Knuth both very high and very low Chi² values are critical:

• If the Chi² value is below the 1% quantile or above the 99% quantile, the tested sample can be considered "non-random".
• If the Chi² value is between the 1% and 5% quantiles or between the 95% and 99% quantiles, the sample should be considered "suspect".
• If the Chi² value is between the 5% and 10% quantiles or between the 90% and 95% quantiles, the sample is considered "almost suspect".
  

In a way this sounds logical. Consider a die that is rolled 60 times and it comes out to exactly ten ones, ten twos etc., This would yield a Chi² value of zero and thus qualify the die as "non-random".

If a run with 10.000 random numbers is grouped into ten categories (0...0,1; 0,1...0,2; 0,2...0,3; ..., 0,9...1) the critical Chi² values for 9 degrees of freedom are:

Code:

 1%: 2,088   99%: 21,666  ("non-random")
 5%: 3,325   95%: 16,919  ("suspect")
10%: 4,168   90%: 14,684  ("almost suspect")

So let's say that for ten categories the Chi² value should fall somewhere between 4 and 15 (that's 9%/91%), or at least between 3 and 18 (3,5%/96,4%).

I did a few tests with the 9821 generator and a number of runs with 10.000 random numbers each. In most cases the Chi² value was well within the mentioned limits. Only sometimes especially the lower limit was not quite met, i.e. the numbers were a bit "too evenly distributed". The (x+pi)^3 and (x+pi)^5 RNGs were no match in this regard.

Additional tests with 20 and finally 100 categories showed essentially the same the results.

What do the math experts say?

Dieter

[Paul Dale]

There are a number of programs available for testing randomness. Knuth is a decent start but there has been a lot more work done over the intervening years.

NIST describes their current process in SP 800-90B draft 2. They also have a Python script that does the analysis. The description is a draft and the examples contain errors but most of the algorithms are properly described, I believe the Python script has correct implementations although there is one case where two sided normal distribution critical values are used and they probably should be single sided. These tests make a reasonable effort, do not require a huge amount of data and attempt to provide an estimate of the entropy present in the samples. However, it isn't difficult to trick them into passing.

The ent program is another testing utility that is very fast and doesn't require a lot of data. In my experience it overestimates the entropy present and shouldn't be given too much weight. A failure here is a definite indication of non-randomness.

The dieharder suite is the defacto standard for testing randomness. It requires a substantial amount of data and can be very slow. It provides a decent battery of tests. Due to the number of tests and the critical values used, some tests usually produce a weak result even for true random numbers.

The current state of the art is the TestU01 suite. It requires a lot of data and is quite discriminating. It usually runs a little faster than the Dieharder suite.

My \( x_{n+1} = frac ( e^{2 x_n + 1} ) \) passes the NIST tests and ent. Dieharder hints at there being a problem with a couple too many weak results and TestU01 fails it outright.


Pauli

[Dieter]

(11-14-2016 10:22 PM)Paul Dale Wrote:  There are a number of programs available for testing randomness. Knuth is a decent start but there has been a lot more work done over the intervening years.
(...)

Thank you for all these links. Will be a lot to read, but I fear much of this is beyond my scope. #-)

(11-14-2016 10:22 PM)Paul Dale Wrote:  My \( x_{n+1} = frac ( e^{2 x_n + 1} ) \) passes the NIST tests and ent. Dieharder hints at there being a problem with a couple too many weak results and TestU01 fails it outright.

If it's good enough for NIST, why not use it for the 67 Games ROM? It also can handle zero end even negative seeds. And a 41 implementation can be done in one stack level as well:

Code:

LBL "R"
RCL 00
FRC
ST+ X
SIGN
ST+ L
X<> L
e^x
FRC
STO 00
END

The only question is speed – the exponential is not exactly fast. I did a test on V41 set to approx. "hardware speed" and it came out about 15% slower than the pi^3 generator. This is acceptable (compare this to the 9821 generator which is >60% slower). Do you think we may use your RNG for the 67 Games ROM?

Edit: I did a short test and found very high Chi² values. Looking at the number distribution there were several cases where the range 0,2...0,3 had extremely low counts while 0,9+ was very high. Here is an example for 100.000 RNs starting with seed=0,181932:

Code:

0,0...    4091
0,1...   12881
0,2...    1534 (!)
0,3...    9148
0,4...    7835
0,5...   10334
0,6...    9083
0,7...   12908
0,8...   10417
0,9...   21769 (!)

The test was done in Excel VBA with rounding to calculator precision (round 2x+1 to 9 digits, calculate the exponential, round this to 8 resp. 9 digits, return fractional part). Running the test with native 15-digit precision returned less extreme results, but still very high Chi² values (around 300 for n=100.000). What am I doing wrong?

Dieter

[Paul Dale]

I'll check some more test results for this function tomorrow if I remember -- as I mentioned it fails TestU01 and is marginal for the dieharder suite.

I never claimed this generator was unbiased, just that it passed the current (at the time) NIST randomness tests with a high entropy result (as close to 8 bits per byte as was possible). It was designed to be presented to the CMVP conference earlier this year as a simple deterministic function that passed the tests but showed they weren't great. It was based on the first draft tests, but the function passes the second draft also, and at the same level of apparent entropy

I'll also try to come up with something similarly simple but less biassed. Entropy is my day job currently...


Pauli

[Namir]

Pauli,

I believe that all (most?) of the PRNG discussed in this thread are, shall we say, light weight PRNG that should work well/satisfactorily with calculator programs, simulations, and games. I don't expect ANY of them to pass the TESTU01 benchmark. Two years ago I tinkered with PRNGs using first Excel and then Matlab (was faster). I was able to (accidentally) find a method that generated calculation noise (in the fractional part of the results) with good random numbers. I checked the method with Excel VBA and realized that the random calculation noise (generated when you take the fractional part of large numbers) was particular to Matlab! I did use several tests (not the full sweet) for randomness and calculated an overall factor that gauged the general performance. You can see the results on my web site. The lesson I learned from this experience is to use integers as the core of PRNG algorithms and calculated the uniform random number only as the last step. Using mainly-floating-point methods puts you at the mercy of the programming language (and version) and its floating-point precision.

Namir

[Dieter]

(11-13-2016 10:46 PM)I Wrote:  According to my results the ^5 version is not better than the ^3 version, as far as Chi² and autocorrelation are concerned. I may be wrong, but just to be sure...

I did another test of the (x+pi)^5 generator and got a somewhat interesting result.
The test was done in Excel VBA, simulating the 41's 10-digit precision with this code:

Code:

    x = Round(r + 3.141592654, 9)
    r = x ^ 5
    d = 9 - Int(Log(r) * lge)  ' remaining digits right of the decimal point, lge = 0,43429...
    r = Round(r - Int(r), d)

And here are the results from 10 runs of 1000 random numbers each:

Code:

 Chi² =   10,40     33,44     33,44     33,44     33,44     33,44     33,44     33,44     33,44     33,44
---------------------------------------------------------------------------------------------------------
0,0..0,1:   106       104       104       104       104       104       104       104       104       104
0,1..0,2:   118       130       130       130       130       130       130       130       130       130
0,2..0,3:   106       120       120       120       120       120       120       120       120       120
0,3..0,4:    87        76        76        76        76        76        76        76        76        76
0,4..0,5:   106        98        98        98        98        98        98        98        98        98
0,5..0,6:    83        72        72        72        72        72        72        72        72        72
0,6..0,7:    98       110       110       110       110       110       110       110       110       110
0,7..0,8:   104       104       104       104       104       104       104       104       104       104
0,8..0,9:    89       108       108       108       108       108       108       108       108       108
0,9..1,0:   103        78        78        78        78        78        78        78        78        78

After the first run (starting with r=0,12345) and Chi²=10,4 all subsequent runs with 1000 numbers yielded the same results. I assume this is due to a rather short period. This behaviour occurs with different seeds.

Now, we have learned that the autocorrelation of the ^5 generator is a bit better than with the ^3 version. But there seems to be another property that makes the ^5 generator ..."not the first choice". Gene and Namir, what do you think?

BTW I could not find something similar with the (x+pi)^3 generator.
At least not yet. ;-)

Dieter

[Paul Dale]

(11-15-2016 07:12 AM)Dieter Wrote:  Edit: I did a short test and found very high Chi² values. Looking at the number distribution there were several cases where the range 0,2...0,3 had extremely low counts while 0,9+ was very high.

You didn't include the \( frac ( 2^{20} x_n ) \) part. The leading bits don't pass the randomness tests but going down the mantissa a bit does.

- Pauli