| 1 | // MersenneTwister.h | 
| 2 | // Mersenne Twister random number generator -- a C++ class MTRand | 
| 3 | // Based on code by Makoto Matsumoto, Takuji Nishimura, and Shawn Cokus | 
| 4 | // Richard J. Wagner  v1.0  15 May 2003  rjwagner@writeme.com | 
| 5 |  | 
| 6 | // The Mersenne Twister is an algorithm for generating random numbers.  It | 
| 7 | // was designed with consideration of the flaws in various other generators. | 
| 8 | // The period, 2^19937-1, and the order of equidistribution, 623 dimensions, | 
| 9 | // are far greater.  The generator is also fast; it avoids multiplication and | 
| 10 | // division, and it benefits from caches and pipelines.  For more information | 
| 11 | // see the inventors' web page at http://www.math.keio.ac.jp/~matumoto/emt.html | 
| 12 |  | 
| 13 | // Reference | 
| 14 | // M. Matsumoto and T. Nishimura, "Mersenne Twister: A 623-Dimensionally | 
| 15 | // Equidistributed Uniform Pseudo-Random Number Generator", ACM Transactions on | 
| 16 | // Modeling and Computer Simulation, Vol. 8, No. 1, January 1998, pp 3-30. | 
| 17 |  | 
| 18 | // Copyright (C) 1997 - 2002, Makoto Matsumoto and Takuji Nishimura, | 
| 19 | // Copyright (C) 2000 - 2003, Richard J. Wagner | 
| 20 | // All rights reserved. | 
| 21 | // | 
| 22 | // Redistribution and use in source and binary forms, with or without | 
| 23 | // modification, are permitted provided that the following conditions | 
| 24 | // are met: | 
| 25 | // | 
| 26 | //   1. Redistributions of source code must retain the above copyright | 
| 27 | //      notice, this list of conditions and the following disclaimer. | 
| 28 | // | 
| 29 | //   2. Redistributions in binary form must reproduce the above copyright | 
| 30 | //      notice, this list of conditions and the following disclaimer in the | 
| 31 | //      documentation and/or other materials provided with the distribution. | 
| 32 | // | 
| 33 | //   3. The names of its contributors may not be used to endorse or promote | 
| 34 | //      products derived from this software without specific prior written | 
| 35 | //      permission. | 
| 36 | // | 
| 37 | // THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS | 
| 38 | // "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT | 
| 39 | // LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR | 
| 40 | // A PARTICULAR PURPOSE ARE DISCLAIMED.  IN NO EVENT SHALL THE COPYRIGHT OWNER OR | 
| 41 | // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, | 
| 42 | // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, | 
| 43 | // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR | 
| 44 | // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF | 
| 45 | // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING | 
| 46 | // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS | 
| 47 | // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. | 
| 48 |  | 
| 49 | // The original code included the following notice: | 
| 50 | // | 
| 51 | //     When you use this, send an email to: matumoto@math.keio.ac.jp | 
| 52 | //     with an appropriate reference to your work. | 
| 53 | // | 
| 54 | // It would be nice to CC: rjwagner@writeme.com and Cokus@math.washington.edu | 
| 55 | // when you write. | 
| 56 |  | 
| 57 | #ifndef MERSENNETWISTER_H | 
| 58 | #define MERSENNETWISTER_H | 
| 59 |  | 
| 60 | // Not thread safe (unless auto-initialization is avoided and each thread has | 
| 61 | // its own MTRand object) | 
| 62 |  | 
| 63 | #include <cassert> | 
| 64 | #include <iostream> | 
| 65 | #include <limits.h> | 
| 66 | #include <stdio.h> | 
| 67 | #include <time.h> | 
| 68 | #include <math.h> | 
| 69 | #include <vector> | 
| 70 | namespace oopse { | 
| 71 |  | 
| 72 | class MTRand { | 
| 73 | // Data | 
| 74 | public: | 
| 75 | typedef unsigned long uint32;  // unsigned integer type, at least 32 bits | 
| 76 |  | 
| 77 | enum { N = 624 };       // length of state vector | 
| 78 | enum { SAVE = N + 1 };  // length of array for save() | 
| 79 |  | 
| 80 | private: | 
| 81 | enum { M = 397 };  // period parameter | 
| 82 |  | 
| 83 | uint32 state[N];   // internal state | 
| 84 | uint32 *pNext;     // next value to get from state | 
| 85 | int left;          // number of values left before reload needed | 
| 86 | int nstrides_; | 
| 87 | int stride_; | 
| 88 |  | 
| 89 | //Methods | 
| 90 | public: | 
| 91 | MTRand( const uint32& oneSeed, int nstrides = 1, int stride = 0);  // initialize with a simple uint32 | 
| 92 | MTRand( uint32 *const bigSeed, uint32 const seedLength = N, int nstrides = 1, int stride = 0);  // or an array | 
| 93 | MTRand(int nstrides = 1, int stride = 0);  // auto-initialize with /dev/urandom or time() and clock() | 
| 94 |  | 
| 95 | // Do NOT use for CRYPTOGRAPHY without securely hashing several returned | 
| 96 | // values together, otherwise the generator state can be learned after | 
| 97 | // reading 624 consecutive values. | 
| 98 |  | 
| 99 | // Access to 32-bit random numbers | 
| 100 | double rand();                          // real number in [0,1] | 
| 101 | double rand( const double& n );         // real number in [0,n] | 
| 102 | double randExc();                       // real number in [0,1) | 
| 103 | double randExc( const double& n );      // real number in [0,n) | 
| 104 | double randDblExc();                    // real number in (0,1) | 
| 105 | double randDblExc( const double& n );   // real number in (0,n) | 
| 106 | uint32 randInt();                       // integer in [0,2^32-1] (modified for striding) | 
| 107 | uint32 rawRandInt();                    // original randInt | 
| 108 | uint32 randInt( const uint32& n );      // integer in [0,n] for n < 2^32 | 
| 109 | double operator()() { return rand(); }  // same as rand() | 
| 110 |  | 
| 111 | // Access to 53-bit random numbers (capacity of IEEE double precision) | 
| 112 | double rand53();  // real number in [0,1) | 
| 113 |  | 
| 114 | // Access to nonuniform random number distributions | 
| 115 | double randNorm( const double& mean = 0.0, const double& variance = 0.0 ); | 
| 116 |  | 
| 117 | // Re-seeding functions with same behavior as initializers | 
| 118 | void seed( const uint32 oneSeed ); | 
| 119 | void seed( uint32 *const bigSeed, const uint32 seedLength = N ); | 
| 120 | void seed(); | 
| 121 |  | 
| 122 | std::vector<uint32>generateSeeds(); | 
| 123 |  | 
| 124 | // Saving and loading generator state | 
| 125 | void save( uint32* saveArray ) const;  // to array of size SAVE | 
| 126 | void load( uint32 *const loadArray );  // from such array | 
| 127 | friend std::ostream& operator<<( std::ostream& os, const MTRand& mtrand ); | 
| 128 | friend std::istream& operator>>( std::istream& is, MTRand& mtrand ); | 
| 129 |  | 
| 130 | protected: | 
| 131 | void initialize( const uint32 oneSeed ); | 
| 132 | void reload(); | 
| 133 | uint32 hiBit( const uint32& u ) const { return u & 0x80000000UL; } | 
| 134 | uint32 loBit( const uint32& u ) const { return u & 0x00000001UL; } | 
| 135 | uint32 loBits( const uint32& u ) const { return u & 0x7fffffffUL; } | 
| 136 | uint32 mixBits( const uint32& u, const uint32& v ) const | 
| 137 | { return hiBit(u) | loBits(v); } | 
| 138 | uint32 twist( const uint32& m, const uint32& s0, const uint32& s1 ) const | 
| 139 | { return m ^ (mixBits(s0,s1)>>1) ^ (-loBit(s1) & 0x9908b0dfUL); } | 
| 140 | static uint32 hash( time_t t, clock_t c ); | 
| 141 | }; | 
| 142 |  | 
| 143 |  | 
| 144 | inline MTRand::MTRand( const uint32& oneSeed, int nstrides, int stride) : nstrides_(nstrides), stride_(stride) { | 
| 145 | assert(stride_ < nstrides_ && stride_ >= 0); | 
| 146 | seed(oneSeed); | 
| 147 | } | 
| 148 |  | 
| 149 | inline MTRand::MTRand( uint32 *const bigSeed, const uint32 seedLength, int nstrides, int stride) : nstrides_(nstrides), stride_(stride) { | 
| 150 | assert(stride_ < nstrides_ && stride_ >= 0); | 
| 151 | seed(bigSeed,seedLength); | 
| 152 | } | 
| 153 |  | 
| 154 | inline MTRand::MTRand(int nstrides, int stride) : nstrides_(nstrides), stride_(stride){ | 
| 155 | assert(stride_ < nstrides_ && stride_ >= 0); | 
| 156 | seed(); | 
| 157 | } | 
| 158 |  | 
| 159 | inline double MTRand::rand() | 
| 160 | { return double(randInt()) * (1.0/4294967295.0); } | 
| 161 |  | 
| 162 | inline double MTRand::rand( const double& n ) | 
| 163 | { return rand() * n; } | 
| 164 |  | 
| 165 | inline double MTRand::randExc() | 
| 166 | { return double(randInt()) * (1.0/4294967296.0); } | 
| 167 |  | 
| 168 | inline double MTRand::randExc( const double& n ) | 
| 169 | { return randExc() * n; } | 
| 170 |  | 
| 171 | inline double MTRand::randDblExc() | 
| 172 | { return ( double(randInt()) + 0.5 ) * (1.0/4294967296.0); } | 
| 173 |  | 
| 174 | inline double MTRand::randDblExc( const double& n ) | 
| 175 | { return randDblExc() * n; } | 
| 176 |  | 
| 177 | inline double MTRand::rand53() | 
| 178 | { | 
| 179 | uint32 a = randInt() >> 5, b = randInt() >> 6; | 
| 180 | return ( a * 67108864.0 + b ) * (1.0/9007199254740992.0);  // by Isaku Wada | 
| 181 | } | 
| 182 |  | 
| 183 | inline double MTRand::randNorm( const double& mean, const double& variance ) | 
| 184 | { | 
| 185 | // Return a real number from a normal (Gaussian) distribution with given | 
| 186 | // mean and variance by Box-Muller method | 
| 187 | double r = sqrt( -2.0 * log( 1.0-randDblExc()) ) * variance; | 
| 188 | double phi = 2.0 * 3.14159265358979323846264338328 * randExc(); | 
| 189 | return mean + r * cos(phi); | 
| 190 | } | 
| 191 |  | 
| 192 | /** | 
| 193 | * This function is modified from the original to allow for random | 
| 194 | * streams on parallel jobs.  It now takes numbers from by striding | 
| 195 | * through the random stream and picking up only one of the random | 
| 196 | * numbers per nstrides_.  The number it picks is the stride_'th | 
| 197 | * number in the stride sequence. | 
| 198 | */ | 
| 199 | inline MTRand::uint32 MTRand::randInt() { | 
| 200 |  | 
| 201 | std::vector<uint32> ranNums(nstrides_); | 
| 202 |  | 
| 203 | for (int i = 0; i < nstrides_; ++i) { | 
| 204 | ranNums[i] = rawRandInt(); | 
| 205 | } | 
| 206 |  | 
| 207 | return ranNums[stride_]; | 
| 208 | } | 
| 209 |  | 
| 210 | /** | 
| 211 | * This is the original randInt function which implements the mersenne | 
| 212 | * twister. | 
| 213 | */ | 
| 214 | inline MTRand::uint32 MTRand::rawRandInt() | 
| 215 | { | 
| 216 | // Pull a 32-bit integer from the generator state | 
| 217 | // Every other access function simply transforms the numbers extracted here | 
| 218 |  | 
| 219 | if( left == 0 ) reload(); | 
| 220 | --left; | 
| 221 |  | 
| 222 | register uint32 s1; | 
| 223 | s1 = *pNext++; | 
| 224 | s1 ^= (s1 >> 11); | 
| 225 | s1 ^= (s1 <<  7) & 0x9d2c5680UL; | 
| 226 | s1 ^= (s1 << 15) & 0xefc60000UL; | 
| 227 | return ( s1 ^ (s1 >> 18) ); | 
| 228 | } | 
| 229 |  | 
| 230 | inline MTRand::uint32 MTRand::randInt( const uint32& n ) | 
| 231 | { | 
| 232 | // Find which bits are used in n | 
| 233 | // Optimized by Magnus Jonsson (magnus@smartelectronix.com) | 
| 234 | uint32 used = n; | 
| 235 | used |= used >> 1; | 
| 236 | used |= used >> 2; | 
| 237 | used |= used >> 4; | 
| 238 | used |= used >> 8; | 
| 239 | used |= used >> 16; | 
| 240 |  | 
| 241 | // Draw numbers until one is found in [0,n] | 
| 242 | uint32 i; | 
| 243 | do | 
| 244 | i = randInt() & used;  // toss unused bits to shorten search | 
| 245 | while( i > n ); | 
| 246 | return i; | 
| 247 | } | 
| 248 |  | 
| 249 |  | 
| 250 | inline void MTRand::seed( const uint32 oneSeed ) | 
| 251 | { | 
| 252 | // Seed the generator with a simple uint32 | 
| 253 | initialize(oneSeed); | 
| 254 | reload(); | 
| 255 | } | 
| 256 |  | 
| 257 |  | 
| 258 | inline void MTRand::seed( uint32 *const bigSeed, const uint32 seedLength ) | 
| 259 | { | 
| 260 | // Seed the generator with an array of uint32's | 
| 261 | // There are 2^19937-1 possible initial states.  This function allows | 
| 262 | // all of those to be accessed by providing at least 19937 bits (with a | 
| 263 | // default seed length of N = 624 uint32's).  Any bits above the lower 32 | 
| 264 | // in each element are discarded. | 
| 265 | // Just call seed() if you want to get array from /dev/urandom | 
| 266 | initialize(19650218UL); | 
| 267 | register int i = 1; | 
| 268 | register uint32 j = 0; | 
| 269 | register int k = ( N > seedLength ? N : seedLength ); | 
| 270 | for( ; k; --k ) | 
| 271 | { | 
| 272 | state[i] = | 
| 273 | state[i] ^ ( (state[i-1] ^ (state[i-1] >> 30)) * 1664525UL ); | 
| 274 | state[i] += ( bigSeed[j] & 0xffffffffUL ) + j; | 
| 275 | state[i] &= 0xffffffffUL; | 
| 276 | ++i;  ++j; | 
| 277 | if( i >= N ) { state[0] = state[N-1];  i = 1; } | 
| 278 | if( j >= seedLength ) j = 0; | 
| 279 | } | 
| 280 | for( k = N - 1; k; --k ) | 
| 281 | { | 
| 282 | state[i] = | 
| 283 | state[i] ^ ( (state[i-1] ^ (state[i-1] >> 30)) * 1566083941UL ); | 
| 284 | state[i] -= i; | 
| 285 | state[i] &= 0xffffffffUL; | 
| 286 | ++i; | 
| 287 | if( i >= N ) { state[0] = state[N-1];  i = 1; } | 
| 288 | } | 
| 289 | state[0] = 0x80000000UL;  // MSB is 1, assuring non-zero initial array | 
| 290 | reload(); | 
| 291 | } | 
| 292 |  | 
| 293 |  | 
| 294 | inline void MTRand::seed() | 
| 295 | { | 
| 296 | std::vector<uint32> seeds; | 
| 297 |  | 
| 298 | seeds = generateSeeds(); | 
| 299 |  | 
| 300 | if (seeds.size() == 1) { | 
| 301 | seed( seeds[0] ); | 
| 302 | } else { | 
| 303 | seed( &seeds[0], seeds.size() ); | 
| 304 | } | 
| 305 | } | 
| 306 |  | 
| 307 |  | 
| 308 | inline std::vector<MTRand::uint32> MTRand::generateSeeds() { | 
| 309 | // Seed the generator with an array from /dev/urandom if available | 
| 310 | // Otherwise use a hash of time() and clock() values | 
| 311 |  | 
| 312 | std::vector<uint32> bigSeed; | 
| 313 |  | 
| 314 | // First try getting an array from /dev/urandom | 
| 315 | FILE* urandom = fopen( "/dev/urandom", "rb" ); | 
| 316 | if( urandom ) | 
| 317 | { | 
| 318 | bigSeed.resize(N); | 
| 319 | register uint32 *s = &bigSeed[0]; | 
| 320 | register int i = N; | 
| 321 | register bool success = true; | 
| 322 | while( success && i-- ) | 
| 323 | success = fread( s++, sizeof(uint32), 1, urandom ); | 
| 324 | fclose(urandom); | 
| 325 | if( success ) { return bigSeed; } | 
| 326 | } | 
| 327 |  | 
| 328 | // Was not successful, so use time() and clock() instead | 
| 329 |  | 
| 330 | bigSeed.push_back(hash( time(NULL), clock())); | 
| 331 | return bigSeed; | 
| 332 | } | 
| 333 |  | 
| 334 |  | 
| 335 | inline void MTRand::initialize( const uint32 seed ) | 
| 336 | { | 
| 337 | // Initialize generator state with seed | 
| 338 | // See Knuth TAOCP Vol 2, 3rd Ed, p.106 for multiplier. | 
| 339 | // In previous versions, most significant bits (MSBs) of the seed affect | 
| 340 | // only MSBs of the state array.  Modified 9 Jan 2002 by Makoto Matsumoto. | 
| 341 | register uint32 *s = state; | 
| 342 | register uint32 *r = state; | 
| 343 | register int i = 1; | 
| 344 | *s++ = seed & 0xffffffffUL; | 
| 345 | for( ; i < N; ++i ) | 
| 346 | { | 
| 347 | *s++ = ( 1812433253UL * ( *r ^ (*r >> 30) ) + i ) & 0xffffffffUL; | 
| 348 | r++; | 
| 349 | } | 
| 350 | } | 
| 351 |  | 
| 352 |  | 
| 353 | inline void MTRand::reload() | 
| 354 | { | 
| 355 | // Generate N new values in state | 
| 356 | // Made clearer and faster by Matthew Bellew (matthew.bellew@home.com) | 
| 357 | register uint32 *p = state; | 
| 358 | register int i; | 
| 359 | for( i = N - M; i--; ++p ) | 
| 360 | *p = twist( p[M], p[0], p[1] ); | 
| 361 | for( i = M; --i; ++p ) | 
| 362 | *p = twist( p[M-N], p[0], p[1] ); | 
| 363 | *p = twist( p[M-N], p[0], state[0] ); | 
| 364 |  | 
| 365 | left = N, pNext = state; | 
| 366 | } | 
| 367 |  | 
| 368 |  | 
| 369 | inline MTRand::uint32 MTRand::hash( time_t t, clock_t c ) | 
| 370 | { | 
| 371 | // Get a uint32 from t and c | 
| 372 | // Better than uint32(x) in case x is floating point in [0,1] | 
| 373 | // Based on code by Lawrence Kirby (fred@genesis.demon.co.uk) | 
| 374 |  | 
| 375 | static uint32 differ = 0;  // guarantee time-based seeds will change | 
| 376 |  | 
| 377 | uint32 h1 = 0; | 
| 378 | unsigned char *p = (unsigned char *) &t; | 
| 379 | for( size_t i = 0; i < sizeof(t); ++i ) | 
| 380 | { | 
| 381 | h1 *= UCHAR_MAX + 2U; | 
| 382 | h1 += p[i]; | 
| 383 | } | 
| 384 | uint32 h2 = 0; | 
| 385 | p = (unsigned char *) &c; | 
| 386 | for( size_t j = 0; j < sizeof(c); ++j ) | 
| 387 | { | 
| 388 | h2 *= UCHAR_MAX + 2U; | 
| 389 | h2 += p[j]; | 
| 390 | } | 
| 391 | return ( h1 + differ++ ) ^ h2; | 
| 392 | } | 
| 393 |  | 
| 394 |  | 
| 395 | inline void MTRand::save( uint32* saveArray ) const | 
| 396 | { | 
| 397 | register uint32 *sa = saveArray; | 
| 398 | register const uint32 *s = state; | 
| 399 | register int i = N; | 
| 400 | for( ; i--; *sa++ = *s++ ) {} | 
| 401 | *sa = left; | 
| 402 | } | 
| 403 |  | 
| 404 |  | 
| 405 | inline void MTRand::load( uint32 *const loadArray ) | 
| 406 | { | 
| 407 | register uint32 *s = state; | 
| 408 | register uint32 *la = loadArray; | 
| 409 | register int i = N; | 
| 410 | for( ; i--; *s++ = *la++ ) {} | 
| 411 | left = *la; | 
| 412 | pNext = &state[N-left]; | 
| 413 | } | 
| 414 |  | 
| 415 |  | 
| 416 | inline std::ostream& operator<<( std::ostream& os, const MTRand& mtrand ) | 
| 417 | { | 
| 418 | register const MTRand::uint32 *s = mtrand.state; | 
| 419 | register int i = mtrand.N; | 
| 420 | for( ; i--; os << *s++ << "\t" ) {} | 
| 421 | return os << mtrand.left; | 
| 422 | } | 
| 423 |  | 
| 424 |  | 
| 425 | inline std::istream& operator>>( std::istream& is, MTRand& mtrand ) | 
| 426 | { | 
| 427 | register MTRand::uint32 *s = mtrand.state; | 
| 428 | register int i = mtrand.N; | 
| 429 | for( ; i--; is >> *s++ ) {} | 
| 430 | is >> mtrand.left; | 
| 431 | mtrand.pNext = &mtrand.state[mtrand.N-mtrand.left]; | 
| 432 | return is; | 
| 433 | } | 
| 434 |  | 
| 435 | } | 
| 436 | #endif  // MERSENNETWISTER_H | 
| 437 |  | 
| 438 | // Change log: | 
| 439 | // | 
| 440 | // v0.1 - First release on 15 May 2000 | 
| 441 | //      - Based on code by Makoto Matsumoto, Takuji Nishimura, and Shawn Cokus | 
| 442 | //      - Translated from C to C++ | 
| 443 | //      - Made completely ANSI compliant | 
| 444 | //      - Designed convenient interface for initialization, seeding, and | 
| 445 | //        obtaining numbers in default or user-defined ranges | 
| 446 | //      - Added automatic seeding from /dev/urandom or time() and clock() | 
| 447 | //      - Provided functions for saving and loading generator state | 
| 448 | // | 
| 449 | // v0.2 - Fixed bug which reloaded generator one step too late | 
| 450 | // | 
| 451 | // v0.3 - Switched to clearer, faster reload() code from Matthew Bellew | 
| 452 | // | 
| 453 | // v0.4 - Removed trailing newline in saved generator format to be consistent | 
| 454 | //        with output format of built-in types | 
| 455 | // | 
| 456 | // v0.5 - Improved portability by replacing static const int's with enum's and | 
| 457 | //        clarifying return values in seed(); suggested by Eric Heimburg | 
| 458 | //      - Removed MAXINT constant; use 0xffffffffUL instead | 
| 459 | // | 
| 460 | // v0.6 - Eliminated seed overflow when uint32 is larger than 32 bits | 
| 461 | //      - Changed integer [0,n] generator to give better uniformity | 
| 462 | // | 
| 463 | // v0.7 - Fixed operator precedence ambiguity in reload() | 
| 464 | //      - Added access for real numbers in (0,1) and (0,n) | 
| 465 | // | 
| 466 | // v0.8 - Included time.h header to properly support time_t and clock_t | 
| 467 | // | 
| 468 | // v1.0 - Revised seeding to match 26 Jan 2002 update of Nishimura and Matsumoto | 
| 469 | //      - Allowed for seeding with arrays of any length | 
| 470 | //      - Added access for real numbers in [0,1) with 53-bit resolution | 
| 471 | //      - Added access for real numbers from normal (Gaussian) distributions | 
| 472 | //      - Increased overall speed by optimizing twist() | 
| 473 | //      - Doubled speed of integer [0,n] generation | 
| 474 | //      - Fixed out-of-range number generation on 64-bit machines | 
| 475 | //      - Improved portability by substituting literal constants for long enum's | 
| 476 | //      - Changed license from GNU LGPL to BSD |