package edu.cornell.lassp.houle.RngPack; import java.util.*; import java.io.Serializable; // // RngPack 1.1a by Paul Houle // http://www.honeylocust.com/RngPack/ // /** * * RANLUX is an advanced pseudo-random number generator based on the * RCARRY algorithm proposed in 1991 by Marsaglia and Zaman. * RCARRY * used a subtract-and-borrow algorithm with a period on the order of * 10¹⁷¹ but still had detectable correlations between * numbers. Martin Luescher proposed the RANLUX * algorithm in 1993; RANLUX generates pseudo-random numbers using * RCARRY but throws away numbers to destroy correlations. * Thus RANLUX trades execution speed for quality; by choosing a * larger luxury setting one gets better random numbers slower. * By the tests availible at the time it was proposed, * RANLUX at the default luxury setting appears to be a * significant advance quality over previous * generators. * *
*
* * * * * * * * * * * * * * * * * *

* LUXURY LEVELS *
level	p
0	24	equivalent to the original `RCARRY` of Marsaglia * and Zaman, very long period, but fails many tests.
1	48	considerable improvement in quality over level 0, * now passes the gap test, but still fails spectral test.
2	97	passes all known tests, but theoretically still * defective.
3	223	* DEFAULT VALUE. Any theoretically possible * correlations have very small chance of being observed.
4	389	highest possible luxury, all 24 bits chaotic.

* *
* * VALIDATION * * * The Java version of RANLUX has been verified against published * values of numbers 1-5 and 101-105 produced by the reference implementation * of RANLUX for the following initial conditions: * *

Default initialization: Ranlux() *
Initialization with: Ranlux(0,0) *
Initialization with: Ranlux(389,1) *
Initialization with: Ranlux(75,0) *

* References: *

* M. Luscher, Computer Physics Communications 79 (1994) 100 *
* F. James, Computer Physics Communications 79 (1994) 111 *
About RANLUX random number generator: Excerpts from discussion in the Usenet news groups *
Miller's FORTRAN 90 implementation of RANLUX with test code *

* * *

* * Source code is available. * * @author Paul Houle (E-mail: paul@honeylocust.com) * @version 1.1a */ public class Ranlux extends RandomSeedable implements Serializable { /** * Maximum luxury level: maxlev=4 */ public static final int maxlev=4; /** * Default luxury level: lxdflt=3 */ public static final int lxdflt=3; static final int igiga = 1000000000, jsdflt = 314159265; static final int twop12 = 4096,itwo24 = 1 << 24,icons = 2147483563; static final int ndskip[]={0,24,73,199,365}; int iseeds[],isdext[],next[]; int luxlev=lxdflt,nskip,inseed,jseed; int in24 = 0, kount = 0, mkount = 0, i24 = 24, j24 = 10; float seeds[],carry= (float) 0.0,twom24,twom12; boolean diagOn=false; /** * Default initialization of RANLUX. Uses default seed * jsdflt=314159265 and luxury level 3. */ public Ranlux() { init_arrays(); rluxdef(); }; /** * Initialize RANLUX with specified luxury level * and seed. * * @param lux luxury level from 0-4. * @param ins seed, a positive integer. * */ public Ranlux(int lux,int ins) { init_arrays(); rluxgo(lux,Math.abs(ins)); }; /** * Initialize RANLUX with specified luxury level * and seed. * * @param lux luxury level from 0-4. * @param ins seed, a positive long. * */ public Ranlux(int lux,long ins) { init_arrays(); rluxgo(lux,Math.abs((int) (ins%Integer.MAX_VALUE))); }; /** * * Initialize RANLUX with default luxury level * and a specified seed. * * @param ins seed, a positive integer */ public Ranlux(int ins) { init_arrays(); rluxgo(lxdflt,Math.abs(ins)); }; /** * * Initialize RANLUX with default luxury level * and a specified seed. * * @param ins seed, a positive integer */ public Ranlux(long ins) { init_arrays(); rluxgo(lxdflt,Math.abs((int) (ins%Integer.MAX_VALUE))); }; /** * * Initialize RANLUX with specified luxury level * and a Date object. Can be used to conveniently initialize RANLUX * from the clock, * *

* RandomElement e = Ranlux(4,new Date());
*

* * @param lux luxury level from 0-4. * @param d date used to generate seed * */ public Ranlux(int lux,Date d) { init_arrays(); rluxgo(lux,(int) (ClockSeed(d)%Integer.MAX_VALUE) ); }; /** * * Initialize RANLUX with default luxury level * and a Date object. Can be used to conveniently initialize RANLUX * from the clock, * *

* RandomElement e = Ranlux(new Date());
*

* * @param d date used to generate seed * */ public Ranlux(Date d) { init_arrays(); rluxgo(lxdflt,(int) (ClockSeed(d)%Integer.MAX_VALUE) ); }; /** * Turns diagnostic messages on and off. If setDiag(true) is * called, RANLUX will print diagnostic information to * System.err * * @param b diagnostic message status */ public void setDiag(boolean b) { diagOn=b; }; /** * * The random number generator. * * @returns a pseudo-random double in the range (0,1) */ public final double raw() { int i,k,lp; float uni,out; uni=seeds[j24]-seeds[i24]-carry; if (uni < (float) 0.0) { uni=uni+ (float) 1.0; carry = twom24; } else carry = (float) 0.0; seeds[i24]=uni; i24=next[i24]; j24=next[j24]; out=uni; if (uni=igiga) { mkount++; kount -= igiga; }; return out; }; private void init_arrays() { /* * * converted from fortran: fortran arrays start at 1, java arrays start at * 0. Here we take the low road to compatibility... We declare arrays that * are one bigger than the fortran code. This wastes three ints and a double; * If you're porting this to a Commodore 64 you might need the space. * */ iseeds = new int[24+1]; isdext = new int[25+1]; next = new int[24+1]; seeds = new float[24+1]; }; private void rluxdef() { int lp,i,k; jseed = jsdflt; inseed = jseed; diag("RANLUX DEFAULT INITIALIZATION: "+jseed); luxlev = lxdflt; nskip = ndskip[luxlev]; lp = nskip + 24; in24 = 0; kount = 0; mkount = 0; diag("RANLUX DEFAULT LUXURY LEVEL = "+luxlev+" p = "+lp); twom24 = (float) 1.0; for(i = 1;i <= 24;i++) { twom24 = twom24 * (float) 0.5; k = jseed / 53668; jseed = 40014 * (jseed-k*53668) - k * 12211; if (jseed<0) jseed = jseed + icons; iseeds[i] = jseed % itwo24; }; twom12 = twom24 * (float) 4096.0; for(i = 1 ; i<=24;i++) { seeds[i] = iseeds[i] * twom24; next[i] = i - 1; }; next[1] = 24; i24 = 24; j24 = 10; carry = (float) 0.0; if (seeds[24] == 0.0) carry = twom24; }; private final void rluxgo(int lux,int ins) { int ilx, i, iouter, isk, k, inner, izip, izip2; float uni; if (lux<00) { luxlev = lxdflt; } else if (lux<=maxlev) { luxlev = lux; } else if (lux<24 || lux>2000) { luxlev = maxlev; diag("RANLUX ILLEGAL LUXURY RLUXGO: "+lux); } else { luxlev = lux; for(ilx=0;ilx<=maxlev;ilx++) if (lux == ndskip[ilx]+24) luxlev=ilx; }; if (luxlev <= maxlev) { nskip = ndskip[luxlev]; diag("RANLUX LUXURY LEVEL SET BY RLUXGO : "+luxlev+" P= "+(nskip+24)); } else { nskip = luxlev-24; diag("RANLUX P-VALUE SET BY RLUXGO TO: "+luxlev); }; in24=0; if (ins<0) diag("Illegal initialization by RLUXGO, negative input seed"); if (ins>0) { jseed = ins; diag("RANLUX INITIALIZED BY RLUXGO FROM SEED "+jseed); } else { jseed=jsdflt; diag("RANLUX INITIALIZED BY RLUXGO FROM DEFAULT SEED"); }; inseed = jseed; twom24 = (float) 1.0; for(i=1;i<=24;i++) { twom24 = twom24 * (float) 0.5; k = jseed / 53668; jseed = 40014 * (jseed-k*53668) - k * 12211; if (jseed < 0) jseed = jseed + icons; iseeds[i] = jseed % itwo24; }; twom12=twom24*4096; for(i=1;i<=24;i++) { seeds[i] = iseeds[i] * twom24; next[i] = i - 1; }; next[1]=24; i24=24; j24=10; carry = (float) 0.0; if (seeds[24] == 0.0) carry = twom24; kount = 0; mkount = 0; }; private void diag(String s) { if (diagOn) System.err.println(s); }; };