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 */ 017 package org.apache.commons.math.random; 018 019 import java.io.Serializable; 020 021 022 /** This class implements a powerful pseudo-random number generator 023 * developed by Makoto Matsumoto and Takuji Nishimura during 024 * 1996-1997. 025 026 * <p>This generator features an extremely long period 027 * (2<sup>19937</sup>-1) and 623-dimensional equidistribution up to 32 028 * bits accuracy. The home page for this generator is located at <a 029 * href="http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/emt.html"> 030 * http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/emt.html</a>.</p> 031 032 * <p>This generator is described in a paper by Makoto Matsumoto and 033 * Takuji Nishimura in 1998: <a 034 * href="http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/ARTICLES/mt.pdf">Mersenne 035 * Twister: A 623-Dimensionally Equidistributed Uniform Pseudo-Random 036 * Number Generator</a>, ACM Transactions on Modeling and Computer 037 * Simulation, Vol. 8, No. 1, January 1998, pp 3--30</p> 038 039 * <p>This class is mainly a Java port of the 2002-01-26 version of 040 * the generator written in C by Makoto Matsumoto and Takuji 041 * Nishimura. Here is their original copyright:</p> 042 043 * <table border="0" width="80%" cellpadding="10" align="center" bgcolor="#E0E0E0"> 044 * <tr><td>Copyright (C) 1997 - 2002, Makoto Matsumoto and Takuji Nishimura, 045 * All rights reserved.</td></tr> 046 047 * <tr><td>Redistribution and use in source and binary forms, with or without 048 * modification, are permitted provided that the following conditions 049 * are met: 050 * <ol> 051 * <li>Redistributions of source code must retain the above copyright 052 * notice, this list of conditions and the following disclaimer.</li> 053 * <li>Redistributions in binary form must reproduce the above copyright 054 * notice, this list of conditions and the following disclaimer in the 055 * documentation and/or other materials provided with the distribution.</li> 056 * <li>The names of its contributors may not be used to endorse or promote 057 * products derived from this software without specific prior written 058 * permission.</li> 059 * </ol></td></tr> 060 061 * <tr><td><strong>THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND 062 * CONTRIBUTORS "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, 063 * INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF 064 * MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE 065 * DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS 066 * BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, 067 * OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, 068 * PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR 069 * PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY 070 * OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT 071 * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE 072 * USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH 073 * DAMAGE.</strong></td></tr> 074 * </table> 075 076 * @version $Revision: 902203 $ $Date: 2010-01-22 13:27:41 -0500 (Fri, 22 Jan 2010) $ 077 * @since 2.0 078 079 */ 080 public class MersenneTwister extends BitsStreamGenerator implements Serializable { 081 082 /** Serializable version identifier. */ 083 private static final long serialVersionUID = 8661194735290153518L; 084 085 /** Size of the bytes pool. */ 086 private static final int N = 624; 087 088 /** Period second parameter. */ 089 private static final int M = 397; 090 091 /** X * MATRIX_A for X = {0, 1}. */ 092 private static final int[] MAG01 = { 0x0, 0x9908b0df }; 093 094 /** Bytes pool. */ 095 private int[] mt; 096 097 /** Current index in the bytes pool. */ 098 private int mti; 099 100 /** Creates a new random number generator. 101 * <p>The instance is initialized using the current time as the 102 * seed.</p> 103 */ 104 public MersenneTwister() { 105 mt = new int[N]; 106 setSeed(System.currentTimeMillis()); 107 } 108 109 /** Creates a new random number generator using a single int seed. 110 * @param seed the initial seed (32 bits integer) 111 */ 112 public MersenneTwister(int seed) { 113 mt = new int[N]; 114 setSeed(seed); 115 } 116 117 /** Creates a new random number generator using an int array seed. 118 * @param seed the initial seed (32 bits integers array), if null 119 * the seed of the generator will be related to the current time 120 */ 121 public MersenneTwister(int[] seed) { 122 mt = new int[N]; 123 setSeed(seed); 124 } 125 126 /** Creates a new random number generator using a single long seed. 127 * @param seed the initial seed (64 bits integer) 128 */ 129 public MersenneTwister(long seed) { 130 mt = new int[N]; 131 setSeed(seed); 132 } 133 134 /** Reinitialize the generator as if just built with the given int seed. 135 * <p>The state of the generator is exactly the same as a new 136 * generator built with the same seed.</p> 137 * @param seed the initial seed (32 bits integer) 138 */ 139 @Override 140 public void setSeed(int seed) { 141 // we use a long masked by 0xffffffffL as a poor man unsigned int 142 long longMT = seed; 143 mt[0]= (int) longMT; 144 for (mti = 1; mti < N; ++mti) { 145 // See Knuth TAOCP Vol2. 3rd Ed. P.106 for multiplier. 146 // initializer from the 2002-01-09 C version by Makoto Matsumoto 147 longMT = (1812433253l * (longMT ^ (longMT >> 30)) + mti) & 0xffffffffL; 148 mt[mti]= (int) longMT; 149 } 150 } 151 152 /** Reinitialize the generator as if just built with the given int array seed. 153 * <p>The state of the generator is exactly the same as a new 154 * generator built with the same seed.</p> 155 * @param seed the initial seed (32 bits integers array), if null 156 * the seed of the generator will be related to the current time 157 */ 158 @Override 159 public void setSeed(int[] seed) { 160 161 if (seed == null) { 162 setSeed(System.currentTimeMillis()); 163 return; 164 } 165 166 setSeed(19650218); 167 int i = 1; 168 int j = 0; 169 170 for (int k = Math.max(N, seed.length); k != 0; k--) { 171 long l0 = (mt[i] & 0x7fffffffl) | ((mt[i] < 0) ? 0x80000000l : 0x0l); 172 long l1 = (mt[i-1] & 0x7fffffffl) | ((mt[i-1] < 0) ? 0x80000000l : 0x0l); 173 long l = (l0 ^ ((l1 ^ (l1 >> 30)) * 1664525l)) + seed[j] + j; // non linear 174 mt[i] = (int) (l & 0xffffffffl); 175 i++; j++; 176 if (i >= N) { 177 mt[0] = mt[N - 1]; 178 i = 1; 179 } 180 if (j >= seed.length) { 181 j = 0; 182 } 183 } 184 185 for (int k = N - 1; k != 0; k--) { 186 long l0 = (mt[i] & 0x7fffffffl) | ((mt[i] < 0) ? 0x80000000l : 0x0l); 187 long l1 = (mt[i-1] & 0x7fffffffl) | ((mt[i-1] < 0) ? 0x80000000l : 0x0l); 188 long l = (l0 ^ ((l1 ^ (l1 >> 30)) * 1566083941l)) - i; // non linear 189 mt[i] = (int) (l & 0xffffffffL); 190 i++; 191 if (i >= N) { 192 mt[0] = mt[N - 1]; 193 i = 1; 194 } 195 } 196 197 mt[0] = 0x80000000; // MSB is 1; assuring non-zero initial array 198 199 } 200 201 /** Reinitialize the generator as if just built with the given long seed. 202 * <p>The state of the generator is exactly the same as a new 203 * generator built with the same seed.</p> 204 * @param seed the initial seed (64 bits integer) 205 */ 206 @Override 207 public void setSeed(long seed) { 208 setSeed(new int[] { (int) (seed >>> 32), (int) (seed & 0xffffffffl) }); 209 } 210 211 /** Generate next pseudorandom number. 212 * <p>This method is the core generation algorithm. It is used by all the 213 * public generation methods for the various primitive types {@link 214 * #nextBoolean()}, {@link #nextBytes(byte[])}, {@link #nextDouble()}, 215 * {@link #nextFloat()}, {@link #nextGaussian()}, {@link #nextInt()}, 216 * {@link #next(int)} and {@link #nextLong()}.</p> 217 * @param bits number of random bits to produce 218 * @return random bits generated 219 */ 220 @Override 221 protected int next(int bits) { 222 223 int y; 224 225 if (mti >= N) { // generate N words at one time 226 int mtNext = mt[0]; 227 for (int k = 0; k < N - M; ++k) { 228 int mtCurr = mtNext; 229 mtNext = mt[k + 1]; 230 y = (mtCurr & 0x80000000) | (mtNext & 0x7fffffff); 231 mt[k] = mt[k + M] ^ (y >>> 1) ^ MAG01[y & 0x1]; 232 } 233 for (int k = N - M; k < N - 1; ++k) { 234 int mtCurr = mtNext; 235 mtNext = mt[k + 1]; 236 y = (mtCurr & 0x80000000) | (mtNext & 0x7fffffff); 237 mt[k] = mt[k + (M - N)] ^ (y >>> 1) ^ MAG01[y & 0x1]; 238 } 239 y = (mtNext & 0x80000000) | (mt[0] & 0x7fffffff); 240 mt[N - 1] = mt[M - 1] ^ (y >>> 1) ^ MAG01[y & 0x1]; 241 242 mti = 0; 243 } 244 245 y = mt[mti++]; 246 247 // tempering 248 y ^= y >>> 11; 249 y ^= (y << 7) & 0x9d2c5680; 250 y ^= (y << 15) & 0xefc60000; 251 y ^= y >>> 18; 252 253 return y >>> (32 - bits); 254 255 } 256 257 }