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ç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 }