001/* 002 * Licensed to the Apache Software Foundation (ASF) under one or more 003 * contributor license agreements. See the NOTICE file distributed with 004 * this work for additional information regarding copyright ownership. 005 * The ASF licenses this file to You under the Apache License, Version 2.0 006 * (the "License"); you may not use this file except in compliance with 007 * the License. You may obtain a copy of the License at 008 * 009 * http://www.apache.org/licenses/LICENSE-2.0 010 * 011 * Unless required by applicable law or agreed to in writing, software 012 * distributed under the License is distributed on an "AS IS" BASIS, 013 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. 014 * See the License for the specific language governing permissions and 015 * limitations under the License. 016 */ 017package org.apache.commons.rng.core.source32; 018 019import java.util.Arrays; 020import org.apache.commons.rng.core.util.NumberFactory; 021 022/** 023 * This abstract class implements the WELL class of pseudo-random number 024 * generator from François Panneton, Pierre L'Ecuyer and Makoto 025 * Matsumoto. 026 * <p> 027 * This generator is described in a paper by François Panneton, 028 * Pierre L'Ecuyer and Makoto Matsumoto 029 * <a href="http://www.iro.umontreal.ca/~lecuyer/myftp/papers/wellrng.pdf"> 030 * Improved Long-Period Generators Based on Linear Recurrences Modulo 2</a> 031 * ACM Transactions on Mathematical Software, 32, 1 (2006). 032 * The errata for the paper are in 033 * <a href="http://www.iro.umontreal.ca/~lecuyer/myftp/papers/wellrng-errata.txt">wellrng-errata.txt</a>. 034 * </p> 035 * 036 * @see <a href="http://www.iro.umontreal.ca/~panneton/WELLRNG.html">WELL Random number generator</a> 037 * 038 * @since 1.0 039 */ 040public abstract class AbstractWell extends IntProvider { 041 /** Block size. */ 042 private static final int BLOCK_SIZE = 32; 043 /** Current index in the bytes pool. */ 044 protected int index; 045 /** Bytes pool. */ 046 protected final int[] v; 047 048 /** 049 * Creates an instance with the given {@code seed}. 050 * 051 * @param k Number of bits in the pool (not necessarily a multiple of 32). 052 * @param seed Initial seed. 053 */ 054 protected AbstractWell(final int k, 055 final int[] seed) { 056 final int r = calculateBlockCount(k); 057 v = new int[r]; 058 index = 0; 059 060 // Initialize the pool content. 061 setSeedInternal(seed); 062 } 063 064 /** {@inheritDoc} */ 065 @Override 066 protected byte[] getStateInternal() { 067 final int[] s = Arrays.copyOf(v, v.length + 1); 068 s[v.length] = index; 069 070 return NumberFactory.makeByteArray(s); 071 } 072 073 /** {@inheritDoc} */ 074 @Override 075 protected void setStateInternal(byte[] s) { 076 checkStateSize(s, (v.length + 1) * 4); 077 078 final int[] tmp = NumberFactory.makeIntArray(s); 079 System.arraycopy(tmp, 0, v, 0, v.length); 080 index = tmp[v.length]; 081 } 082 083 /** 084 * Initializes the generator with the given {@code seed}. 085 * 086 * @param seed Seed. Cannot be null. 087 */ 088 private void setSeedInternal(final int[] seed) { 089 System.arraycopy(seed, 0, v, 0, Math.min(seed.length, v.length)); 090 091 if (seed.length < v.length) { 092 for (int i = seed.length; i < v.length; ++i) { 093 final long current = v[i - seed.length]; 094 v[i] = (int) ((1812433253L * (current ^ (current >> 30)) + i) & 0xffffffffL); 095 } 096 } 097 098 index = 0; 099 } 100 101 /** 102 * Calculate the number of 32-bits blocks. 103 * 104 * @param k Number of bits in the pool (not necessarily a multiple of 32). 105 * @return the number of 32-bits blocks. 106 */ 107 private static int calculateBlockCount(final int k) { 108 // The bits pool contains k bits, k = r w - p where r is the number 109 // of w bits blocks, w is the block size (always 32 in the original paper) 110 // and p is the number of unused bits in the last block. 111 return (k + BLOCK_SIZE - 1) / BLOCK_SIZE; 112 } 113 114 /** 115 * Inner class used to store the indirection index table which is fixed for a given 116 * type of WELL class of pseudo-random number generator. 117 */ 118 protected static final class IndexTable { 119 /** Index indirection table giving for each index its predecessor taking table size into account. */ 120 private final int[] iRm1; 121 /** Index indirection table giving for each index its second predecessor taking table size into account. */ 122 private final int[] iRm2; 123 /** Index indirection table giving for each index the value index + m1 taking table size into account. */ 124 private final int[] i1; 125 /** Index indirection table giving for each index the value index + m2 taking table size into account. */ 126 private final int[] i2; 127 /** Index indirection table giving for each index the value index + m3 taking table size into account. */ 128 private final int[] i3; 129 130 /** Creates a new pre-calculated indirection index table. 131 * @param k number of bits in the pool (not necessarily a multiple of 32) 132 * @param m1 first parameter of the algorithm 133 * @param m2 second parameter of the algorithm 134 * @param m3 third parameter of the algorithm 135 */ 136 public IndexTable(final int k, final int m1, final int m2, final int m3) { 137 138 final int r = calculateBlockCount(k); 139 140 // precompute indirection index tables. These tables are used for optimizing access 141 // they allow saving computations like "(j + r - 2) % r" with costly modulo operations 142 iRm1 = new int[r]; 143 iRm2 = new int[r]; 144 i1 = new int[r]; 145 i2 = new int[r]; 146 i3 = new int[r]; 147 for (int j = 0; j < r; ++j) { 148 iRm1[j] = (j + r - 1) % r; 149 iRm2[j] = (j + r - 2) % r; 150 i1[j] = (j + m1) % r; 151 i2[j] = (j + m2) % r; 152 i3[j] = (j + m3) % r; 153 } 154 } 155 156 /** 157 * Returns the predecessor of the given index modulo the table size. 158 * @param index the index to look at 159 * @return (index - 1) % table size 160 */ 161 public int getIndexPred(final int index) { 162 return iRm1[index]; 163 } 164 165 /** 166 * Returns the second predecessor of the given index modulo the table size. 167 * @param index the index to look at 168 * @return (index - 2) % table size 169 */ 170 public int getIndexPred2(final int index) { 171 return iRm2[index]; 172 } 173 174 /** 175 * Returns index + M1 modulo the table size. 176 * @param index the index to look at 177 * @return (index + M1) % table size 178 */ 179 public int getIndexM1(final int index) { 180 return i1[index]; 181 } 182 183 /** 184 * Returns index + M2 modulo the table size. 185 * @param index the index to look at 186 * @return (index + M2) % table size 187 */ 188 public int getIndexM2(final int index) { 189 return i2[index]; 190 } 191 192 /** 193 * Returns index + M3 modulo the table size. 194 * @param index the index to look at 195 * @return (index + M3) % table size 196 */ 197 public int getIndexM3(final int index) { 198 return i3[index]; 199 } 200 } 201}