#ifndef MTWIST_H
#define MTWIST_H
/*
* $Id: mtwist.h,v 1.24 2012-12-31 22:22:03-08 geoff Exp $
*
* Header file for C/C++ use of the Mersenne-Twist pseudo-RNG. See
* http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/emt.html for full
* information.
*
* Author of this header file: Geoff Kuenning, March 18, 2001.
*
* IMPORTANT NOTE: this implementation assumes a modern compiler. In
* particular, it assumes that the "inline" keyword is available, and
* that the "stdint.h" header file is present.
*
* The variables above are defined in an inverted sense because I
* expect that most modern compilers will support these features. By
* inverting the sense, this common case will require no special
* compiler flags.
*
* IMPORTANT NOTE: this software requires access to a 32-bit type.
* The Mersenne Twist algorithms are not guaranteed to produce correct
* results with a 64-bit type.
*
* The executable part of this software is based on LGPL-ed code by
* Takuji Nishimura. The header file is therefore also distributed
* under the LGPL:
*
* This library is free software; you can redistribute it and/or
* modify it under the terms of the GNU Library General Public License
* as published by the Free Software Foundation; either version 2 of
* the License, or (at your option) any later version.
*
* This library is distributed in the hope that it will be useful, but
* WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
* Library General Public License for more details. You should have
* received a copy of the GNU Library General Public License along
* with this library; if not, write to the Free Foundation, Inc., 59
* Temple Place, Suite 330, Boston, MA 02111-1307 USA
*
* $Log: mtwist.h,v $
* Revision 1.24 2012-12-31 22:22:03-08 geoff
* Fix the out-of-bounds bug in mt_llrand and mt_ldrand that were
* overlooked because I assumed they used mts_*.
*
* Revision 1.23 2013-01-01 01:18:52-08 geoff
* Fix a lot of compiler warnings.
*
* Revision 1.22 2012-12-30 16:24:49-08 geoff
* Declare the new versions of the /dev/random and urandom seeding
* functions, which now return the seed chosen.
*
* Revision 1.21 2012-09-23 23:15:40-07 geoff
* Fix an array index violation found by valgrind and reported by David
* Chapman; under some circumstances statevec[-1] could be accessed and
* used in random-number generation. The bug only affects the *_llrand
* and *_ldrand functions.
*
* Revision 1.20 2010-12-10 03:28:18-08 geoff
* Add support for GENERATE_CODE_IN_HEADER. Fix the URL for the original
* Web page.
*
* Revision 1.19 2010-06-24 20:53:59+12 geoff
* Switch to using types from stdint.h. Get rid of all compilation
* options.
*
* Revision 1.18 2010-06-24 00:29:38-07 geoff
* Do a better job of auto-determining MT_MACHINE_BITS.
*
* Revision 1.17 2007-10-26 00:21:06-07 geoff
* Introduce, document, and use the new mt_u32bit_t type so that the code
* will compile and run on 64-bit platforms (although it does not
* currently use the 64-bit Mersenne Twist algorithm).
*
* Revision 1.16 2005/11/11 08:21:39 geoff
* If possible, try to infer MT_MACHINE_BITS from limits.h.
*
* Revision 1.15 2003/09/11 23:56:20 geoff
* Allow stdio references in C++ files; it turns out that ANSI has
* blessed it. Declare the various functions as external even if they're
* inlined or being compiled directly (in mtwist.c). Get rid of a #ifdef
* that can't ever be true.
*
* Revision 1.14 2003/09/11 05:50:53 geoff
* Don't allow stdio references from C++, since they're not guaranteed to
* work on all compilers. Disable inlining using the MT_INLINE keyword
* rather than #defining inline, since doing the latter can affect other
* files and functions than our own.
*
* Revision 1.13 2003/07/01 23:29:29 geoff
* Refer to streams from the standard library using the correct namespace.
*
* Revision 1.12 2002/10/30 07:39:54 geoff
* Declare the new seeding functions.
*
* Revision 1.11 2001/06/19 00:41:16 geoff
* For consistency with other C++ types, don't put out a newline after
* the saved data.
*
* Revision 1.10 2001/06/18 10:09:24 geoff
* Fix some places where I forgot to set one of the result values. Make
* the C++ state vector protected so the random-distributions package can
* pass it to the C functions.
*
* Revision 1.9 2001/06/18 05:40:12 geoff
* Prefix the compile options with MT_.
*
* Revision 1.8 2001/06/14 10:26:59 geoff
* Invert the sense of the #define flags so that the default is the
* normal case (if gcc is normal!). Also default MT_MACHINE_BITS to 32.
*
* Revision 1.7 2001/06/14 10:10:38 geoff
* Move the critical-path PRNG code into the header file so that it can
* be inlined. Add saving/loading of state. Add functions to seed based
* on /dev/random or the time. Add the function-call operator in the C++
* code.
*
* Revision 1.6 2001/06/11 10:00:04 geoff
* Add declarations of the refresh and /dev/random seeding functions.
* Change getstate to return a complete state pointer, since knowing the
* position in the state vector is critical to restoring the state.
*
* Revision 1.5 2001/04/23 08:36:03 geoff
* Remember to zero the state pointer when constructing, since otherwise
* proper initialization won't happen.
*
* Revision 1.4 2001/04/14 01:33:32 geoff
* Clarify the license
*
* Revision 1.3 2001/04/14 01:04:54 geoff
* Add a C++ class, mt_prng, that makes usage more convenient for C++
* programmers.
*
* Revision 1.2 2001/04/09 08:45:00 geoff
* Fix the name in the #ifndef wrapper, and clean up some outdated comments.
*
* Revision 1.1 2001/04/07 09:43:41 geoff
* Initial revision
*
*/
#include <stdio.h>
#ifdef __cplusplus
#include <iostream>
#endif /* __cplusplus */
#define __STDC_LIMIT_MACROS
#include <stdint.h>
/*
* The following value is a fundamental parameter of the algorithm.
* It was found experimentally using methods described in Matsumoto
* and Nishimura's paper. It is exceedingly magic; don't change it.
*/
#define MT_STATE_SIZE 624 /* Size of the MT state vector */
/*
* Internal state for an MT RNG. The user can keep multiple mt_state
* structures around as a way of generating multiple streams of random
* numbers.
*
* In Matsumoto and Nishimura's original paper, the state vector was
* processed in a forward direction. I have reversed the state vector
* in this implementation. The reason for the reversal is that it
* allows the critical path to use a test against zero instead of a
* test against 624 to detect the need to refresh the state. on most
* machines, testing against zero is slightly faster. It also means
* that a state that has been set to all zeros will be correctly
* detected as needing initialization; this means that setting a state
* vector to zero (either with memset or by statically allocating it)
* will cause the RNG to operate properly.
*/
typedef struct
{
uint32_t statevec[MT_STATE_SIZE];
/* Vector holding current state */
int stateptr; /* Next state entry to be used */
int initialized; /* NZ if state was initialized */
}
mt_state;
#ifdef __cplusplus
extern "C"
{
#endif
/*
* Functions for manipulating any generator (given a state pointer).
*/
extern void mts_mark_initialized(mt_state* state);
/* Mark a PRNG state as initialized */
extern void mts_seed32(mt_state* state, uint32_t seed);
/* Set random seed for any generator */
extern void mts_seed32new(mt_state* state, uint32_t seed);
/* Set random seed for any generator */
extern void mts_seedfull(mt_state* state,
uint32_t seeds[MT_STATE_SIZE]);
/* Set complicated seed for any gen. */
extern uint32_t mts_seed(mt_state* state);
/* Choose seed from random input. */
/* ..Prefers /dev/urandom; uses time */
/* ..if /dev/urandom unavailable. */
/* ..Only gives 32 bits of entropy. */
/* ..Returns seed usable with seed32 */
extern uint32_t mts_goodseed(mt_state* state);
/* Choose seed from more random */
/* ..input than mts_seed. Prefers */
/* ../dev/random; uses time if that */
/* ..is unavailable. Only gives 32 */
/* ..bits of entropy. */
/* ..Returns seed usable with seed32 */
extern void mts_bestseed(mt_state* state);
/* Choose seed from extremely random */
/* ..input (can be *very* slow). */
/* ..Prefers /dev/random and reads */
/* ..the entire state from there. */
/* ..If /dev/random is unavailable, */
/* ..falls back to mt_goodseed(). */
/* ..Not usually worth the cost. */
extern void mts_refresh(mt_state* state);
/* Generate 624 more random values */
extern int mts_savestate(FILE* statefile, mt_state* state);
/* Save state to a file (ASCII). */
/* ..Returns NZ if succeeded. */
extern int mts_loadstate(FILE* statefile, mt_state* state);
/* Load state from a file (ASCII). */
/* ..Returns NZ if succeeded. */
/*
* Functions for manipulating the default generator.
*/
extern void mt_seed32(uint32_t seed);
/* Set random seed for default gen. */
extern void mt_seed32new(uint32_t seed);
/* Set random seed for default gen. */
extern void mt_seedfull(uint32_t seeds[MT_STATE_SIZE]);
/* Set complicated seed for default */
extern uint32_t mt_seed(void); /* Choose seed from random input. */
/* ..Prefers /dev/urandom; uses time */
/* ..if /dev/urandom unavailable. */
/* ..Only gives 32 bits of entropy. */
extern uint32_t mt_goodseed(void);
/* Choose seed from more random */
/* ..input than mts_seed. Prefers */
/* ../dev/random; uses time if that */
/* ..is unavailable. Only gives 32 */
/* ..bits of entropy. */
extern void mt_bestseed(void);
/* Choose seed from extremely random */
/* ..input (can be *very* slow). */
/* ..Prefers /dev/random and reads */
/* ..the entire state from there. */
/* ..If /dev/random is unavailable, */
/* ..falls back to mt_goodseed(). */
/* ..Not usually worth the cost. */
extern mt_state* mt_getstate(void);
/* Get current state of default */
/* ..generator */
extern int mt_savestate(FILE* statefile);
/* Save state to a file (ASCII) */
/* ..Returns NZ if succeeded. */
extern int mt_loadstate(FILE* statefile);
/* Load state from a file (ASCII) */
/* ..Returns NZ if succeeded. */
#ifdef __cplusplus
}
#endif
/*
* Functions for generating random numbers. The actual code of the
* functions is given in this file so that it can be declared inline.
* For compilers that don't have the inline feature, mtwist.c will
* incorporate this file with some clever #defining so that the code
* actually gets compiled. In that case, however, "extern"
* definitions will be needed here, so we give them.
*/
#ifdef __cplusplus
#endif /* __cplusplus */
extern uint32_t mts_lrand(mt_state* state);
/* Generate 32-bit value, any gen. */
#ifdef UINT64_MAX
extern uint64_t mts_llrand(mt_state* state);
/* Generate 64-bit value, any gen. */
#endif /* UINT64_MAX */
extern double mts_drand(mt_state* state);
/* Generate floating value, any gen. */
/* Fast, with only 32-bit precision */
extern double mts_ldrand(mt_state* state);
/* Generate floating value, any gen. */
/* Slower, with 64-bit precision */
extern uint32_t mt_lrand(void); /* Generate 32-bit random value */
#ifdef UINT64_MAX
extern uint64_t mt_llrand(void);
/* Generate 64-bit random value */
#endif /* UINT64_MAX */
extern double mt_drand(void);
/* Generate floating value */
/* Fast, with only 32-bit precision */
extern double mt_ldrand(void);
/* Generate floating value */
/* Slower, with 64-bit precision */
/*
* Tempering parameters. These are perhaps the most magic of all the magic
* values in the algorithm. The values are again experimentally determined.
* The values generated by the recurrence relation (constants above) are not
* equidistributed in 623-space. For some reason, the tempering process
* produces that effect. Don't ask me why. Read the paper if you can
* understand the math. Or just trust these magic numbers.
*/
#define MT_TEMPERING_MASK_B 0x9d2c5680
#define MT_TEMPERING_MASK_C 0xefc60000
#define MT_TEMPERING_SHIFT_U(y) \
(y >> 11)
#define MT_TEMPERING_SHIFT_S(y) \
(y << 7)
#define MT_TEMPERING_SHIFT_T(y) \
(y << 15)
#define MT_TEMPERING_SHIFT_L(y) \
(y >> 18)
/*
* Macros to do the tempering. MT_PRE_TEMPER does all but the last step;
* it's useful for situations where the final step can be incorporated
* into a return statement. MT_FINAL_TEMPER does that final step (not as
* an assignment). MT_TEMPER does the entire process. Note that
* MT_PRE_TEMPER and MT_TEMPER both modify their arguments.
*/
#define MT_PRE_TEMPER(value) \
do \
{ \
value ^= MT_TEMPERING_SHIFT_U(value); \
value ^= MT_TEMPERING_SHIFT_S(value) & MT_TEMPERING_MASK_B; \
value ^= MT_TEMPERING_SHIFT_T(value) & MT_TEMPERING_MASK_C; \
} \
while (0)
#define MT_FINAL_TEMPER(value) \
((value) ^ MT_TEMPERING_SHIFT_L(value))
#define MT_TEMPER(value) \
do \
{ \
value ^= MT_TEMPERING_SHIFT_U(value); \
value ^= MT_TEMPERING_SHIFT_S(value) & MT_TEMPERING_MASK_B; \
value ^= MT_TEMPERING_SHIFT_T(value) & MT_TEMPERING_MASK_C; \
value ^= MT_TEMPERING_SHIFT_L(value); \
} \
while (0)
/*
* The Mersenne Twist PRNG makes it default state available as an
* external variable. This feature is undocumented, but is useful to
* use because it allows us to avoid implementing every randistr function
* twice. (In fact, the feature was added to enable randistrs.c to be
* written. It would be better to write in C++, where I could control
* the access to the state.)
*/
extern mt_state mt_default_state;
/* State of the default generator */
extern double mt_32_to_double;
/* Multiplier to convert long to dbl */
extern double mt_64_to_double;
/* Mult'r to cvt long long to dbl */
/*
* In gcc, inline functions must be declared extern or they'll produce
* assembly code (and thus linking errors). We have to work around
* that difficulty with the MT_EXTERN define.
*/
#ifndef MT_EXTERN
#ifdef __cplusplus
#define MT_EXTERN /* C++ doesn't need static */
#else /* __cplusplus */
#define MT_EXTERN extern /* C (at least gcc) needs extern */
#endif /* __cplusplus */
#endif /* MT_EXTERN */
/*
* Make it possible for mtwist.c to disable the inline keyword. We
* use our own keyword so that we don't interfere with inlining in
* C/C++ header files, above.
*/
#ifndef MT_INLINE
#define MT_INLINE inline /* Compiler has inlining */
#endif /* MT_INLINE */
/*
* Try to guess whether the compiler is one (like gcc) that requires
* inline code to be available in the header file, or a smarter one
* that gets inlines directly from object files. But if we've been
* given the information, trust it.
*/
#ifndef MT_GENERATE_CODE_IN_HEADER
#ifdef __GNUC__
#define MT_GENERATE_CODE_IN_HEADER 1
#endif /* __GNUC__ */
#if defined(__INTEL_COMPILER) || defined(_MSC_VER)
#define MT_GENERATE_CODE_IN_HEADER 0
#endif /* __INTEL_COMPILER || _MSC_VER */
#endif /* MT_GENERATE_CODE_IN_HEADER */
#if MT_GENERATE_CODE_IN_HEADER
/*
* Generate a random number in the range 0 to 2^32-1, inclusive, working
* from a given state vector.
*
* The generator is optimized for speed. The primary optimization is that
* the pseudorandom numbers are generated in batches of MT_STATE_SIZE. This
* saves the cost of a modulus operation in the critical path.
*/
MT_EXTERN MT_INLINE uint32_t mts_lrand(
register mt_state* state) /* State for the PRNG */
{
register uint32_t random_value; /* Pseudorandom value generated */
if (state->stateptr <= 0)
mts_refresh(state);
random_value = state->statevec[--state->stateptr];
MT_PRE_TEMPER(random_value);
return MT_FINAL_TEMPER(random_value);
}
#ifdef UINT64_MAX
/*
* Generate a random number in the range 0 to 2^64-1, inclusive, working
* from a given state vector.
*
* According to Matsumoto and Nishimura, such a number can be generated by
* simply concatenating two 32-bit pseudorandom numbers. Who am I to argue?
*
* Note that there is a slight inefficiency here: if the 624-entry state is
* recycled on the second call to mts_lrand, there will be an unnecessary
* check to see if the state has been initialized. The cost of that check
* seems small (since it happens only once every 624 random numbers, and
* never if only 64-bit numbers are being generated), so I didn't bother to
* optimize it out. Doing so would be messy, since it would require two
* nearly-identical internal implementations of mts_lrand.
*/
MT_EXTERN MT_INLINE uint64_t mts_llrand(
register mt_state* state) /* State for the PRNG */
{
register uint32_t random_value_1; /* 1st pseudorandom value generated */
register uint32_t random_value_2; /* 2nd pseudorandom value generated */
/*
* For maximum speed, we'll handle the two overflow cases
* together. That will save us one test in the common case, at
* the expense of an extra one in the overflow case.
*/
if (--state->stateptr <= 0)
{
if (state->stateptr < 0)
{
mts_refresh(state);
random_value_1 = state->statevec[--state->stateptr];
}
else
{
random_value_1 = state->statevec[state->stateptr];
mts_refresh(state);
}
}
else
random_value_1 = state->statevec[state->stateptr];
MT_TEMPER(random_value_1);
random_value_2 = state->statevec[--state->stateptr];
MT_PRE_TEMPER(random_value_2);
return ((uint64_t) random_value_1 << 32)
| (uint64_t) MT_FINAL_TEMPER(random_value_2);
}
#endif /* UINT64_MAX */
/*
* Generate a double-precision random number between 0 (inclusive) and 1.0
* (exclusive). This function is optimized for speed, but it only generates
* 32 bits of precision. Use mts_ldrand to get 64 bits of precision.
*/
MT_EXTERN MT_INLINE double mts_drand(
register mt_state* state) /* State for the PRNG */
{
register uint32_t random_value; /* Pseudorandom value generated */
if (state->stateptr <= 0)
mts_refresh(state);
random_value = state->statevec[--state->stateptr];
MT_TEMPER(random_value);
return random_value * mt_32_to_double;
}
/*
* Generate a double-precision random number between 0 (inclusive) and 1.0
* (exclusive). This function generates 64 bits of precision. Use
* mts_drand for more speed but less precision.
*/
MT_EXTERN MT_INLINE double mts_ldrand(
register mt_state* state) /* State for the PRNG */
{
#ifdef UINT64_MAX
uint64_t final_value; /* Final (integer) value */
#endif /* UINT64_MAX */
register uint32_t random_value_1; /* 1st pseudorandom value generated */
register uint32_t random_value_2; /* 2nd pseudorandom value generated */
/*
* For maximum speed, we'll handle the two overflow cases
* together. That will save us one test in the common case, at
* the expense of an extra one in the overflow case.
*/
if (--state->stateptr <= 0)
{
if (state->stateptr < 0)
{
mts_refresh(state);
random_value_1 = state->statevec[--state->stateptr];
}
else
{
random_value_1 = state->statevec[state->stateptr];
mts_refresh(state);
}
}
else
random_value_1 = state->statevec[state->stateptr];
MT_TEMPER(random_value_1);
random_value_2 = state->statevec[--state->stateptr];
MT_TEMPER(random_value_2);
#ifdef UINT64_MAX
final_value = ((uint64_t) random_value_1 << 32) | (uint64_t) random_value_2;
return final_value * mt_64_to_double;
#else /* UINT64_MAX */
return random_value_1 * mt_32_to_double + random_value_2 * mt_64_to_double;
#endif /* UINT64_MAX */
}
/*
* Generate a random number in the range 0 to 2^32-1, inclusive, working
* from the default state vector.
*
* See mts_lrand for full commentary.
*/
MT_EXTERN MT_INLINE uint32_t mt_lrand(void)
{
register uint32_t random_value; /* Pseudorandom value generated */
if (mt_default_state.stateptr <= 0)
mts_refresh(&mt_default_state);
random_value = mt_default_state.statevec[--mt_default_state.stateptr];
MT_PRE_TEMPER(random_value);
return MT_FINAL_TEMPER(random_value);
}
#ifdef UINT64_MAX
/*
* Generate a random number in the range 0 to 2^64-1, inclusive, working
* from the default state vector.
*
* See mts_llrand for full commentary.
*/
MT_EXTERN MT_INLINE uint64_t mt_llrand(void)
{
register uint32_t random_value_1; /* 1st pseudorandom value generated */
register uint32_t random_value_2; /* 2nd pseudorandom value generated */
/*
* For maximum speed, we'll handle the two overflow cases
* together. That will save us one test in the common case, at
* the expense of an extra one in the overflow case.
*/
if (--mt_default_state.stateptr <= 0)
{
if (mt_default_state.stateptr < 0)
{
mts_refresh(&mt_default_state);
random_value_1 =
mt_default_state.statevec[--mt_default_state.stateptr];
}
else
{
random_value_1 =
mt_default_state.statevec[mt_default_state.stateptr];
mts_refresh(&mt_default_state);
}
}
else
random_value_1 = mt_default_state.statevec[mt_default_state.stateptr];
MT_TEMPER(random_value_1);
random_value_2 = mt_default_state.statevec[--mt_default_state.stateptr];
MT_PRE_TEMPER(random_value_2);
return ((uint64_t) random_value_1 << 32)
| (uint64_t) MT_FINAL_TEMPER(random_value_2);
}
#endif /* UINT64_MAX */
/*
* Generate a double-precision random number between 0 (inclusive) and 1.0
* (exclusive). This function is optimized for speed, but it only generates
* 32 bits of precision. Use mt_ldrand to get 64 bits of precision.
*/
MT_EXTERN MT_INLINE double mt_drand(void)
{
register uint32_t random_value; /* Pseudorandom value generated */
if (mt_default_state.stateptr <= 0)
mts_refresh(&mt_default_state);
random_value = mt_default_state.statevec[--mt_default_state.stateptr];
MT_TEMPER(random_value);
return random_value * mt_32_to_double;
}
/*
* Generate a double-precision random number between 0 (inclusive) and 1.0
* (exclusive). This function generates 64 bits of precision. Use
* mts_drand for more speed but less precision.
*/
MT_EXTERN MT_INLINE double mt_ldrand(void)
{
#ifdef UINT64_MAX
uint64_t final_value; /* Final (integer) value */
#endif /* UINT64_MAX */
register uint32_t random_value_1; /* 1st pseudorandom value generated */
register uint32_t random_value_2; /* 2nd pseudorandom value generated */
/*
* For maximum speed, we'll handle the two overflow cases
* together. That will save us one test in the common case, at
* the expense of an extra one in the overflow case.
*/
if (--mt_default_state.stateptr <= 0)
{
if (mt_default_state.stateptr < 0)
{
mts_refresh(&mt_default_state);
random_value_1 =
mt_default_state.statevec[--mt_default_state.stateptr];
}
else
{
random_value_1 =
mt_default_state.statevec[mt_default_state.stateptr];
mts_refresh(&mt_default_state);
}
}
else
random_value_1 = mt_default_state.statevec[mt_default_state.stateptr];
MT_TEMPER(random_value_1);
random_value_2 = mt_default_state.statevec[--mt_default_state.stateptr];
MT_TEMPER(random_value_2);
#ifdef UINT64_MAX
final_value = ((uint64_t) random_value_1 << 32) | (uint64_t) random_value_2;
return final_value * mt_64_to_double;
#else /* UINT64_MAX */
return random_value_1 * mt_32_to_double + random_value_2 * mt_64_to_double;
#endif /* UINT64_MAX */
}
#endif /* MT_GENERATE_CODE_IN_HEADER */
#ifdef __cplusplus
/*
* C++ interface to the Mersenne Twist PRNG. This class simply
* provides a more C++-ish way to access the PRNG. Only state-based
* functions are provided. All functions are inlined, both for speed
* and so that the same implementation code can be used in C and C++.
*/
class mt_prng
{
friend class mt_empirical_distribution;
public:
/*
* Constructors and destructors. The default constructor
* leaves initialization (seeding) for later unless pickSeed
* is true, in which case the seed is chosen based on either
* /dev/urandom (if available) or the system time. The other
* constructors accept either a 32-bit seed, or a full
* 624-integer seed.
*/
mt_prng( // Default constructor
bool pickSeed = false)
// True to get seed from /dev/urandom
// ..or time
{
state.stateptr = 0;
state.initialized = 0;
if (pickSeed)
(void)mts_seed(&state);
}
mt_prng(uint32_t newseed)
// Construct with 32-bit seeding
{
state.stateptr = 0;
state.initialized = 0;
mts_seed32(&state, newseed);
}
mt_prng(uint32_t seeds[MT_STATE_SIZE])
// Construct with full seeding
{
state.stateptr = 0;
state.initialized = 0;
mts_seedfull(&state, seeds);
}
~mt_prng() { }
/*
* Copy and assignment are best left defaulted.
*/
/*
* PRNG seeding functions.
*/
void seed32(uint32_t newseed)
// Set 32-bit random seed
{
mts_seed32(&state, newseed);
}
void seed32new(uint32_t newseed)
// Set 32-bit random seed
{
mts_seed32new(&state, newseed);
}
void seedfull(uint32_t seeds[MT_STATE_SIZE])
// Set complicated random seed
{
mts_seedfull(&state, seeds);
}
uint32_t seed() // Choose seed from random input
{
return mts_seed(&state);
}
uint32_t goodseed() // Choose better seed from random input
{
return mts_goodseed(&state);
}
void bestseed() // Choose best seed from random input
{
mts_bestseed(&state);
}
friend std::ostream&
operator<<(std::ostream& stream, const mt_prng& rng);
friend std::istream&
operator>>(std::istream& stream, mt_prng& rng);
/*
* PRNG generation functions
*/
uint32_t lrand() // Generate 32-bit pseudo-random value
{
return mts_lrand(&state);
}
#ifdef UINT64_MAX
uint64_t llrand() // Generate 64-bit pseudo-random value
{
return mts_llrand(&state);
}
#endif /* UINT64_MAX */
double drand() // Generate fast 32-bit floating value
{
return mts_drand(&state);
}
double ldrand() // Generate slow 64-bit floating value
{
return mts_ldrand(&state);
}
/*
* Following Richard J. Wagner's example, we overload the
* function-call operator to return a 64-bit floating value.
* That allows the common use of the PRNG to be simplified as
* in the following example:
*
* mt_prng ranno(true);
* // ...
* coinFlip = ranno() >= 0.5 ? heads : tails;
*/
double operator()()
{
return mts_drand(&state);
}
protected:
/*
* Protected data
*/
mt_state state; // Current state of the PRNG
};
#if MT_GENERATE_CODE_IN_HEADER
/*
* Save state to a stream. See mts_savestate.
*/
MT_INLINE std::ostream& operator<<(
std::ostream& stream, // Stream to save to
const mt_prng& rng) // PRNG to save
{
for (int i = MT_STATE_SIZE; --i >= 0; )
{
if (!(stream << rng.state.statevec[i] << ' '))
return stream;
}
return stream << rng.state.stateptr;
}
/*
* Restore state from a stream. See mts_loadstate.
*/
MT_INLINE std::istream& operator>>(
std::istream& stream, // Stream to laod from
mt_prng& rng) // PRNG to load
{
rng.state.initialized = rng.state.stateptr = 0;
for (int i = MT_STATE_SIZE; --i >= 0; )
{
if (!(stream >> rng.state.statevec[i]))
return stream;
}
if (!(stream >> rng.state.stateptr))
{
rng.state.stateptr = 0;
return stream;
}
/*
* If the state is invalid, all we can do is to make it uninitialized.
*/
if (rng.state.stateptr < 0 || rng.state.stateptr > MT_STATE_SIZE)
{
rng.state.stateptr = 0;
return stream;
}
mts_mark_initialized(&rng.state);
return stream;
}
#endif /* MT_GENERATE_CODE_IN_HEADER */
#endif /* __cplusplus */
#endif /* MTWIST_H */