| 101 |
|
double randExc( const double& n ); // real number in [0,n) |
| 102 |
|
double randDblExc(); // real number in (0,1) |
| 103 |
|
double randDblExc( const double& n ); // real number in (0,n) |
| 104 |
< |
uint32 randInt(); // integer in [0,2^32-1] |
| 104 |
> |
uint32 randInt(); // integer in [0,2^32-1] (modified for striding) |
| 105 |
> |
uint32 rawRandInt(); // original randInt |
| 106 |
|
uint32 randInt( const uint32& n ); // integer in [0,n] for n < 2^32 |
| 107 |
|
double operator()() { return rand(); } // same as rand() |
| 108 |
|
|
| 185 |
|
return mean + r * cos(phi); |
| 186 |
|
} |
| 187 |
|
|
| 188 |
< |
inline MTRand::uint32 MTRand::randInt() |
| 188 |
> |
/** |
| 189 |
> |
* This function is modified from the original to allow for random |
| 190 |
> |
* streams on parallel jobs. It now takes numbers from by striding |
| 191 |
> |
* through the random stream and picking up only one of the random |
| 192 |
> |
* numbers per nstrides_. The number it picks is the stride_'th |
| 193 |
> |
* number in the stride sequence. |
| 194 |
> |
*/ |
| 195 |
> |
inline MTRand::uint32 MTRand::randInt() { |
| 196 |
> |
|
| 197 |
> |
uint32 ranNums[nstrides_]; |
| 198 |
> |
|
| 199 |
> |
for (int i = 0; i < nstrides_; ++i) { |
| 200 |
> |
ranNums[i] = rawRandInt(); |
| 201 |
> |
} |
| 202 |
> |
|
| 203 |
> |
return ranNums[stride_]; |
| 204 |
> |
} |
| 205 |
> |
|
| 206 |
> |
/** |
| 207 |
> |
* This is the original randInt function which implements the mersenne |
| 208 |
> |
* twister. |
| 209 |
> |
*/ |
| 210 |
> |
inline MTRand::uint32 MTRand::rawRandInt() |
| 211 |
|
{ |
| 212 |
|
// Pull a 32-bit integer from the generator state |
| 213 |
|
// Every other access function simply transforms the numbers extracted here |
| 214 |
< |
|
| 215 |
< |
uint32 ranNums[nstrides_]; |
| 216 |
< |
|
| 217 |
< |
for (int i = 0; i < nstrides_; ++i) { |
| 218 |
< |
if( left == 0 ) { |
| 219 |
< |
reload(); |
| 220 |
< |
} |
| 221 |
< |
|
| 222 |
< |
--left; |
| 223 |
< |
|
| 201 |
< |
register uint32 s1; |
| 202 |
< |
s1 = *pNext++; |
| 203 |
< |
s1 ^= (s1 >> 11); |
| 204 |
< |
s1 ^= (s1 << 7) & 0x9d2c5680UL; |
| 205 |
< |
s1 ^= (s1 << 15) & 0xefc60000UL; |
| 206 |
< |
ranNums[i] = s1 ^ (s1 >> 18); |
| 207 |
< |
} |
| 208 |
< |
|
| 209 |
< |
return ranNums[stride_]; |
| 214 |
> |
|
| 215 |
> |
if( left == 0 ) reload(); |
| 216 |
> |
--left; |
| 217 |
> |
|
| 218 |
> |
register uint32 s1; |
| 219 |
> |
s1 = *pNext++; |
| 220 |
> |
s1 ^= (s1 >> 11); |
| 221 |
> |
s1 ^= (s1 << 7) & 0x9d2c5680UL; |
| 222 |
> |
s1 ^= (s1 << 15) & 0xefc60000UL; |
| 223 |
> |
return ( s1 ^ (s1 >> 18) ); |
| 224 |
|
} |
| 225 |
|
|
| 226 |
|
inline MTRand::uint32 MTRand::randInt( const uint32& n ) |
| 289 |
|
|
| 290 |
|
inline void MTRand::seed() |
| 291 |
|
{ |
| 292 |
< |
// Seed the generator with an array from /dev/urandom if available |
| 293 |
< |
// Otherwise use a hash of time() and clock() values |
| 294 |
< |
|
| 295 |
< |
// First try getting an array from /dev/urandom |
| 296 |
< |
FILE* urandom = fopen( "/dev/urandom", "rb" ); |
| 297 |
< |
if( urandom ) |
| 298 |
< |
{ |
| 299 |
< |
uint32 bigSeed[N]; |
| 300 |
< |
register uint32 *s = bigSeed; |
| 287 |
< |
register int i = N; |
| 288 |
< |
register bool success = true; |
| 289 |
< |
while( success && i-- ) |
| 290 |
< |
success = fread( s++, sizeof(uint32), 1, urandom ); |
| 291 |
< |
fclose(urandom); |
| 292 |
< |
if( success ) { seed( bigSeed, N ); return; } |
| 293 |
< |
} |
| 294 |
< |
|
| 295 |
< |
// Was not successful, so use time() and clock() instead |
| 296 |
< |
seed( hash( time(NULL), clock() ) ); |
| 292 |
> |
vector<uint32> seeds; |
| 293 |
> |
|
| 294 |
> |
seeds = generateSeeds(); |
| 295 |
> |
|
| 296 |
> |
if (seeds.size() == 1) { |
| 297 |
> |
seed( seeds[0] ); |
| 298 |
> |
} else { |
| 299 |
> |
seed( &seeds[0], seeds.size() ); |
| 300 |
> |
} |
| 301 |
|
} |
| 302 |
|
|
| 303 |
|
|
| 304 |
+ |
inline vector<uint32> MTRand::generateSeeds() { |
| 305 |
+ |
// Seed the generator with an array from /dev/urandom if available |
| 306 |
+ |
// Otherwise use a hash of time() and clock() values |
| 307 |
+ |
|
| 308 |
+ |
vector<uint32> bigSeed; |
| 309 |
+ |
|
| 310 |
+ |
// First try getting an array from /dev/urandom |
| 311 |
+ |
FILE* urandom = fopen( "/dev/urandom", "rb" ); |
| 312 |
+ |
if( urandom ) |
| 313 |
+ |
{ |
| 314 |
+ |
bigSeed.resize(N); |
| 315 |
+ |
register uint32 *s = &bigSeed[0]; |
| 316 |
+ |
register int i = N; |
| 317 |
+ |
register bool success = true; |
| 318 |
+ |
while( success && i-- ) |
| 319 |
+ |
success = fread( s++, sizeof(uint32), 1, urandom ); |
| 320 |
+ |
fclose(urandom); |
| 321 |
+ |
if( success ) { return bigSeed; } |
| 322 |
+ |
} |
| 323 |
+ |
|
| 324 |
+ |
// Was not successful, so use time() and clock() instead |
| 325 |
+ |
|
| 326 |
+ |
bigSeed.push_back(hash( time(NULL), clock())); |
| 327 |
+ |
return bigSeed; |
| 328 |
+ |
} |
| 329 |
+ |
|
| 330 |
+ |
|
| 331 |
|
inline void MTRand::initialize( const uint32 seed ) |
| 332 |
|
{ |
| 333 |
|
// Initialize generator state with seed |
