Search All of the Math Forum:

Views expressed in these public forums are not endorsed by NCTM or The Math Forum.

Notice: We are no longer accepting new posts, but the forums will continue to be readable.

Topic: RNGs: A Super KISS
Replies: 14   Last Post: Feb 26, 2013 10:53 AM

 Messages: [ Previous | Next ]
 rossum Posts: 731 Registered: 12/13/04
Re: RNGs: A Super KISS
Posted: Feb 26, 2013 8:06 AM

On Mon, 25 Feb 2013 20:10:54 -0800 (PST), latexpst@gmail.com wrote:

>On Wednesday, November 4, 2009 2:46:27 AM UTC+11, geo wrote:
>> /*
>> For those mesmerized (or Mersenne-ized?) by a RNG
>> with period 2^19937-1, I offer one here with period
>> 54767*2^1337279---over 10^396564 times as long.
>> It is one of my CMWC (Complimentary-Multiply-With-Carry) RNGs,
>> and is suggested here as one of the components of a
>> super-long-period KISS (Keep-It-Simple-Stupid) RNG.
>>
>> With b=2^32 and a=7010176, and given a 32-bit x, and a 32-bit c, this
>> generator produces a new x,c by forming 64-bit t=a*x+c then replacing:
>> c=top 32 bits of t and x=(b-1)-(bottom 32 bits of t). In C: c=t>>32;
>> x=~t;
>>
>> For many years, CPUs have had machine instructions to form such a
>> 64-bit t and extract the top and bottom halves, but unfortunately
>> only recent Fortran versions have means to easily invoke them.
>>
>> Ability to do those extractions leads to implementations that are
>> simple
>> and extremely fast---some 140 million per second on my desktop PC.
>>
>> Used alone, this generator passes all the Diehard Battery of Tests,
>> but
>> its simplicity makes it well-suited to serve as one of the three
>> components
>> of a KISS RNG, based on the Keep-It-Simple-Stupid principle, and the
>> idea,
>> supported by both theory and practice, that the combination of RNGs
>> based on
>> different mathematical models can be no worse---and is usually
>> better---than
>> any of the components.
>>
>> So here is a complete C version of what might be called a SUPER KISS
>> RNG,
>> combining, by addition mod 2^32, a Congruential RNG, a Xorshift RNG
>> and the super-long-period CMWC RNG:
>> */
>>
>> #include <stdio.h>
>> static unsigned long Q
>> [41790],indx=41790,carry=362436,xcng=1236789,xs=521288629;
>>
>> #define CNG ( xcng=69609*xcng+123 ) /*Congruential*/
>> #define XS ( xs^=xs<<13, xs^=(unsigned)xs>>17, xs^=xs>>5 ) /
>> *Xorshift*/
>> #define SUPR ( indx<41790 ? Q[indx++] : refill() )
>> #define KISS SUPR+CNG+XS
>>
>> int refill( )
>> { int i; unsigned long long t;
>> for(i=0;i<41790;i++) { t=7010176LL*Q[i]+carry; carry=(t>>32); Q[i]=~
>> (t);}
>> indx=1; return (Q[0]);
>> }
>>
>> int main()
>> {unsigned long i,x;
>> for(i=0;i<41790;i++) Q[i]=CNG+XS;
>> for(i=0;i<1000000000;i++) x=KISS;
>> printf(" x=%d.\nDoes x=-872412446?\n",x);
>> }
>>
>> /*
>> Running this program should produce 10^9 KISSes in some 7-15 seconds.
>> You are invited to cut, paste, compile and run for yourself, checking
>> to
>> see if the last value is as designated, (formatted as a signed integer
>> for
>> potential comparisons with systems using signed integers).
>> You may want to report or comment on implementations for other
>> languages.
>>
>> The arithmetic operations are suited for either signed or unsigned
>> integers.
>> Thus, with (64-bit)t=a*x+c, x=t%b in C or x=mod(t,b) in Fortran, and
>> c=c/b in either C or Fortran, but with ways to avoid integer
>> divisions,
>> and subsequent replacement of x by its base-b complement, ~x in C.
>>
>> With b=2^32 and p=54767*2^1337287+1, the SUPR part of this Super KISS
>> uses my CMWC method to produce, in reverse order, the base-b expansion
>> of k/p for some k determined by the values used to seed the Q array.
>> The period is the order of b for that prime p:
>> 54767*2^1337279, about 2^1337294 or 10^402566.
>> (It took a continuous run of 24+ days on an earlier PC to
>> establish that order. My thanks to the wizards behind PFGW
>> and to Phil Carmody for some suggested code.)
>>
>> Even the Q's all zero, should seeding be overlooked in main(),
>> will still produce a sequence of the required period, but will
>> put the user in a strange and exceedingly rare place in the entire
>> sequence. Users should choose a reasonable number of the 1337280
>> random bits that a fully-seeded Q array requires.
>>
>> Using your own choices of merely 87 seed bits, 32 each for xcng,xs
>> and 23 for carry<7010176, then initializing the Q array with
>> for(i=0;i<41790;i++) Q[i]=CNG+XS;
>> should serve well for many applications, but others, such as in
>> Law or Gaming, where a minimum number of possible outcomes may be
>> required, might need more of the 1337280 seed bits for the Q array.
>>
>> As might applications in cryptography: With an unknown but fully-
>> seeded Q array, a particular string of, say, 41000 successive SUPR
>> values will appear at more than 2^20000 locations in the full
>> sequence,
>> making it virtually impossible to get the location of that particular
>> string in the full loop, and thus predict coming or earlier values,
>> even if able to undo the CNG+XS operations.
>> */
>>
>> /*
>> So I again invite you to cut, paste, compile and run the above C
>> program.
>> 1000 million KISSes should be generated, and the specified result
>> appear,
>> by the time you count slowly to fifteen.
>> (Without an optimizing compiler, you may have to count more slowly.)
>> */
>>
>> /* George Marsaglia */

>
>What's the best way to generate double precision numbers between 0 and 1, excluding 0, from this generator? Thanks.

You may have to wait some time for a reply. Unfortunately George died
in 2011: http://en.wikipedia.org/wiki/George_Marsaglia

rossum

Date Subject Author
11/3/09 geo
11/3/09 Dann Corbit
11/3/09 tom@iahu.ca
11/3/09 Dann Corbit