Berack/upo-valpre
0
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases 3 Wiki Activity

Rngs

Browse Source 
- imoprted and fixed for java
This commit is contained in:
Giacomo Bertolazzi
2025-01-17 12:28:11 +01:00
parent 49ac43b1fa
commit f0b6cbfba3
3 changed files with 403 additions and 0 deletions
Download Patch File Download Diff File Expand all files Collapse all files

99
src/main/java/net/berack/upo/valpre/rand/Rng.java Normal file
View File

@@ -0,0 +1,99 @@
package net.berack.upo.valpre.rand;
/**
* This is an Java library for random number generation. The use of this
* library is recommended as a replacement for the Java class Random,
* particularly in simulation applications where the statistical
* 'goodness' of the random number generator is important.
*
* The generator used in this library is a so-called 'Lehmer random number
* generator' which returns a pseudo-random number uniformly distributed
* between 0.0 and 1.0. The period is (m - 1) where m = 2,147,483,647 and
* the smallest and largest possible values are (1 / m) and 1 - (1 / m)
* respectively. For more details see:
*
* "Random Number Generators: Good Ones Are Hard To Find"
* Steve Park and Keith Miller
* Communications of the ACM, October 1988
*
* Note that as of 7-11-90 the multiplier used in this library has changed
* from the previous "minimal standard" 16807 to a new value of 48271. To
* use this library in its old (16807) form change the constants MULTIPLIER
* and CHECK as indicated in the comments.
*
* Name : Rng.java (Random Number Generation - Single Stream)
* Authors : Steve Park & Dave Geyer
* Translated by : Jun Wang & Richard Dutton
* Language : Java
* Latest Revision : 6-10-04
*
* Program rng : Section 2.2
*/
public class Rng {
public final static long CHECK = 399268537L; /* use 1043616065 for the "minimal standard" */
public final static long DEFAULT = 123456789L; /* initial seed, use 0 < DEFAULT < MODULUS */
public final static long MODULUS = 2147483647; /* DON'T CHANGE THIS VALUE */
public final static long MULTIPLIER = 48271; /* use 16807 for the "minimal standard" */
private long seed = DEFAULT; /* seed is the state of the generator */
public Rng() {
}
public Rng(long seed) {
this.putSeed(seed);
}
/**
* Random is a Lehmer generator that returns a pseudo-random real number
* uniformly distributed between 0.0 and 1.0. The period is (m - 1)
* where m = 2,147,483,647 amd the smallest and largest possible values
* are (1 / m) and 1 - (1 / m) respectively.
*/
public double random() {
this.seed = Rng.newSeed(MODULUS, MULTIPLIER, this.seed);
return ((double) this.seed / MODULUS);
}
/**
* Use this (optional) procedure to initialize or reset the state of
* the random number generator according to the following conventions:
* if x > 0 then x is the initial seed (unless too large)
* if x <= 0 then the initial seed is obtained from the system clock
*/
public void putSeed(long x) {
if (x > 0L) {
x = x % MODULUS; /* correct if x is too large */
} else {
x = System.currentTimeMillis();
// x = ((unsigned long) time((time_t *) NULL)) % MODULUS;
}
this.seed = x;
}
/**
* Use this (optional) procedure to get the current state of the random
* number generator.
*/
public long getSeed() {
return this.seed;
}
/**
* Use this (optional) procedure to test for a correct implementation.
*/
public static boolean testRandom() {
Rng rng = new Rng(1); /* set initial state to 1 */
for (var i = 0; i < 10000; i++)
rng.random();
return rng.getSeed() == CHECK;
}
public static long newSeed(long modulus, long multiplier, long seed) {
var Q = modulus / multiplier;
var R = modulus % multiplier;
var t = multiplier * (seed % Q) - R * (seed / Q);
return t > 0 ? t : (t + modulus);
}
}

95
src/main/java/net/berack/upo/valpre/rand/Rngs.java Normal file
View File

@@ -0,0 +1,95 @@
package net.berack.upo.valpre.rand;
/**
* This is an Java library for multi-stream random number generation.
* The use of this library is recommended as a replacement for the Java
* class Random, particularly in simulation applications where the
* statistical 'goodness' of the random number generator is important.
* The library supplies 256 streams of random numbers; use
* selectStream(s) to switch between streams indexed s = 0,1,...,255.
*
* The streams must be initialized. The recommended way to do this is by
* using the function plantSeeds(x) with the value of x used to initialize
* the default stream and all other streams initialized automatically with
* values dependent on the value of x. The following convention is used
* to initialize the default stream:
* if x > 0 then x is the state
* if x < 0 then the state is obtained from the system clock
* if x = 0 then the state is to be supplied interactively.
*
* The generator used in this library is a so-called 'Lehmer random number
* generator' which returns a pseudo-random number uniformly distributed
* 0.0 and 1.0. The period is (m - 1) where m = 2,147,483,647 and the
* smallest and largest possible values are (1 / m) and 1 - (1 / m)
* respectively. For more details see:
*
* "Random Number Generators: Good Ones Are Hard To Find"
* Steve Park and Keith Miller
* Communications of the ACM, October 1988
*
* Name : Rngs.java (Random Number Generation - Multiple Streams)
* Authors : Steve Park & Dave Geyer
* Translated by : Jun Wang & Richard Dutton
* Language : Java
* Latest Revision : 6-10-04
*/
class Rngs {
private final static int STREAMS = 256; /* # of streams, DON'T CHANGE THIS VALUE */
private final static long A256 = 22925; /* jump multiplier, DON'T CHANGE THIS VALUE */
private Rng[] rngs;
private int current = 0;
public Rngs(long seed) {
this.rngs = new Rng[STREAMS];
this.plantSeeds(seed);
}
public Rng getRng() {
return this.rngs[this.current];
}
public Rng getRng(int stream) {
return this.rngs[stream % STREAMS];
}
/**
* Use this function to set the state of all the random number generator
* streams by "planting" a sequence of states (seeds), one per stream,
* with all states dictated by the state of the default stream.
* The sequence of planted states is separated one from the next by
* 8,367,782 calls to Random().
*/
public void plantSeeds(long seed0) {
this.rngs[0] = new Rng(seed0);
for (int j = 1; j < STREAMS; j++) {
seed0 = Rng.newSeed(Rng.MODULUS, A256, seed0);
this.rngs[j] = new Rng(seed0);
}
}
/**
* Use this function to set the current random number generator
* stream -- that stream from which the next random number will come.
*/
public void selectStream(int index) {
this.current = index % STREAMS;
}
/**
* Use this (optional) function to test for a correct implementation.
*/
public static boolean testRandom() {
var rngs = new Rngs(1);
var first = rngs.getRng();
for (int i = 0; i < 10000; i++)
first.random();
var x = first.getSeed(); /* get the new state value */
var ok = (x == Rng.CHECK); /* and check for correctness */
x = rngs.getRng(1).getSeed(); /* get the state of stream 1 */
return ok && (x == A256); /* x should be the jump multiplier */
}
}

209
src/main/java/net/berack/upo/valpre/rand/Rvgs.java Normal file
View File

@@ -0,0 +1,209 @@
package net.berack.upo.valpre.rand;
/* --------------------------------------------------------------------------
* This is a Java library for generating random variates from six discrete
* distributions
*
* Generator Range (x) Mean Variance
*
* bernoulli(p) x = 0,1 p p*(1-p)
* binomial(n, p) x = 0,...,n n*p n*p*(1-p)
* equilikely(a, b) x = a,...,b (a+b)/2 ((b-a+1)*(b-a+1)-1)/12
* Geometric(p) x = 0,... p/(1-p) p/((1-p)*(1-p))
* pascal(n, p) x = 0,... n*p/(1-p) n*p/((1-p)*(1-p))
* poisson(m) x = 0,... m m
*
* and seven continuous distributions
*
* uniform(a, b) a < x < b (a + b)/2 (b - a)*(b - a)/12
* exponential(m) x > 0 m m*m
* erlang(n, b) x > 0 n*b n*b*b
* normal(m, s) all x m s*s
* logNormal(a, b) x > 0 see below
* chiSquare(n) x > 0 n 2*n
* student(n) all x 0 (n > 1) n/(n - 2) (n > 2)
*
* For the a Lognormal(a, b) random variable, the mean and variance are
*
* mean = exp(a + 0.5*b*b)
* variance = (exp(b*b) - 1) * exp(2*a + b*b)
*
* Name : Rvgs.java (Random Variate GeneratorS)
* Authors : Steve Park & Dave Geyer
* Translated by : Richard Dutton & Jun Wang
* Language : Java
* Latest Revision : 7-1-04
* --------------------------------------------------------------------------
*/
public class Rvgs {
private final Rng rng;
// public Rvgs() {
// this.rngs = new Rngs(Rng.DEFAULT);
// }
public Rvgs(Rng rng) {
if (rng == null)
throw new NullPointerException();
this.rng = rng;
}
/**
* Returns 1 with probability p or 0 with probability 1 - p.
* NOTE: use 0.0 < p < 1.0
*/
public long bernoulli(double p) {
return ((this.rng.random() < (1.0 - p)) ? 0 : 1);
}
/**
* Returns a binomial distributed integer between 0 and n inclusive.
* NOTE: use n > 0 and 0.0 < p < 1.0
*/
public long binomial(long n, double p) {
long i, x = 0;
for (i = 0; i < n; i++)
x += bernoulli(p);
return (x);
}
/**
* Returns an equilikely distributed integer between a and b inclusive.
* NOTE: use a < b
*/
public long equilikely(long a, long b) {
return (a + (long) ((b - a + 1) * this.rng.random()));
}
/**
* Returns a geometric distributed non-negative integer.
* NOTE: use 0.0 < p < 1.0
*/
public long geometric(double p) {
return ((long) (Math.log(1.0 - this.rng.random()) / Math.log(p)));
}
/**
* Returns a Pascal distributed non-negative integer.
* NOTE: use n > 0 and 0.0 < p < 1.0
*/
public long pascal(long n, double p) {
long i, x = 0;
for (i = 0; i < n; i++)
x += geometric(p);
return (x);
}
/**
* Returns a Poisson distributed non-negative integer.
* NOTE: use m > 0
*/
public long poisson(double m) {
double t = 0.0;
long x = 0;
while (t < m) {
t += exponential(1.0);
x++;
}
return (x - 1);
}
/**
* Returns a uniformly distributed real number between a and b.
* NOTE: use a < b
*/
public double uniform(double a, double b) {
return (a + (b - a) * this.rng.random());
}
/**
* Returns an exponentially distributed positive real number.
* NOTE: use m > 0.0
*/
public double exponential(double m) {
return (-m * Math.log(1.0 - this.rng.random()));
}
/**
* Returns an Erlang distributed positive real number.
* NOTE: use n > 0 and b > 0.0
*/
public double erlang(long n, double b) {
long i;
double x = 0.0;
for (i = 0; i < n; i++)
x += exponential(b);
return (x);
}
/**
* Returns a normal (Gaussian) distributed real number.
* NOTE: use s > 0.0
*
* Uses a very accurate approximation of the normal idf due to Odeh & Evans,
* J. Applied Statistics, 1974, vol 23, pp 96-97.
*/
public double normal(double m, double s) {
final double p0 = 0.322232431088;
final double q0 = 0.099348462606;
final double p1 = 1.0;
final double q1 = 0.588581570495;
final double p2 = 0.342242088547;
final double q2 = 0.531103462366;
final double p3 = 0.204231210245e-1;
final double q3 = 0.103537752850;
final double p4 = 0.453642210148e-4;
final double q4 = 0.385607006340e-2;
double u, t, p, q, z;
u = this.rng.random();
if (u < 0.5)
t = Math.sqrt(-2.0 * Math.log(u));
else
t = Math.sqrt(-2.0 * Math.log(1.0 - u));
p = p0 + t * (p1 + t * (p2 + t * (p3 + t * p4)));
q = q0 + t * (q1 + t * (q2 + t * (q3 + t * q4)));
if (u < 0.5)
z = (p / q) - t;
else
z = t - (p / q);
return (m + s * z);
}
/**
* Returns a lognormal distributed positive real number.
* NOTE: use b > 0.0
*/
public double logNormal(double a, double b) {
return (Math.exp(a + b * normal(0.0, 1.0)));
}
/**
* Returns a chi-square distributed positive real number.
* NOTE: use n > 0
*/
public double chiSquare(long n) {
long i;
double z, x = 0.0;
for (i = 0; i < n; i++) {
z = normal(0.0, 1.0);
x += z * z;
}
return (x);
}
/**
* Returns a student-t distributed real number.
* NOTE: use n > 0
*/
public double student(long n) {
return (normal(0.0, 1.0) / Math.sqrt(chiSquare(n) / n));
}
}