Rev 10 | Go to most recent revision | Show entire file | Ignore whitespace | Details | Blame | Last modification | View Log | RSS feed
| Rev 10 | Rev 8160 | ||
|---|---|---|---|
| Line 1... | Line 1... | ||
| 1 | /* |
1 | /* |
| 2 | A C-program for MT19937, with initialization improved 2002/1/26. |
2 | A C-program for MT19937, with initialization improved 2002/1/26. |
| 3 | Coded by Takuji Nishimura and Makoto Matsumoto. |
3 | Coded by Takuji Nishimura and Makoto Matsumoto. |
| 4 | 4 | ||
| 5 | Before using, initialize the state by using init_genrand(seed) |
5 | Before using, initialize the state by using init_genrand(seed) |
| 6 | or init_by_array(init_key, key_length). |
6 | or init_by_array(init_key, key_length). |
| 7 | 7 | ||
| 8 | Copyright (C) 1997 - 2002, Makoto Matsumoto and Takuji Nishimura, |
8 | Copyright (C) 1997 - 2002, Makoto Matsumoto and Takuji Nishimura, |
| 9 | All rights reserved. |
9 | All rights reserved. |
| 10 | 10 | ||
| 11 | Redistribution and use in source and binary forms, with or without |
11 | Redistribution and use in source and binary forms, with or without |
| 12 | modification, are permitted provided that the following conditions |
12 | modification, are permitted provided that the following conditions |
| 13 | are met: |
13 | are met: |
| 14 | 14 | ||
| Line 17... | Line 17... | ||
| 17 | 17 | ||
| 18 | 2. Redistributions in binary form must reproduce the above copyright |
18 | 2. Redistributions in binary form must reproduce the above copyright |
| 19 | notice, this list of conditions and the following disclaimer in the |
19 | notice, this list of conditions and the following disclaimer in the |
| 20 | documentation and/or other materials provided with the distribution. |
20 | documentation and/or other materials provided with the distribution. |
| 21 | 21 | ||
| 22 | 3. The names of its contributors may not be used to endorse or promote |
22 | 3. The names of its contributors may not be used to endorse or promote |
| 23 | products derived from this software without specific prior written |
23 | products derived from this software without specific prior written |
| 24 | permission. |
24 | permission. |
| 25 | 25 | ||
| 26 | THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS |
26 | THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS |
| 27 | "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT |
27 | "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT |
| 28 | LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR |
28 | LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR |
| Line 39... | Line 39... | ||
| 39 | Any feedback is very welcome. |
39 | Any feedback is very welcome. |
| 40 | http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/emt.html |
40 | http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/emt.html |
| 41 | email: m-mat @ math.sci.hiroshima-u.ac.jp (remove space) |
41 | email: m-mat @ math.sci.hiroshima-u.ac.jp (remove space) |
| 42 | */ |
42 | */ |
| 43 | 43 | ||
| 44 | /* Period parameters */ |
44 | /* Period parameters */ |
| 45 | #define N 624 |
45 | #define N 624 |
| 46 | #define M 397 |
46 | #define M 397 |
| 47 | #define MATRIX_A 0x9908b0dfUL /* constant vector a */ |
47 | #define MATRIX_A 0x9908b0dfUL /* constant vector a */ |
| 48 | #define UPPER_MASK 0x80000000UL /* most significant w-r bits */ |
48 | #define UPPER_MASK 0x80000000UL /* most significant w-r bits */ |
| 49 | #define LOWER_MASK 0x7fffffffUL /* least significant r bits */ |
49 | #define LOWER_MASK 0x7fffffffUL /* least significant r bits */ |
| Line 54... | Line 54... | ||
| 54 | /* initializes mt[N] with a seed */ |
54 | /* initializes mt[N] with a seed */ |
| 55 | void init_genrand(unsigned long s) |
55 | void init_genrand(unsigned long s) |
| 56 | { |
56 | { |
| 57 | mt[0]= s & 0xffffffffUL; |
57 | mt[0]= s & 0xffffffffUL; |
| 58 | for (mti=1; mti<N; mti++) { |
58 | for (mti=1; mti<N; mti++) { |
| 59 | mt[mti] = |
59 | mt[mti] = |
| 60 | (1812433253UL * (mt[mti-1] ^ (mt[mti-1] >> 30)) + mti); |
60 | (1812433253UL * (mt[mti-1] ^ (mt[mti-1] >> 30)) + mti); |
| 61 | /* See Knuth TAOCP Vol2. 3rd Ed. P.106 for multiplier. */ |
61 | /* See Knuth TAOCP Vol2. 3rd Ed. P.106 for multiplier. */ |
| 62 | /* In the previous versions, MSBs of the seed affect */ |
62 | /* In the previous versions, MSBs of the seed affect */ |
| 63 | /* only MSBs of the array mt[]. */ |
63 | /* only MSBs of the array mt[]. */ |
| 64 | /* 2002/01/09 modified by Makoto Matsumoto */ |
64 | /* 2002/01/09 modified by Makoto Matsumoto */ |
| 65 | mt[mti] &= 0xffffffffUL; |
65 | mt[mti] &= 0xffffffffUL; |
| Line 91... | Line 91... | ||
| 91 | mt[i] &= 0xffffffffUL; /* for WORDSIZE > 32 machines */ |
91 | mt[i] &= 0xffffffffUL; /* for WORDSIZE > 32 machines */ |
| 92 | i++; |
92 | i++; |
| 93 | if (i>=N) { mt[0] = mt[N-1]; i=1; } |
93 | if (i>=N) { mt[0] = mt[N-1]; i=1; } |
| 94 | } |
94 | } |
| 95 | 95 | ||
| 96 | mt[0] = 0x80000000UL; /* MSB is 1; assuring non-zero initial array */ |
96 | mt[0] = 0x80000000UL; /* MSB is 1; assuring non-zero initial array */ |
| 97 | } |
97 | } |
| 98 | 98 | ||
| 99 | /* generates a random number on [0,0xffffffff]-interval */ |
99 | /* generates a random number on [0,0xffffffff]-interval */ |
| 100 | unsigned long genrand_int32(void) |
100 | unsigned long genrand_int32(void) |
| 101 | { |
101 | { |
| Line 120... | Line 120... | ||
| 120 | y = (mt[N-1]&UPPER_MASK)|(mt[0]&LOWER_MASK); |
120 | y = (mt[N-1]&UPPER_MASK)|(mt[0]&LOWER_MASK); |
| 121 | mt[N-1] = mt[M-1] ^ (y >> 1) ^ mag01[y & 0x1UL]; |
121 | mt[N-1] = mt[M-1] ^ (y >> 1) ^ mag01[y & 0x1UL]; |
| 122 | 122 | ||
| 123 | mti = 0; |
123 | mti = 0; |
| 124 | } |
124 | } |
| 125 | 125 | ||
| 126 | y = mt[mti++]; |
126 | y = mt[mti++]; |
| 127 | 127 | ||
| 128 | /* Tempering */ |
128 | /* Tempering */ |
| 129 | y ^= (y >> 11); |
129 | y ^= (y >> 11); |
| 130 | y ^= (y << 7) & 0x9d2c5680UL; |
130 | y ^= (y << 7) & 0x9d2c5680UL; |
| Line 134... | Line 134... | ||
| 134 | return y; |
134 | return y; |
| 135 | } |
135 | } |
| 136 | 136 | ||
| 137 | /* generates a random number on [0,0x7fffffff]-interval */ |
137 | /* generates a random number on [0,0x7fffffff]-interval */ |
| 138 | long genrand_int31(void) |
138 | long genrand_int31(void) |
| 139 | { |
139 | { |
| 140 | return (long)(genrand_int32()>>1); |
140 | return (long)(genrand_int32()>>1); |
| 141 | } |
141 | } |
| 142 | 142 | ||
| 143 | /* generates a random number on [0,1]-real-interval */ |
143 | /* generates a random number on [0,1]-real-interval */ |
| 144 | double genrand_real1(void) |
144 | double genrand_real1(void) |
| 145 | { |
145 | { |
| 146 | return genrand_int32()*(1.0/4294967295.0); |
146 | return genrand_int32()*(1.0/4294967295.0); |
| 147 | /* divided by 2^32-1 */ |
147 | /* divided by 2^32-1 */ |
| 148 | } |
148 | } |
| 149 | 149 | ||
| 150 | /* generates a random number on [0,1)-real-interval */ |
150 | /* generates a random number on [0,1)-real-interval */ |
| 151 | double genrand_real2(void) |
151 | double genrand_real2(void) |
| 152 | { |
152 | { |
| 153 | return genrand_int32()*(1.0/4294967296.0); |
153 | return genrand_int32()*(1.0/4294967296.0); |
| 154 | /* divided by 2^32 */ |
154 | /* divided by 2^32 */ |
| 155 | } |
155 | } |
| 156 | 156 | ||
| 157 | /* generates a random number on (0,1)-real-interval */ |
157 | /* generates a random number on (0,1)-real-interval */ |
| 158 | double genrand_real3(void) |
158 | double genrand_real3(void) |
| 159 | { |
159 | { |
| 160 | return (((double)genrand_int32()) + 0.5)*(1.0/4294967296.0); |
160 | return (((double)genrand_int32()) + 0.5)*(1.0/4294967296.0); |
| 161 | /* divided by 2^32 */ |
161 | /* divided by 2^32 */ |
| 162 | } |
162 | } |
| 163 | 163 | ||
| 164 | /* generates a random number on [0,1) with 53-bit resolution*/ |
164 | /* generates a random number on [0,1) with 53-bit resolution*/ |
| 165 | double genrand_res53(void) |
165 | double genrand_res53(void) |
| 166 | { |
166 | { |
| 167 | unsigned long a=genrand_int32()>>5, b=genrand_int32()>>6; |
167 | unsigned long a=genrand_int32()>>5, b=genrand_int32()>>6; |
| 168 | return(a*67108864.0+b)*(1.0/9007199254740992.0); |
168 | return(a*67108864.0+b)*(1.0/9007199254740992.0); |
| 169 | } |
169 | } |
| 170 | /* These real versions are due to Isaku Wada, 2002/01/09 added */ |
170 | /* These real versions are due to Isaku Wada, 2002/01/09 added */ |
| 171 | 171 | ||