A Universal Random Number Generator.

[B. Narasimhan]

As I have seen quite a clamour for good random number generators I am enclosing 	
the code for the UNIVERSAL RANDOM NUMBER GENERATOR coded in Microsoft C. It is
not the best implementation, but I hope it is clear enough for the user to 
improve upon. For a reference see:

         "Toward a Universal Random Number Generator" 
                           by
          G. Marsaglia, Arif Zaman and Wai Wan Tsang 

to appear in Letters in Statistics and Probability sometime this year.
In my opinion, this is the best RNG existing today and has passed ALL
stringent tests put to it besides being completely portable from ETA10
down to CP/M machines.

B. Narasimhan       e-mail: naras at stat.fsu.edu on Internet.


        /*******************************************************
        * This program tests the Universal Random Number       *
        * generator.                                           *
        *                                                      *
        * SOFTWARE PROVIDED BY                                 *
        *  Dr. G. Marsaglia and B. Narasimhan                  *
        *      Dept. of Statistics, Florida State University,  *
        *      Tallahassee, FL 32306.                          *
        *      Ph. (904)-644-3218.                             *
        *******************************************************/

        /*----------------------------------------------------*/
        /* Application Notes:                                 */
        /* This program has been compiled and tested with     */
        /* Microsoft C 5.0. Hopefully, it is portable. Note   */
        /* the code can be made more efficient in easy-to-see */
        /* spots. We have avoided these obvious changes in    */
        /* interest of readability.                           */
        /*----------------------------------------------------*/

#include <stdio.h>
#include <stdlib.h>
#include <math.h>

#define TAB_SIZE 97
#define NO_OF_BITS 24
#define MOD_1 179
#define MOD_2 169
#define CONST_1 53
#define CONST_2 362436
#define CONST_3 7654321
#define CONST_4 16777213
#define CONST_5 96
#define CONST_6 32
#define CONST_7 16777216.0
#define MASK 0xffffff

unsigned long UNI_TABL [TAB_SIZE];
signed long C_VAL;
short I_IND, J_IND;

SET_UNI(a, b, c, d)
int a, b, c, d;
{
    int i, j;
    ldiv_t div_result;
    unsigned long sum;

    for (i = 0; i < TAB_SIZE ; i++) {
	sum = 0;
	for (j = NO_OF_BITS; j > 0; j--) {
	    div_result = ldiv(a*b,MOD_1);
	    div_result = ldiv(div_result.rem * c,MOD_1);
	    a = b;
	    b = c;
	    c = div_result.rem;
	    div_result = ldiv(CONST_1 * d + 1, MOD_2);
	    d = div_result.rem;
	    sum *= 2;
	    div_result = ldiv(d*c,64);
	    if (div_result.rem >= 32)
		sum += 1;
	}
	UNI_TABL[i] = sum;
    }
    C_VAL = CONST_2;
    I_IND = CONST_5;
    J_IND = CONST_6;
}
float UNI(void)
{
    unsigned long uni_val;

    uni_val = UNI_TABL[I_IND] - UNI_TABL[J_IND];
    uni_val &= MASK;
    UNI_TABL[I_IND] = uni_val;

    if (I_IND == 0)
	I_IND = CONST_5;
    else
	I_IND--;

    if (J_IND == 0)
	J_IND = CONST_5;
    else
	J_IND--;

    C_VAL -= CONST_3;
    if (C_VAL < 0)
	C_VAL += CONST_4;
    uni_val -= C_VAL;
    uni_val &= MASK;
    return (float) uni_val / CONST_7;
}

main()
{
    int i;
    float temp;
    unsigned long stemp;

    SET_UNI(12,34,56,78);

    for (i = 20000; i > 0; i--)
	temp = UNI();            /* Generate 20000 unis. */

    for (i = 6; i > 0; i--) {
	stemp = UNI() * CONST_7;
	printf("\t%lu\n",stemp);     /* Print next 6 of them as integers. */
   }
}