diff options
| author | Eric Andersen <andersen@codepoet.org> | 2002-03-19 12:19:07 +0000 | 
|---|---|---|
| committer | Eric Andersen <andersen@codepoet.org> | 2002-03-19 12:19:07 +0000 | 
| commit | e2d6080d4d663c4c8bee9df1d1b6a87fa1944a22 (patch) | |
| tree | 2c2f42427ef6f02241a5a4ea442d8fbea50a9325 /libc/stdlib | |
| parent | 54d956c541aa6ea5a8e39d3db8bb3d4f3c9f4bb2 (diff) | |
Merge glibc random, which gets us a much better RNG and a
reentrant one as well.  It is not much bigger than what we
had, so...
 -Erik
Diffstat (limited to 'libc/stdlib')
| -rw-r--r-- | libc/stdlib/Makefile | 2 | ||||
| -rw-r--r-- | libc/stdlib/random.c | 264 | ||||
| -rw-r--r-- | libc/stdlib/random_r.c | 365 | 
3 files changed, 607 insertions, 24 deletions
diff --git a/libc/stdlib/Makefile b/libc/stdlib/Makefile index 8d7f949df..71136c2fa 100644 --- a/libc/stdlib/Makefile +++ b/libc/stdlib/Makefile @@ -38,7 +38,7 @@ MSRC2=atexit.c  MOBJ2=atexit.o on_exit.o __exit_handler.o exit.o  CSRC =	abort.c getenv.c mktemp.c qsort.c realpath.c bsearch.c \ -	mkstemp.c putenv.c rand.c random.c setenv.c system.c div.c \ +	mkstemp.c putenv.c rand.c random.c random_r.c setenv.c system.c div.c \  	ldiv.c getpt.c ptsname.c grantpt.c unlockpt.c gcvt.c  CSRC+=  drand48.c drand48-iter.c drand48_r.c erand48.c erand48_r.c \  	jrand48.c jrand48_r.c lrand48.c lrand48_r.c mrand48.c mrand48_r.c \ diff --git a/libc/stdlib/random.c b/libc/stdlib/random.c index cbd4206ae..bc20d1e1b 100644 --- a/libc/stdlib/random.c +++ b/libc/stdlib/random.c @@ -1,37 +1,255 @@ -#include <stdlib.h> +/* + * Copyright (c) 1983 Regents of the University of California. + * All rights reserved. + * + * Redistribution and use in source and binary forms are permitted + * provided that the above copyright notice and this paragraph are + * duplicated in all such forms and that any documentation, + * advertising materials, and other materials related to such + * distribution and use acknowledge that the software was developed + * by the University of California, Berkeley.  The name of the + * University may not be used to endorse or promote products derived + * from this software without specific prior written permission. + * THIS SOFTWARE IS PROVIDED ``AS IS'' AND WITHOUT ANY EXPRESS OR + * IMPLIED WARRANTIES, INCLUDING, WITHOUT LIMITATION, THE IMPLIED + * WARRANTIES OF MERCHANTIBILITY AND FITNESS FOR A PARTICULAR PURPOSE. + */  /* - * This generator is a combination of three linear congruential generators - * with periods or 2^15-405, 2^15-1041 and 2^15-1111. It has a period that - * is the product of these three numbers. + * This is derived from the Berkeley source: + *	@(#)random.c	5.5 (Berkeley) 7/6/88 + * It was reworked for the GNU C Library by Roland McGrath. + * Rewritten to use reentrant functions by Ulrich Drepper, 1995.   */ -static long int seed1 = 1; -static long int seed2 = 1; -static long int seed3 = 1; +#define _GNU_SOURCE +#include <features.h> +#include <limits.h> +#include <stddef.h> +#include <stdlib.h> +#ifdef __UCLIBC_HAS_THREADS__ +#include <pthread.h> +/* POSIX.1c requires that there is mutual exclusion for the `rand' and +   `srand' functions to prevent concurrent calls from modifying common +   data.  */ +static pthread_mutex_t lock = PTHREAD_RECURSIVE_MUTEX_INITIALIZER_NP; +#else +#define pthread_mutex_lock(x) +#define pthread_mutex_unlock(x) +#endif + +/* An improved random number generation package.  In addition to the standard +   rand()/srand() like interface, this package also has a special state info +   interface.  The initstate() routine is called with a seed, an array of +   bytes, and a count of how many bytes are being passed in; this array is +   then initialized to contain information for random number generation with +   that much state information.  Good sizes for the amount of state +   information are 32, 64, 128, and 256 bytes.  The state can be switched by +   calling the setstate() function with the same array as was initialized +   with initstate().  By default, the package runs with 128 bytes of state +   information and generates far better random numbers than a linear +   congruential generator.  If the amount of state information is less than +   32 bytes, a simple linear congruential R.N.G. is used.  Internally, the +   state information is treated as an array of longs; the zeroth element of +   the array is the type of R.N.G. being used (small integer); the remainder +   of the array is the state information for the R.N.G.  Thus, 32 bytes of +   state information will give 7 longs worth of state information, which will +   allow a degree seven polynomial.  (Note: The zeroth word of state +   information also has some other information stored in it; see setstate +   for details).  The random number generation technique is a linear feedback +   shift register approach, employing trinomials (since there are fewer terms +   to sum up that way).  In this approach, the least significant bit of all +   the numbers in the state table will act as a linear feedback shift register, +   and will have period 2^deg - 1 (where deg is the degree of the polynomial +   being used, assuming that the polynomial is irreducible and primitive). +   The higher order bits will have longer periods, since their values are +   also influenced by pseudo-random carries out of the lower bits.  The +   total period of the generator is approximately deg*(2**deg - 1); thus +   doubling the amount of state information has a vast influence on the +   period of the generator.  Note: The deg*(2**deg - 1) is an approximation +   only good for large deg, when the period of the shift register is the +   dominant factor.  With deg equal to seven, the period is actually much +   longer than the 7*(2**7 - 1) predicted by this formula.  */ + + + +/* For each of the currently supported random number generators, we have a +   break value on the amount of state information (you need at least this many +   bytes of state info to support this random number generator), a degree for +   the polynomial (actually a trinomial) that the R.N.G. is based on, and +   separation between the two lower order coefficients of the trinomial.  */ + +/* Linear congruential.  */ +#define	TYPE_0		0 +#define	BREAK_0		8 +#define	DEG_0		0 +#define	SEP_0		0 + +/* x**7 + x**3 + 1.  */ +#define	TYPE_1		1 +#define	BREAK_1		32 +#define	DEG_1		7 +#define	SEP_1		3 -#define CRANK(a,b,c,m,s) 	\ -	q = s/a;		\ -	s = b*(s-a*q) - c*q;	\ -	if(s<0) s+=m; +/* x**15 + x + 1.  */ +#define	TYPE_2		2 +#define	BREAK_2		64 +#define	DEG_2		15 +#define	SEP_2		1 -long int random() +/* x**31 + x**3 + 1.  */ +#define	TYPE_3		3 +#define	BREAK_3		128 +#define	DEG_3		31 +#define	SEP_3		3 + +/* x**63 + x + 1.  */ +#define	TYPE_4		4 +#define	BREAK_4		256 +#define	DEG_4		63 +#define	SEP_4		1 + + +/* Array versions of the above information to make code run faster. +   Relies on fact that TYPE_i == i.  */ + +#define	MAX_TYPES	5	/* Max number of types above.  */ + + +/* Initially, everything is set up as if from: +	initstate(1, randtbl, 128); +   Note that this initialization takes advantage of the fact that srandom +   advances the front and rear pointers 10*rand_deg times, and hence the +   rear pointer which starts at 0 will also end up at zero; thus the zeroth +   element of the state information, which contains info about the current +   position of the rear pointer is just +	(MAX_TYPES * (rptr - state)) + TYPE_3 == TYPE_3.  */ + +static int32_t randtbl[DEG_3 + 1] =  { -	register long int q; +    TYPE_3, -	CRANK(206, 157, 31, 32363, seed1); -	CRANK(217, 146, 45, 31727, seed2); -	CRANK(222, 142, 133, 31657, seed3); +    -1726662223, 379960547, 1735697613, 1040273694, 1313901226, +    1627687941, -179304937, -2073333483, 1780058412, -1989503057, +    -615974602, 344556628, 939512070, -1249116260, 1507946756, +    -812545463, 154635395, 1388815473, -1926676823, 525320961, +    -1009028674, 968117788, -123449607, 1284210865, 435012392, +    -2017506339, -911064859, -370259173, 1132637927, 1398500161, +    -205601318, +}; -	return seed1 ^ seed2 ^ seed3; + +static struct random_data unsafe_state = +{ +    /* FPTR and RPTR are two pointers into the state info, a front and a rear +       pointer.  These two pointers are always rand_sep places aparts, as they +       cycle through the state information.  (Yes, this does mean we could get +       away with just one pointer, but the code for random is more efficient +       this way).  The pointers are left positioned as they would be from the call: +       initstate(1, randtbl, 128); +       (The position of the rear pointer, rptr, is really 0 (as explained above +       in the initialization of randtbl) because the state table pointer is set +       to point to randtbl[1] (as explained below).)  */ + +    fptr : &randtbl[SEP_3 + 1], +    rptr : &randtbl[1], + +    /* The following things are the pointer to the state information table, +       the type of the current generator, the degree of the current polynomial +       being used, and the separation between the two pointers. +       Note that for efficiency of random, we remember the first location of +       the state information, not the zeroth.  Hence it is valid to access +       state[-1], which is used to store the type of the R.N.G. +       Also, we remember the last location, since this is more efficient than +       indexing every time to find the address of the last element to see if +       the front and rear pointers have wrapped.  */ + +    state : &randtbl[1], + +    rand_type : TYPE_3, +    rand_deg : DEG_3, +    rand_sep : SEP_3, + +    end_ptr : &randtbl[sizeof (randtbl) / sizeof (randtbl[0])] +}; + + +/* Initialize the random number generator based on the given seed.  If the +   type is the trivial no-state-information type, just remember the seed. +   Otherwise, initializes state[] based on the given "seed" via a linear +   congruential generator.  Then, the pointers are set to known locations +   that are exactly rand_sep places apart.  Lastly, it cycles the state +   information a given number of times to get rid of any initial dependencies +   introduced by the L.C.R.N.G.  Note that the initialization of randtbl[] +   for default usage relies on values produced by this routine.  */ +void srandom (unsigned int x) +{ +    pthread_mutex_lock(&lock); +    srandom_r (x, &unsafe_state); +    pthread_mutex_unlock(&lock); +} +weak_alias (srandom, srand) + +/* Initialize the state information in the given array of N bytes for +   future random number generation.  Based on the number of bytes we +   are given, and the break values for the different R.N.G.'s, we choose +   the best (largest) one we can and set things up for it.  srandom is +   then called to initialize the state information.  Note that on return +   from srandom, we set state[-1] to be the type multiplexed with the current +   value of the rear pointer; this is so successive calls to initstate won't +   lose this information and will be able to restart with setstate. +   Note: The first thing we do is save the current state, if any, just like +   setstate so that it doesn't matter when initstate is called. +   Returns a pointer to the old state.  */ +char * initstate (unsigned int seed, char *arg_state, size_t n) +{ +    int32_t *ostate; + +    pthread_mutex_lock(&lock); +    ostate = &unsafe_state.state[-1]; +    initstate_r (seed, arg_state, n, &unsafe_state); +    pthread_mutex_unlock(&lock); +    return (char *) ostate;  } -void srandom(unsigned int seed) +/* Restore the state from the given state array. +   Note: It is important that we also remember the locations of the pointers +   in the current state information, and restore the locations of the pointers +   from the old state information.  This is done by multiplexing the pointer +   location into the zeroth word of the state information. Note that due +   to the order in which things are done, it is OK to call setstate with the +   same state as the current state +   Returns a pointer to the old state information.  */ +char * setstate (char *arg_state)  { -	seed &= RAND_MAX; -	seed1 = seed % 32362 + 1; -	seed2 = seed % 31726 + 1; -	seed3 = seed % 31656 + 1; +    int32_t *ostate; + +    pthread_mutex_lock(&lock); +    ostate = &unsafe_state.state[-1]; +    if (setstate_r (arg_state, &unsafe_state) < 0) +	ostate = NULL; +    pthread_mutex_unlock(&lock); +    return (char *) ostate; +} + +/* If we are using the trivial TYPE_0 R.N.G., just do the old linear +   congruential bit.  Otherwise, we do our fancy trinomial stuff, which is the +   same in all the other cases due to all the global variables that have been +   set up.  The basic operation is to add the number at the rear pointer into +   the one at the front pointer.  Then both pointers are advanced to the next +   location cyclically in the table.  The value returned is the sum generated, +   reduced to 31 bits by throwing away the "least random" low bit. +   Note: The code takes advantage of the fact that both the front and +   rear pointers can't wrap on the same call by not testing the rear +   pointer if the front one has wrapped.  Returns a 31-bit random number.  */ + +long int random () +{ +  int32_t retval; + +  pthread_mutex_lock(&lock); +  random_r (&unsafe_state, &retval); +  pthread_mutex_unlock(&lock); +  return retval;  } -weak_alias(srandom, srand); diff --git a/libc/stdlib/random_r.c b/libc/stdlib/random_r.c new file mode 100644 index 000000000..fce3a98ab --- /dev/null +++ b/libc/stdlib/random_r.c @@ -0,0 +1,365 @@ +/* + * Copyright (c) 1983 Regents of the University of California. + * All rights reserved. + * + * Redistribution and use in source and binary forms are permitted + * provided that the above copyright notice and this paragraph are + * duplicated in all such forms and that any documentation, + * advertising materials, and other materials related to such + * distribution and use acknowledge that the software was developed + * by the University of California, Berkeley.  The name of the + * University may not be used to endorse or promote products derived + * from this software without specific prior written permission. + * THIS SOFTWARE IS PROVIDED ``AS IS'' AND WITHOUT ANY EXPRESS OR + * IMPLIED WARRANTIES, INCLUDING, WITHOUT LIMITATION, THE IMPLIED + * WARRANTIES OF MERCHANTIBILITY AND FITNESS FOR A PARTICULAR PURPOSE. + */ + +/* + * This is derived from the Berkeley source: + *	@(#)random.c	5.5 (Berkeley) 7/6/88 + * It was reworked for the GNU C Library by Roland McGrath. + * Rewritten to be reentrant by Ulrich Drepper, 1995 + */ + +#define _GNU_SOURCE +#include <features.h> +#include <errno.h> +#include <limits.h> +#include <stddef.h> +#include <stdlib.h> + + + +/* An improved random number generation package.  In addition to the standard +   rand()/srand() like interface, this package also has a special state info +   interface.  The initstate() routine is called with a seed, an array of +   bytes, and a count of how many bytes are being passed in; this array is +   then initialized to contain information for random number generation with +   that much state information.  Good sizes for the amount of state +   information are 32, 64, 128, and 256 bytes.  The state can be switched by +   calling the setstate() function with the same array as was initialized +   with initstate().  By default, the package runs with 128 bytes of state +   information and generates far better random numbers than a linear +   congruential generator.  If the amount of state information is less than +   32 bytes, a simple linear congruential R.N.G. is used.  Internally, the +   state information is treated as an array of longs; the zeroth element of +   the array is the type of R.N.G. being used (small integer); the remainder +   of the array is the state information for the R.N.G.  Thus, 32 bytes of +   state information will give 7 longs worth of state information, which will +   allow a degree seven polynomial.  (Note: The zeroth word of state +   information also has some other information stored in it; see setstate +   for details).  The random number generation technique is a linear feedback +   shift register approach, employing trinomials (since there are fewer terms +   to sum up that way).  In this approach, the least significant bit of all +   the numbers in the state table will act as a linear feedback shift register, +   and will have period 2^deg - 1 (where deg is the degree of the polynomial +   being used, assuming that the polynomial is irreducible and primitive). +   The higher order bits will have longer periods, since their values are +   also influenced by pseudo-random carries out of the lower bits.  The +   total period of the generator is approximately deg*(2**deg - 1); thus +   doubling the amount of state information has a vast influence on the +   period of the generator.  Note: The deg*(2**deg - 1) is an approximation +   only good for large deg, when the period of the shift register is the +   dominant factor.  With deg equal to seven, the period is actually much +   longer than the 7*(2**7 - 1) predicted by this formula.  */ + + + +/* For each of the currently supported random number generators, we have a +   break value on the amount of state information (you need at least this many +   bytes of state info to support this random number generator), a degree for +   the polynomial (actually a trinomial) that the R.N.G. is based on, and +   separation between the two lower order coefficients of the trinomial.  */ + +/* Linear congruential.  */ +#define	TYPE_0		0 +#define	BREAK_0		8 +#define	DEG_0		0 +#define	SEP_0		0 + +/* x**7 + x**3 + 1.  */ +#define	TYPE_1		1 +#define	BREAK_1		32 +#define	DEG_1		7 +#define	SEP_1		3 + +/* x**15 + x + 1.  */ +#define	TYPE_2		2 +#define	BREAK_2		64 +#define	DEG_2		15 +#define	SEP_2		1 + +/* x**31 + x**3 + 1.  */ +#define	TYPE_3		3 +#define	BREAK_3		128 +#define	DEG_3		31 +#define	SEP_3		3 + +/* x**63 + x + 1.  */ +#define	TYPE_4		4 +#define	BREAK_4		256 +#define	DEG_4		63 +#define	SEP_4		1 + + +/* Array versions of the above information to make code run faster. +   Relies on fact that TYPE_i == i.  */ + +#define	MAX_TYPES	5	/* Max number of types above.  */ + +struct random_poly_info +{ +    int seps[MAX_TYPES]; +    int degrees[MAX_TYPES]; +}; + +static const struct random_poly_info random_poly_info = +{ +    { SEP_0, SEP_1, SEP_2, SEP_3, SEP_4 }, +    { DEG_0, DEG_1, DEG_2, DEG_3, DEG_4 } +}; + + + + +/* Initialize the random number generator based on the given seed.  If the +   type is the trivial no-state-information type, just remember the seed. +   Otherwise, initializes state[] based on the given "seed" via a linear +   congruential generator.  Then, the pointers are set to known locations +   that are exactly rand_sep places apart.  Lastly, it cycles the state +   information a given number of times to get rid of any initial dependencies +   introduced by the L.C.R.N.G.  Note that the initialization of randtbl[] +   for default usage relies on values produced by this routine.  */ +int srandom_r (unsigned int seed, struct random_data *buf) +{ +    int type; +    int32_t *state; +    long int i; +    long int word; +    int32_t *dst; +    int kc; + +    if (buf == NULL) +	goto fail; +    type = buf->rand_type; +    if ((unsigned int) type >= MAX_TYPES) +	goto fail; + +    state = buf->state; +    /* We must make sure the seed is not 0.  Take arbitrarily 1 in this case.  */ +    if (seed == 0) +	seed = 1; +    state[0] = seed; +    if (type == TYPE_0) +	goto done; + +    dst = state; +    word = seed; +    kc = buf->rand_deg; +    for (i = 1; i < kc; ++i) +    { +	/* This does: +	   state[i] = (16807 * state[i - 1]) % 2147483647; +	   but avoids overflowing 31 bits.  */ +	long int hi = word / 127773; +	long int lo = word % 127773; +	word = 16807 * lo - 2836 * hi; +	if (word < 0) +	    word += 2147483647; +	*++dst = word; +    } + +    buf->fptr = &state[buf->rand_sep]; +    buf->rptr = &state[0]; +    kc *= 10; +    while (--kc >= 0) +    { +	int32_t discard; +	(void) random_r (buf, &discard); +    } + +done: +    return 0; + +fail: +    return -1; +} + +/* Initialize the state information in the given array of N bytes for +   future random number generation.  Based on the number of bytes we +   are given, and the break values for the different R.N.G.'s, we choose +   the best (largest) one we can and set things up for it.  srandom is +   then called to initialize the state information.  Note that on return +   from srandom, we set state[-1] to be the type multiplexed with the current +   value of the rear pointer; this is so successive calls to initstate won't +   lose this information and will be able to restart with setstate. +   Note: The first thing we do is save the current state, if any, just like +   setstate so that it doesn't matter when initstate is called. +   Returns a pointer to the old state.  */ +int initstate_r (seed, arg_state, n, buf) +     unsigned int seed; +     char *arg_state; +     size_t n; +     struct random_data *buf; +{ +    int type; +    int degree; +    int separation; +    int32_t *state; + +    if (buf == NULL) +	goto fail; + +    if (n >= BREAK_3) +	type = n < BREAK_4 ? TYPE_3 : TYPE_4; +    else if (n < BREAK_1) +    { +	if (n < BREAK_0) +	{ +	    __set_errno (EINVAL); +	    goto fail; +	} +	type = TYPE_0; +    } +    else +	type = n < BREAK_2 ? TYPE_1 : TYPE_2; + +    degree = random_poly_info.degrees[type]; +    separation = random_poly_info.seps[type]; + +    buf->rand_type = type; +    buf->rand_sep = separation; +    buf->rand_deg = degree; +    state = &((int32_t *) arg_state)[1];	/* First location.  */ +    /* Must set END_PTR before srandom.  */ +    buf->end_ptr = &state[degree]; + +    buf->state = state; + +    srandom_r (seed, buf); + +    state[-1] = TYPE_0; +    if (type != TYPE_0) +	state[-1] = (buf->rptr - state) * MAX_TYPES + type; + +    return 0; + +fail: +    __set_errno (EINVAL); +    return -1; +} + +/* Restore the state from the given state array. +   Note: It is important that we also remember the locations of the pointers +   in the current state information, and restore the locations of the pointers +   from the old state information.  This is done by multiplexing the pointer +   location into the zeroth word of the state information. Note that due +   to the order in which things are done, it is OK to call setstate with the +   same state as the current state +   Returns a pointer to the old state information.  */ +int setstate_r (char *arg_state, struct random_data *buf) +{ +    int32_t *new_state = 1 + (int32_t *) arg_state; +    int type; +    int old_type; +    int32_t *old_state; +    int degree; +    int separation; + +    if (arg_state == NULL || buf == NULL) +	goto fail; + +    old_type = buf->rand_type; +    old_state = buf->state; +    if (old_type == TYPE_0) +	old_state[-1] = TYPE_0; +    else +	old_state[-1] = (MAX_TYPES * (buf->rptr - old_state)) + old_type; + +    type = new_state[-1] % MAX_TYPES; +    if (type < TYPE_0 || type > TYPE_4) +	goto fail; + +    buf->rand_deg = degree = random_poly_info.degrees[type]; +    buf->rand_sep = separation = random_poly_info.seps[type]; +    buf->rand_type = type; + +    if (type != TYPE_0) +    { +	int rear = new_state[-1] / MAX_TYPES; +	buf->rptr = &new_state[rear]; +	buf->fptr = &new_state[(rear + separation) % degree]; +    } +    buf->state = new_state; +    /* Set end_ptr too.  */ +    buf->end_ptr = &new_state[degree]; + +    return 0; + +fail: +    __set_errno (EINVAL); +    return -1; +} + +/* If we are using the trivial TYPE_0 R.N.G., just do the old linear +   congruential bit.  Otherwise, we do our fancy trinomial stuff, which is the +   same in all the other cases due to all the global variables that have been +   set up.  The basic operation is to add the number at the rear pointer into +   the one at the front pointer.  Then both pointers are advanced to the next +   location cyclically in the table.  The value returned is the sum generated, +   reduced to 31 bits by throwing away the "least random" low bit. +   Note: The code takes advantage of the fact that both the front and +   rear pointers can't wrap on the same call by not testing the rear +   pointer if the front one has wrapped.  Returns a 31-bit random number.  */ + +int random_r (buf, result) +     struct random_data *buf; +     int32_t *result; +{ +    int32_t *state; + +    if (buf == NULL || result == NULL) +	goto fail; + +    state = buf->state; + +    if (buf->rand_type == TYPE_0) +    { +	int32_t val = state[0]; +	val = ((state[0] * 1103515245) + 12345) & 0x7fffffff; +	state[0] = val; +	*result = val; +    } +    else +    { +	int32_t *fptr = buf->fptr; +	int32_t *rptr = buf->rptr; +	int32_t *end_ptr = buf->end_ptr; +	int32_t val; + +	val = *fptr += *rptr; +	/* Chucking least random bit.  */ +	*result = (val >> 1) & 0x7fffffff; +	++fptr; +	if (fptr >= end_ptr) +	{ +	    fptr = state; +	    ++rptr; +	} +	else +	{ +	    ++rptr; +	    if (rptr >= end_ptr) +		rptr = state; +	} +	buf->fptr = fptr; +	buf->rptr = rptr; +    } +    return 0; + +fail: +    __set_errno (EINVAL); +    return -1; +} +  | 
