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 * https://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
18 /*
19 * This is not the original file distributed by the Apache Software Foundation
20 * It has been modified by the Hipparchus project
21 */
22 package org.hipparchus.random;
23
24 import java.io.Serializable;
25
26 import org.hipparchus.util.FastMath;
27
28
29 /**
30 * This abstract class implements the WELL class of pseudo-random number generator
31 * from François Panneton, Pierre L'Ecuyer and Makoto Matsumoto.
32 * <p>
33 * This generator is described in a paper by François Panneton,
34 * Pierre L'Ecuyer and Makoto Matsumoto
35 * <a href="http://www.iro.umontreal.ca/~lecuyer/myftp/papers/wellrng.pdf">
36 * Improved Long-Period Generators Based on Linear Recurrences Modulo 2</a>
37 * ACM Transactions on Mathematical Software, 32, 1 (2006). The errata for the paper
38 * are in <a href="http://www.iro.umontreal.ca/~lecuyer/myftp/papers/wellrng-errata.txt">
39 * wellrng-errata.txt</a>.
40 *
41 * @see <a href="http://www.iro.umontreal.ca/~panneton/WELLRNG.html">WELL Random number generator</a>
42 */
43 public abstract class AbstractWell extends IntRandomGenerator implements Serializable {
44
45 /** Serializable version identifier. */
46 private static final long serialVersionUID = 20150223L;
47
48 /** Current index in the bytes pool. */
49 protected int index;
50
51 /** Bytes pool. */
52 protected final int[] v;
53
54 /** Creates a new random number generator.
55 * <p>The instance is initialized using the current time plus the
56 * system identity hash code of this instance as the seed.</p>
57 * @param k number of bits in the pool (not necessarily a multiple of 32)
58 */
59 protected AbstractWell(final int k) {
60 this(k, null);
61 }
62
63 /** Creates a new random number generator using a single int seed.
64 * @param k number of bits in the pool (not necessarily a multiple of 32)
65 * @param seed the initial seed (32 bits integer)
66 */
67 protected AbstractWell(final int k, final int seed) {
68 this(k, new int[] { seed });
69 }
70
71 /**
72 * Creates a new random number generator using an int array seed.
73 * @param k number of bits in the pool (not necessarily a multiple of 32)
74 * @param seed the initial seed (32 bits integers array), if null
75 * the seed of the generator will be related to the current time
76 */
77 protected AbstractWell(final int k, final int[] seed) {
78
79 final int r = calculateBlockCount(k);
80 this.v = new int[r];
81 this.index = 0;
82
83 // initialize the pool content
84 setSeed(seed);
85 }
86
87 /**
88 * Creates a new random number generator using a single long seed.
89 * @param k number of bits in the pool (not necessarily a multiple of 32)
90 * @param seed the initial seed (64 bits integer)
91 */
92 protected AbstractWell(final int k, final long seed) {
93 this(k, new int[] { (int) (seed >>> 32), (int) (seed & 0xffffffffl) });
94 }
95
96 /**
97 * Reinitialize the generator as if just built with the given int array seed.
98 * <p>
99 * The state of the generator is exactly the same as a new
100 * generator built with the same seed.
101 *
102 * @param seed the initial seed (32 bits integers array). If null
103 * the seed of the generator will be the system time plus the system identity
104 * hash code of the instance.
105 */
106 @Override
107 public void setSeed(final int[] seed) {
108 if (seed == null) {
109 setSeed(System.currentTimeMillis() + System.identityHashCode(this));
110 return;
111 }
112
113 System.arraycopy(seed, 0, v, 0, FastMath.min(seed.length, v.length));
114
115 if (seed.length < v.length) {
116 for (int i = seed.length; i < v.length; ++i) {
117 final long l = v[i - seed.length];
118 v[i] = (int) ((1812433253l * (l ^ (l >> 30)) + i) & 0xffffffffL);
119 }
120 }
121
122 index = 0;
123 clearCache(); // Clear normal deviate cache
124 }
125
126 /**
127 * Calculate the number of 32-bits blocks.
128 * @param k number of bits in the pool (not necessarily a multiple of 32)
129 * @return the number of 32-bits blocks
130 */
131 private static int calculateBlockCount(final int k) {
132 // the bits pool contains k bits, k = r w - p where r is the number
133 // of w bits blocks, w is the block size (always 32 in the original paper)
134 // and p is the number of unused bits in the last block
135 final int w = 32;
136 return (k + w - 1) / w;
137 }
138
139 /**
140 * Inner class used to store the indirection index table which is fixed
141 * for a given type of WELL class of pseudo-random number generator.
142 */
143 protected static final class IndexTable {
144 /**
145 * Index indirection table giving for each index its predecessor
146 * taking table size into account.
147 */
148 private final int[] iRm1;
149
150 /**
151 * Index indirection table giving for each index its second predecessor
152 * taking table size into account.
153 */
154 private final int[] iRm2;
155
156 /**
157 * Index indirection table giving for each index the value index + m1
158 * taking table size into account.
159 */
160 private final int[] i1;
161
162 /**
163 * Index indirection table giving for each index the value index + m2
164 * taking table size into account.
165 */
166 private final int[] i2;
167
168 /**
169 * Index indirection table giving for each index the value index + m3
170 * taking table size into account.
171 */
172 private final int[] i3;
173
174 /**
175 * Creates a new pre-calculated indirection index table.
176 * @param k number of bits in the pool (not necessarily a multiple of 32)
177 * @param m1 first parameter of the algorithm
178 * @param m2 second parameter of the algorithm
179 * @param m3 third parameter of the algorithm
180 */
181 public IndexTable(final int k, final int m1, final int m2, final int m3) {
182
183 final int r = calculateBlockCount(k);
184
185 // precompute indirection index tables. These tables are used for optimizing access
186 // they allow saving computations like "(j + r - 2) % r" with costly modulo operations
187 iRm1 = new int[r];
188 iRm2 = new int[r];
189 i1 = new int[r];
190 i2 = new int[r];
191 i3 = new int[r];
192 for (int j = 0; j < r; ++j) {
193 iRm1[j] = (j + r - 1) % r;
194 iRm2[j] = (j + r - 2) % r;
195 i1[j] = (j + m1) % r;
196 i2[j] = (j + m2) % r;
197 i3[j] = (j + m3) % r;
198 }
199 }
200
201 /**
202 * Returns the predecessor of the given index modulo the table size.
203 * @param index the index to look at
204 * @return (index - 1) % table size
205 */
206 public int getIndexPred(final int index) {
207 return iRm1[index];
208 }
209
210 /**
211 * Returns the second predecessor of the given index modulo the table size.
212 * @param index the index to look at
213 * @return (index - 2) % table size
214 */
215 public int getIndexPred2(final int index) {
216 return iRm2[index];
217 }
218
219 /**
220 * Returns index + M1 modulo the table size.
221 * @param index the index to look at
222 * @return (index + M1) % table size
223 */
224 public int getIndexM1(final int index) {
225 return i1[index];
226 }
227
228 /**
229 * Returns index + M2 modulo the table size.
230 * @param index the index to look at
231 * @return (index + M2) % table size
232 */
233 public int getIndexM2(final int index) {
234 return i2[index];
235 }
236
237 /**
238 * Returns index + M3 modulo the table size.
239 * @param index the index to look at
240 * @return (index + M3) % table size
241 */
242 public int getIndexM3(final int index) {
243 return i3[index];
244 }
245 }
246 }