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.source64; 018 019import java.util.Arrays; 020import org.apache.commons.rng.core.util.NumberFactory; 021 022/** 023 * This class provides the 64-bits version of the originally 32-bits 024 * {@link org.apache.commons.rng.core.source32.MersenneTwister 025 * Mersenne Twister}. 026 * 027 * <p> 028 * This class is mainly a Java port of 029 * <a href="http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/emt64.html"> 030 * the 2014/2/23 version of the generator 031 * </a> written in C by Takuji Nishimura and Makoto Matsumoto. 032 * </p> 033 * 034 * <p> 035 * Here is their original copyright: 036 * </p> 037 * 038 * <table style="background-color: #E0E0E0; width: 80%"> 039 * <caption>Mersenne Twister licence</caption> 040 * <tr><td style="padding: 10px">Copyright (C) 2004, Makoto Matsumoto and Takuji Nishimura, 041 * All rights reserved.</td></tr> 042 * 043 * <tr><td style="padding: 10px">Redistribution and use in source and binary forms, with or without 044 * modification, are permitted provided that the following conditions 045 * are met: 046 * <ol> 047 * <li>Redistributions of source code must retain the above copyright 048 * notice, this list of conditions and the following disclaimer.</li> 049 * <li>Redistributions in binary form must reproduce the above copyright 050 * notice, this list of conditions and the following disclaimer in the 051 * documentation and/or other materials provided with the distribution.</li> 052 * <li>The names of its contributors may not be used to endorse or promote 053 * products derived from this software without specific prior written 054 * permission.</li> 055 * </ol></td></tr> 056 * 057 * <tr><td style="padding: 10px"><strong>THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND 058 * CONTRIBUTORS "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, 059 * INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF 060 * MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE 061 * DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS 062 * BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, 063 * OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, 064 * PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR 065 * PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY 066 * OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT 067 * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE 068 * USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH 069 * DAMAGE.</strong></td></tr> 070 * </table> 071 * 072 * @see <a href="https://en.wikipedia.org/wiki/Mersenne_Twister">Mersenne Twister (Wikipedia)</a> 073 * @since 1.0 074 */ 075public class MersenneTwister64 extends LongProvider { 076 /** Size of the bytes pool. */ 077 private static final int NN = 312; 078 /** Period second parameter. */ 079 private static final int MM = 156; 080 /** X * MATRIX_A for X = {0, 1}. */ 081 private static final long[] MAG01 = {0x0L, 0xb5026f5aa96619e9L}; 082 /** Most significant 33 bits. */ 083 private static final long UM = 0xffffffff80000000L; 084 /** Least significant 31 bits. */ 085 private static final long LM = 0x7fffffffL; 086 /** Bytes pool. */ 087 private final long[] mt = new long[NN]; 088 /** Current index in the bytes pool. */ 089 private int mti; 090 091 /** 092 * Creates a new random number generator. 093 * 094 * @param seed Initial seed. 095 */ 096 public MersenneTwister64(long[] seed) { 097 setSeedInternal(seed); 098 } 099 100 /** {@inheritDoc} */ 101 @Override 102 protected byte[] getStateInternal() { 103 final long[] s = Arrays.copyOf(mt, NN + 1); 104 s[NN] = mti; 105 106 return composeStateInternal(NumberFactory.makeByteArray(s), 107 super.getStateInternal()); 108 } 109 110 /** {@inheritDoc} */ 111 @Override 112 protected void setStateInternal(byte[] s) { 113 final byte[][] c = splitStateInternal(s, (NN + 1) * 8); 114 115 final long[] tmp = NumberFactory.makeLongArray(c[0]); 116 System.arraycopy(tmp, 0, mt, 0, NN); 117 mti = (int) tmp[NN]; 118 119 super.setStateInternal(c[1]); 120 } 121 122 /** 123 * Initializes the generator with the given seed. 124 * 125 * @param inputSeed Initial seed. 126 */ 127 private void setSeedInternal(long[] inputSeed) { 128 // Accept empty seed. 129 final long[] seed = (inputSeed.length == 0) ? new long[1] : inputSeed; 130 131 initState(19650218L); 132 int i = 1; 133 int j = 0; 134 135 for (int k = Math.max(NN, seed.length); k != 0; k--) { 136 final long mm1 = mt[i - 1]; 137 mt[i] = (mt[i] ^ ((mm1 ^ (mm1 >>> 62)) * 0x369dea0f31a53f85L)) + seed[j] + j; // non linear 138 i++; 139 j++; 140 if (i >= NN) { 141 mt[0] = mt[NN - 1]; 142 i = 1; 143 } 144 if (j >= seed.length) { 145 j = 0; 146 } 147 } 148 for (int k = NN - 1; k != 0; k--) { 149 final long mm1 = mt[i - 1]; 150 mt[i] = (mt[i] ^ ((mm1 ^ (mm1 >>> 62)) * 0x27bb2ee687b0b0fdL)) - i; // non linear 151 i++; 152 if (i >= NN) { 153 mt[0] = mt[NN - 1]; 154 i = 1; 155 } 156 } 157 158 mt[0] = 0x8000000000000000L; // MSB is 1; assuring non-zero initial array 159 } 160 161 /** 162 * Initialize the internal state of this instance. 163 * 164 * @param seed Seed. 165 */ 166 private void initState(long seed) { 167 mt[0] = seed; 168 for (mti = 1; mti < NN; mti++) { 169 final long mm1 = mt[mti - 1]; 170 mt[mti] = 0x5851f42d4c957f2dL * (mm1 ^ (mm1 >>> 62)) + mti; 171 } 172 } 173 174 /** {@inheritDoc} */ 175 @Override 176 public long next() { 177 long x; 178 179 if (mti >= NN) { // generate NN words at one time 180 for (int i = 0; i < NN - MM; i++) { 181 x = (mt[i] & UM) | (mt[i + 1] & LM); 182 mt[i] = mt[i + MM] ^ (x >>> 1) ^ MAG01[(int)(x & 0x1L)]; 183 } 184 for (int i = NN - MM; i < NN - 1; i++) { 185 x = (mt[i] & UM) | (mt[i + 1] & LM); 186 mt[i] = mt[ i + (MM - NN)] ^ (x >>> 1) ^ MAG01[(int)(x & 0x1L)]; 187 } 188 189 x = (mt[NN - 1] & UM) | (mt[0] & LM); 190 mt[NN - 1] = mt[MM - 1] ^ (x >>> 1) ^ MAG01[(int)(x & 0x1L)]; 191 192 mti = 0; 193 } 194 195 x = mt[mti++]; 196 197 x ^= (x >>> 29) & 0x5555555555555555L; 198 x ^= (x << 17) & 0x71d67fffeda60000L; 199 x ^= (x << 37) & 0xfff7eee000000000L; 200 x ^= x >>> 43; 201 202 return x; 203 } 204}