diff options
Diffstat (limited to 'src/lib/libc/stdlib/random.c')
| -rw-r--r-- | src/lib/libc/stdlib/random.c | 417 |
1 files changed, 417 insertions, 0 deletions
diff --git a/src/lib/libc/stdlib/random.c b/src/lib/libc/stdlib/random.c new file mode 100644 index 0000000000..5ce7c90ee9 --- /dev/null +++ b/src/lib/libc/stdlib/random.c | |||
| @@ -0,0 +1,417 @@ | |||
| 1 | /* | ||
| 2 | * Copyright (c) 1983 Regents of the University of California. | ||
| 3 | * All rights reserved. | ||
| 4 | * | ||
| 5 | * Redistribution and use in source and binary forms, with or without | ||
| 6 | * modification, are permitted provided that the following conditions | ||
| 7 | * are met: | ||
| 8 | * 1. Redistributions of source code must retain the above copyright | ||
| 9 | * notice, this list of conditions and the following disclaimer. | ||
| 10 | * 2. Redistributions in binary form must reproduce the above copyright | ||
| 11 | * notice, this list of conditions and the following disclaimer in the | ||
| 12 | * documentation and/or other materials provided with the distribution. | ||
| 13 | * 3. All advertising materials mentioning features or use of this software | ||
| 14 | * must display the following acknowledgement: | ||
| 15 | * This product includes software developed by the University of | ||
| 16 | * California, Berkeley and its contributors. | ||
| 17 | * 4. Neither the name of the University nor the names of its contributors | ||
| 18 | * may be used to endorse or promote products derived from this software | ||
| 19 | * without specific prior written permission. | ||
| 20 | * | ||
| 21 | * THIS SOFTWARE IS PROVIDED BY THE REGENTS AND CONTRIBUTORS ``AS IS'' AND | ||
| 22 | * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE | ||
| 23 | * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE | ||
| 24 | * ARE DISCLAIMED. IN NO EVENT SHALL THE REGENTS OR CONTRIBUTORS BE LIABLE | ||
| 25 | * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL | ||
| 26 | * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS | ||
| 27 | * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) | ||
| 28 | * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT | ||
| 29 | * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY | ||
| 30 | * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF | ||
| 31 | * SUCH DAMAGE. | ||
| 32 | */ | ||
| 33 | |||
| 34 | #if defined(LIBC_SCCS) && !defined(lint) | ||
| 35 | static char *rcsid = "$OpenBSD: random.c,v 1.9 2000/04/04 14:27:00 millert Exp $"; | ||
| 36 | #endif /* LIBC_SCCS and not lint */ | ||
| 37 | |||
| 38 | #include <sys/types.h> | ||
| 39 | #include <sys/time.h> | ||
| 40 | #include <fcntl.h> | ||
| 41 | #include <stdio.h> | ||
| 42 | #include <stdlib.h> | ||
| 43 | #include <unistd.h> | ||
| 44 | |||
| 45 | /* | ||
| 46 | * random.c: | ||
| 47 | * | ||
| 48 | * An improved random number generation package. In addition to the standard | ||
| 49 | * rand()/srand() like interface, this package also has a special state info | ||
| 50 | * interface. The initstate() routine is called with a seed, an array of | ||
| 51 | * bytes, and a count of how many bytes are being passed in; this array is | ||
| 52 | * then initialized to contain information for random number generation with | ||
| 53 | * that much state information. Good sizes for the amount of state | ||
| 54 | * information are 32, 64, 128, and 256 bytes. The state can be switched by | ||
| 55 | * calling the setstate() routine with the same array as was initiallized | ||
| 56 | * with initstate(). By default, the package runs with 128 bytes of state | ||
| 57 | * information and generates far better random numbers than a linear | ||
| 58 | * congruential generator. If the amount of state information is less than | ||
| 59 | * 32 bytes, a simple linear congruential R.N.G. is used. | ||
| 60 | * | ||
| 61 | * Internally, the state information is treated as an array of longs; the | ||
| 62 | * zeroeth element of the array is the type of R.N.G. being used (small | ||
| 63 | * integer); the remainder of the array is the state information for the | ||
| 64 | * R.N.G. Thus, 32 bytes of state information will give 7 longs worth of | ||
| 65 | * state information, which will allow a degree seven polynomial. (Note: | ||
| 66 | * the zeroeth word of state information also has some other information | ||
| 67 | * stored in it -- see setstate() for details). | ||
| 68 | * | ||
| 69 | * The random number generation technique is a linear feedback shift register | ||
| 70 | * approach, employing trinomials (since there are fewer terms to sum up that | ||
| 71 | * way). In this approach, the least significant bit of all the numbers in | ||
| 72 | * the state table will act as a linear feedback shift register, and will | ||
| 73 | * have period 2^deg - 1 (where deg is the degree of the polynomial being | ||
| 74 | * used, assuming that the polynomial is irreducible and primitive). The | ||
| 75 | * higher order bits will have longer periods, since their values are also | ||
| 76 | * influenced by pseudo-random carries out of the lower bits. The total | ||
| 77 | * period of the generator is approximately deg*(2**deg - 1); thus doubling | ||
| 78 | * the amount of state information has a vast influence on the period of the | ||
| 79 | * generator. Note: the deg*(2**deg - 1) is an approximation only good for | ||
| 80 | * large deg, when the period of the shift register is the dominant factor. | ||
| 81 | * With deg equal to seven, the period is actually much longer than the | ||
| 82 | * 7*(2**7 - 1) predicted by this formula. | ||
| 83 | */ | ||
| 84 | |||
| 85 | /* | ||
| 86 | * For each of the currently supported random number generators, we have a | ||
| 87 | * break value on the amount of state information (you need at least this | ||
| 88 | * many bytes of state info to support this random number generator), a degree | ||
| 89 | * for the polynomial (actually a trinomial) that the R.N.G. is based on, and | ||
| 90 | * the separation between the two lower order coefficients of the trinomial. | ||
| 91 | */ | ||
| 92 | #define TYPE_0 0 /* linear congruential */ | ||
| 93 | #define BREAK_0 8 | ||
| 94 | #define DEG_0 0 | ||
| 95 | #define SEP_0 0 | ||
| 96 | |||
| 97 | #define TYPE_1 1 /* x**7 + x**3 + 1 */ | ||
| 98 | #define BREAK_1 32 | ||
| 99 | #define DEG_1 7 | ||
| 100 | #define SEP_1 3 | ||
| 101 | |||
| 102 | #define TYPE_2 2 /* x**15 + x + 1 */ | ||
| 103 | #define BREAK_2 64 | ||
| 104 | #define DEG_2 15 | ||
| 105 | #define SEP_2 1 | ||
| 106 | |||
| 107 | #define TYPE_3 3 /* x**31 + x**3 + 1 */ | ||
| 108 | #define BREAK_3 128 | ||
| 109 | #define DEG_3 31 | ||
| 110 | #define SEP_3 3 | ||
| 111 | |||
| 112 | #define TYPE_4 4 /* x**63 + x + 1 */ | ||
| 113 | #define BREAK_4 256 | ||
| 114 | #define DEG_4 63 | ||
| 115 | #define SEP_4 1 | ||
| 116 | |||
| 117 | /* | ||
| 118 | * Array versions of the above information to make code run faster -- | ||
| 119 | * relies on fact that TYPE_i == i. | ||
| 120 | */ | ||
| 121 | #define MAX_TYPES 5 /* max number of types above */ | ||
| 122 | |||
| 123 | static int degrees[MAX_TYPES] = { DEG_0, DEG_1, DEG_2, DEG_3, DEG_4 }; | ||
| 124 | static int seps [MAX_TYPES] = { SEP_0, SEP_1, SEP_2, SEP_3, SEP_4 }; | ||
| 125 | |||
| 126 | /* | ||
| 127 | * Initially, everything is set up as if from: | ||
| 128 | * | ||
| 129 | * initstate(1, &randtbl, 128); | ||
| 130 | * | ||
| 131 | * Note that this initialization takes advantage of the fact that srandom() | ||
| 132 | * advances the front and rear pointers 10*rand_deg times, and hence the | ||
| 133 | * rear pointer which starts at 0 will also end up at zero; thus the zeroeth | ||
| 134 | * element of the state information, which contains info about the current | ||
| 135 | * position of the rear pointer is just | ||
| 136 | * | ||
| 137 | * MAX_TYPES * (rptr - state) + TYPE_3 == TYPE_3. | ||
| 138 | */ | ||
| 139 | |||
| 140 | static long randtbl[DEG_3 + 1] = { | ||
| 141 | TYPE_3, | ||
| 142 | 0x991539b1, 0x16a5bce3, 0x6774a4cd, 0x3e01511e, 0x4e508aaa, 0x61048c05, | ||
| 143 | 0xf5500617, 0x846b7115, 0x6a19892c, 0x896a97af, 0xdb48f936, 0x14898454, | ||
| 144 | 0x37ffd106, 0xb58bff9c, 0x59e17104, 0xcf918a49, 0x09378c83, 0x52c7a471, | ||
| 145 | 0x8d293ea9, 0x1f4fc301, 0xc3db71be, 0x39b44e1c, 0xf8a44ef9, 0x4c8b80b1, | ||
| 146 | 0x19edc328, 0x87bf4bdd, 0xc9b240e5, 0xe9ee4b1b, 0x4382aee7, 0x535b6b41, | ||
| 147 | 0xf3bec5da, | ||
| 148 | }; | ||
| 149 | |||
| 150 | /* | ||
| 151 | * fptr and rptr are two pointers into the state info, a front and a rear | ||
| 152 | * pointer. These two pointers are always rand_sep places aparts, as they | ||
| 153 | * cycle cyclically through the state information. (Yes, this does mean we | ||
| 154 | * could get away with just one pointer, but the code for random() is more | ||
| 155 | * efficient this way). The pointers are left positioned as they would be | ||
| 156 | * from the call | ||
| 157 | * | ||
| 158 | * initstate(1, randtbl, 128); | ||
| 159 | * | ||
| 160 | * (The position of the rear pointer, rptr, is really 0 (as explained above | ||
| 161 | * in the initialization of randtbl) because the state table pointer is set | ||
| 162 | * to point to randtbl[1] (as explained below). | ||
| 163 | */ | ||
| 164 | static long *fptr = &randtbl[SEP_3 + 1]; | ||
| 165 | static long *rptr = &randtbl[1]; | ||
| 166 | |||
| 167 | /* | ||
| 168 | * The following things are the pointer to the state information table, the | ||
| 169 | * type of the current generator, the degree of the current polynomial being | ||
| 170 | * used, and the separation between the two pointers. Note that for efficiency | ||
| 171 | * of random(), we remember the first location of the state information, not | ||
| 172 | * the zeroeth. Hence it is valid to access state[-1], which is used to | ||
| 173 | * store the type of the R.N.G. Also, we remember the last location, since | ||
| 174 | * this is more efficient than indexing every time to find the address of | ||
| 175 | * the last element to see if the front and rear pointers have wrapped. | ||
| 176 | */ | ||
| 177 | static long *state = &randtbl[1]; | ||
| 178 | static int rand_type = TYPE_3; | ||
| 179 | static int rand_deg = DEG_3; | ||
| 180 | static int rand_sep = SEP_3; | ||
| 181 | static long *end_ptr = &randtbl[DEG_3 + 1]; | ||
| 182 | |||
| 183 | /* | ||
| 184 | * srandom: | ||
| 185 | * | ||
| 186 | * Initialize the random number generator based on the given seed. If the | ||
| 187 | * type is the trivial no-state-information type, just remember the seed. | ||
| 188 | * Otherwise, initializes state[] based on the given "seed" via a linear | ||
| 189 | * congruential generator. Then, the pointers are set to known locations | ||
| 190 | * that are exactly rand_sep places apart. Lastly, it cycles the state | ||
| 191 | * information a given number of times to get rid of any initial dependencies | ||
| 192 | * introduced by the L.C.R.N.G. Note that the initialization of randtbl[] | ||
| 193 | * for default usage relies on values produced by this routine. | ||
| 194 | */ | ||
| 195 | void | ||
| 196 | srandom(x) | ||
| 197 | u_int x; | ||
| 198 | { | ||
| 199 | register long int test; | ||
| 200 | register int i; | ||
| 201 | ldiv_t val; | ||
| 202 | |||
| 203 | if (rand_type == TYPE_0) | ||
| 204 | state[0] = x; | ||
| 205 | else { | ||
| 206 | state[0] = x; | ||
| 207 | for (i = 1; i < rand_deg; i++) { | ||
| 208 | /* | ||
| 209 | * Implement the following, without overflowing 31 bits: | ||
| 210 | * | ||
| 211 | * state[i] = (16807 * state[i - 1]) % 2147483647; | ||
| 212 | * | ||
| 213 | * 2^31-1 (prime) = 2147483647 = 127773*16807+2836 | ||
| 214 | */ | ||
| 215 | val = ldiv(state[i-1], 127773); | ||
| 216 | test = 16807 * val.rem - 2836 * val.quot; | ||
| 217 | state[i] = test + (test < 0 ? 2147483647 : 0); | ||
| 218 | } | ||
| 219 | fptr = &state[rand_sep]; | ||
| 220 | rptr = &state[0]; | ||
| 221 | for (i = 0; i < 10 * rand_deg; i++) | ||
| 222 | (void)random(); | ||
| 223 | } | ||
| 224 | } | ||
| 225 | |||
| 226 | /* | ||
| 227 | * srandomdev: | ||
| 228 | * | ||
| 229 | * Many programs choose the seed value in a totally predictable manner. | ||
| 230 | * This often causes problems. We seed the generator using the much more | ||
| 231 | * secure arandom(4) interface. Note that this particular seeding | ||
| 232 | * procedure can generate states which are impossible to reproduce by | ||
| 233 | * calling srandom() with any value, since the succeeding terms in the | ||
| 234 | * state buffer are no longer derived from the LC algorithm applied to | ||
| 235 | * a fixed seed. | ||
| 236 | */ | ||
| 237 | void | ||
| 238 | srandomdev() | ||
| 239 | { | ||
| 240 | int fd; | ||
| 241 | size_t len; | ||
| 242 | |||
| 243 | if (rand_type == TYPE_0) | ||
| 244 | len = sizeof(state[0]); | ||
| 245 | else | ||
| 246 | len = rand_deg * sizeof(state[0]); | ||
| 247 | |||
| 248 | if ((fd = open("/dev/arandom", O_RDONLY, 0)) != -1 && | ||
| 249 | read(fd, (void *) state, len) == (ssize_t) len) { | ||
| 250 | close(fd); | ||
| 251 | } else { | ||
| 252 | struct timeval tv; | ||
| 253 | u_int junk; | ||
| 254 | |||
| 255 | /* XXX - this could be better */ | ||
| 256 | gettimeofday(&tv, NULL); | ||
| 257 | srandom(getpid() ^ tv.tv_sec ^ tv.tv_usec ^ junk); | ||
| 258 | if (fd != -1) | ||
| 259 | close(fd); | ||
| 260 | return; | ||
| 261 | } | ||
| 262 | |||
| 263 | if (rand_type != TYPE_0) { | ||
| 264 | fptr = &state[rand_sep]; | ||
| 265 | rptr = &state[0]; | ||
| 266 | } | ||
| 267 | } | ||
| 268 | |||
| 269 | /* | ||
| 270 | * initstate: | ||
| 271 | * | ||
| 272 | * Initialize the state information in the given array of n bytes for future | ||
| 273 | * random number generation. Based on the number of bytes we are given, and | ||
| 274 | * the break values for the different R.N.G.'s, we choose the best (largest) | ||
| 275 | * one we can and set things up for it. srandom() is then called to | ||
| 276 | * initialize the state information. | ||
| 277 | * | ||
| 278 | * Note that on return from srandom(), we set state[-1] to be the type | ||
| 279 | * multiplexed with the current value of the rear pointer; this is so | ||
| 280 | * successive calls to initstate() won't lose this information and will be | ||
| 281 | * able to restart with setstate(). | ||
| 282 | * | ||
| 283 | * Note: the first thing we do is save the current state, if any, just like | ||
| 284 | * setstate() so that it doesn't matter when initstate is called. | ||
| 285 | * | ||
| 286 | * Returns a pointer to the old state. | ||
| 287 | */ | ||
| 288 | char * | ||
| 289 | initstate(seed, arg_state, n) | ||
| 290 | u_int seed; /* seed for R.N.G. */ | ||
| 291 | char *arg_state; /* pointer to state array */ | ||
| 292 | size_t n; /* # bytes of state info */ | ||
| 293 | { | ||
| 294 | register char *ostate = (char *)(&state[-1]); | ||
| 295 | |||
| 296 | if (rand_type == TYPE_0) | ||
| 297 | state[-1] = rand_type; | ||
| 298 | else | ||
| 299 | state[-1] = MAX_TYPES * (rptr - state) + rand_type; | ||
| 300 | if (n < BREAK_0) | ||
| 301 | return(NULL); | ||
| 302 | if (n < BREAK_1) { | ||
| 303 | rand_type = TYPE_0; | ||
| 304 | rand_deg = DEG_0; | ||
| 305 | rand_sep = SEP_0; | ||
| 306 | } else if (n < BREAK_2) { | ||
| 307 | rand_type = TYPE_1; | ||
| 308 | rand_deg = DEG_1; | ||
| 309 | rand_sep = SEP_1; | ||
| 310 | } else if (n < BREAK_3) { | ||
| 311 | rand_type = TYPE_2; | ||
| 312 | rand_deg = DEG_2; | ||
| 313 | rand_sep = SEP_2; | ||
| 314 | } else if (n < BREAK_4) { | ||
| 315 | rand_type = TYPE_3; | ||
| 316 | rand_deg = DEG_3; | ||
| 317 | rand_sep = SEP_3; | ||
| 318 | } else { | ||
| 319 | rand_type = TYPE_4; | ||
| 320 | rand_deg = DEG_4; | ||
| 321 | rand_sep = SEP_4; | ||
| 322 | } | ||
| 323 | state = &(((long *)arg_state)[1]); /* first location */ | ||
| 324 | end_ptr = &state[rand_deg]; /* must set end_ptr before srandom */ | ||
| 325 | srandom(seed); | ||
| 326 | if (rand_type == TYPE_0) | ||
| 327 | state[-1] = rand_type; | ||
| 328 | else | ||
| 329 | state[-1] = MAX_TYPES*(rptr - state) + rand_type; | ||
| 330 | return(ostate); | ||
| 331 | } | ||
| 332 | |||
| 333 | /* | ||
| 334 | * setstate: | ||
| 335 | * | ||
| 336 | * Restore the state from the given state array. | ||
| 337 | * | ||
| 338 | * Note: it is important that we also remember the locations of the pointers | ||
| 339 | * in the current state information, and restore the locations of the pointers | ||
| 340 | * from the old state information. This is done by multiplexing the pointer | ||
| 341 | * location into the zeroeth word of the state information. | ||
| 342 | * | ||
| 343 | * Note that due to the order in which things are done, it is OK to call | ||
| 344 | * setstate() with the same state as the current state. | ||
| 345 | * | ||
| 346 | * Returns a pointer to the old state information. | ||
| 347 | */ | ||
| 348 | char * | ||
| 349 | setstate(arg_state) | ||
| 350 | const char *arg_state; | ||
| 351 | { | ||
| 352 | register long *new_state = (long *)arg_state; | ||
| 353 | register int type = new_state[0] % MAX_TYPES; | ||
| 354 | register int rear = new_state[0] / MAX_TYPES; | ||
| 355 | char *ostate = (char *)(&state[-1]); | ||
| 356 | |||
| 357 | if (rand_type == TYPE_0) | ||
| 358 | state[-1] = rand_type; | ||
| 359 | else | ||
| 360 | state[-1] = MAX_TYPES * (rptr - state) + rand_type; | ||
| 361 | switch(type) { | ||
| 362 | case TYPE_0: | ||
| 363 | case TYPE_1: | ||
| 364 | case TYPE_2: | ||
| 365 | case TYPE_3: | ||
| 366 | case TYPE_4: | ||
| 367 | rand_type = type; | ||
| 368 | rand_deg = degrees[type]; | ||
| 369 | rand_sep = seps[type]; | ||
| 370 | break; | ||
| 371 | default: | ||
| 372 | return(NULL); | ||
| 373 | } | ||
| 374 | state = &new_state[1]; | ||
| 375 | if (rand_type != TYPE_0) { | ||
| 376 | rptr = &state[rear]; | ||
| 377 | fptr = &state[(rear + rand_sep) % rand_deg]; | ||
| 378 | } | ||
| 379 | end_ptr = &state[rand_deg]; /* set end_ptr too */ | ||
| 380 | return(ostate); | ||
| 381 | } | ||
| 382 | |||
| 383 | /* | ||
| 384 | * random: | ||
| 385 | * | ||
| 386 | * If we are using the trivial TYPE_0 R.N.G., just do the old linear | ||
| 387 | * congruential bit. Otherwise, we do our fancy trinomial stuff, which is | ||
| 388 | * the same in all the other cases due to all the global variables that have | ||
| 389 | * been set up. The basic operation is to add the number at the rear pointer | ||
| 390 | * into the one at the front pointer. Then both pointers are advanced to | ||
| 391 | * the next location cyclically in the table. The value returned is the sum | ||
| 392 | * generated, reduced to 31 bits by throwing away the "least random" low bit. | ||
| 393 | * | ||
| 394 | * Note: the code takes advantage of the fact that both the front and | ||
| 395 | * rear pointers can't wrap on the same call by not testing the rear | ||
| 396 | * pointer if the front one has wrapped. | ||
| 397 | * | ||
| 398 | * Returns a 31-bit random number. | ||
| 399 | */ | ||
| 400 | long | ||
| 401 | random() | ||
| 402 | { | ||
| 403 | long i; | ||
| 404 | |||
| 405 | if (rand_type == TYPE_0) | ||
| 406 | i = state[0] = (state[0] * 1103515245 + 12345) & 0x7fffffff; | ||
| 407 | else { | ||
| 408 | *fptr += *rptr; | ||
| 409 | i = (*fptr >> 1) & 0x7fffffff; /* chucking least random bit */ | ||
| 410 | if (++fptr >= end_ptr) { | ||
| 411 | fptr = state; | ||
| 412 | ++rptr; | ||
| 413 | } else if (++rptr >= end_ptr) | ||
| 414 | rptr = state; | ||
| 415 | } | ||
| 416 | return(i); | ||
| 417 | } | ||
