commons-commits mailing list archives

Site index · List index
Message view « Date » · « Thread »
Top « Date » · « Thread »
From l..@apache.org
Subject svn commit: r999135 - in /commons/proper/math/trunk/src/main/java/org/apache/commons/math/random: AbstractWell.java Well1024a.java Well19937a.java Well19937c.java Well44497a.java Well44497b.java Well512a.java
Date Mon, 20 Sep 2010 22:06:07 GMT
Author: luc
Date: Mon Sep 20 22:06:07 2010
New Revision: 999135

URL: http://svn.apache.org/viewvc?rev=999135&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/trunk/src/main/java/org/apache/commons/math/random/AbstractWell.java
    commons/proper/math/trunk/src/main/java/org/apache/commons/math/random/Well1024a.java
    commons/proper/math/trunk/src/main/java/org/apache/commons/math/random/Well19937a.java
    commons/proper/math/trunk/src/main/java/org/apache/commons/math/random/Well19937c.java
    commons/proper/math/trunk/src/main/java/org/apache/commons/math/random/Well44497a.java
    commons/proper/math/trunk/src/main/java/org/apache/commons/math/random/Well44497b.java
    commons/proper/math/trunk/src/main/java/org/apache/commons/math/random/Well512a.java

Modified: commons/proper/math/trunk/src/main/java/org/apache/commons/math/random/AbstractWell.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math/random/AbstractWell.java?rev=999135&r1=999134&r2=999135&view=diff
==============================================================================
--- commons/proper/math/trunk/src/main/java/org/apache/commons/math/random/AbstractWell.java
(original)
+++ commons/proper/math/trunk/src/main/java/org/apache/commons/math/random/AbstractWell.java
Mon Sep 20 22:06:07 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&ccedil;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/trunk/src/main/java/org/apache/commons/math/random/Well1024a.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math/random/Well1024a.java?rev=999135&r1=999134&r2=999135&view=diff
==============================================================================
--- commons/proper/math/trunk/src/main/java/org/apache/commons/math/random/Well1024a.java
(original)
+++ commons/proper/math/trunk/src/main/java/org/apache/commons/math/random/Well1024a.java
Mon Sep 20 22:06:07 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/trunk/src/main/java/org/apache/commons/math/random/Well19937a.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math/random/Well19937a.java?rev=999135&r1=999134&r2=999135&view=diff
==============================================================================
--- commons/proper/math/trunk/src/main/java/org/apache/commons/math/random/Well19937a.java
(original)
+++ commons/proper/math/trunk/src/main/java/org/apache/commons/math/random/Well19937a.java
Mon Sep 20 22:06:07 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/trunk/src/main/java/org/apache/commons/math/random/Well19937c.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math/random/Well19937c.java?rev=999135&r1=999134&r2=999135&view=diff
==============================================================================
--- commons/proper/math/trunk/src/main/java/org/apache/commons/math/random/Well19937c.java
(original)
+++ commons/proper/math/trunk/src/main/java/org/apache/commons/math/random/Well19937c.java
Mon Sep 20 22:06:07 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/trunk/src/main/java/org/apache/commons/math/random/Well44497a.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math/random/Well44497a.java?rev=999135&r1=999134&r2=999135&view=diff
==============================================================================
--- commons/proper/math/trunk/src/main/java/org/apache/commons/math/random/Well44497a.java
(original)
+++ commons/proper/math/trunk/src/main/java/org/apache/commons/math/random/Well44497a.java
Mon Sep 20 22:06:07 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/trunk/src/main/java/org/apache/commons/math/random/Well44497b.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math/random/Well44497b.java?rev=999135&r1=999134&r2=999135&view=diff
==============================================================================
--- commons/proper/math/trunk/src/main/java/org/apache/commons/math/random/Well44497b.java
(original)
+++ commons/proper/math/trunk/src/main/java/org/apache/commons/math/random/Well44497b.java
Mon Sep 20 22:06:07 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/trunk/src/main/java/org/apache/commons/math/random/Well512a.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math/random/Well512a.java?rev=999135&r1=999134&r2=999135&view=diff
==============================================================================
--- commons/proper/math/trunk/src/main/java/org/apache/commons/math/random/Well512a.java (original)
+++ commons/proper/math/trunk/src/main/java/org/apache/commons/math/random/Well512a.java Mon
Sep 20 22:06:07 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);
     }
 
 }



Mime
View raw message