Date: Fri, 06 Jun 1997 08:45:42 -0700
From: Bill Currie <billc AT blackmagic DOT tait DOT co DOT nz>
Subject: One method for small random nuimbers
To: djgpp AT delorie DOT com
Reply-to: billc AT blackmagic DOT tait DOT co DOT nz
Message-id: <339830A6.672@blackmagic.tait.co.nz>
Organization: Tait Electronics NZ
MIME-version: 1.0
Content-type: multipart/mixed; boundary="------------2D4C1B034A5E"
References: <2 DOT 2 DOT 32 DOT 19970604162434 DOT 0069d50c AT gate> <3395B097 DOT 5077 AT cornell DOT edu>
Precedence: bulk

This is a multi-part message in MIME format.

--------------2D4C1B034A5E
Content-Type: text/plain; charset=us-ascii
Content-Transfer-Encoding: 7bit

I noticed in the thread about random numbers (Random numbers/George)
that at one stage people were discussing how to get smaller ranges of
random numbers.  Now, I don't lay claims to my method being any better
than anybody elses, but it DOES work for me (it seemed to significantly
improve the randomness of my generated numbers).

I've attached it for those interested in it (as well as the random
number generator which I got from DDJ back in 91-93 (can't remeber) and
modified to produce 32 bit random numbers.

My method works by extracing the required number of bits to cover the
desired range (eg 3 bits for a range of 0-5) from a random number and
then taking the modulo of the desired range (6 in the example). Nothing
special so far, but what I do when it comes to the next time a number is
desired (any range, it probably only improves the randomness anyway) is
extract the needed bits from the SAME number, but ignoring the bits used
last time.  Basicaly, I treat the number received from the main RNG as a
bit stream and pull off the bits as needed.

Of course it isn't perfect, but it did improve things for me and I hope
that somebody else can find it useful.

Bill
-- 
Leave others their otherness.

--------------2D4C1B034A5E
Content-Type: text/plain; charset=us-ascii; name="RNUM.C"
Content-Transfer-Encoding: 7bit
Content-Disposition: inline; filename="RNUM.C"

#include <time.h>
#include <stdlib.h>
#include "r250.h"

unsigned long rnum(unsigned long max)
{
	static unsigned long bit_bank=0;
	static char bit_balance=-1;
	unsigned long mask=0,n=max;
	char bits=1,ob=0;

	if (bit_balance<0) {
		unsigned long l;
		time(&l);
		srand(l);
		r250_init((rand()<<16)|rand());
		bit_balance=0;
	}
	do {
		if (bits!=ob) {
			for (bits=mask=0,n=max; n; n>>=1) mask|=1<<bits++;
			ob=bits;
		}
		n=0;
		if (bit_balance<bits) {
			if (bit_balance) {
				char t=bits-bit_balance;
				n=(bit_bank&(mask>>t))<<t;
				mask>>=bit_balance;
				bits-=bit_balance;
			}
			bit_bank=r250();
			bit_balance=32;
		}
		n|=bit_bank&mask;
		bit_bank>>=bits;
		bit_balance-=bits;
	} while (n>max);
	return n;
}

--------------2D4C1B034A5E
Content-Type: text/plain; charset=us-ascii; name="R250.C"
Content-Transfer-Encoding: 7bit
Content-Disposition: inline; filename="R250.C"

/******************************************************************************
*  Module:  r250.c   Description: implements R250 random number generator,
*  from S. Kirkpatrick and E. Stoll, Journal of Computational Physics, 40,
*  p. 517 (1981). Written by:    W. L. Maier
******************************************************************************/

#include <stdlib.h>

/**** Satic variables ****/
static unsigned int r250_buffer[250];
static int r250_index;

/**** Function prototypes  ****/
void r250_init(int seed);
unsigned int r250();
double dr250();

/**** Function: r250_init  Description: initializes r250 random number
																generator. ****/

void r250_init(int seed)
{
/*---------------------------------------------------------------------------*/
	int		   j, k;
	unsigned int mask;
	unsigned int msb;
/*---------------------------------------------------------------------------*/
	srand(seed);
	r250_index=0;
	for (j = 0; j < 250; j++)	  /* Fill the r250 buffer with 32-bit values */
		r250_buffer[j] = (rand()<<16)|rand();
	msb = 0x80000000;		/* To turn on the diagonal bit	 */
	mask = 0xffffffff;		/* To turn off the leftmost bits */
	for (j = 0; j < 32; j++)
		{
		k = 11 * j + 3;				/* select a word to operate on		  */
		r250_buffer[k%250] &=mask;	/* Turn off bits left of the diagonal */
		r250_buffer[k%250] |= msb;	/* Turn on the diagonal bit			  */
		mask >>=1;
		msb >>=1;
		}
}
/**** Function:  r250  Description: returns a random unsigned integer. ****/
unsigned int r250()
{
/*---------------------------------------------------------------------------*/
	register int	j;
	register unsigned int new_rand;
/*---------------------------------------------------------------------------*/
	if (r250_index >= 146)
		j = r250_index - 147;	  /* Wrap pointer around */
	else
		j = r250_index +103;

	new_rand = r250_buffer[r250_index] ^ r250_buffer[j];
	r250_buffer[r250_index] = new_rand;
	if (r250_index >= 249)		/* Increment pointer for next time */
		r250_index=0;
	else
		r250_index++;

	return new_rand;
}
/**** Function:  r250  Description: returns a random double in range 0-1. ****/
double dr250()
{
/*---------------------------------------------------------------------------*/
	register int	j;
	register unsigned int new_rand;
/*---------------------------------------------------------------------------*/
	if (r250_index >= 146)
		j = r250_index - 147;	  /* Wrap pointer around */
	else
		j = r250_index +103;

	new_rand = r250_buffer[r250_index] ^ r250_buffer[j];
	r250_buffer[r250_index] = new_rand;
	if (r250_index >= 249)		/* Increment pointer for next time */
		r250_index=0;
	else
		r250_index++;

	return new_rand / (double)0xffffffff;/* Return a number in 0.0 to 1.0 */
}

--------------2D4C1B034A5E--