View Javadoc
1   /*
2    * Licensed to the Apache Software Foundation (ASF) under one or more
3    * contributor license agreements.  See the NOTICE file distributed with
4    * this work for additional information regarding copyright ownership.
5    * The ASF licenses this file to You under the Apache License, Version 2.0
6    * (the "License"); you may not use this file except in compliance with
7    * the License.  You may obtain a copy of the License at
8    *
9    *      http://www.apache.org/licenses/LICENSE-2.0
10   *
11   * Unless required by applicable law or agreed to in writing, software
12   * distributed under the License is distributed on an "AS IS" BASIS,
13   * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14   * See the License for the specific language governing permissions and
15   * limitations under the License.
16   */
17  package org.apache.commons.rng.core.source32;
18  
19  import java.util.Arrays;
20  import org.apache.commons.rng.core.util.NumberFactory;
21  
22  /**
23   * This abstract class implements the WELL class of pseudo-random number
24   * generator from François Panneton, Pierre L'Ecuyer and Makoto
25   * Matsumoto.
26   * <p>
27   * This generator is described in a paper by Fran&ccedil;ois Panneton,
28   * Pierre L'Ecuyer and Makoto Matsumoto
29   * <a href="http://www.iro.umontreal.ca/~lecuyer/myftp/papers/wellrng.pdf">
30   * Improved Long-Period Generators Based on Linear Recurrences Modulo 2</a>
31   * ACM Transactions on Mathematical Software, 32, 1 (2006).
32   * The errata for the paper are in
33   * <a href="http://www.iro.umontreal.ca/~lecuyer/myftp/papers/wellrng-errata.txt">wellrng-errata.txt</a>.
34   * </p>
35   *
36   * @see <a href="http://www.iro.umontreal.ca/~panneton/WELLRNG.html">WELL Random number generator</a>
37   *
38   * @since 1.0
39   */
40  public abstract class AbstractWell extends IntProvider {
41      /** Block size. */
42      private static final int BLOCK_SIZE = 32;
43      /** Current index in the bytes pool. */
44      protected int index;
45      /** Bytes pool. */
46      protected final int[] v;
47  
48      /**
49       * Creates an instance with the given {@code seed}.
50       *
51       * @param k Number of bits in the pool (not necessarily a multiple of 32).
52       * @param seed Initial seed.
53       */
54      protected AbstractWell(final int k,
55                             final int[] seed) {
56          final int r = calculateBlockCount(k);
57          v = new int[r];
58          index = 0;
59  
60          // Initialize the pool content.
61          setSeedInternal(seed);
62      }
63  
64      /** {@inheritDoc} */
65      @Override
66      protected byte[] getStateInternal() {
67          final int[] s = Arrays.copyOf(v, v.length + 1);
68          s[v.length] = index;
69  
70          return composeStateInternal(NumberFactory.makeByteArray(s),
71                                      super.getStateInternal());
72      }
73  
74      /** {@inheritDoc} */
75      @Override
76      protected void setStateInternal(byte[] s) {
77          final byte[][] c = splitStateInternal(s, (v.length + 1) * 4);
78  
79          final int[] tmp = NumberFactory.makeIntArray(c[0]);
80          System.arraycopy(tmp, 0, v, 0, v.length);
81          index = tmp[v.length];
82  
83          super.setStateInternal(c[1]);
84      }
85  
86      /**
87       * Initializes the generator with the given {@code seed}.
88       *
89       * @param seed Seed. Cannot be null.
90       */
91      private void setSeedInternal(final int[] seed) {
92          System.arraycopy(seed, 0, v, 0, Math.min(seed.length, v.length));
93  
94          if (seed.length < v.length) {
95              for (int i = seed.length; i < v.length; ++i) {
96                  final long current = v[i - seed.length];
97                  v[i] = (int) ((1812433253L * (current ^ (current >> 30)) + i) & 0xffffffffL);
98              }
99          }
100 
101         index = 0;
102     }
103 
104     /**
105      * Calculate the number of 32-bits blocks.
106      *
107      * @param k Number of bits in the pool (not necessarily a multiple of 32).
108      * @return the number of 32-bits blocks.
109      */
110     private static int calculateBlockCount(final int k) {
111         // The bits pool contains k bits, k = r w - p where r is the number
112         // of w bits blocks, w is the block size (always 32 in the original paper)
113         // and p is the number of unused bits in the last block.
114         return (k + BLOCK_SIZE - 1) / BLOCK_SIZE;
115     }
116 
117     /**
118      * Inner class used to store the indirection index table which is fixed for a given
119      * type of WELL class of pseudo-random number generator.
120      */
121     protected static final class IndexTable {
122         /** Index indirection table giving for each index its predecessor taking table size into account. */
123         private final int[] iRm1;
124         /** Index indirection table giving for each index its second predecessor taking table size into account. */
125         private final int[] iRm2;
126         /** Index indirection table giving for each index the value index + m1 taking table size into account. */
127         private final int[] i1;
128         /** Index indirection table giving for each index the value index + m2 taking table size into account. */
129         private final int[] i2;
130         /** Index indirection table giving for each index the value index + m3 taking table size into account. */
131         private final int[] i3;
132 
133         /** Creates a new pre-calculated indirection index table.
134          * @param k number of bits in the pool (not necessarily a multiple of 32)
135          * @param m1 first parameter of the algorithm
136          * @param m2 second parameter of the algorithm
137          * @param m3 third parameter of the algorithm
138          */
139         public IndexTable(final int k, final int m1, final int m2, final int m3) {
140 
141             final int r = calculateBlockCount(k);
142 
143             // precompute indirection index tables. These tables are used for optimizing access
144             // they allow saving computations like "(j + r - 2) % r" with costly modulo operations
145             iRm1 = new int[r];
146             iRm2 = new int[r];
147             i1 = new int[r];
148             i2 = new int[r];
149             i3 = new int[r];
150             for (int j = 0; j < r; ++j) {
151                 iRm1[j] = (j + r - 1) % r;
152                 iRm2[j] = (j + r - 2) % r;
153                 i1[j] = (j + m1) % r;
154                 i2[j] = (j + m2) % r;
155                 i3[j] = (j + m3) % r;
156             }
157         }
158 
159         /**
160          * Returns the predecessor of the given index modulo the table size.
161          * @param index the index to look at
162          * @return (index - 1) % table size
163          */
164         public int getIndexPred(final int index) {
165             return iRm1[index];
166         }
167 
168         /**
169          * Returns the second predecessor of the given index modulo the table size.
170          * @param index the index to look at
171          * @return (index - 2) % table size
172          */
173         public int getIndexPred2(final int index) {
174             return iRm2[index];
175         }
176 
177         /**
178          * Returns index + M1 modulo the table size.
179          * @param index the index to look at
180          * @return (index + M1) % table size
181          */
182         public int getIndexM1(final int index) {
183             return i1[index];
184         }
185 
186         /**
187          * Returns index + M2 modulo the table size.
188          * @param index the index to look at
189          * @return (index + M2) % table size
190          */
191         public int getIndexM2(final int index) {
192             return i2[index];
193         }
194 
195         /**
196          * Returns index + M3 modulo the table size.
197          * @param index the index to look at
198          * @return (index + M3) % table size
199          */
200         public int getIndexM3(final int index) {
201             return i3[index];
202         }
203     }
204 }