1 /**********************************************************************
2
3 random.c -
4
5 $Author: michal $
6 $Date: 2002/08/21 15:47:54 $
7 created at: Fri Dec 24 16:39:21 JST 1993
8
9 Copyright (C) 1993-2002 Yukihiro Matsumoto
10
11 **********************************************************************/
12
13 /*
14 This is based on trimmed version of MT19937. To get the original version,
15 contact <http://www.math.keio.ac.jp/~matumoto/emt.html>.
16
17 The original copyright notice follows.
18
19 A C-program for MT19937, with initialization improved 2002/2/10.
20 Coded by Takuji Nishimura and Makoto Matsumoto.
21 This is a faster version by taking Shawn Cokus's optimization,
22 Matthe Bellew's simplification, Isaku Wada's real version.
23
24 Before using, initialize the state by using init_genrand(seed)
25 or init_by_array(init_key, key_length).
26
27 Copyright (C) 1997 - 2002, Makoto Matsumoto and Takuji Nishimura,
28 All rights reserved.
29
30 Redistribution and use in source and binary forms, with or without
31 modification, are permitted provided that the following conditions
32 are met:
33
34 1. Redistributions of source code must retain the above copyright
35 notice, this list of conditions and the following disclaimer.
36
37 2. Redistributions in binary form must reproduce the above copyright
38 notice, this list of conditions and the following disclaimer in the
39 documentation and/or other materials provided with the distribution.
40
41 3. The names of its contributors may not be used to endorse or promote
42 products derived from this software without specific prior written
43 permission.
44
45 THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
46 "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
47 LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
48 A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR
49 CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
50 EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
51 PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
52 PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
53 LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
54 NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
55 SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
56
57
58 Any feedback is very welcome.
59 http://www.math.keio.ac.jp/matumoto/emt.html
60 email: matumoto@math.keio.ac.jp
61 */
62
63 /* Period parameters */
64 #define N 624
65 #define M 397
66 #define MATRIX_A 0x9908b0dfUL /* constant vector a */
67 #define UMASK 0x80000000UL /* most significant w-r bits */
68 #define LMASK 0x7fffffffUL /* least significant r bits */
69 #define MIXBITS(u,v) ( ((u) & UMASK) | ((v) & LMASK) )
70 #define TWIST(u,v) ((MIXBITS(u,v) >> 1) ^ ((v)&1UL ? MATRIX_A : 0UL))
71
72 static unsigned long state[N]; /* the array for the state vector */
73 static int left = 1;
74 static int initf = 0;
75 static unsigned long *next;
76
77 /* initializes state[N] with a seed */
78 static void
79 init_genrand(s)
80 unsigned long s;
81 {
82 int j;
83 state[0]= s & 0xffffffffUL;
84 for (j=1; j<N; j++) {
85 state[j] = (1812433253UL * (state[j-1] ^ (state[j-1] >> 30)) + j);
86 /* See Knuth TAOCP Vol2. 3rd Ed. P.106 for multiplier. */
87 /* In the previous versions, MSBs of the seed affect */
88 /* only MSBs of the array state[]. */
89 /* 2002/01/09 modified by Makoto Matsumoto */
90 state[j] &= 0xffffffffUL; /* for >32 bit machines */
91 }
92 left = 1; initf = 1;
93 }
94
95 static void
96 next_state()
97 {
98 unsigned long *p=state;
99 int j;
100
101 /* if init_genrand() has not been called, */
102 /* a default initial seed is used */
103 if (initf==0) init_genrand(5489UL);
104
105 left = N;
106 next = state;
107
108 for (j=N-M+1; --j; p++)
109 *p = p[M] ^ TWIST(p[0], p[1]);
110
111 for (j=M; --j; p++)
112 *p = p[M-N] ^ TWIST(p[0], p[1]);
113
114 *p = p[M-N] ^ TWIST(p[0], state[0]);
115 }
116
117 /* generates a random number on [0,1)-real-interval */
118 static double
119 genrand_real()
120 {
121 unsigned long y;
122
123 if (--left == 0) next_state();
124 y = *next++;
125
126 /* Tempering */
127 y ^= (y >> 11);
128 y ^= (y << 7) & 0x9d2c5680UL;
129 y ^= (y << 15) & 0xefc60000UL;
130 y ^= (y >> 18);
131
132 return (double)y * (1.0/4294967296.0);
133 /* divided by 2^32 */
134 }
135
136 #undef N
137 #undef M
138
139 /* These real versions are due to Isaku Wada, 2002/01/09 added */
140
141 #include "ruby.h"
142
143 #ifdef HAVE_UNISTD_H
144 #include <unistd.h>
145 #endif
146 #include <time.h>
147 #ifndef NT
148 #ifdef HAVE_SYS_TIME_H
149 # include <sys/time.h>
150 #else
151 struct timeval {
152 long tv_sec; /* seconds */
153 long tv_usec; /* and microseconds */
154 };
155 #endif
156 #endif /* NT */
157
158 static int first = 1;
159
160 static int
161 rand_init(seed)
162 unsigned long seed;
163 {
164 static unsigned long saved_seed;
165 unsigned long old;
166
167 first = 0;
168 init_genrand(seed);
169 old = saved_seed;
170 saved_seed = seed;
171
172 return old;
173 }
174
175 static unsigned long
176 random_seed()
177 {
178 static int n = 0;
179 struct timeval tv;
180
181 gettimeofday(&tv, 0);
182 return tv.tv_sec ^ tv.tv_usec ^ getpid() ^ n++;
183 }
184
185 static VALUE
186 rb_f_srand(argc, argv, obj)
187 int argc;
188 VALUE *argv;
189 VALUE obj;
190 {
191 VALUE sd;
192 unsigned long seed, old;
193
194 rb_secure(4);
195 if (rb_scan_args(argc, argv, "01", &sd) == 0) {
196 seed = random_seed();
197 }
198 else {
199 seed = NUM2ULONG(sd);
200 }
201 old = rand_init(seed);
202
203 return ULONG2NUM(old);
204 }
205
206 static VALUE
207 rb_f_rand(argc, argv, obj)
208 int argc;
209 VALUE *argv;
210 VALUE obj;
211 {
212 VALUE vmax;
213 long val, max;
214
215 rb_scan_args(argc, argv, "01", &vmax);
216 if (first) {
217 rand_init(random_seed());
218 }
219 switch (TYPE(vmax)) {
220 case T_FLOAT:
221 if (RFLOAT(vmax)->value <= LONG_MAX && RFLOAT(vmax)->value >= LONG_MIN) {
222 max = (long)RFLOAT(vmax)->value;
223 break;
224 }
225 vmax = rb_dbl2big(RFLOAT(vmax)->value);
226 /* fall through */
227 {
228 long len = RBIGNUM(vmax)->len;
229 double *buf = ALLOCA_N(double, len);
230
231 while (len--) {
232 buf[len] = genrand_real();
233 }
234 return rb_big_rand(vmax, buf);
235 }
236 case T_NIL:
237 max = 0;
238 break;
239 default:
240 max = NUM2LONG(vmax);
241 break;
242 }
243
244 if (max == 0) {
245 return rb_float_new(genrand_real());
246 }
247 if (max < 0) max = -max;
248 val = max*genrand_real();
249
250 return LONG2NUM(val);
251 }
252
253 void
254 Init_Random()
255 {
256 rb_define_global_function("srand", rb_f_srand, -1);
257 rb_define_global_function("rand", rb_f_rand, -1);
258 }