G

#### geo

diverse comments over various newsgroups that I

thought it better to use a separate posting to

bring suggested changes to interested potential

users via sci.math, comp.lang.c, comp.lang.fortran.

At least two readers have pointed out special cases that

might arise from the code I suggested for the KISS4691.

christain.bau pointed out that (x<<13)+c

can exceed 2^32 when the rightmost 19 bits of x are 1's

and c is its maximum possible, a-1=8192=2^13.

Changing

c=(t<x)+(x>>19); to c=(t<=x)+(x>>19); fixes that.

But io_x points out that this overlooks the case x=c=0,

for which c=(t<=x)+(x>>19) would incorrectly make c=1

rather than 0.

This is perhaps best handled by merely, when x=0,

setting Q[j]=c; c=0; return Q[j];

Below is a revised listing of the C version;

changes should be transparent for versions

kindly provided by interested viewers---in Fortran,

PL/1 and others.

While neither correction is likely to have an

appreciable effect on use of KISS4691 as an RNG,

both corrections are necessary to ensure that

the period of MWC() will be the value (> 10^45191)

called for by the underlying mathematics of MWC,

(assuming the two period-1 seed sets,

c=0 and each Q

*=0*

c=8192 and each Q

c=8192 and each Q

*=2^32-1*

are not used).

And the results after 10^9 calls to MWC() then

10^9 call to KISS as MWC()+CNG+XS,

confirmed in various versions, are not changed.

George Marsaglia

----------------------------------------------------------

static unsigned long xs=521288629,xcng=362436069,Q[4691];

unsigned long MWC(void) /*238 million/second*/

{unsigned long t,x; static c=0,j=4691;

j=(j<4690)? j+1:0;

x=Q[j]; if(x==0) {Q[j]=c; c=0; return Q[j];}

t=(x<<13)+c+x; c=(t<=x)+(x>>19);

return (Q[j]=t);

}

#define CNG ( xcng=69069*xcng+123 )

#define XS ( xs^=(xs<<13), xs^=(xs>>17), xs^=(xs<<5) )

#define KISS ( MWC()+CNG+XS ) /*138 million/sec*/

#include <stdio.h>

int main()

{unsigned long i,x;

for(i=0;i<4691;i++) Qare not used).

And the results after 10^9 calls to MWC() then

10^9 call to KISS as MWC()+CNG+XS,

confirmed in various versions, are not changed.

George Marsaglia

----------------------------------------------------------

static unsigned long xs=521288629,xcng=362436069,Q[4691];

unsigned long MWC(void) /*238 million/second*/

{unsigned long t,x; static c=0,j=4691;

j=(j<4690)? j+1:0;

x=Q[j]; if(x==0) {Q[j]=c; c=0; return Q[j];}

t=(x<<13)+c+x; c=(t<=x)+(x>>19);

return (Q[j]=t);

}

#define CNG ( xcng=69069*xcng+123 )

#define XS ( xs^=(xs<<13), xs^=(xs>>17), xs^=(xs<<5) )

#define KISS ( MWC()+CNG+XS ) /*138 million/sec*/

#include <stdio.h>

int main()

{unsigned long i,x;

for(i=0;i<4691;i++) Q

*=CNG+XS;*

for(i=0;i<1000000000;i++) x=MWC();

printf(" MWC result=3740121002 ?\n%22u\n",x);

for(i=0;i<1000000000;i++) x=KISS;

printf("KISS result=2224631993 ?\n%22u\n",x);

}

---------------------------------------------------------for(i=0;i<1000000000;i++) x=MWC();

printf(" MWC result=3740121002 ?\n%22u\n",x);

for(i=0;i<1000000000;i++) x=KISS;

printf("KISS result=2224631993 ?\n%22u\n",x);

}

---------------------------------------------------------