Go to most recent revision | Details | Last modification | View Log | RSS feed
Rev | Author | Line No. | Line |
---|---|---|---|
1882 | clevermous | 1 | /* Copyright (C) 1995 DJ Delorie, see COPYING.DJ for details */ |
2 | /* This is file RANDOM.C */ |
||
3 | /* This file may have been modified by DJ Delorie (Jan 1995). If so, |
||
4 | ** these modifications are Coyright (C) 1993 DJ Delorie, 24 Kirsten Ave, |
||
5 | ** Rochester NH, 03867-2954, USA. |
||
6 | */ |
||
7 | |||
8 | /* |
||
9 | * Copyright (c) 1983 Regents of the University of California. |
||
10 | * All rights reserved. |
||
11 | * |
||
12 | * Redistribution and use in source and binary forms are permitted |
||
13 | * provided that: (1) source distributions retain this entire copyright |
||
14 | * notice and comment, and (2) distributions including binaries display |
||
15 | * the following acknowledgement: ``This product includes software |
||
16 | * developed by the University of California, Berkeley and its contributors'' |
||
17 | * in the documentation or other materials provided with the distribution |
||
18 | * and in all advertising materials mentioning features or use of this |
||
19 | * software. Neither the name of the University nor the names of its |
||
20 | * contributors may be used to endorse or promote products derived |
||
21 | * from this software without specific prior written permission. |
||
22 | * THIS SOFTWARE IS PROVIDED ``AS IS'' AND WITHOUT ANY EXPRESS OR |
||
23 | * IMPLIED WARRANTIES, INCLUDING, WITHOUT LIMITATION, THE IMPLIED |
||
24 | * WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE. |
||
25 | */ |
||
26 | |||
27 | #include |
||
28 | |||
29 | /* |
||
30 | * random.c: |
||
31 | * An improved random number generation package. In addition to the standard |
||
32 | * rand()/srand() like interface, this package also has a special state info |
||
33 | * interface. The initstate() routine is called with a seed, an array of |
||
34 | * bytes, and a count of how many bytes are being passed in; this array is then |
||
35 | * initialized to contain information for random number generation with that |
||
36 | * much state information. Good sizes for the amount of state information are |
||
37 | * 32, 64, 128, and 256 bytes. The state can be switched by calling the |
||
38 | * setstate() routine with the same array as was initiallized with initstate(). |
||
39 | * By default, the package runs with 128 bytes of state information and |
||
40 | * generates far better random numbers than a linear congruential generator. |
||
41 | * If the amount of state information is less than 32 bytes, a simple linear |
||
42 | * congruential R.N.G. is used. |
||
43 | * Internally, the state information is treated as an array of longs; the |
||
44 | * zeroeth element of the array is the type of R.N.G. being used (small |
||
45 | * integer); the remainder of the array is the state information for the |
||
46 | * R.N.G. Thus, 32 bytes of state information will give 7 longs worth of |
||
47 | * state information, which will allow a degree seven polynomial. (Note: the |
||
48 | * zeroeth word of state information also has some other information stored |
||
49 | * in it -- see setstate() for details). |
||
50 | * The random number generation technique is a linear feedback shift register |
||
51 | * approach, employing trinomials (since there are fewer terms to sum up that |
||
52 | * way). In this approach, the least significant bit of all the numbers in |
||
53 | * the state table will act as a linear feedback shift register, and will have |
||
54 | * period 2^deg - 1 (where deg is the degree of the polynomial being used, |
||
55 | * assuming that the polynomial is irreducible and primitive). The higher |
||
56 | * order bits will have longer periods, since their values are also influenced |
||
57 | * by pseudo-random carries out of the lower bits. The total period of the |
||
58 | * generator is approximately deg*(2**deg - 1); thus doubling the amount of |
||
59 | * state information has a vast influence on the period of the generator. |
||
60 | * Note: the deg*(2**deg - 1) is an approximation only good for large deg, |
||
61 | * when the period of the shift register is the dominant factor. With deg |
||
62 | * equal to seven, the period is actually much longer than the 7*(2**7 - 1) |
||
63 | * predicted by this formula. |
||
64 | */ |
||
65 | |||
66 | |||
67 | |||
68 | /* |
||
69 | * For each of the currently supported random number generators, we have a |
||
70 | * break value on the amount of state information (you need at least this |
||
71 | * many bytes of state info to support this random number generator), a degree |
||
72 | * for the polynomial (actually a trinomial) that the R.N.G. is based on, and |
||
73 | * the separation between the two lower order coefficients of the trinomial. |
||
74 | */ |
||
75 | |||
76 | #define TYPE_0 0 /* linear congruential */ |
||
77 | #define BREAK_0 8 |
||
78 | #define DEG_0 0 |
||
79 | #define SEP_0 0 |
||
80 | |||
81 | #define TYPE_1 1 /* x**7 + x**3 + 1 */ |
||
82 | #define BREAK_1 32 |
||
83 | #define DEG_1 7 |
||
84 | #define SEP_1 3 |
||
85 | |||
86 | #define TYPE_2 2 /* x**15 + x + 1 */ |
||
87 | #define BREAK_2 64 |
||
88 | #define DEG_2 15 |
||
89 | #define SEP_2 1 |
||
90 | |||
91 | #define TYPE_3 3 /* x**31 + x**3 + 1 */ |
||
92 | #define BREAK_3 128 |
||
93 | #define DEG_3 31 |
||
94 | #define SEP_3 3 |
||
95 | |||
96 | #define TYPE_4 4 /* x**63 + x + 1 */ |
||
97 | #define BREAK_4 256 |
||
98 | #define DEG_4 63 |
||
99 | #define SEP_4 1 |
||
100 | |||
101 | |||
102 | /* |
||
103 | * Array versions of the above information to make code run faster -- relies |
||
104 | * on fact that TYPE_i == i. |
||
105 | */ |
||
106 | |||
107 | #define MAX_TYPES 5 /* max number of types above */ |
||
108 | static int degrees[MAX_TYPES] = { DEG_0, DEG_1, DEG_2, DEG_3, DEG_4 }; |
||
109 | static int seps[MAX_TYPES] = { SEP_0, SEP_1, SEP_2, SEP_3, SEP_4 }; |
||
110 | |||
111 | /* |
||
112 | * Initially, everything is set up as if from : |
||
113 | * initstate( 1, &randtbl, 128 ); |
||
114 | * Note that this initialization takes advantage of the fact that srandom() |
||
115 | * advances the front and rear pointers 10*rand_deg times, and hence the |
||
116 | * rear pointer which starts at 0 will also end up at zero; thus the zeroeth |
||
117 | * element of the state information, which contains info about the current |
||
118 | * position of the rear pointer is just |
||
119 | * MAX_TYPES*(rptr - state) + TYPE_3 == TYPE_3. |
||
120 | */ |
||
121 | |||
122 | static unsigned long randtbl[DEG_3 + 1] = { TYPE_3, |
||
123 | 0x9a319039U, 0x32d9c024U, 0x9b663182U, 0x5da1f342U, |
||
124 | 0xde3b81e0U, 0xdf0a6fb5U, 0xf103bc02U, 0x48f340fbU, |
||
125 | 0x7449e56bU, 0xbeb1dbb0U, 0xab5c5918U, 0x946554fdU, |
||
126 | 0x8c2e680fU, 0xeb3d799fU, 0xb11ee0b7U, 0x2d436b86U, |
||
127 | 0xda672e2aU, 0x1588ca88U, 0xe369735dU, 0x904f35f7U, |
||
128 | 0xd7158fd6U, 0x6fa6f051U, 0x616e6b96U, 0xac94efdcU, |
||
129 | 0x36413f93U, 0xc622c298U, 0xf5a42ab8U, 0x8a88d77bU, |
||
130 | 0xf5ad9d0eU, 0x8999220bU, 0x27fb47b9U }; |
||
131 | |||
132 | /* |
||
133 | * fptr and rptr are two pointers into the state info, a front and a rear |
||
134 | * pointer. These two pointers are always rand_sep places aparts, as they cycle |
||
135 | * cyclically through the state information. (Yes, this does mean we could get |
||
136 | * away with just one pointer, but the code for random() is more efficient this |
||
137 | * way). The pointers are left positioned as they would be from the call |
||
138 | * initstate( 1, randtbl, 128 ) |
||
139 | * (The position of the rear pointer, rptr, is really 0 (as explained above |
||
140 | * in the initialization of randtbl) because the state table pointer is set |
||
141 | * to point to randtbl[1] (as explained below). |
||
142 | */ |
||
143 | |||
144 | static long *fptr = &randtbl[ SEP_3 + 1 ]; |
||
145 | static long *rptr = &randtbl[ 1 ]; |
||
146 | |||
147 | /* |
||
148 | * The following things are the pointer to the state information table, |
||
149 | * the type of the current generator, the degree of the current polynomial |
||
150 | * being used, and the separation between the two pointers. |
||
151 | * Note that for efficiency of random(), we remember the first location of |
||
152 | * the state information, not the zeroeth. Hence it is valid to access |
||
153 | * state[-1], which is used to store the type of the R.N.G. |
||
154 | * Also, we remember the last location, since this is more efficient than |
||
155 | * indexing every time to find the address of the last element to see if |
||
156 | * the front and rear pointers have wrapped. |
||
157 | */ |
||
158 | |||
159 | static long *state = &randtbl[ 1 ]; |
||
160 | static int rand_type = TYPE_3; |
||
161 | static int rand_deg = DEG_3; |
||
162 | static int rand_sep = SEP_3; |
||
163 | static long *end_ptr = &randtbl[ DEG_3 + 1 ]; |
||
164 | |||
165 | /* |
||
166 | * srandom: |
||
167 | * Initialize the random number generator based on the given seed. If the |
||
168 | * type is the trivial no-state-information type, just remember the seed. |
||
169 | * Otherwise, initializes state[] based on the given "seed" via a linear |
||
170 | * congruential generator. Then, the pointers are set to known locations |
||
171 | * that are exactly rand_sep places apart. Lastly, it cycles the state |
||
172 | * information a given number of times to get rid of any initial dependencies |
||
173 | * introduced by the L.C.R.N.G. |
||
174 | * Note that the initialization of randtbl[] for default usage relies on |
||
175 | * values produced by this routine. |
||
176 | */ |
||
177 | |||
178 | int |
||
179 | srandom(int x) |
||
180 | { |
||
181 | int i, j; |
||
182 | |||
183 | if (rand_type == TYPE_0) |
||
184 | { |
||
185 | state[ 0 ] = x; |
||
186 | } |
||
187 | else |
||
188 | { |
||
189 | j = 1; |
||
190 | state[ 0 ] = x; |
||
191 | for (i = 1; i < rand_deg; i++) |
||
192 | { |
||
193 | state[i] = 1103515245*state[i - 1] + 12345; |
||
194 | } |
||
195 | fptr = &state[rand_sep]; |
||
196 | rptr = &state[0]; |
||
197 | for( i = 0; i < 10*rand_deg; i++ ) |
||
198 | random(); |
||
199 | } |
||
200 | return 0; |
||
201 | } |
||
202 | |||
203 | /* |
||
204 | * initstate: |
||
205 | * Initialize the state information in the given array of n bytes for |
||
206 | * future random number generation. Based on the number of bytes we |
||
207 | * are given, and the break values for the different R.N.G.'s, we choose |
||
208 | * the best (largest) one we can and set things up for it. srandom() is |
||
209 | * then called to initialize the state information. |
||
210 | * Note that on return from srandom(), we set state[-1] to be the type |
||
211 | * multiplexed with the current value of the rear pointer; this is so |
||
212 | * successive calls to initstate() won't lose this information and will |
||
213 | * be able to restart with setstate(). |
||
214 | * Note: the first thing we do is save the current state, if any, just like |
||
215 | * setstate() so that it doesn't matter when initstate is called. |
||
216 | * Returns a pointer to the old state. |
||
217 | */ |
||
218 | |||
219 | char * |
||
220 | initstate (unsigned seed, char *arg_state, int n) |
||
221 | { |
||
222 | char *ostate = (char *)(&state[ -1 ]); |
||
223 | |||
224 | if (rand_type == TYPE_0) |
||
225 | state[-1] = rand_type; |
||
226 | else |
||
227 | state[-1] = MAX_TYPES * (rptr - state) + rand_type; |
||
228 | if (n < BREAK_1) |
||
229 | { |
||
230 | if (n < BREAK_0) |
||
231 | return 0; |
||
232 | rand_type = TYPE_0; |
||
233 | rand_deg = DEG_0; |
||
234 | rand_sep = SEP_0; |
||
235 | } |
||
236 | else |
||
237 | { |
||
238 | if (n < BREAK_2) |
||
239 | { |
||
240 | rand_type = TYPE_1; |
||
241 | rand_deg = DEG_1; |
||
242 | rand_sep = SEP_1; |
||
243 | } |
||
244 | else |
||
245 | { |
||
246 | if (n < BREAK_3) |
||
247 | { |
||
248 | rand_type = TYPE_2; |
||
249 | rand_deg = DEG_2; |
||
250 | rand_sep = SEP_2; |
||
251 | } |
||
252 | else |
||
253 | { |
||
254 | if (n < BREAK_4) |
||
255 | { |
||
256 | rand_type = TYPE_3; |
||
257 | rand_deg = DEG_3; |
||
258 | rand_sep = SEP_3; |
||
259 | } |
||
260 | else |
||
261 | { |
||
262 | rand_type = TYPE_4; |
||
263 | rand_deg = DEG_4; |
||
264 | rand_sep = SEP_4; |
||
265 | } |
||
266 | } |
||
267 | } |
||
268 | } |
||
269 | state = &(((long *)arg_state)[1]); /* first location */ |
||
270 | end_ptr = &state[rand_deg]; /* must set end_ptr before srandom */ |
||
271 | srandom(seed); |
||
272 | if (rand_type == TYPE_0) |
||
273 | state[-1] = rand_type; |
||
274 | else |
||
275 | state[-1] = MAX_TYPES * (rptr - state) + rand_type; |
||
276 | return ostate; |
||
277 | } |
||
278 | |||
279 | /* |
||
280 | * setstate: |
||
281 | * Restore the state from the given state array. |
||
282 | * Note: it is important that we also remember the locations of the pointers |
||
283 | * in the current state information, and restore the locations of the pointers |
||
284 | * from the old state information. This is done by multiplexing the pointer |
||
285 | * location into the zeroeth word of the state information. |
||
286 | * Note that due to the order in which things are done, it is OK to call |
||
287 | * setstate() with the same state as the current state. |
||
288 | * Returns a pointer to the old state information. |
||
289 | */ |
||
290 | |||
291 | char * |
||
292 | setstate(char *arg_state) |
||
293 | { |
||
294 | long *new_state = (long *)arg_state; |
||
295 | int type = new_state[0]%MAX_TYPES; |
||
296 | int rear = new_state[0]/MAX_TYPES; |
||
297 | char *ostate = (char *)( &state[ -1 ] ); |
||
298 | |||
299 | if (rand_type == TYPE_0) |
||
300 | state[-1] = rand_type; |
||
301 | else |
||
302 | state[-1] = MAX_TYPES * (rptr - state) + rand_type; |
||
303 | switch (type) |
||
304 | { |
||
305 | case TYPE_0: |
||
306 | case TYPE_1: |
||
307 | case TYPE_2: |
||
308 | case TYPE_3: |
||
309 | case TYPE_4: |
||
310 | rand_type = type; |
||
311 | rand_deg = degrees[ type ]; |
||
312 | rand_sep = seps[ type ]; |
||
313 | break; |
||
314 | } |
||
315 | state = &new_state[ 1 ]; |
||
316 | if (rand_type != TYPE_0) |
||
317 | { |
||
318 | rptr = &state[rear]; |
||
319 | fptr = &state[(rear + rand_sep)%rand_deg]; |
||
320 | } |
||
321 | end_ptr = &state[rand_deg]; /* set end_ptr too */ |
||
322 | return ostate; |
||
323 | } |
||
324 | |||
325 | /* |
||
326 | * random: |
||
327 | * If we are using the trivial TYPE_0 R.N.G., just do the old linear |
||
328 | * congruential bit. Otherwise, we do our fancy trinomial stuff, which is the |
||
329 | * same in all ther other cases due to all the global variables that have been |
||
330 | * set up. The basic operation is to add the number at the rear pointer into |
||
331 | * the one at the front pointer. Then both pointers are advanced to the next |
||
332 | * location cyclically in the table. The value returned is the sum generated, |
||
333 | * reduced to 31 bits by throwing away the "least random" low bit. |
||
334 | * Note: the code takes advantage of the fact that both the front and |
||
335 | * rear pointers can't wrap on the same call by not testing the rear |
||
336 | * pointer if the front one has wrapped. |
||
337 | * Returns a 31-bit random number. |
||
338 | */ |
||
339 | |||
340 | long |
||
341 | random(void) |
||
342 | { |
||
343 | long i; |
||
344 | |||
345 | if (rand_type == TYPE_0) |
||
346 | { |
||
347 | i = state[0] = ( state[0]*1103515245 + 12345 )&0x7fffffff; |
||
348 | } |
||
349 | else |
||
350 | { |
||
351 | *fptr += *rptr; |
||
352 | i = (*fptr >> 1)&0x7fffffff; /* chucking least random bit */ |
||
353 | if (++fptr >= end_ptr ) |
||
354 | { |
||
355 | fptr = state; |
||
356 | ++rptr; |
||
357 | } |
||
358 | else |
||
359 | { |
||
360 | if (++rptr >= end_ptr) |
||
361 | rptr = state; |
||
362 | } |
||
363 | } |
||
364 | return i; |
||
365 | }>>>>>>> |