You are viewing a plain text version of this content. The canonical link for it is here.
Posted to commits@commons.apache.org by lu...@apache.org on 2010/09/21 00:06:40 UTC
svn commit: r999136 - in
/commons/proper/math/branches/MATH_2_X/src/main/java/org/apache/commons/math/random:
AbstractWell.java Well1024a.java Well19937a.java Well19937c.java
Well44497a.java Well44497b.java Well512a.java
Author: luc
Date: Mon Sep 20 22:06:40 2010
New Revision: 999136
URL: http://svn.apache.org/viewvc?rev=999136&view=rev
Log:
improved Well pseudo random numbers generators performances by inlining transforms and using indirection arrays to avoid index computation
Modified:
commons/proper/math/branches/MATH_2_X/src/main/java/org/apache/commons/math/random/AbstractWell.java
commons/proper/math/branches/MATH_2_X/src/main/java/org/apache/commons/math/random/Well1024a.java
commons/proper/math/branches/MATH_2_X/src/main/java/org/apache/commons/math/random/Well19937a.java
commons/proper/math/branches/MATH_2_X/src/main/java/org/apache/commons/math/random/Well19937c.java
commons/proper/math/branches/MATH_2_X/src/main/java/org/apache/commons/math/random/Well44497a.java
commons/proper/math/branches/MATH_2_X/src/main/java/org/apache/commons/math/random/Well44497b.java
commons/proper/math/branches/MATH_2_X/src/main/java/org/apache/commons/math/random/Well512a.java
Modified: commons/proper/math/branches/MATH_2_X/src/main/java/org/apache/commons/math/random/AbstractWell.java
URL: http://svn.apache.org/viewvc/commons/proper/math/branches/MATH_2_X/src/main/java/org/apache/commons/math/random/AbstractWell.java?rev=999136&r1=999135&r2=999136&view=diff
==============================================================================
--- commons/proper/math/branches/MATH_2_X/src/main/java/org/apache/commons/math/random/AbstractWell.java (original)
+++ commons/proper/math/branches/MATH_2_X/src/main/java/org/apache/commons/math/random/AbstractWell.java Mon Sep 20 22:06:40 2010
@@ -38,29 +38,22 @@ public abstract class AbstractWell exten
/** Serializable version identifier. */
private static final long serialVersionUID = -8068371019303673353L;
- /** Number of bits blocks in the pool. */
- private final int r;
-
/** Bit mask preserving the first w - p bits in a w bits block. */
- private final int mp;
+ protected final int mp;
/** Bit mask preserving the last p bits in a w bits block. */
- private final int mpTilde;
-
- /** First parameter of the algorithm. */
- private final int m1;
-
- /** Second parameter of the algorithm. */
- private final int m2;
-
- /** Third parameter of the algorithm. */
- private final int m3;
+ protected final int mpTilde;
/** Current index in the bytes pool. */
- private int index;
+ protected int index;
/** Bytes pool. */
- private final int[] v;
+ protected final int[] v;
+ protected final int[] iRm1;
+ protected final int[] iRm2;
+ protected final int[] i1;
+ protected final int[] i2;
+ protected final int[] i3;
/** Creates a new random number generator.
* <p>The instance is initialized using the current time as the
@@ -99,18 +92,29 @@ public abstract class AbstractWell exten
// of w bits blocks, w is the block size (always 32 in the original paper)
// and p is the number of unused bits in the last block
final int w = 32;
- this.r = (k + w - 1) / w;
+ final int r = (k + w - 1) / w;
final int p = r * w - k;
// set up generator parameters
this.mp = (-1) << p;
this.mpTilde = ~mp;
- this.m1 = m1;
- this.m2 = m2;
- this.m3 = m3;
this.v = new int[r];
this.index = 0;
+ // set up indirection indices
+ iRm1 = new int[r];
+ iRm2 = new int[r];
+ i1 = new int[r];
+ i2 = new int[r];
+ i3 = new int[r];
+ for (int j = 0; j < r; ++j) {
+ iRm1[j] = (j + r - 1) % r;
+ iRm2[j] = (j + r - 2) % r;
+ i1[j] = (j + m1) % r;
+ i2[j] = (j + m2) % r;
+ i3[j] = (j + m3) % r;
+ }
+
// initialize the pool content
setSeed(seed);
@@ -171,164 +175,7 @@ public abstract class AbstractWell exten
setSeed(new int[] { (int) (seed >>> 32), (int) (seed & 0xffffffffl) });
}
- /** Generate next pseudorandom number.
- * <p>This method is the core generation algorithm. It is used by all the
- * public generation methods for the various primitive types {@link
- * #nextBoolean()}, {@link #nextBytes(byte[])}, {@link #nextDouble()},
- * {@link #nextFloat()}, {@link #nextGaussian()}, {@link #nextInt()},
- * {@link #next(int)} and {@link #nextLong()}.</p>
- * <p>This implementation is the general WELL algorithm, described in
- * a paper by François Panneton, Pierre L'Ecuyer and Makoto Matsumoto
- * <a href="http://www.iro.umontreal.ca/~lecuyer/myftp/papers/wellrng.pdf">Improved
- * Long-Period Generators Based on Linear Recurrences Modulo 2</a> ACM
- * Transactions on Mathematical Software, 32, 1 (2006).</p>
- * @param bits number of random bits to produce
- * @return random bits generated
- */
- protected int next(final int bits) {
-
- final int iRm1 = (index + r - 1) % r;
- final int iRm2 = (index + r - 2) % r;
- final int i1 = (index + m1) % r;
- final int i2 = (index + m2) % r;
- final int i3 = (index + m3) % r;
-
- final int z0 = (mp & v[iRm1]) ^ (mpTilde & v[iRm2]);
- final int z1 = t0(v[index]) ^ t1(v[i1]);
- final int z2 = t2(v[i2]) ^ t3(v[i3]);
- final int z3 = z1 ^ z2;
- final int z4 = t4(z0) ^ t5(z1) ^ t6(z2) ^ t7(z3);
-
- v[index] = z3;
- v[iRm1] = z4;
- v[iRm2] &= mp;
- index = iRm1;
-
- return z4 >>> (32 - bits);
-
- }
-
- /** Apply transform M<sub>0</sub> to a bits block.
- * @param x bits block to apply transform to
- * @return M<sub>0</sub>(x)
- */
- protected int m0(final int x) {
- return 0;
- }
-
- /** Apply transform M<sub>1</sub> to a bits block.
- * @param x bits block to apply transform to
- * @return M<sub>1</sub>(x)
- */
- protected int m1(final int x) {
- return x;
- }
-
- /** Apply transform M<sub>2</sub> to a bits block.
- * @param t parameter of the transform
- * @param x bits block to apply transform to
- * @return M<sub>2, t</sub>(x)
- */
- protected int m2(final int t, final int x) {
- return (t >= 0) ? (x >>> t) : (x << -t);
- }
-
- /** Apply transform M<sub>3</sub> to a bits block.
- * @param t parameter of the transform
- * @param x bits block to apply transform to
- * @return M<sub>3, t</sub>(x)
- */
- protected int m3(final int t, final int x) {
- return x ^ ((t >= 0) ? (x >>> t) : (x << -t));
- }
-
- /** Apply transform M<sub>4</sub> to a bits block.
- * @param a parameter of the transform
- * @param x bits block to apply transform to
- * @return M<sub>4, a</sub>(x)
- */
- protected int m4(final int a, final int x) {
- final int shiftedX = x >>> 1;
- return ((x & 0x80000000) != 0) ? (shiftedX ^ a) : shiftedX;
- }
-
- /** Apply transform M<sub>5</sub> to a bits block.
- * @param t first parameter of the transform
- * @param b second parameter of the transform
- * @param x bits block to apply transform to
- * @return M<sub>5, t, b</sub>(x)
- */
- protected int m5(final int t, final int b, final int x) {
- // table I of the paper specifies that a left shift for positive t and
- // a right shift for negative t, however, reference implementation does
- // the opposite (and in fact this transform is used only by Well512a
- // with t = -28). Here, we follow the reference implementation with a
- // left shift for NEGATIVE t
- final int shiftedX = (t >= 0) ? (x >>> t) : (x << -t);
- return x ^ (shiftedX & b);
- }
-
- /** Apply transform M<sub>6</sub> to a bits block.
- * @param q first parameter of the transform
- * @param dsMask second parameter of the transform as a bit mask
- * @param tMask third parameter of the transform as a bit mask
- * @param a fourth parameter of the transform
- * @param x bits block to apply transform to
- * @return M<sub>6, q, s, t, a</sub>(x)
- */
- protected int m6(final int q, final int dsMask, final int tMask, final int a, final int x) {
- final int lShiftedX = x << q;
- final int rShiftedX = x >>> (32 - q);
- final int z = (lShiftedX ^ rShiftedX) & dsMask;
- return ((x & tMask) != 0) ? (z ^ a) : z;
- }
-
- /** Apply transform T<sub>0</sub> to a bits block.
- * @param vi0 bits block to apply transform to
- * @return T<sub>0</sub> v<sub>i,0</sub>
- */
- protected abstract int t0(int vi0);
-
- /** Apply transform T<sub>1</sub> to a bits block.
- * @param vim1 bits block to apply transform to
- * @return T<sub>1</sub> v<sub>i,m1</sub>
- */
- protected abstract int t1(int vim1);
-
- /** Apply transform T<sub>2</sub> to a bits block.
- * @param vim2 bits block to apply transform to
- * @return T<sub>2</sub> v<sub>i,m2</sub>
- */
- protected abstract int t2(int vim2);
-
- /** Apply transform T<sub>3</sub> to a bits block.
- * @param vim3 bits block to apply transform to
- * @return T<sub>3</sub> v<sub>i,m3</sub>
- */
- protected abstract int t3(int vim3);
-
- /** Apply transform T<sub>4</sub> to a bits block.
- * @param z0 bits block to apply transform to
- * @return T<sub>4</sub> z<sub>0</sub>
- */
- protected abstract int t4(int z0);
-
- /** Apply transform T<sub>5</sub> to a bits block.
- * @param z1 bits block to apply transform to
- * @return T<sub>5</sub> z<sub>1</sub>
- */
- protected abstract int t5(int z1);
-
- /** Apply transform T<sub>6</sub> to a bits block.
- * @param z2 bits block to apply transform to
- * @return T<sub>6</sub> z<sub>2</sub>
- */
- protected abstract int t6(int z2);
-
- /** Apply transform T<sub>7</sub> to a bits block.
- * @param z3 bits block to apply transform to
- * @return T<sub>7</sub> z<sub>3</sub>
- */
- protected abstract int t7(int z3);
+ /** {@inheritDoc} */
+ protected abstract int next(final int bits);
}
Modified: commons/proper/math/branches/MATH_2_X/src/main/java/org/apache/commons/math/random/Well1024a.java
URL: http://svn.apache.org/viewvc/commons/proper/math/branches/MATH_2_X/src/main/java/org/apache/commons/math/random/Well1024a.java?rev=999136&r1=999135&r2=999136&view=diff
==============================================================================
--- commons/proper/math/branches/MATH_2_X/src/main/java/org/apache/commons/math/random/Well1024a.java (original)
+++ commons/proper/math/branches/MATH_2_X/src/main/java/org/apache/commons/math/random/Well1024a.java Mon Sep 20 22:06:40 2010
@@ -34,7 +34,7 @@ package org.apache.commons.math.random;
public class Well1024a extends AbstractWell {
/** Serializable version identifier. */
- private static final long serialVersionUID = -5403981908127539981L;
+ private static final long serialVersionUID = 5680173464174485492L;
/** Number of bits in the pool. */
private static final int K = 1024;
@@ -79,43 +79,28 @@ public class Well1024a extends AbstractW
}
/** {@inheritDoc} */
- protected int t0(final int vi0) {
- return m1(vi0);
- }
-
- /** {@inheritDoc} */
- protected int t1(final int vim1) {
- return m3(8, vim1);
- }
-
- /** {@inheritDoc} */
- protected int t2(final int vim2) {
- return m3(-19, vim2);
- }
-
- /** {@inheritDoc} */
- protected int t3(final int vim3) {
- return m3(-14, vim3);
- }
+ protected int next(final int bits) {
- /** {@inheritDoc} */
- protected int t4(final int z0) {
- return m3(-11, z0);
- }
+ final int indexRm1 = iRm1[index];
+ final int indexRm2 = iRm2[index];
- /** {@inheritDoc} */
- protected int t5(final int z1) {
- return m3(-7, z1);
- }
+ final int v0 = v[index];
+ final int vM1 = v[i1[index]];
+ final int vM2 = v[i2[index]];
+ final int vM3 = v[i3[index]];
+
+ final int z0 = v[indexRm1];
+ final int z1 = v0 ^ (vM1 ^ (vM1 >>> 8));
+ final int z2 = (vM2 ^ (vM2 << 19)) ^ (vM3 ^ (vM3 << 14));
+ final int z3 = z1 ^ z2;
+ final int z4 = (z0 ^ (z0 << 11)) ^ (z1 ^ (z1 << 7)) ^ (z2 ^ (z2 << 13));
+
+ v[index] = z3;
+ v[indexRm1] = z4;
+ v[indexRm2] &= mp;
+ index = indexRm1;
- /** {@inheritDoc} */
- protected int t6(final int z2) {
- return m3(-13, z2);
- }
+ return z4 >>> (32 - bits);
- /** {@inheritDoc} */
- protected int t7(final int z3) {
- return m0(z3);
}
-
}
Modified: commons/proper/math/branches/MATH_2_X/src/main/java/org/apache/commons/math/random/Well19937a.java
URL: http://svn.apache.org/viewvc/commons/proper/math/branches/MATH_2_X/src/main/java/org/apache/commons/math/random/Well19937a.java?rev=999136&r1=999135&r2=999136&view=diff
==============================================================================
--- commons/proper/math/branches/MATH_2_X/src/main/java/org/apache/commons/math/random/Well19937a.java (original)
+++ commons/proper/math/branches/MATH_2_X/src/main/java/org/apache/commons/math/random/Well19937a.java Mon Sep 20 22:06:40 2010
@@ -34,7 +34,7 @@ package org.apache.commons.math.random;
public class Well19937a extends AbstractWell {
/** Serializable version identifier. */
- private static final long serialVersionUID = -8052371714518610855L;
+ private static final long serialVersionUID = -7462102162223815419L;
/** Number of bits in the pool. */
private static final int K = 19937;
@@ -79,43 +79,28 @@ public class Well19937a extends Abstract
}
/** {@inheritDoc} */
- protected int t0(final int vi0) {
- return m3(-25, vi0);
- }
-
- /** {@inheritDoc} */
- protected int t1(final int vim1) {
- return m3(27, vim1);
- }
-
- /** {@inheritDoc} */
- protected int t2(final int vim2) {
- return m2(9, vim2);
- }
-
- /** {@inheritDoc} */
- protected int t3(final int vim3) {
- return m3(1, vim3);
- }
+ protected int next(final int bits) {
- /** {@inheritDoc} */
- protected int t4(final int z0) {
- return m1(z0);
- }
+ final int indexRm1 = iRm1[index];
+ final int indexRm2 = iRm2[index];
- /** {@inheritDoc} */
- protected int t5(final int z1) {
- return m3(-9, z1);
- }
+ final int v0 = v[index];
+ final int vM1 = v[i1[index]];
+ final int vM2 = v[i2[index]];
+ final int vM3 = v[i3[index]];
+
+ final int z0 = (0x80000000 & v[indexRm1]) ^ (0x7FFFFFFF & v[indexRm2]);
+ final int z1 = (v0 ^ (v0 << 25)) ^ (vM1 ^ (vM1 >>> 27));
+ final int z2 = (vM2 >>> 9) ^ (vM3 ^ (vM3 >>> 1));
+ final int z3 = z1 ^ z2;
+ final int z4 = z0 ^ (z1 ^ (z1 << 9)) ^ (z2 ^ (z2 << 21)) ^ (z3 ^ (z3 >>> 21));
+
+ v[index] = z3;
+ v[indexRm1] = z4;
+ v[indexRm2] &= mp;
+ index = indexRm1;
- /** {@inheritDoc} */
- protected int t6(final int z2) {
- return m3(-21, z2);
- }
+ return z4 >>> (32 - bits);
- /** {@inheritDoc} */
- protected int t7(final int z3) {
- return m3(21, z3);
}
-
}
Modified: commons/proper/math/branches/MATH_2_X/src/main/java/org/apache/commons/math/random/Well19937c.java
URL: http://svn.apache.org/viewvc/commons/proper/math/branches/MATH_2_X/src/main/java/org/apache/commons/math/random/Well19937c.java?rev=999136&r1=999135&r2=999136&view=diff
==============================================================================
--- commons/proper/math/branches/MATH_2_X/src/main/java/org/apache/commons/math/random/Well19937c.java (original)
+++ commons/proper/math/branches/MATH_2_X/src/main/java/org/apache/commons/math/random/Well19937c.java Mon Sep 20 22:06:40 2010
@@ -31,23 +31,36 @@ package org.apache.commons.math.random;
* @since 2.2
*/
-public class Well19937c extends Well19937a {
+public class Well19937c extends AbstractWell {
/** Serializable version identifier. */
private static final long serialVersionUID = -7203498180754925124L;
+ /** Number of bits in the pool. */
+ private static final int K = 19937;
+
+ /** First parameter of the algorithm. */
+ private static final int M1 = 70;
+
+ /** Second parameter of the algorithm. */
+ private static final int M2 = 179;
+
+ /** Third parameter of the algorithm. */
+ private static final int M3 = 449;
+
/** Creates a new random number generator.
* <p>The instance is initialized using the current time as the
* seed.</p>
*/
public Well19937c() {
+ super(K, M1, M2, M3);
}
/** Creates a new random number generator using a single int seed.
* @param seed the initial seed (32 bits integer)
*/
public Well19937c(int seed) {
- super(seed);
+ super(K, M1, M2, M3, seed);
}
/** Creates a new random number generator using an int array seed.
@@ -55,29 +68,45 @@ public class Well19937c extends Well1993
* the seed of the generator will be related to the current time
*/
public Well19937c(int[] seed) {
- super(seed);
+ super(K, M1, M2, M3, seed);
}
/** Creates a new random number generator using a single long seed.
* @param seed the initial seed (64 bits integer)
*/
public Well19937c(long seed) {
- super(seed);
+ super(K, M1, M2, M3, seed);
}
/** {@inheritDoc} */
protected int next(final int bits) {
- // compute raw value given by WELL19937a generator
- // which is NOT maximally-equidistributed
- int z = super.next(32);
+ final int indexRm1 = iRm1[index];
+ final int indexRm2 = iRm2[index];
+
+ final int v0 = v[index];
+ final int vM1 = v[i1[index]];
+ final int vM2 = v[i2[index]];
+ final int vM3 = v[i3[index]];
+
+ final int z0 = (0x80000000 & v[indexRm1]) ^ (0x7FFFFFFF & v[indexRm2]);
+ final int z1 = (v0 ^ (v0 << 25)) ^ (vM1 ^ (vM1 >>> 27));
+ final int z2 = (vM2 >>> 9) ^ (vM3 ^ (vM3 >>> 1));
+ final int z3 = z1 ^ z2;
+ int z4 = z0 ^ (z1 ^ (z1 << 9)) ^ (z2 ^ (z2 << 21)) ^ (z3 ^ (z3 >>> 21));
+
+ v[index] = z3;
+ v[indexRm1] = z4;
+ v[indexRm2] &= mp;
+ index = indexRm1;
+
// add Matsumoto-Kurita tempering
// to get a maximally-equidistributed generator
- z = z ^ ((z << 7) & 0xe46e1700);
- z = z ^ ((z << 15) & 0x9b868000);
+ z4 = z4 ^ ((z4 << 7) & 0xe46e1700);
+ z4 = z4 ^ ((z4 << 15) & 0x9b868000);
- return z >>> (32 - bits);
+ return z4 >>> (32 - bits);
}
Modified: commons/proper/math/branches/MATH_2_X/src/main/java/org/apache/commons/math/random/Well44497a.java
URL: http://svn.apache.org/viewvc/commons/proper/math/branches/MATH_2_X/src/main/java/org/apache/commons/math/random/Well44497a.java?rev=999136&r1=999135&r2=999136&view=diff
==============================================================================
--- commons/proper/math/branches/MATH_2_X/src/main/java/org/apache/commons/math/random/Well44497a.java (original)
+++ commons/proper/math/branches/MATH_2_X/src/main/java/org/apache/commons/math/random/Well44497a.java Mon Sep 20 22:06:40 2010
@@ -34,7 +34,7 @@ package org.apache.commons.math.random;
public class Well44497a extends AbstractWell {
/** Serializable version identifier. */
- private static final long serialVersionUID = 5154222742730470272L;
+ private static final long serialVersionUID = -3859207588353972099L;
/** Number of bits in the pool. */
private static final int K = 44497;
@@ -79,46 +79,30 @@ public class Well44497a extends Abstract
}
/** {@inheritDoc} */
- protected int t0(final int vi0) {
- return m3(-24, vi0);
- }
-
- /** {@inheritDoc} */
- protected int t1(final int vim1) {
- return m3(30, vim1);
- }
-
- /** {@inheritDoc} */
- protected int t2(final int vim2) {
- return m3(-10, vim2);
- }
-
- /** {@inheritDoc} */
- protected int t3(final int vim3) {
- return m2(-26, vim3);
- }
+ protected int next(final int bits) {
- /** {@inheritDoc} */
- protected int t4(final int z0) {
- return m1(z0);
- }
+ final int indexRm1 = iRm1[index];
+ final int indexRm2 = iRm2[index];
- /** {@inheritDoc} */
- protected int t5(final int z1) {
- return m3(20, z1);
- }
+ final int v0 = v[index];
+ final int vM1 = v[i1[index]];
+ final int vM2 = v[i2[index]];
+ final int vM3 = v[i3[index]];
+
+ final int z0 = (0xFFFF8000 & v[indexRm1]) ^ (0x00007FFF & v[indexRm2]);
+ final int z1 = (v0 ^ (v0 << 24)) ^ (vM1 ^ (vM1 >>> 30));
+ final int z2 = (vM2 ^ (vM2 << 10)) ^ (vM3 << 26);
+ final int z3 = z1 ^ z2;
+ final int z2Prime = ((z2 << 9) ^ (z2 >>> 23)) & 0xfbffffff;
+ final int z2Second = ((z2 & 0x00020000) != 0) ? (z2Prime ^ 0xb729fcec) : z2Prime;
+ final int z4 = z0 ^ (z1 ^ (z1 >>> 20)) ^ z2Second ^ z3;
+
+ v[index] = z3;
+ v[indexRm1] = z4;
+ v[indexRm2] &= mp;
+ index = indexRm1;
- /** {@inheritDoc} */
- protected int t6(final int z2) {
- // table II of the paper specifies t6 to be m6(9, d14, t5, 0xb729fcec, z2)
- // however, the reference implementation uses m6(9, d26, t17, 0xb729fcec, z2).
- // Here, we follow the reference implementation
- return m6(9, (-1) ^ (0x1 << 26), 0x1 << 17, 0xb729fcec, z2);
- }
+ return z4 >>> (32 - bits);
- /** {@inheritDoc} */
- protected int t7(final int z3) {
- return m1(z3);
}
-
}
Modified: commons/proper/math/branches/MATH_2_X/src/main/java/org/apache/commons/math/random/Well44497b.java
URL: http://svn.apache.org/viewvc/commons/proper/math/branches/MATH_2_X/src/main/java/org/apache/commons/math/random/Well44497b.java?rev=999136&r1=999135&r2=999136&view=diff
==============================================================================
--- commons/proper/math/branches/MATH_2_X/src/main/java/org/apache/commons/math/random/Well44497b.java (original)
+++ commons/proper/math/branches/MATH_2_X/src/main/java/org/apache/commons/math/random/Well44497b.java Mon Sep 20 22:06:40 2010
@@ -31,23 +31,36 @@ package org.apache.commons.math.random;
* @since 2.2
*/
-public class Well44497b extends Well44497a {
+public class Well44497b extends AbstractWell {
/** Serializable version identifier. */
private static final long serialVersionUID = 4032007538246675492L;
+ /** Number of bits in the pool. */
+ private static final int K = 44497;
+
+ /** First parameter of the algorithm. */
+ private static final int M1 = 23;
+
+ /** Second parameter of the algorithm. */
+ private static final int M2 = 481;
+
+ /** Third parameter of the algorithm. */
+ private static final int M3 = 229;
+
/** Creates a new random number generator.
* <p>The instance is initialized using the current time as the
* seed.</p>
*/
public Well44497b() {
+ super(K, M1, M2, M3);
}
/** Creates a new random number generator using a single int seed.
* @param seed the initial seed (32 bits integer)
*/
public Well44497b(int seed) {
- super(seed);
+ super(K, M1, M2, M3, seed);
}
/** Creates a new random number generator using an int array seed.
@@ -55,14 +68,14 @@ public class Well44497b extends Well4449
* the seed of the generator will be related to the current time
*/
public Well44497b(int[] seed) {
- super(seed);
+ super(K, M1, M2, M3, seed);
}
/** Creates a new random number generator using a single long seed.
* @param seed the initial seed (64 bits integer)
*/
public Well44497b(long seed) {
- super(seed);
+ super(K, M1, M2, M3, seed);
}
/** {@inheritDoc} */
@@ -70,14 +83,33 @@ public class Well44497b extends Well4449
// compute raw value given by WELL44497a generator
// which is NOT maximally-equidistributed
- int z = super.next(32);
+ final int indexRm1 = iRm1[index];
+ final int indexRm2 = iRm2[index];
+
+ final int v0 = v[index];
+ final int vM1 = v[i1[index]];
+ final int vM2 = v[i2[index]];
+ final int vM3 = v[i3[index]];
+
+ final int z0 = (0xFFFF8000 & v[indexRm1]) ^ (0x00007FFF & v[indexRm2]);
+ final int z1 = (v0 ^ (v0 << 24)) ^ (vM1 ^ (vM1 >>> 30));
+ final int z2 = (vM2 ^ (vM2 << 10)) ^ (vM3 << 26);
+ final int z3 = z1 ^ z2;
+ final int z2Prime = ((z2 << 9) ^ (z2 >>> 23)) & 0xfbffffff;
+ final int z2Second = ((z2 & 0x00020000) != 0) ? (z2Prime ^ 0xb729fcec) : z2Prime;
+ int z4 = z0 ^ (z1 ^ (z1 >>> 20)) ^ z2Second ^ z3;
+
+ v[index] = z3;
+ v[indexRm1] = z4;
+ v[indexRm2] &= mp;
+ index = indexRm1;
// add Matsumoto-Kurita tempering
// to get a maximally-equidistributed generator
- z = z ^ ((z << 7) & 0x93dd1400);
- z = z ^ ((z << 15) & 0xfa118000);
+ z4 = z4 ^ ((z4 << 7) & 0x93dd1400);
+ z4 = z4 ^ ((z4 << 15) & 0xfa118000);
- return z >>> (32 - bits);
+ return z4 >>> (32 - bits);
}
Modified: commons/proper/math/branches/MATH_2_X/src/main/java/org/apache/commons/math/random/Well512a.java
URL: http://svn.apache.org/viewvc/commons/proper/math/branches/MATH_2_X/src/main/java/org/apache/commons/math/random/Well512a.java?rev=999136&r1=999135&r2=999136&view=diff
==============================================================================
--- commons/proper/math/branches/MATH_2_X/src/main/java/org/apache/commons/math/random/Well512a.java (original)
+++ commons/proper/math/branches/MATH_2_X/src/main/java/org/apache/commons/math/random/Well512a.java Mon Sep 20 22:06:40 2010
@@ -34,7 +34,7 @@ package org.apache.commons.math.random;
public class Well512a extends AbstractWell {
/** Serializable version identifier. */
- private static final long serialVersionUID = 8706771840051210473L;
+ private static final long serialVersionUID = -6104179812103820574L;
/** Number of bits in the pool. */
private static final int K = 512;
@@ -79,46 +79,27 @@ public class Well512a extends AbstractWe
}
/** {@inheritDoc} */
- protected int t0(final int vi0) {
- return m3(-16, vi0);
- }
+ protected int next(final int bits) {
- /** {@inheritDoc} */
- protected int t1(final int vim1) {
- return m3(-15, vim1);
- }
+ final int indexRm1 = iRm1[index];
- /** {@inheritDoc} */
- protected int t2(final int vim2) {
- return m3(11, vim2);
- }
+ final int vi = v[index];
+ final int vi1 = v[i1[index]];
+ final int vi2 = v[i2[index]];
+ final int z0 = v[indexRm1];
+ // m3: x ^ ((t >= 0) ? (x >>> t) : (x << -t));
+
+ final int z1 = (vi ^ (vi << 16)) ^ (vi1 ^ (vi1 << 15));
+ final int z2 = vi2 ^ (vi2 >>> 11);
+ final int z3 = z1 ^ z2;
+ final int z4 = (z0 ^ (z0 << 2)) ^ (z1 ^ (z1 << 18)) ^ (z2 << 28) ^ (z3 ^ ((z3 << 5) & 0xda442d24));
+
+ v[index] = z3;
+ v[indexRm1] = z4;
+ index = indexRm1;
- /** {@inheritDoc} */
- protected int t3(final int vim3) {
- return m0(vim3);
- }
+ return z4 >>> (32 - bits);
- /** {@inheritDoc} */
- protected int t4(final int z0) {
- return m3(-2, z0);
- }
-
- /** {@inheritDoc} */
- protected int t5(final int z1) {
- return m3(-18, z1);
- }
-
- /** {@inheritDoc} */
- protected int t6(final int z2) {
- // table II of the paper specifies t6 to be m3(-28, z2)
- // however, the reference implementation uses m2(-28, z2).
- // Here, we follow the reference implementation
- return m2(-28, z2);
- }
-
- /** {@inheritDoc} */
- protected int t7(final int z3) {
- return m5(-5, 0xda442d24, z3);
}
}