/*-------------------------------------------------------------*/

testKISS32: procedure options (main);

   /* George Marsaglia's Pseudo-Random Number Generator having a period of */
   /* 5*2^1320481*(2^32-1).                                                */
   /* 32-bit PRNGs are generated directly.                                 */
   /* The procedure FILL must be called first before invoking KISS32.      */

KISS32: procedure() returns(fixed binary (32) unsigned) options (reorder);
   declare x fixed binary (32) unsigned;
   declare Q(41265) fixed binary (32) unsigned external,
           indx fixed binary (31) static initial (41266);
   declare (carry initial (362),
           xcng  initial (1236789),
           xs    initial (521288629)) fixed binary (32) unsigned external;

   if indx <= 41265 then do; x=Q(indx); indx=indx+1; end;
   else do; x=refill(); indx = 2; end;
   xcng=xcng*69069+123;
   xs=ieor(xs,isll(xs,13));
   xs=ieor(xs,isrl(xs,17));
   xs=ieor(xs,isll(xs,5));
   return (x+xcng+xs);
end KISS32;

refill: procedure() returns(fixed binary(32) unsigned) options (reorder);
   declare (i,s,z,h) fixed binary (32) unsigned;
   declare Q(41265) fixed binary (32) unsigned external;
   declare (carry initial (362),
           xcng  initial (1236789),
           xs    initial (521288629)) fixed binary (32) unsigned external;

   do i = 1 to 41265;
      h = iand(carry,1);
      z = isrl(isll(Q(i),9),1) + isrl(isll(Q(i),7),1) + isrl(carry,1);
      carry=isrl(Q(i),23)+isrl(Q(i),25)+isrl(z,31);
      Q(i)=inot(isll(z,1)+h);
   end;
   return (Q(1));
end refill;

fill: procedure options (reorder);
   declare Q(41265) fixed binary (32) unsigned external;
   declare (carry initial (362),
           xcng  initial (1236789),
           xs    initial (521288629)) fixed binary (32) unsigned external;
   declare i fixed binary (31);

   do i=1 to 41265;   /* fill Q with Congruential+Xorshift */
      xcng=xcng*69069+123;
      xs=ieor(xs,isll(xs,13));
      xs=ieor(xs,isrl(xs,17));
      xs=ieor(xs,isll(xs,5));
      Q(i)=xcng+xs;
   end;
end fill;

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

   call fill;
   do i=1 to 1000000000;  x=KISS32();  end;
   put skip edit ( 'Does x = 1809478889 ?', 'x =', x)
      (a, skip, x(5), a, f(11));
 end testKISS32;