00001 /* 00002 * Copyright (c) 1983 Regents of the University of California. 00003 * All rights reserved. 00004 * 00005 * Redistribution and use in source and binary forms, with or without 00006 * modification, are permitted provided that the following conditions 00007 * are met: 00008 * 1. Redistributions of source code must retain the above copyright 00009 * notice, this list of conditions and the following disclaimer. 00010 * 2. Redistributions in binary form must reproduce the above copyright 00011 * notice, this list of conditions and the following disclaimer in the 00012 * documentation and/or other materials provided with the distribution. 00013 * 3. [rescinded 22 July 1999] 00014 * 4. Neither the name of the University nor the names of its contributors 00015 * may be used to endorse or promote products derived from this software 00016 * without specific prior written permission. 00017 * 00018 * THIS SOFTWARE IS PROVIDED BY THE REGENTS AND CONTRIBUTORS ``AS IS'' AND 00019 * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE 00020 * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE 00021 * ARE DISCLAIMED. IN NO EVENT SHALL THE REGENTS OR CONTRIBUTORS BE LIABLE 00022 * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL 00023 * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS 00024 * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) 00025 * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT 00026 * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY 00027 * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF 00028 * SUCH DAMAGE. 00029 */ 00030 00031 /* 00032 * This is derived from the Berkeley source: 00033 * @(#)random.c 5.5 (Berkeley) 7/6/88 00034 * It was reworked for the GNU C Library by Roland McGrath. 00035 */ 00036 00037 /* 00038 00039 @deftypefn Supplement {long int} random (void) 00040 @deftypefnx Supplement void srandom (unsigned int @var{seed}) 00041 @deftypefnx Supplement void* initstate (unsigned int @var{seed}, void *@var{arg_state}, unsigned long @var{n}) 00042 @deftypefnx Supplement void* setstate (void *@var{arg_state}) 00043 00044 Random number functions. @code{random} returns a random number in the 00045 range 0 to @code{LONG_MAX}. @code{srandom} initializes the random 00046 number generator to some starting point determined by @var{seed} 00047 (else, the values returned by @code{random} are always the same for each 00048 run of the program). @code{initstate} and @code{setstate} allow fine-grained 00049 control over the state of the random number generator. 00050 00051 @end deftypefn 00052 00053 */ 00054 00055 #include <errno.h> 00056 00057 #if 0 00058 00059 #include <ansidecl.h> 00060 #include <limits.h> 00061 #include <stddef.h> 00062 #include <stdlib.h> 00063 00064 #else 00065 00066 #define ULONG_MAX ((unsigned long)(~0L)) /* 0xFFFFFFFF for 32-bits */ 00067 #define LONG_MAX ((long)(ULONG_MAX >> 1)) /* 0x7FFFFFFF for 32-bits*/ 00068 00069 #ifdef __STDC__ 00070 # define PTR void * 00071 # ifndef NULL 00072 # define NULL (void *) 0 00073 # endif 00074 #else 00075 # define PTR char * 00076 # ifndef NULL 00077 # define NULL (void *) 0 00078 # endif 00079 #endif 00080 00081 #endif 00082 00083 long int random (); 00084 00085 /* An improved random number generation package. In addition to the standard 00086 rand()/srand() like interface, this package also has a special state info 00087 interface. The initstate() routine is called with a seed, an array of 00088 bytes, and a count of how many bytes are being passed in; this array is 00089 then initialized to contain information for random number generation with 00090 that much state information. Good sizes for the amount of state 00091 information are 32, 64, 128, and 256 bytes. The state can be switched by 00092 calling the setstate() function with the same array as was initiallized 00093 with initstate(). By default, the package runs with 128 bytes of state 00094 information and generates far better random numbers than a linear 00095 congruential generator. If the amount of state information is less than 00096 32 bytes, a simple linear congruential R.N.G. is used. Internally, the 00097 state information is treated as an array of longs; the zeroeth element of 00098 the array is the type of R.N.G. being used (small integer); the remainder 00099 of the array is the state information for the R.N.G. Thus, 32 bytes of 00100 state information will give 7 longs worth of state information, which will 00101 allow a degree seven polynomial. (Note: The zeroeth word of state 00102 information also has some other information stored in it; see setstate 00103 for details). The random number generation technique is a linear feedback 00104 shift register approach, employing trinomials (since there are fewer terms 00105 to sum up that way). In this approach, the least significant bit of all 00106 the numbers in the state table will act as a linear feedback shift register, 00107 and will have period 2^deg - 1 (where deg is the degree of the polynomial 00108 being used, assuming that the polynomial is irreducible and primitive). 00109 The higher order bits will have longer periods, since their values are 00110 also influenced by pseudo-random carries out of the lower bits. The 00111 total period of the generator is approximately deg*(2**deg - 1); thus 00112 doubling the amount of state information has a vast influence on the 00113 period of the generator. Note: The deg*(2**deg - 1) is an approximation 00114 only good for large deg, when the period of the shift register is the 00115 dominant factor. With deg equal to seven, the period is actually much 00116 longer than the 7*(2**7 - 1) predicted by this formula. */ 00117 00118 00119 00120 /* For each of the currently supported random number generators, we have a 00121 break value on the amount of state information (you need at least thi 00122 bytes of state info to support this random number generator), a degree for 00123 the polynomial (actually a trinomial) that the R.N.G. is based on, and 00124 separation between the two lower order coefficients of the trinomial. */ 00125 00126 /* Linear congruential. */ 00127 #define TYPE_0 0 00128 #define BREAK_0 8 00129 #define DEG_0 0 00130 #define SEP_0 0 00131 00132 /* x**7 + x**3 + 1. */ 00133 #define TYPE_1 1 00134 #define BREAK_1 32 00135 #define DEG_1 7 00136 #define SEP_1 3 00137 00138 /* x**15 + x + 1. */ 00139 #define TYPE_2 2 00140 #define BREAK_2 64 00141 #define DEG_2 15 00142 #define SEP_2 1 00143 00144 /* x**31 + x**3 + 1. */ 00145 #define TYPE_3 3 00146 #define BREAK_3 128 00147 #define DEG_3 31 00148 #define SEP_3 3 00149 00150 /* x**63 + x + 1. */ 00151 #define TYPE_4 4 00152 #define BREAK_4 256 00153 #define DEG_4 63 00154 #define SEP_4 1 00155 00156 00157 /* Array versions of the above information to make code run faster. 00158 Relies on fact that TYPE_i == i. */ 00159 00160 #define MAX_TYPES 5 /* Max number of types above. */ 00161 00162 static int degrees[MAX_TYPES] = { DEG_0, DEG_1, DEG_2, DEG_3, DEG_4 }; 00163 static int seps[MAX_TYPES] = { SEP_0, SEP_1, SEP_2, SEP_3, SEP_4 }; 00164 00165 00166 00167 /* Initially, everything is set up as if from: 00168 initstate(1, randtbl, 128); 00169 Note that this initialization takes advantage of the fact that srandom 00170 advances the front and rear pointers 10*rand_deg times, and hence the 00171 rear pointer which starts at 0 will also end up at zero; thus the zeroeth 00172 element of the state information, which contains info about the current 00173 position of the rear pointer is just 00174 (MAX_TYPES * (rptr - state)) + TYPE_3 == TYPE_3. */ 00175 00176 static long int randtbl[DEG_3 + 1] = 00177 { TYPE_3, 00178 0x9a319039, 0x32d9c024, 0x9b663182, 0x5da1f342, 00179 0xde3b81e0, 0xdf0a6fb5, 0xf103bc02, 0x48f340fb, 00180 0x7449e56b, 0xbeb1dbb0, 0xab5c5918, 0x946554fd, 00181 0x8c2e680f, 0xeb3d799f, 0xb11ee0b7, 0x2d436b86, 00182 0xda672e2a, 0x1588ca88, 0xe369735d, 0x904f35f7, 00183 0xd7158fd6, 0x6fa6f051, 0x616e6b96, 0xac94efdc, 00184 0x36413f93, 0xc622c298, 0xf5a42ab8, 0x8a88d77b, 00185 0xf5ad9d0e, 0x8999220b, 0x27fb47b9 00186 }; 00187 00188 /* FPTR and RPTR are two pointers into the state info, a front and a rear 00189 pointer. These two pointers are always rand_sep places aparts, as they 00190 cycle through the state information. (Yes, this does mean we could get 00191 away with just one pointer, but the code for random is more efficient 00192 this way). The pointers are left positioned as they would be from the call: 00193 initstate(1, randtbl, 128); 00194 (The position of the rear pointer, rptr, is really 0 (as explained above 00195 in the initialization of randtbl) because the state table pointer is set 00196 to point to randtbl[1] (as explained below).) */ 00197 00198 static long int *fptr = &randtbl[SEP_3 + 1]; 00199 static long int *rptr = &randtbl[1]; 00200 00201 00202 00203 /* The following things are the pointer to the state information table, 00204 the type of the current generator, the degree of the current polynomial 00205 being used, and the separation between the two pointers. 00206 Note that for efficiency of random, we remember the first location of 00207 the state information, not the zeroeth. Hence it is valid to access 00208 state[-1], which is used to store the type of the R.N.G. 00209 Also, we remember the last location, since this is more efficient than 00210 indexing every time to find the address of the last element to see if 00211 the front and rear pointers have wrapped. */ 00212 00213 static long int *state = &randtbl[1]; 00214 00215 static int rand_type = TYPE_3; 00216 static int rand_deg = DEG_3; 00217 static int rand_sep = SEP_3; 00218 00219 static long int *end_ptr = &randtbl[sizeof(randtbl) / sizeof(randtbl[0])]; 00220 00221 /* Initialize the random number generator based on the given seed. If the 00222 type is the trivial no-state-information type, just remember the seed. 00223 Otherwise, initializes state[] based on the given "seed" via a linear 00224 congruential generator. Then, the pointers are set to known locations 00225 that are exactly rand_sep places apart. Lastly, it cycles the state 00226 information a given number of times to get rid of any initial dependencies 00227 introduced by the L.C.R.N.G. Note that the initialization of randtbl[] 00228 for default usage relies on values produced by this routine. */ 00229 void 00230 srandom (x) 00231 unsigned int x; 00232 { 00233 state[0] = x; 00234 if (rand_type != TYPE_0) 00235 { 00236 register long int i; 00237 for (i = 1; i < rand_deg; ++i) 00238 state[i] = (1103515145 * state[i - 1]) + 12345; 00239 fptr = &state[rand_sep]; 00240 rptr = &state[0]; 00241 for (i = 0; i < 10 * rand_deg; ++i) 00242 random(); 00243 } 00244 } 00245 00246 /* Initialize the state information in the given array of N bytes for 00247 future random number generation. Based on the number of bytes we 00248 are given, and the break values for the different R.N.G.'s, we choose 00249 the best (largest) one we can and set things up for it. srandom is 00250 then called to initialize the state information. Note that on return 00251 from srandom, we set state[-1] to be the type multiplexed with the current 00252 value of the rear pointer; this is so successive calls to initstate won't 00253 lose this information and will be able to restart with setstate. 00254 Note: The first thing we do is save the current state, if any, just like 00255 setstate so that it doesn't matter when initstate is called. 00256 Returns a pointer to the old state. */ 00257 PTR 00258 initstate (seed, arg_state, n) 00259 unsigned int seed; 00260 PTR arg_state; 00261 unsigned long n; 00262 { 00263 PTR ostate = (PTR) &state[-1]; 00264 00265 if (rand_type == TYPE_0) 00266 state[-1] = rand_type; 00267 else 00268 state[-1] = (MAX_TYPES * (rptr - state)) + rand_type; 00269 if (n < BREAK_1) 00270 { 00271 if (n < BREAK_0) 00272 { 00273 errno = EINVAL; 00274 return NULL; 00275 } 00276 rand_type = TYPE_0; 00277 rand_deg = DEG_0; 00278 rand_sep = SEP_0; 00279 } 00280 else if (n < BREAK_2) 00281 { 00282 rand_type = TYPE_1; 00283 rand_deg = DEG_1; 00284 rand_sep = SEP_1; 00285 } 00286 else if (n < BREAK_3) 00287 { 00288 rand_type = TYPE_2; 00289 rand_deg = DEG_2; 00290 rand_sep = SEP_2; 00291 } 00292 else if (n < BREAK_4) 00293 { 00294 rand_type = TYPE_3; 00295 rand_deg = DEG_3; 00296 rand_sep = SEP_3; 00297 } 00298 else 00299 { 00300 rand_type = TYPE_4; 00301 rand_deg = DEG_4; 00302 rand_sep = SEP_4; 00303 } 00304 00305 state = &((long int *) arg_state)[1]; /* First location. */ 00306 /* Must set END_PTR before srandom. */ 00307 end_ptr = &state[rand_deg]; 00308 srandom(seed); 00309 if (rand_type == TYPE_0) 00310 state[-1] = rand_type; 00311 else 00312 state[-1] = (MAX_TYPES * (rptr - state)) + rand_type; 00313 00314 return ostate; 00315 } 00316 00317 /* Restore the state from the given state array. 00318 Note: It is important that we also remember the locations of the pointers 00319 in the current state information, and restore the locations of the pointers 00320 from the old state information. This is done by multiplexing the pointer 00321 location into the zeroeth word of the state information. Note that due 00322 to the order in which things are done, it is OK to call setstate with the 00323 same state as the current state 00324 Returns a pointer to the old state information. */ 00325 00326 PTR 00327 setstate (arg_state) 00328 PTR arg_state; 00329 { 00330 register long int *new_state = (long int *) arg_state; 00331 register int type = new_state[0] % MAX_TYPES; 00332 register int rear = new_state[0] / MAX_TYPES; 00333 PTR ostate = (PTR) &state[-1]; 00334 00335 if (rand_type == TYPE_0) 00336 state[-1] = rand_type; 00337 else 00338 state[-1] = (MAX_TYPES * (rptr - state)) + rand_type; 00339 00340 switch (type) 00341 { 00342 case TYPE_0: 00343 case TYPE_1: 00344 case TYPE_2: 00345 case TYPE_3: 00346 case TYPE_4: 00347 rand_type = type; 00348 rand_deg = degrees[type]; 00349 rand_sep = seps[type]; 00350 break; 00351 default: 00352 /* State info munged. */ 00353 errno = EINVAL; 00354 return NULL; 00355 } 00356 00357 state = &new_state[1]; 00358 if (rand_type != TYPE_0) 00359 { 00360 rptr = &state[rear]; 00361 fptr = &state[(rear + rand_sep) % rand_deg]; 00362 } 00363 /* Set end_ptr too. */ 00364 end_ptr = &state[rand_deg]; 00365 00366 return ostate; 00367 } 00368 00369 /* If we are using the trivial TYPE_0 R.N.G., just do the old linear 00370 congruential bit. Otherwise, we do our fancy trinomial stuff, which is the 00371 same in all ther other cases due to all the global variables that have been 00372 set up. The basic operation is to add the number at the rear pointer into 00373 the one at the front pointer. Then both pointers are advanced to the next 00374 location cyclically in the table. The value returned is the sum generated, 00375 reduced to 31 bits by throwing away the "least random" low bit. 00376 Note: The code takes advantage of the fact that both the front and 00377 rear pointers can't wrap on the same call by not testing the rear 00378 pointer if the front one has wrapped. Returns a 31-bit random number. */ 00379 00380 long int 00381 random () 00382 { 00383 if (rand_type == TYPE_0) 00384 { 00385 state[0] = ((state[0] * 1103515245) + 12345) & LONG_MAX; 00386 return state[0]; 00387 } 00388 else 00389 { 00390 long int i; 00391 *fptr += *rptr; 00392 /* Chucking least random bit. */ 00393 i = (*fptr >> 1) & LONG_MAX; 00394 ++fptr; 00395 if (fptr >= end_ptr) 00396 { 00397 fptr = state; 00398 ++rptr; 00399 } 00400 else 00401 { 00402 ++rptr; 00403 if (rptr >= end_ptr) 00404 rptr = state; 00405 } 00406 return i; 00407 } 00408 }
1.5.6