Re: Ensuring the long period of the KISS4691 RNG.

Re: Ensuring the long period of the KISS4691 RNG.

[robin]

"geo" <gmarsaglia@gmail.com> wrote in message news:08db9776-d933-44a5-a20f-7a137e81cce5@5g2000yqz.googlegroups.com...
| My posting for KISS4691 has become so mixed with
| 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.

Below is my PL/I version, using unsigned arithmetic,
which I confirm delivers the same test results as before.

(NOSIZE, NOFOFL):
RNG: PROCEDURE OPTIONS (MAIN, REORDER);

   declare (xs initial (521288629), xcng initial (362436069),
            Q(0:4690) ) static fixed binary (32) unsigned;

MWC: procedure () returns (fixed binary (32) unsigned);
   declare (t,x,i) fixed binary (32) unsigned;
   declare (c initial (0), j initial (4691) ) fixed binary (32) unsigned static;

   if j < hbound(Q,1) then j = j + 1; else j = 0;
   x = Q(j); if x = 0 then do; Q(j) = c; c = 0; return(Q(j)); end;
   t = isll(x,13)+c+x; c = (t<=x) + isrl(x,19);
   Q(j)=t;
   return (t);
end MWC;

CNG: procedure returns (fixed binary (32) unsigned);
  xcng=bin(69069)*xcng+bin(123);
  return (xcng);
end CNG;

XXS: procedure returns (fixed binary (32) unsigned);
   xs = ieor (xs, isll(xs, 13) );
   xs = ieor (xs, isrl(xs, 17) );
   xs = ieor (xs, isll(xs,  5) );
   return (xs);
end XXS;

KISS: procedure returns (fixed binary (32) unsigned);
   return ( MWC()+CNG+XXS );
end KISS;

   declare (i,x) fixed binary (32) unsigned;

   Q = CNG + XXS; /* Initialize: */
   do i = 1 to 1000000000; x=MWC(); end;
   put skip edit (" MWC result=3740121002",x) (a, skip, f(22));
   do i = 1 to 1000000000; x=KISS; end;
   put skip edit ("KISS result=2224631993",x) (a, skip, f(22));

end RNG;

Re: Ensuring the long period of the KISS4691 RNG.

[robin]

"geo" <gmarsaglia@gmail.com> wrote in message news:08db9776-d933-44a5-a20f-7a137e81cce5@5g2000yqz.googlegroups.com...

| 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.

Below is my PL/I version, using signed arithmetic,
which I confirm delivers the same test results as previously,
and which are the same results as the signed version.

(NOSIZE, NOFOFL):
RNG: PROCEDURE OPTIONS (MAIN, REORDER);

   declare (xs initial (521288629), xcng initial (362436069),
            Q(0:4690) ) static fixed binary (31);

MWC: procedure () returns (fixed binary (31));
   declare (t,x,i) fixed binary (31);
   declare (c initial (0), j initial (4691) ) fixed binary (31) static;
   declare t1 fixed binary (31);

   if j < hbound(Q,1) then j = j + 1; else j = 0;
   x = Q(j); if x = 0 then do; Q(j) = c; c = 0; return(Q(j)); end;
   t = isll(x,13)+c+x;
   /* Tests x >= t by computing x - t in two parts. */
   t1 = isrl(x, 2) - isrl(t, 2);
   if t1  = 0 then t1 = iand(x, 3) - iand(t, 3);
   if t1 >= 0 then t1 = 1; else t1 = 0;
   c = t1 + isrl(x, 19);
   Q(j)=t;
   return (t);
end MWC;

CNG: procedure returns (fixed binary (31));
  xcng=bin(69069)*xcng+bin(123);
  return (xcng);
end CNG;

XXS: procedure returns (fixed binary (31));
   xs = ieor (xs, isll(xs, 13) );
   xs = ieor (xs, isrl(xs, 17) );
   xs = ieor (xs, isll(xs,  5) );
   return (xs);
end XXS;

KISS: procedure returns (fixed binary (31));
   return ( MWC()+CNG+XXS );
end KISS;

   declare (i,x) fixed binary (31);
   declare y fixed decimal (11);

   Q = CNG+XXS; /* Initialize. */
   do i = 1 to 1000000000; x=MWC(); end;
   put skip edit (" Expected MWC result = 3740121002", 'computed =', x)
      (a, skip, x(12), a, f(11));
   y = iand(x, 2147483647);
   if x < 0 then y = y + 2147483648;
   put skip edit (y) (x(11), f(22)); put skip;
   do i = 1 to 1000000000; x=KISS; end;
   put skip edit ("Expected KISS result = 2224631993", 'computed =', x)
      (a, skip, x(12), a, f(11));
   y = iand(x, 2147483647);
   if x < 0 then y = y + 2147483648;
   put skip edit (y) (x(11), f(22));

end RNG; 

Re: Ensuring the long period of the KISS4691 RNG.

[robin]

"robin" <robin51@dodo.com.au> wrote in message news:4c722ae5$0$1261$c30e37c6@exi-reader.telstra.net...
| "geo" <gmarsaglia@gmail.com> wrote in message news:08db9776-d933-44a5-a20f-7a137e81cce5@5g2000yqz.googlegroups.com...
|
|| 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.
|
| Below is my PL/I version, using signed arithmetic,
| which I confirm delivers the same test results as previously,
| and which are the same results as the signed version.

um...unsigned version. 

Re: Ensuring the long period of the KISS4691 RNG.

[io_x]

"robin" <robin51@dodo.com.au> ha scritto nel messaggio 
news:4c722e84$0$1261$c30e37c6@exi-reader.telstra.net...
> "robin" <robin51@dodo.com.au> wrote in message 
> news:4c722ae5$0$1261$c30e37c6@exi-reader.telstra.net...
> | "geo" <gmarsaglia@gmail.com> wrote in message 
> news:08db9776-d933-44a5-a20f-7a137e81cce5@5g2000yqz.googlegroups.com...
> |
> || 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.
> |
> | Below is my PL/I version, using signed arithmetic,
> | which I confirm delivers the same test results as previously,
> | and which are the same results as the signed version.
>
> um...unsigned version.

there is someone here can do the signed one in C?

even the "*" change unsinged-signed
i suspect there are many iussue on that

NEVER USE SIGNED!