looking for uniformly distributed (0..1] random number generator

[Eliot Moss]

... and here's our version, based on the same source (Note: a previous posting
of mine had some of the citation information wrong, which I have now corrected
below) ....							Eliot
------------------------------------------------------------------------
/*
      $RCSfile: random.h,v $
      $Revision: 1.3 $
      $Date: 90/01/21 14:08:34 $
      $Author: moss $
      $State: Exp $
      $Locker: moss $
*/
/* header file for random.c */
# ifndef __RANDOM__
# define __RANDOM__
# ifdef __STDC__
extern double Random(void);		/* Random Generator */
extern int RandomInt(int modulus);	/* Returns integer 0 <= r < modulus */
extern int TimeRandomize(void);		/* Assign seeds using the system time */
extern int SetRandomSeeds(int s1,int s2);	/* Assign seeds. */
# else
extern double Random();		/* Random Generator */
extern int RandomInt();		/* Returns integer 0 <= r < modulus */
extern int TimeRandomize();	/* Assign seeds using the system time */
extern int SetRandomSeeds();	/* Assign seeds. */
# endif
# endif
------------------------------------------------------------------------
/*
      $RCSfile: random.c,v $
      $Revision: 1.3 $
      $Date: 90/01/19 06:42:14 $
      $Author: nayeri $
      $State: Exp $
      $Locker: moss $
*/
/* 
* 
* random.h 	Random Generator routines for 32-bit machines.
* Author 	Farshad Nayeri
* Date 		July 15, 1989
* Description	These are a set of library routines for 
*		random generation.  With the current constants
*		the period is about 2E18.
* Disclaimer
*		The following routines are correct to the best
*		of my knowledge, however I will not responsible 
*		for direct or indirect implications of using these 
*		routines in case of any error.
* See Also
*		"Efficient and Portable Combined Random Number 
*		Generators" by Pierre L'Ecuyer, Comm. of ACM, 
*		Vol 31, Number 6, June 1988.
*
*/

#ifdef VMS
#  include <time.h>
#else
#  include <sys/types.h>
#  include <sys/timeb.h>
#endif

#define NORMAL 0	 
#define ERROR (-1)	/* TimeRandomize depeneds on this value */
#define MAX(s1,s2) ((s1) > (s2) ? (s1) : (s2))
#define MIN(s1,s2) ((s1) < (s2) ? (s1) : (s2))

#define MINSEED1 	1
#define MAXSEED1	2147483562
#define DIVISOR1 	53668			/* constants for random */
#define SEED1KMULT 	40014			/* generator  */
#define K1MULT 		12211

#define MINSEED2 	1
#define DIVISOR2 	52274
#define MAXSEED2 	2147483398
#define SEED2KMULT 	40692
#define K2MULT 		3791

#define TIMEFACTOR 	10000

double normalizer = (1.0 / ((double)MAX(MAXSEED1,MAXSEED2)));
						/* normalization factor */
double maxSeed1and2 = MAX(MAXSEED1,MAXSEED2);

static int seed1 = MAXSEED1 / 2 ,		/* first and second seed*/
         seed2 = MAXSEED2 / 2;


/*
*
* DESCRIPTION
*	Random will return a random double precision real value
*	between 0.0 and 1.0 (See Efficient and Portable Combined
*	Random Number Generators by Pierre L'Ecuyer, Comm. of ACM,
*	Vol 31, Number 6, June 1988. )
* NOTES
*	This is a 32-bit generator.
*/

# ifdef __STDC__
double Random(void)
# else
double Random()
# endif

  {
  int z,k;

  k = seed1 / DIVISOR1;
  seed1 = SEED1KMULT * ( seed1 - k * DIVISOR1) - k * K1MULT;
  if (seed1 < 0) seed1 += MAXSEED1 + 1;
              /* seed1 = seed1 + MAXSEED1 + 1 */


  k = seed2 / DIVISOR2;
  seed2 = SEED2KMULT * ( seed2 - k * DIVISOR2) - k * K2MULT;
  if (seed2 < 0) seed2 += MAXSEED2 + 1;
	      /* seed2 = seed2 + MAXSEED2 + 1 */

  z = seed1 - seed2;
  if (z < 1) z += maxSeed1and2;
  return ((double) z  * normalizer);

  } 


/*
*
* DESCRIPTION
*	RandomInt returns a random integer between 0 and the modulus-1
*	Since Random only returns values r > 0, 1-r < 1, so 
*	this function will always return values between 0 and modulus-1
*
*/

# ifdef __STDC__
int RandomInt(int modulus)
# else
int RandomInt(modulus)
  int modulus;
# endif
  {
  return((int) ((1.0 - Random()) * modulus));
  }



/*
*
* DESCRIPTION 
*	If s1 and s2 are valid, set the seeds, otherwise returns 
*	error and keep the old seeds.
* RETURN VALUE
*	0 if everything is O.K.
*	-1 if there is an error.
*
*/

# ifdef __STDC__
int SetRandomSeeds(int s1,int s2)
# else
int SetRandomSeeds(s1,s2)
  int s1;
  int s2;
# endif
  {

  if ((s1>=MINSEED1) && (s1<=MAXSEED1)) 
    seed1 = s1;
  else {
    return(ERROR);
    }
 
  if ((s2>=MINSEED2) && (s2<=MAXSEED2)) 
    seed2 = s2;
  else {
    return(ERROR);
    }

  return(NORMAL);
  } 


/*
* 
* DESCRIPTION
*	Use the system time to set the seeds. 
* RETURN VALUE
*	0 if everything is O.K.
*	-1 if there is a problem with time routine.
*
*/
# ifdef __STDC__
int TimeRandomize(void)
# else
int TimeRandomize()
# endif
  {
  unsigned short theTime = 0;
  struct timeb tp;
 
  ftime(&tp);
  theTime = tp.millitm / 10;	/* millitm is always multiple of 10 */

  if ((theTime & 001) == 0)
     SetRandomSeeds((int)(MIN(seed1,seed2) +  theTime * TIMEFACTOR) 
							% MAXSEED1 + MINSEED1,
                    (int)(MAX(seed1,seed2)             	% MAXSEED2 + MINSEED2));
  else
    SetRandomSeeds((int)(MAX(seed1,seed2) +  theTime * TIMEFACTOR) 
							% MAXSEED1 + MINSEED1,
                   (int)(MIN(seed1,seed2)             	% MAXSEED2 + MINSEED2));
  return(NORMAL);
  }
------------------------------------------------------------------------
--

		J. Eliot B. Moss, Assistant Professor
		Department of Computer and Information Science
		Lederle Graduate Research Center
		University of Massachusetts
		Amherst, MA  01003
		(413) 545-4206, 545-1249 (fax); Moss at cs.umass.edu