Date: Fri, 23 Sep 1994 11:24:33 -0700 (PDT)
From: Gordon Hogenson <ghogenso AT u DOT washington DOT edu>
Subject: Re: Problem with rand() (fwd)
To: ZRS <roalz AT varano DOT ing DOT unico DOT it>
Cc: djgpp <djgpp AT sun DOT soe DOT clarkson DOT edu>

On Fri, 23 Sep 1994, ZRS wrote:

> 
> It seems that rand() returns alternately an even and an odd number
> (I.E. 0 141 2900 12717 ...)
> 
> 1.) Is this the recommended (ANSI) behaviour ??
> 
> BUT... I can get only EVEN numbers for 'result'.
> I do know I could use RNGs from libg++, but then I must have libg++ 
> anywhere the Die is ported!
> 

Never, never trust a system's rand() function!  Most implementations
of rand() are hopelessly flawed.  ANSI wasn't much help either--
in their specification, they saw fit to include an "example"
implementation, which was a very poor generator indeed.  Some
library implementors further decided to add 'extra flourishes' of
randomness, by putting in masks here and there, or some such
nonsense, most of which make things worse.   IBM, once contacted about 
the poor quality of thier rand() function, was reported to have 
made the statement, "We guarantee that each number is random
individually, but we do not guarantee that more than one of them
is random." (!??!?)

There is a good discussion of random numbers in "Numerical Recipes"
(2nd edition or later).  The following bit of code is the
random number generator we use for our Monte Carlo simulations,
which is basically taken from Numerical Recipes.  It returns
double values from 0 to 1, but you can easily modify it to return
integers if you like.

The authors of Numerical Recipes are so confident in this generator
that they offer $1000 to the first person to write a better one or
come up with a test that it fails.

---------------------------------random.c---------------------
/* The following codes implements a useful random number generator.
This verson was obtained from the parasoft parallel one, used on 
tranputers.

'pranset' seeds the generator, 'pranf' carries out the recursion,
(pranf is not called from client code).  Clients use the 'ran0()' function 
which does a "shuffle"  */
  


static int PR_AAA;
static int PR_BBB;
static int PR_randx;

#define	PR_MULT	1103515245
#define	PR_ADD 	12345
#define	PR_MASK	( 0x7fffffff )
#define	PR_TWOTO31	2147483648.0

void pranset(long seed) 
{
    PR_AAA = 1;
    PR_BBB = 0;
   
    
	PR_AAA = (PR_MULT * PR_AAA) & PR_MASK;
	PR_BBB = (PR_MULT * PR_BBB + PR_ADD) & PR_MASK;
	
	    PR_randx = (PR_AAA*seed + PR_BBB) & PR_MASK;
}

double pranf(void)
/*
	Return a random double in [0, 1.0)
*/
{
	double retvalue;

	retvalue=PR_randx / PR_TWOTO31;
	PR_randx = (PR_AAA*PR_randx + PR_BBB)& PR_MASK;
	return( retvalue );
}


double ran0(void)
{
 static double z,v[98],dum;
 static int iff=0;
 int j;
  if (iff == 0)
    { iff = 1;
      for (j=1; j<=197;j++) dum= pranf();
      for (j=1;j<=97;j++) v[j] = pranf();
      z = pranf();
    }
  j=1+(int)(97.0*z);
  z=v[j];
  v[j]=pranf();
  return z;
}