Re: [Dev-C++] random numbers with double precision

[Alfred P. R. <al...@ya...>]

Pardon all, in advance, for the length of this reply!
 
Don't waste your time Anton, with the rand() function.  I spent a week last month going through RNGs, developing an implementation of random number generators that pass the DieHard tests.  The system rand() function is severely deficient, and fails ALL Diehard tests!
 
I wound up using The Mother (integer arithemetic) and an implementation of the Wichman and Hill (double precision) procedure RNG.  The original implementation of both, as available on the Internet, may return both integer and double precesion values concurrently.  An excellent page on RNGs is Random Numbers and Monte Carlo Methods, Random Number Generators.
 
As both are public domain, and therefore not copy-rightable by me, I'll give you the coding to both.  These PASS ALL of the DieHard tests in an excellent fashion, with very good period.  I've statistics on both of these, BTW, created using a 300M unsigned 16 bit integer test file using the constant seed -1 (I've also the code available for generating the test file, running DieHard, and pattern repetition checking.)
 
George Marsaglia's "The Mother"  RNG:

// Mother RNG related defines.
#define m16Long 65536L    /* 2^16 */
#define m16Mask 0xFFFF    /* Mask for lower 16 bits */
#define m15Mask 0x7FFF    /* Mask for lower 15 bits */
#define m31Mask 0x7FFFFFFF   /* Mask for 31 bits */
#define m32Double 4294967295.0  /* 2^32-1 */
// Mother related variables.
bool topword;   // Return to Mother RNG word if true.
short mother1[10];  // Arrays of values for Mother RNG.
short mother2[10];
long current_mother; // Current 32 bit Mother RNG random number.

The Mother Itself:
unsigned int mrand ( void ) {
/*
  Name: Modified Mother Random Number Generator
  Copyright: 2005, Happy Cat Technologies
  Author: Alfred P. Reaud
  Date: 01/02/05 22:34
  Description: This member is an unsigned integer implementation of George Marsaglia's
 The Mother of All Random Number Generators.  The seed is set at the time of
 this class's construction from the system clock. Verification has occurred
 using the Diehard battery of tests, with no significant differences encountered
 between this implementation using a timer generated seed and the test
 implementation using a constant seed.
  RETURNS: An unsigned integer is returned representing a 16 bit random number.
  VERIFIED: 02/02/05 22:39  Passes the DIEHARD RNG TESTS! 200 nsec/number.
*/
unsigned long number1, number2;
 // Check for which word to return.
 if (topword) {
 // Return the top 16 bits.
 topword = false;
 return (unsigned int)((current_mother >> 16) & m16Mask);
 }
 else {
  // Generate a new random number, save, and return the bottom 16 bits.
  topword = true;
  /* Move elements 1 to 8 to 2 to 9 */
  memmove( mother1 + 2, mother1 + 1, 8 * sizeof(short) );
  memmove( mother2 + 2, mother2 + 1, 8 * sizeof(short) );
  /* Put the carry values in numberi */
  number1 = mother1[0];
  number2 = mother2[0];
  /* Form the linear combinations */
  number1 += 1941 * mother1[2] + 1860 * mother1[3] + 1812 * mother1[4]
    + 1776 * mother1[5] + 1492 * mother1[6] + 1215 * mother1[7]
    + 1066 * mother1[8] + 12013 * mother1[9];
  number2 += 1111 * mother2[2] + 2222 * mother2[3] + 3333 * mother2[4]
    + 4444 * mother2[5] + 5555 * mother2[6] + 6666 * mother2[7]
    + 7777 * mother2[8] + 9272 * mother2[9];
  /* Save the high bits of numberi as the new carry */
  mother1[0] = number1 / m16Long;
  mother2[0] = number2 / m16Long;
  /* Put the low bits of numberi into motheri[1] */
  mother1[1] = m16Mask & number1;
  mother2[1] = m16Mask & number2;
  current_mother = (((long)mother1[1]) << 16) + (long)mother2[1];
  return (unsigned int)(current_mother & m16Mask);
 }
}
 
Initializing "The Mother":
{
short n = 0;
short *p = NULL;
unsigned short sNumber = 0;
unsigned long seed_dWord = 0;
unsigned long number = 0;
 // Initialize the Mother RNG.
 seed_dWord = (rand() << 16) + rand();
 sNumber = seed_dWord & m16Mask; /* The low 16 bits */
 number = seed_dWord & m31Mask; /* Only want 31 bits */
 p = mother1;
 for ( n = 18; n--; ) {
  /* One line multiply-with-carry */
  number = 30903 * sNumber + (number >> 16);
  *p++ = sNumber = number & m16Mask;
  if (n == 9) p = mother2;
 }
 /* Make cary 15 bits */
 mother1[0] &= m15Mask;
 mother2[0] &= m15Mask;
}

 
Wichman and Hil RNG:
#define TWORAISED16 65536 // 2^16.  Used to extract a random unsigned integer.
unsigned int secure_random ( void ) {
/*
  Name: Secure Random Member Function
  Copyright: 2004-2005, Happy Cat Technologies
  Author: Alfred P. Reaud
  Date: 28/01/05 21:26
  Description: This member implements a more secure random number generator using
 the Wichman and Hill procedure.
  RETURNS: An unsigned integer is returned representing a 16 bit random number.
  VERIFIED: 28/01/05 21:26  Passes the DIEHARD RNG TESTS! 975 nsec/number.
*/
int ix, iy, iz;
double result, rand_f, i_part;
 // Create 3 random numbers.
 ix = (171 * rand()) % 30269;
 iy = (172 * rand()) % 30307;
 iz = (170 * rand()) % 30323;
 // Sum them as a double.
 result = ((double)ix / 30269.0) + ((double)iy / 30307.0) + ((double)iz / 30323.0);
 // Make a random number that will be 16 bits wide.
 rand_f = modf(result, &i_part) * TWORAISED16;
 result = modf(rand_f, &i_part);
 return ((unsigned int)i_part) % TWORAISED16;
}

Anton Zilman <zi...@ma...> wrote:
Hi,
 
  I need to produce random numbes uniformly distributed between 0 and 1. I can use (rand()%MAX_INT)/MAX_INT, but this gives me only 30 thousand something different numbers, because rand() is an integer.
  So, if I cannot simulate something probability of which is less than 1/MAX_INT. Is there an analog of drand() function in Dec C/C++?
  If anyone can help, it is apprecated!
 
  Anton




Alfred P. Reaud