Account Options

  1. Sign in
The old Google Groups will be going away soon, but your browser is incompatible with the new version.
Google Groups Home
« Groups Home
Message from discussion 64-bit KISS RNGs
The group you are posting to is a Usenet group. Messages posted to this group will make your email address visible to anyone on the Internet.
Your reply message has not been sent.
Your post was successful
 
From:
To:
Cc:
Followup To:
Add Cc | Add Followup-to | Edit Subject
Subject:
Validation:
For verification purposes please type the characters you see in the picture below or the numbers you hear by clicking the accessibility icon. Listen and type the numbers you hear
 
geo  
View profile  
 More options Feb 28 2009, 7:30 am
Newsgroups: sci.math, comp.lang.c, comp.lang.fortran
From: geo <gmarsag...@gmail.com>
Date: Sat, 28 Feb 2009 04:30:48 -0800 (PST)
Local: Sat, Feb 28 2009 7:30 am
Subject: 64-bit KISS RNGs
64-bit KISS RNGs

Consistent with the Keep It Simple Stupid (KISS) principle,
I have previously suggested 32-bit KISS Random Number
Generators (RNGs) that seem to have been frequently adopted.

Having had requests for 64-bit KISSes, and now that
64-bit integers are becoming more available, I will
describe here a 64-bit KISS RNG, with comments on
implementation for various languages, speed, periods
and performance after extensive tests of randomness.

This 64-bit KISS RNG has three components, each nearly
good enough to serve alone.    The components are:
Multiply-With-Carry (MWC), period (2^121+2^63-1)
Xorshift (XSH), period 2^64-1
Congruential (CNG), period 2^64

Compact C and Fortran listings are given below.  They
can be cut, pasted, compiled and run to see if, after
100 million calls, results agree with that provided
by theory, assuming the default seeds.

Users may want to put the content in other forms, and,
for general use, provide means to set the 250 seed bits
required in the variables x,y,z (64 bits) and c (58 bits)
that have been given default values in the test versions.

The C version uses #define macros to enumerate the few
instructions that MWC, XSH and CNG require.   The KISS
macro adds MWC+XSH+CNG mod 2^64, so that KISS can be
inserted at any place in a C program where a random 64-bit
integer is required.
Fortran's requirement that integers be signed makes the
necessary code more complicated, hence a function KISS().

C version; test by invoking macro KISS 100 million times
-----------------------------------------------------------------
#include <stdio.h>

static unsigned long long
x=1234567890987654321ULL,c=123456123456123456ULL,
y=362436362436362436ULL,z=1066149217761810ULL,t;

#define MWC (t=(x<<58)+c, c=(x>>6), x+=t, c+=(x<t), x)
#define XSH ( y^=(y<<13), y^=(y>>17), y^=(y<<43) )
#define CNG ( z=6906969069LL*z+1234567 )
#define KISS (MWC+XSH+CNG)

 int main(void)
 {int i;
  for(i=0;i<100000000;i++) t=KISS;
  (t==1666297717051644203ULL) ?
  printf("100 million uses of KISS OK"):
  printf("Fail");
 }

---------------------------------------------------------------
Fortran version; test by calling KISS() 100 million times
---------------------------------------------------------------
 program testkiss
 implicit integer*8(a-z)
 do i=1,100000000; t=KISS(); end do
 if(t.eq.1666297717051644203_8) then
 print*,"100 million calls to KISS() OK"
 else; print*,"Fail"
 end if; end

 function KISS()
 implicit integer*8(a-z)
 data x,y,z,c /1234567890987654321_8, 362436362436362436_8,&
 1066149217761810_8, 123456123456123456_8/
 save x,y,z,c
  m(x,k)=ieor(x,ishft(x,k))  !statement function
  s(x)=ishft(x,-63)             !statement function
 t=ishft(x,58)+c
 if(s(x).eq.s(t)) then; c=ishft(x,-6)+s(x)
      else; c=ishft(x,-6)+1-s(x+t); endif
 x=t+x
 y=m(m(m(y,13_8),-17_8),43_8)
 z=6906969069_8*z+1234567
 KISS=x+y+z
 return; end
---------------------------------------------------------------

Output from using the macro KISS or the function KISS() is
MWC+XSH+CNG mod 2^64.

CNG is easily implemented on machines with 64-bit integers,
as arithmetic is automatically mod 2^64, whether integers
are considered signed or unsigned.   The CNG statement is
    z=6906969069*z+1234567.
When I established the lattice structure of congruential
generators in the 60's, a search produced 69069 as an easy-
to-remember multiplier with nearly cubic lattices in 2,3,4,5-
space, so I tried concatenating, using 6906969069 as
my first test multiplier. Remarkably---a seemingly one in many
hundreds chance---it turned out to also have excellent lattice
structure in 2,3,4,5-space, so that's the one chosen.
(I doubt if lattice structure of CNG has much influence on the
composite 64-bit KISS produced via MWC+XSH+CNG mod 2^64.)

XSH, the Xorshift component, described in
     www.jstatsoft.org/v08/i14/paper
uses three invocations of an integer "xor"ed with a shifted
version of itself.
The XSH component used for this KISS is, in C notation:
   y^=(y<<13); y^=(y>>17); y^=(y<<43)
with Fortran equivalents y=ieor(y,ishft(y,13)), etc., although
this can be effected by a Fortran statement function:
     f(y,k)=ieor(y,ishft(y,k))
     y=f(f(f(y,13),-17),43)
As with lattice structure, choice of the triple 13,-17,43 is
probably of no particular importance; any one of the 275 full-
period triples listed in the above article is likely to provide
a  satisfactory component XSH for the composite MWC+XSH+CNG.

The choice of multiplier 'a' for the multiply-with-carry (MWC)
component of KISS is not so easily made.  In effect, a multiply-
with-carry sequence has a current value x and current "carry" c,
and from each given x,c a new x,c pair is constructed by forming
  t=a*x+c, then x=t mod b=2^64 and c=floor(t/b).
This is easily implemented for 32-bit computers that permit
forming a*t+c in 64 bits, from which the new x is the bottom and
the new c the top 32-bits.

When a,x and c are 64-bits, not many computers seem to have an easy
way to form t=a*x+c in 128 bits, then extract the top and bottom
64-bit segments.  For that reason, special choices for 'a' are
needed among those that satisfy the general requirement that
p=a*b-1 is a prime for which b=2^64 has order (p-1)/2.

My choice---and the only one of this form---is a=2^58+1. Then the
top 64 bits of an imagined 128-bit t=a*x+c may be obtained as
(using C notation)  (x>>6)+ 1 or 0, depending
on whether the 64-bit parts of (x<<58)+c+x cause an overflow.
Since (x<<58)+c cannot itself cause overflow (c will always be <a),
we get the carry as c=(x>>6) plus overflow from (x<<58)+x.

This is easily done in C with unsigned integers, using  a different
kind of 't':   t=(x<<58)+c; c=(x>>6); x=t+x; c=c+(x<t);
For Fortran and others that permit only signed integers, more work
is needed.
Equivalent mod 2^64 versions of t=(x<<58)+c and c=(x>>6) are easy,
and if s(x) represents (x>>63) in C or ishft(x,-66) in Fortran,
then for signed integers, the new carry c comes from the rule
 if s(x) equals s(t) then c=(x>>6)+s(x) else c=(x>>6)+1-s(x+t)

Speed:
A C version of this KISS RNG takes 18 nanosecs for each
64-bit random number on my desktop (Vista) PC, thus
producing KISSes at a rate exceeding 55 million per second.
Fortran or other integers-must-be-signed compilers might get
"only" around 40 million per second.

Setting seeds:
Use of KISS or KISS() as a general 64-bit RNG requires specifying
3*64+58=250 bits for seeds, 64 bits each for x,y,z and 58 for c,
resulting in a composite sequence with period around 2^250.
The actual period is
  (2^250+2^192+2^64-2^186-2^129)/6 ~= 2^(247.42)  or 10^(74.48).
We "lose" 1+1.58=2.58 bits from maximum possible period, one bit
because b=2^64, a square, cannot be a primitive root of p=ab-1,
so the best possible order for b is (p-1)/2.
The periods of MWC and  XSH have gcd 3=2^1.58, so another  1.58
bits  are "lost" from the best possible period we could expect
from 250 seed bits.

Some users may think 250 seed bits are an unreasonable requirement.
A good seeding  procedure might be to assume the default seed
values then let the user choose none, one, two,..., or all
of x,y,z, and c to be reseeded.

Tests:
Latest tests in The Diehard Battery, available at
   http://i.cs.hku.hk/~diehard/
were applied extensively. Those tests that specifically required
32-bit integers were applied to the leftmost 32 bits
(e,g, KISS>>32;), then to the middle 32-bits ((KISS<<16)>>32;)
then to the rightmost 32 bits, ( (KISS<<32)>>32).
There were no extremes in the more than 700 p-values returned
by the tests, nor indeed for similar tests applied to just two of the
KISS components: MWC+XSH, then MWC+CNG, then XSH+CNG.

The simplicity, speed, period around 2^250 and performance on
tests of randomness---as well as ability to produce exactly
the same 64-bit patterns, whether considered signed or unsigned
integers---make this 64-bit KISS well worth considering for
adoption or adaption to languages other than C or Fortran,
as has been done for 32-bit KISSes.

George Marsaglia


 
You must Sign in before you can post messages.
To post a message you must first join this group.
Please update your nickname on the subscription settings page before posting.
You do not have the permission required to post.