• Advertisement
Sign in to follow this  

"R250'' algorithm of Kirkpatrick-Stoll (random number)

This topic is 3891 days old which is more than the 365 day threshold we allow for new replies. Please post a new topic.

If you intended to correct an error in the post then please contact us.

Recommended Posts

R250 algorithm of Kirkpatrick-Stoll (random number generator) http://homepage.univie.ac.at/Franz.Vesely/cp0102/dx/node37.html Could someone explain the above, (dumbed * 2 + " " + down) with simple short code.
// net code, this code works fine, and i believe its based of the 
// same algorithm
#include "stdafx.h"
#include <iostream>
#include <stdlib.h>



unsigned int r250();
unsigned int r250n(unsigned n);
double dr250();


/**** Static variables ****/
static int r250_index = 0;
static unsigned int r250_buffer[250] = {
	15301,57764,10921,56345,19316,43154,54727,49252,32360,49582,
	26124,25833,34404,11030,26232,13965,16051,63635,55860,5184,
	15931,39782,16845,11371,38624,10328,9139,1684,48668,59388,
	13297,1364,56028,15687,63279,27771,5277,44628,31973,46977,
	16327,23408,36065,52272,33610,61549,58364,3472,21367,56357,
	56345,54035,7712,55884,39774,10241,50164,47995,1718,46887,
	47892,6010,29575,54972,30458,21966,54449,10387,4492,644,
	57031,41607,61820,54588,40849,54052,59875,43128,50370,44691,
	286,12071,3574,61384,15592,45677,9711,23022,35256,45493,
	48913,146,9053,5881,36635,43280,53464,8529,34344,64955,
	38266,12730,101,16208,12607,58921,22036,8221,31337,11984,
	20290,26734,19552,48,31940,43448,34762,53344,60664,12809,
	57318,17436,44730,19375,30,17425,14117,5416,23853,55783,
	57995,32074,26526,2192,11447,11,53446,35152,64610,64883,
	26899,25357,7667,3577,39414,51161,4,58427,57342,58557,
	53233,1066,29237,36808,19370,17493,37568,3,61468,38876,
	17586,64937,21716,56472,58160,44955,55221,63880,1,32200,
	62066,22911,24090,10438,40783,36364,14999,2489,43284,9898,
	39612,9245,593,34857,41054,30162,65497,53340,27209,45417,
	37497,4612,58397,52910,56313,62716,22377,40310,15190,34471,
	64005,18090,11326,50839,62901,59284,5580,15231,9467,13161,
	58500,7259,317,50968,2962,23006,32280,6994,18751,5148,
	52739,49370,51892,18552,52264,54031,2804,17360,1919,19639,
	2323,9448,43821,11022,45500,31509,49180,35598,38883,19754,
	987,11521,55494,38056,20664,2629,50986,31009,54043,59743
	};

static unsigned myrand();
static void mysrand(unsigned newseed);

/**** Function: r250_init  
	Description: initializes r250 random number generator. ****/
void r250_init(int seed)
{
/*---------------------------------------------------------------------------*/
    int        j, k;
    unsigned int mask;
    unsigned int msb;
/*---------------------------------------------------------------------------*/
    mysrand(seed);
    r250_index = 0;
    for (j = 0; j < 250; j++)     /* Fill the r250 buffer with 15-bit values */
        r250_buffer[j] = myrand();
    for (j = 0; j < 250; j++)     /* Set some of the MS bits to 1 */
        if (myrand() > 16384)
            r250_buffer[j] |= 0x8000;
    msb = 0x8000;       /* To turn on the diagonal bit   */
    mask = 0xffff;      /* To turn off the leftmost bits */
    for (j = 0; j < 16; j++)
        {
        k = 11 * j + 3;             /* Select a word to operate on        */
        r250_buffer[k] &= mask;     /* Turn off bits left of the diagonal */
        r250_buffer[k] |= msb;      /* Turn on the diagonal bit           */
        mask >>= 1;
        msb >>= 1;
        }
}

/**** Function: r250 Description: returns a random unsigned integer k
				uniformly distributed in the interval 0 <= k < 65536.  ****/
unsigned int r250()
{
/*---------------------------------------------------------------------------*/
    register int    j;
    register unsigned int new_rand;
/*---------------------------------------------------------------------------*/
    if (r250_index >= 147)
        j = r250_index - 147;      /* Wrap pointer around */
    else
        j = r250_index + 103;

    new_rand = r250_buffer[r250_index] ^= r250_buffer[j];

    if (r250_index >= 249)      /* Increment pointer for next time */
        r250_index = 0;
    else
        r250_index++;

    return new_rand;
}

/**** Function: r250n Description: returns a random unsigned integer k
					uniformly distributed in the interval 0 <= k < n ****/
unsigned int r250n(unsigned n)
{
/*---------------------------------------------------------------------------*/
    register int    j;
    register unsigned int new_rand, limit;
/*---------------------------------------------------------------------------*/
	limit = (65535U/n)*n;
	do 
        {new_rand = r250();
        if (r250_index >= 147)
            j = r250_index - 147;      /* Wrap pointer around */
        else
            j = r250_index + 103;

        new_rand = r250_buffer[r250_index] ^= r250_buffer[j];

        if (r250_index >= 249)      /* Increment pointer for next time */
            r250_index = 0;
        else
            r250_index++;
        } while(new_rand >= limit);
    return new_rand%n;
}

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

    new_rand = r250_buffer[r250_index] ^= r250_buffer[j];

    if (r250_index >= 249)      /* Increment pointer for next time */
        r250_index = 0;
    else
        r250_index++;
    return new_rand / 65536.;   /* Return a number in [0.0 to 1.0) */
}


/***  linear congruent pseudorandom number generator for initialization ***/

static unsigned long seed=1;

	/*	Return a pseudorandom number in the interval 0 <= n < 32768.
		This produces the following sequence of pseudorandom 
		numbers:
		346, 130, 10982, 1090...  (9996 numbers skipped) ...23369,
		2020, 5703, 12762, 10828, 16252, 28648, 27041, 23444, 6604... 
	*/ 
static unsigned myrand()
{
	seed = seed*0x015a4e35L + 1;
	return (seed>>16)&0x7fff;
}

	/*	Initialize above linear congruent pseudorandom number generator */
static void mysrand(unsigned newseed)
{	seed = newseed;
}


int _tmain(int argc, _TCHAR* argv[])
{
r250_init(256);
unsigned int X;
X = myrand();
int i = 0;
while(i++ < 10){
		std::cout << "Double "<< dr250() << std::endl;
		std::cout << "Int "<< r250() << std::endl;
		std::cout << "X " << X << std::endl;
		std::cout << "Unsigned " << r250n(X) << std::endl;
}

	return 0;
}


Share this post


Link to post
Share on other sites
Advertisement
What exactly about that material are you not comprehending? As is, it is written as simple as it can get without being repetitive and over simplified. The only difference between the code you posted and the algorithm as presented by the link is the code is generating multiple bits per iteration.

Share this post


Link to post
Share on other sites
Sign in to follow this  

  • Advertisement