I made a new pseudo-random number generator

[p_m21987]

PLEASE NOTE: This project is purely for research and educational purposes, it must NOT be used for any practical purposes or applications. I'm posting it here for discussion, as a potential example of what NOT to do, and just as an interesting curiosity. I developed this function in the hopes that it will NEVER be used.
For practical purposes, you should NEVER attempt to write your own pseudo-random number generator, even if you are a highly experienced mathematician. It's very easy to make something with serious biases, and the consequences could be very serious.

Hello,

Recently I've been experimenting with 'dieharder' and 'practrand' statistical testing suite software, and researching different methods that are used to generate pseudo-random numbers, and I independently developed a technique that I've been calling "nutting".

'Nutting' is to take the result of a multiplication fmod 1.0.
This design uses two state variables which increase at different rates, and eventually wrap back around. (The wrap-around event is known as the 'sideways nutting action'.)
I tested it with a C implementation that writes a stream of 32bit integers to stdout, and fed that into 'practrand'. I also ran in parallel a separate version of it with the output bits in reverse order, and tested that with practrand too.

Both instances passed the full 32TB of testing from practrand. This is the point where you should be especially uncomfortable and suspicious, because it doesn't mean there aren't serious problems, it just means that practrand was not able to detect them.

The source code for the C implementation is here:

Code: Select all

double fnf( double v ){
 return v-(unsigned int)v;
}

double nn1 = 0.0;
double nn2 = 0.0;
// SuperMegaNutter, version 'b7'.
double SuperMegaNutter_b7(){
 nn1 += 0.499743730631690299468026764666361;
 if( nn1 > 2147483647.0 ){
  nn1 = fnf( nn1 );
 }
 nn2 -= 0.619712099029093592809216927916784+0.01*fnf(nn1);
 if( nn2 < 0.0 ){
  nn2 += 2147483647.0;
 }
 double a = fnf( nn1 * fnf( nn1 * 0.70635556640556476940272083033647 ) * fnf( nn1 * 0.96241840600902081782173280676402 ) );
 double b = fnf( nn2 * fnf( nn2 * 0.41059769134948776074925592402104 ) * fnf( nn2 * 0.83598420020319692166434447425947 ) );
 nn1 -= fnf(a+b)*0.249871861;
 return fnf( a * 2.79218049951371579168124556907984 - a
              + 
             b
              +
             fnf(nn1+nn2)*1.1235804235832785923589346329  );
}

#include <stdio.h>

int main(){
 unsigned int o;
 while(1){
  o = (unsigned int)(SuperMegaNutter_b7()*(double)0x100000000);
  #ifdef FLIP
  unsigned int flip = 0;
  for (int i = 0; i < 32; i++) {
   flip = (flip << 1) | ((o >> i) & 1);
  }
  fwrite(&flip, 4, 1, stdout);
  #else
  fwrite( &o, 4, 1, stdout );
  #endif
 }
}

Here is a BBC BASIC adaptation:

Code: Select all

      GCOL 0, -16: CLG: OFF
      REPEAT
        GCOL 0, 16 * FNSuperMegaNutter_b7
        PLOT 69, 1280 * FNSuperMegaNutter_b7, 1024 * FNSuperMegaNutter_b7
      UNTIL FALSE

      DEF FNf(v)=v-INTv

      REM SuperMegaNutter, version 'b7'.
      DEF FNSuperMegaNutter_b7
      nn1 += 0.499743730631690299468026764666361
      IF nn1 > 2147483647.0 THEN nn1 = FNf( nn1 )
      nn2 -= 0.619712099029093592809216927916784 + 0.01 * FNf( nn1 )
      IF nn2 < 0.0 THEN nn2 += 2147483647.0
      LOCAL a,b
      a = FNf( nn1 * FNf( nn1 * 0.70635556640556476940272083033647 ) * FNf( nn1 * 0.96241840600902081782173280676402 ) )
      b = FNf( nn2 * FNf( nn2 * 0.41059769134948776074925592402104 ) * FNf( nn2 * 0.83598420020319692166434447425947 ) )
      nn1 -= FNf( a + b ) * 0.249871861
      = FNf( a * 2.79218049951371579168124556907984 - a + b + FNf( nn1 + nn2 ) * 1.1235804235832785923589346329 )

To anyone who's interested, I'd like to welcome your comments and thoughts. I hope we can have a positive discussion.

Regards,
pm

[Hated Moron]

To anyone who's interested, I'd like to welcome your comments and thoughts. I hope we can have a positive discussion.

Given that it uses floating-point calculations, have you checked that it works just as well with 64-bit floats (as used by the Raspberry Pi, Android, iOS and in-browser editions of BBC BASIC) as it does with 80-bit floats (as used by the Windows, MacOS and Linux-86 editions of BBC BASIC)?

Assuming that its performance doesn't suffer with 64-it floats, the follow-up question is: does it generate the same sequence of numbers (from a given seed) in both 64-bit and 80-bit editions? Doing so is a desirable feature of a PRNG for BBC BASIC.

I well remember the complication that arose from the RNG in David Williams's Forces of Darkness which, because of tiny initial differences magnifying after many iterations, meant that a different Demo Mode data file was needed for 64-bit floats and 80-bit floats!