Close

a POSIXy rand()

A project log for PEAC Pisano with End-Around Carry algorithm

Add X to Y and Y to X, says the song. And carry on.

yann-guidon-ygdesYann Guidon / YGDES 2 days ago0 Comments

After the last log, here's a first draft of a POSIX-compliant PRNG to replace the dumb congruential generator. The combined period should exceed 2^64 yet no heavy operation (like multiplies or moduli) are used. Of course I should test the quality of the output...

#include <stdint.h>

#define RAND_MAX 32767

/* From the spec: If rand() is called before any calls
   to srand() are made, the same sequence shall
   be generated as when srand() is first called
   with a seed value of 1. */
uint32_t rand_LFSR    =  1,
         rand_PEAC_XC =  1;
uint16_t rand_PEAC_Y  =  0;

void srand(unsigned seed) {
  if (seed) rand_LFSR = seed;
      else  rand_LFSR = 0x89ABCDEFUL;
  rand_PEAC_XC  =  seed;
  rand_PEAC_Y   = ~seed & 1;
}

int rand(void) {
  // LFSR PRNG :
  uint32_t t = rand_LFSR & 1;
  rand_LFSR >>= 1;
  if (t)
    rand_LFSR ^= 0x82608EDB; // CRC32 poly

  // PEAC scrambler :
  uint16_t rand_X = rand_PEAC_XC; // 16 LSB
  rand_PEAC_XC = (rand_PEAC_XC >> 16)
               + rand_X + rand_PEAC_Y;
  rand_PEAC_Y  = rand_X + rand_LFSR; // 16 LSB

  return rand_PEAC_XC & RAND_MAX;
}

Let's see how far this one flies.

The first compile and run go smoothly:

[yg@ygrand]$ gcc -Wall -Os -o test_ygrand test_ygrand.c  && ./test_ygrand 
0000000000000001
0000111011011101
0101100010010100
0100110001001101
0110000110011000
0000110001000000
0000111111001111
0110110100001010
0110001101111111
0100001111011101
0001111011001111
0001111001100110
0101000000111101
0011010111111011
0010100111100101
0111000110110111
0110010010000111
0000000011101101
0011101011001100
0010000000101001
0100110100101110
0110011001110011
...

The PRNG is not even initialised, only using the default value of 1, and it immediately reached a high scrambling state. Calling srand() is a must though, and should call rand() a few times itself for "better results".

Initialising with srand(0) also gives a good-looking sequence with the help of some extra calls to rand().

The Hamming distance between consecutive numbers looks good after 10000000000 iterations. Note that slot 16 is cleared because only 15 bits are provided at each step.

  0 : 305938
  1 : 4575646
  2 : 32047057
  3 : 138854426
  4 : 416575520
  5 : 916453074
  6 : 1527396508
  7 : 1963876973
  8 : 1963793753
  9 : 1527393214
 10 : 916423856
 11 : 416545995
 12 : 138835287
 13 : 32036947
 14 : 4580472
 15 : 305334
 16 : 0

The peak is well centered around slots 7 and 8, the histogram is nicely symmetrical.

For comparison, I made a "fake PRNG" that pulls data out of /dev/random and here is the histogram:

  0 : 305586d
  1 : 4579839d
  2 : 32047565d
  3 : 138857968d
  4 : 416545970d
  5 : 916470025d
  6 : 1527349897d
  7 : 1963799989d
  8 : 1963842436d
  9 : 1527391622d
 10 : 916495473d
 11 : 416539621d
 12 : 138850583d
 13 : 32042754d
 14 : 4575830d
 15 : 304842d
 16 : 0d

It's pretty spot on but that's only one of many aspects to check.

I've put the files there : ygrand.tgz

Discussions