v05i009: Minimal Standard Random Number Generator

[Fuat C. Baran]

Posting-number: Volume 5, Issue 9
Submitted-by: "Fuat C. Baran" <fuat at cunixc.cc.columbia.edu>
Archive-name: random

[Conflict!!!  BSD already has a random() which is better than rand().  ++bsa]

I don't know if this is really worth posting (I just converted
pseudo-code to C from the latest issue of CACM), but in case other
people need this I thought I'd make it available.  It implements a
"minimal standard random number generator" discussed in CACM.  It is
meant as a replacement for rand(3) that comes with 4.x BSD which is
said to be a poor generator according to the article.  It should work
on any system in which the maximum integer is 2^31-1 (when compiled
with -DTEST it checks to see whether or not it was implemented
correctly).

					--Fuat

#! /bin/sh
# This is a shell archive, meaning:
# 1. Remove everything above the #! /bin/sh line.
# 2. Save the resulting text in a file.
# 3. Execute the file with /bin/sh (not csh) to create:
#	random.c
# This archive created: Mon Oct 24 17:27:26 1988
export PATH; PATH=/bin:/usr/bin:$PATH
if test -f 'random.c'
then
	echo shar: "will not over-write existing file 'random.c'"
else
sed 's/^X//' << \SHAR_EOF > 'random.c'
X/**
X ** random.c: Minimal Standard Pseudo-Random Number Generator
X **
X ** Author: Fuat C. Baran, Columbia University, 1988
X **
X ** Based on code in "Random Number Generators: Good Ones are Hard to Find",
X ** by Stephen K. Park and Keith W. Miller in Communications of the ACM,
X ** 31, 10 (Oct. 1988) pp. 1192-1201.
X **
X ** Requirements: maxint must be 2^31 -1 or larger.
X **
X ** Compile with -DTEST and run to see if it was implemented correctly.
X **/
X
X/* some constants we need */
X#define A 16807
X#define M 2147483647			/* Mersenne prime 2^31 -1 */
X#define Q 127773			/* M div A (M / A) */
X#define R 2836				/* M mod A (M % A) */
X
Xstatic long seed = 1;			/* initialize with rand_init() */
X
X
X/*
X * random:
X * minimal standard random number generator.
X * call rand_init first to initialize seed.
X * returns a random float.
X */
X
Xfloat
Xrandom () {
X    long hi, lo;
X    
X    hi = seed / Q;
X    lo = seed % Q;
X    if ((seed = A * lo - R * hi) <= 0)
X        seed += M;
X    return ((float) seed / M);
X}
X
X
X/*
X * rand_init:
X * initialize seed.
X */
X
Xvoid
Xrand_init (s)
Xlong s;
X{
X    seed = s;
X}
X
X
X#ifdef TEST
X
Xmain () {
X    int i;
X    float r;
X
X    rand_init (1);
X    for (i = 0; i < 10000; i++)
X	r = random();
X    printf ("Seed is %ld (should be 1043618065)\n", seed);
X}
X
X#else
X
Xmain () {
X    int i;
X
X    for (i = 0; i < 1000000; i++)
X	random();
X}
X
X#endif /* TEST */
X
X
SHAR_EOF
fi
exit 0
#	End of shell archive



-- 
INTERNET: fuat at columbia.edu          U.S. MAIL: Columbia University
BITNET:   fuat at cunixc.cc.columbia.edu           Center for Computing Activities
USENET:   ...!rutgers!columbia!cunixc!fuat      712 Watson Labs, 612 W115th St.
PHONE:    (212) 280-5128                        New York, NY 10025