57 #ifndef MERSENNETWISTER_H
58 #define MERSENNETWISTER_H
63 #ifdef _MSC_VER // !!! mb pragmas added 19.02.2007
64 #pragma warning(disable: 4146)
65 #pragma warning(disable: 4996)
102 double rand(
const double& n );
104 double randExc(
const double& n );
115 double randNorm(
const double& mean = 0.0,
const double& variance = 0.0 );
137 {
return m ^ (
mixBits(s0,s1)>>1) ^ (-
loBit(s1) & 0x9908b0dfUL); }
146 {
seed(bigSeed,seedLength); }
152 {
return double(
randInt()) * (1.0/4294967295.0); }
155 {
return rand() * n; }
158 {
return double(
randInt()) * (1.0/4294967296.0); }
164 {
return (
double(
randInt()) + 0.5 ) * (1.0/4294967296.0); }
172 return ( a * 67108864.0 + b ) * (1.0/9007199254740992.0);
179 double r = sqrt( -2.0 * log( 1.0-
randDblExc()) ) * variance;
180 double phi = 2.0 * 3.14159265358979323846264338328 *
randExc();
181 return mean + r * cos(phi);
195 s1 ^= (s1 << 7) & 0x9d2c5680UL;
196 s1 ^= (s1 << 15) & 0xefc60000UL;
197 return ( s1 ^ (s1 >> 18) );
240 if( seedLength > k ) k = seedLength;
245 state[i] += ( bigSeed[j] & 0xffffffffUL ) + j;
246 state[i] &= 0xffffffffUL;
248 if( i >=
N ) { state[0] = state[
N-1]; i = 1; }
249 if( j >= seedLength ) j = 0;
251 for( k =
N - 1; k; --k )
256 state[i] &= 0xffffffffUL;
258 if( i >=
N ) { state[0] = state[
N-1]; i = 1; }
260 state[0] = 0x80000000UL;
271 FILE* urandom = fopen(
"/dev/urandom",
"rb" );
275 register uint32 *s = bigSeed;
277 register bool success =
true;
278 while( success && i-- )
279 success = fread( s++,
sizeof(
uint32), 1, urandom )!=0;
281 if( success ) {
seed( bigSeed,
N );
return; }
285 seed(
hash( time(NULL), clock() ) );
298 *s++ = seed & 0xffffffffUL;
301 *s++ = ( 1812433253UL * ( *r ^ (*r >> 30) ) + i ) & 0xffffffffUL;
313 for( i =
N -
M; i--; ++p )
314 *p =
twist( p[
M], p[0], p[1] );
315 for( i = M; --i; ++p )
316 *p =
twist( p[M-
N], p[0], p[1] );
332 unsigned char *p = (
unsigned char *) &t;
333 for(
size_t i = 0; i <
sizeof(t); ++i )
335 h1 *= UCHAR_MAX + 2U;
339 p = (
unsigned char *) &c;
340 for(
size_t j = 0; j <
sizeof(c); ++j )
342 h2 *= UCHAR_MAX + 2U;
345 return ( h1 + differ++ ) ^ h2;
351 register uint32 *sa = saveArray;
354 for( ; i--; *sa++ = *s++ ) {}
362 register uint32 *la = loadArray;
364 for( ; i--; *s++ = *la++ ) {}
373 register int i = mtrand.
N;
374 for( ; i--; os << *s++ <<
"\t" ) {}
375 return os << mtrand.
left;
382 register int i = mtrand.
N;
383 for( ; i--; is >> *s++ ) {}
389 #endif // MERSENNETWISTER_H