Re: svn commit: r300956 - head/lib/libc/stdlib
On Tue, 31 May 2016, Chris Torek wrote: That was what I was complaining about. div.c is for C90 (misspelled "ANSI"). It wasn't misspelled when I wrote it. :-) The 1989 ANSI C :-) standard was formally ratified in Dec 1989, and the draft was pretty firm by the time I wrote the code (which I am sure was also 1989, despite the 1990 copyright; we added or updated all the copyrights at the last minute, for net-2). ISO's quick adoption, and hence the name C90, post-date all of that. The weren't quick enough to be in the same year. I usually spell the year as 90, but this can be confusing when discussing the compiler c89 or the compiler flag -std=c89. Hmm, old gcc doesn't have -std=c90, but clang does. Correct rounding for a positive divisor is towards minus infinity so that the remainder is not negative. This is modulo arithmetic and has good algebraic properties. C99 requires rounding minus infinity, at least for positive divisors, under the extension "reliable integer division". Did you state this backwards? For integer divison I see: Oops. I seem to have misedited a whole clause. "towards" is also missing. When integers are divided, the result of the / operator is the algebraic quotient with any fractional part discarded.[105] If the quotient a/b is representable, the expression (a/b)*b + a%b shall equal a; otherwise, the behavior of both a/b and a%b is undefined. which (as footnote 105 notes) is "truncation towards zero", so that (-1)/2 is 0 and not -1. In C90, the rounding is implementation-defined, so it may be correct, but it is "unreliable" since it can be anything. In C90, div() is specified as giving "reliable" division, and that is what the fixups implement. This now wastes time to change nothing. Right -- as long as the compiler must meet C99 rules, div.c can just use the / and % operators. ache added the ifdef. I checked that compilers (old gcc and current clang on amd64) don't auto-inline div(). It is very suitable for inlining without the fixup. Bruce ___ svn-src-all@freebsd.org mailing list https://lists.freebsd.org/mailman/listinfo/svn-src-all To unsubscribe, send any mail to "svn-src-all-unsubscr...@freebsd.org"
Re: svn commit: r300956 - head/lib/libc/stdlib
>That was what I was complaining about. div.c is for C90 (misspelled >"ANSI"). It wasn't misspelled when I wrote it. :-) The 1989 ANSI C standard was formally ratified in Dec 1989, and the draft was pretty firm by the time I wrote the code (which I am sure was also 1989, despite the 1990 copyright; we added or updated all the copyrights at the last minute, for net-2). ISO's quick adoption, and hence the name C90, post-date all of that. >... div.c hasn't caught up with C99. C99 broke this by changing >the rules for integer division to disallow correct rounding for and >require consistently incorrect rounding. Specifically, C99 requires truncation towards zero, so that (-1)/2 = 0 (not -1) and 1/(-2) = 0 (not -1). You (and I) might prefer (-1)/2 = -1, but most CPU "integer divide" ops provide the other result (so that "x >> 1" and "x / 2" are different when x is negative: "x >> 1" sign extends so that (-1)>>1 == -1). C90 (and C99) both require div() and ldiv() to behave like most CPUs, this this truncation-towards-zero behavior, which makes the two functions kind of pointless in C99. >Correct rounding for a >positive divisor is towards minus infinity so that the remainder is >not negative. This is modulo arithmetic and has good algebraic >properties. C99 requires rounding minus infinity, at least for positive >divisors, under the extension "reliable integer division". Did you state this backwards? For integer divison I see: When integers are divided, the result of the / operator is the algebraic quotient with any fractional part discarded.[105] If the quotient a/b is representable, the expression (a/b)*b + a%b shall equal a; otherwise, the behavior of both a/b and a%b is undefined. which (as footnote 105 notes) is "truncation towards zero", so that (-1)/2 is 0 and not -1. >In C90, >the rounding is implementation-defined, so it may be correct, but it >is "unreliable" since it can be anything. > >In C90, div() is specified as giving "reliable" division, and that is >what the fixups implement. This now wastes time to change nothing. Right -- as long as the compiler must meet C99 rules, div.c can just use the / and % operators. Chris ___ svn-src-all@freebsd.org mailing list https://lists.freebsd.org/mailman/listinfo/svn-src-all To unsubscribe, send any mail to "svn-src-all-unsubscr...@freebsd.org"
Re: svn commit: r300956 - head/lib/libc/stdlib
On 31.05.2016 17:17, Bruce Evans wrote: > Our rand should use just 1, and it is dangerous to change RAND_MAX again, > but can we even change the sequences? Something might depend on > reproducing the old sequences. This question already arises many times. The sequence must be reproducible during program runtime only, it is impossible to change implementation at all otherwise. If someone needs something always reproducible it needs to implement their own PRNG. ___ svn-src-all@freebsd.org mailing list https://lists.freebsd.org/mailman/listinfo/svn-src-all To unsubscribe, send any mail to "svn-src-all-unsubscr...@freebsd.org"
Re: svn commit: r300956 - head/lib/libc/stdlib
On Tue, 31 May 2016, Andrey Chernov wrote: On 31.05.2016 12:58, Bruce Evans wrote: On Tue, 31 May 2016, Andrey Chernov wrote: ... You can download SPRNG library implementing all of them here: http://www.sprng.org/RNG/ For me it is overcomplicated. The general case is certainly too complicated. It would let the user specify all the parameters for the LCG and generate optimal code for the choice on the fly. Maybe also change the parameters with the seed. But even most implementors don't know how to choose the parameters. That's just for the not so good LCG family of RNGs. All mentioned PRNGs with exact parameters (info.h) can be found in the subdirs under SRC, f.e. lcg64 or pmlcfg. They are still looks complicated and use GMP Bignum library for calculations. Our rand should use just 1, and it is dangerous to change RAND_MAX again, but can we even change the sequences? Something might depend on reproducing the old sequences. For fixups, see full explanation in the comment of stdlib/div.c That was what I was complaining about. div.c is for C90 (misspelled "ANSI"). div.c hasn't caught up with C99. C99 broke this by changing ... In C90, div() is specified as giving "reliable" division, and that is what the fixups implement. This now wastes time to change nothing. If the hardware does correct rounding, then the compiler already had to do the fixup and id should have produced good MD code for this. The C code then adds not so good MI code. It seems POSIX says nothing about any fixup and negative results: "These functions shall compute the quotient and remainder of the division of the numerator numer by the denominator denom. If the division is inexact, the resulting quotient is the long integer (for the ldiv( ) function) or long long integer (for the lldiv( ) function) of lesser magnitude that is the nearest to the algebraic quotient. If the result cannot be represented, the behavior is undefined; otherwise, quot * denom+rem shall equal numer." This uses the same not so good wording as C90. "algebraic" means "in the field of ordinary rational numbers", not in a general algebraic object. Then we use the ordering property of the rationals, so we're closer to analysis than algebra. POSIX doesn't define the division operator, and the difference is only there, except for bugs in the wording. In C90, the rounding can go either way so the fixup is needed, but in C99 and C11 the division operator is intended to be specified to do the same rounding as div() used to be specified to do, but with the unimproved wording "the result [of division] is the algebraic quotient with any fractional part discarded". C99 and C11 remove the special wording for div(). This removes all ordering requirements from the standard (1 was moved to a footnote), so C has "unreliable" division again, depending on the unspecified choice of sign for the fractional part. Does it mean that this fixup should be removed? Just ifdef it for < C99. This corresponds to removing the special wording in >= C99. Bruce ___ svn-src-all@freebsd.org mailing list https://lists.freebsd.org/mailman/listinfo/svn-src-all To unsubscribe, send any mail to "svn-src-all-unsubscr...@freebsd.org"
Re: svn commit: r300956 - head/lib/libc/stdlib
This implementation of shifts 32-bit and 64-bit PRNGs looks more reliable than Mersenne primes and don't use Bignums: http://wwwmaths.anu.edu.au/~brent/ftp/random/xorgens305.tar.gz PDF article inside. On 31.05.2016 14:30, Andrey Chernov wrote: > On 31.05.2016 12:58, Bruce Evans wrote: >> On Tue, 31 May 2016, Andrey Chernov wrote: >> >>> On 31.05.2016 9:48, Bruce Evans wrote: > Perhaps you can find some ideas, answers and PRNG comparison in the > original paper: > http://www.firstpr.com.au/dsp/rand31/p1192-park.pdf The ones with Mersenne primes and tweaked Mersenne primes in the reference (lanl?) given by pfg@ seem OK. >>> >>> It should be again correctly implemented to not overflow 64bit (as any >>> other 64bit generator too). >>> >>> Moreover, it have visible patterns in the low order bits. This one from >>> the same site is better >>> http://nuclear.llnl.gov/CNP/rng/rngman/node7.html >>> but still have some pattern too and 61-bit. Next one does not suffer >>> from the pattern at all >>> http://nuclear.llnl.gov/CNP/rng/rngman/node8.html >>> but is 2^64-2^10+1. >> >> That's the tweaked Mersenne prime one. >> >> Surely there is a pattern in the bits for all of these, and the better >> ones just hide the pattern better from normal uses and tests? > > Excepting bad first one from pfg@ > http://nuclear.llnl.gov/CNP/rng/rngman/node4.html > as I understand, they not hide (i.e. drop) the pattern but mix it with > other part in calculations, it cause whole range reduced. > >>> You can download SPRNG library implementing all of them here: >>> http://www.sprng.org/RNG/ >>> For me it is overcomplicated. >> >> The general case is certainly too complicated. It would let the user >> specify all the parameters for the LCG and generate optimal code for >> the choice on the fly. Maybe also change the parameters with the seed. >> But even most implementors don't know how to choose the parameters. >> That's just for the not so good LCG family of RNGs. > > All mentioned PRNGs with exact parameters (info.h) can be found in the > subdirs under SRC, f.e. lcg64 or pmlcfg. They are still looks > complicated and use GMP Bignum library for calculations. ___ svn-src-all@freebsd.org mailing list https://lists.freebsd.org/mailman/listinfo/svn-src-all To unsubscribe, send any mail to "svn-src-all-unsubscr...@freebsd.org"
Re: svn commit: r300956 - head/lib/libc/stdlib
On 31.05.2016 12:58, Bruce Evans wrote: > On Tue, 31 May 2016, Andrey Chernov wrote: > >> On 31.05.2016 9:48, Bruce Evans wrote: Perhaps you can find some ideas, answers and PRNG comparison in the original paper: http://www.firstpr.com.au/dsp/rand31/p1192-park.pdf >>> >>> The ones with Mersenne primes and tweaked Mersenne primes in the >>> reference >>> (lanl?) given by pfg@ seem OK. >> >> It should be again correctly implemented to not overflow 64bit (as any >> other 64bit generator too). >> >> Moreover, it have visible patterns in the low order bits. This one from >> the same site is better >> http://nuclear.llnl.gov/CNP/rng/rngman/node7.html >> but still have some pattern too and 61-bit. Next one does not suffer >> from the pattern at all >> http://nuclear.llnl.gov/CNP/rng/rngman/node8.html >> but is 2^64-2^10+1. > > That's the tweaked Mersenne prime one. > > Surely there is a pattern in the bits for all of these, and the better > ones just hide the pattern better from normal uses and tests? Excepting bad first one from pfg@ http://nuclear.llnl.gov/CNP/rng/rngman/node4.html as I understand, they not hide (i.e. drop) the pattern but mix it with other part in calculations, it cause whole range reduced. >> You can download SPRNG library implementing all of them here: >> http://www.sprng.org/RNG/ >> For me it is overcomplicated. > > The general case is certainly too complicated. It would let the user > specify all the parameters for the LCG and generate optimal code for > the choice on the fly. Maybe also change the parameters with the seed. > But even most implementors don't know how to choose the parameters. > That's just for the not so good LCG family of RNGs. All mentioned PRNGs with exact parameters (info.h) can be found in the subdirs under SRC, f.e. lcg64 or pmlcfg. They are still looks complicated and use GMP Bignum library for calculations. >> For fixups, see full >> explanation in the comment of stdlib/div.c > > That was what I was complaining about. div.c is for C90 (misspelled > "ANSI"). div.c hasn't caught up with C99. C99 broke this by changing > the rules for integer division to disallow correct rounding for and > require consistently incorrect rounding. Correct rounding for a > positive divisor is towards minus infinity so that the remainder is > not negative. This is modulo arithmetic and has good algebraic > properties. C99 requires rounding minus infinity, at least for positive > divisors, under the extension "reliable integer division". In C90, > the rounding is implementation- defined, so it may be correct, but it > is "unreliable" since it can be anything. > > In C90, div() is specified as giving "reliable" division, and that is > what the fixups implement. This now wastes time to change nothing. > If the hardware does correct rounding, then the compiler already > had to do the fixup and id should have produced good MD code for this. > The C code then adds not so good MI code. It seems POSIX says nothing about any fixup and negative results: "These functions shall compute the quotient and remainder of the division of the numerator numer by the denominator denom. If the division is inexact, the resulting quotient is the long integer (for the ldiv( ) function) or long long integer (for the lldiv( ) function) of lesser magnitude that is the nearest to the algebraic quotient. If the result cannot be represented, the behavior is undefined; otherwise, quot * denom+rem shall equal numer." Does it mean that this fixup should be removed? ___ svn-src-all@freebsd.org mailing list https://lists.freebsd.org/mailman/listinfo/svn-src-all To unsubscribe, send any mail to "svn-src-all-unsubscr...@freebsd.org"
Re: svn commit: r300956 - head/lib/libc/stdlib
On Tue, 31 May 2016, Andrey Chernov wrote: On 31.05.2016 9:48, Bruce Evans wrote: Perhaps you can find some ideas, answers and PRNG comparison in the original paper: http://www.firstpr.com.au/dsp/rand31/p1192-park.pdf The ones with Mersenne primes and tweaked Mersenne primes in the reference (lanl?) given by pfg@ seem OK. It should be again correctly implemented to not overflow 64bit (as any other 64bit generator too). Moreover, it have visible patterns in the low order bits. This one from the same site is better http://nuclear.llnl.gov/CNP/rng/rngman/node7.html but still have some pattern too and 61-bit. Next one does not suffer from the pattern at all http://nuclear.llnl.gov/CNP/rng/rngman/node8.html but is 2^64-2^10+1. That's the tweaked Mersenne prime one. Surely there is a pattern in the bits for all of these, and the better ones just hide the pattern better from normal uses and tests? Isn't it well known that the lower bits have more patterns so shouldn't be used? This depends on whether the implementation exposes all its bits. Since the rand() specification and FreeBSD man pages don't say anything about this, it is a bug in the implementation to expose all its bits. The implementation is handicapped by using only 32 bits internally but wanting RAND_MAX to be almost 32 bits. You can download SPRNG library implementing all of them here: http://www.sprng.org/RNG/ For me it is overcomplicated. The general case is certainly too complicated. It would let the user specify all the parameters for the LCG and generate optimal code for the choice on the fly. Maybe also change the parameters with the seed. But even most implementors don't know how to choose the parameters. That's just for the not so good LCG family of RNGs. clang -O uses single "idivl" calculating both quotient and reminder for current code, but for ldiv(3) case call/storage/additional calcs overhead will be added. ldiv(3) does not reduce anything, see stdlib/ldiv.c The extern functions give lots of overhead. Sometimes they get replaced automatically by builtins, so it is unclear when the extern functions are actually used. ldiv.c compiles to idivl for the division part but has lots of extras for the fixups. The fixups do nothing except waste time on most hardware/cstd combinations. IIRC, C99 requires direct division to have the same poor rounding rules as idiv(), so idiv() is not really needed in C99 and the fixups in it are negatively needed. The builtins know what is needed so they don't do the fixups in the usual case that the hardware follows the poor rounding rules. We don't have ldiv() builtin in any cae. Apparently integer division is too fundamental to be a builtin. strings -a `which clang` | grep builtin.*div on amd64 gives 47 builtins for division, but none seem to be for integer division. gcc-3.3 has just 4, all for SSE FP division. For fixups, see full explanation in the comment of stdlib/div.c That was what I was complaining about. div.c is for C90 (misspelled "ANSI"). div.c hasn't caught up with C99. C99 broke this by changing the rules for integer division to disallow correct rounding for and require consistently incorrect rounding. Correct rounding for a positive divisor is towards minus infinity so that the remainder is not negative. This is modulo arithmetic and has good algebraic properties. C99 requires rounding minus infinity, at least for positive divisors, under the extension "reliable integer division". In C90, the rounding is implementation- defined, so it may be correct, but it is "unreliable" since it can be anything. In C90, div() is specified as giving "reliable" division, and that is what the fixups implement. This now wastes time to change nothing. If the hardware does correct rounding, then the compiler already had to do the fixup and id should have produced good MD code for this. The C code then adds not so good MI code. Bruce ___ svn-src-all@freebsd.org mailing list https://lists.freebsd.org/mailman/listinfo/svn-src-all To unsubscribe, send any mail to "svn-src-all-unsubscr...@freebsd.org"
Re: svn commit: r300956 - head/lib/libc/stdlib
On 31.05.2016 9:48, Bruce Evans wrote: >> Perhaps you can find some ideas, answers and PRNG comparison in the >> original paper: >> http://www.firstpr.com.au/dsp/rand31/p1192-park.pdf > > The ones with Mersenne primes and tweaked Mersenne primes in the reference > (lanl?) given by pfg@ seem OK. It should be again correctly implemented to not overflow 64bit (as any other 64bit generator too). Moreover, it have visible patterns in the low order bits. This one from the same site is better http://nuclear.llnl.gov/CNP/rng/rngman/node7.html but still have some pattern too and 61-bit. Next one does not suffer from the pattern at all http://nuclear.llnl.gov/CNP/rng/rngman/node8.html but is 2^64-2^10+1. You can download SPRNG library implementing all of them here: http://www.sprng.org/RNG/ For me it is overcomplicated. >> clang -O uses single "idivl" calculating both quotient and reminder for >> current code, but for ldiv(3) case call/storage/additional calcs >> overhead will be added. ldiv(3) does not reduce anything, see >> stdlib/ldiv.c > > The extern functions give lots of overhead. Sometimes they get replaced > automatically by builtins, so it is unclear when the extern functions are > actually used. ldiv.c compiles to idivl for the division part but has > lots of extras for the fixups. The fixups do nothing except waste time > on most hardware/cstd combinations. IIRC, C99 requires direct division > to have the same poor rounding rules as idiv(), so idiv() is not really > needed in C99 and the fixups in it are negatively needed. The builtins > know what is needed so they don't do the fixups in the usual case that > the hardware follows the poor rounding rules. We don't have ldiv() builtin in any cae. For fixups, see full explanation in the comment of stdlib/div.c ___ svn-src-all@freebsd.org mailing list https://lists.freebsd.org/mailman/listinfo/svn-src-all To unsubscribe, send any mail to "svn-src-all-unsubscr...@freebsd.org"
Re: svn commit: r300956 - head/lib/libc/stdlib
On Mon, 30 May 2016, Andrey Chernov wrote: On 30.05.2016 5:17, Bruce Evans wrote: ... Even 1980's compiler technology was not far from reducing the division to a multiplication. The LCG expression would then reduce to (uintN_t)(A * x + B) where N is either 32 or 64. Perhaps N needs to be 64 even with the small coeefficients, due to the divisor being large and not a power of 2. But if we have 64-bit arithmetic, then we can choose much better coefficients than the C90 32-bit ones or the ACM barely 16-bit ones, and uses A * x + B directly, giving a 64-bit period, and have a chance of our 31-bit RAND_MAX finally working. Perhaps you can find some ideas, answers and PRNG comparison in the original paper: http://www.firstpr.com.au/dsp/rand31/p1192-park.pdf The ones with Mersenne primes and tweaked Mersenne primes in the reference (lanl?) given by pfg@ seem OK. I like a simple power of 2 modulus. I once always used fancy prime pairs for hash tables but found that 1 prime and 1 power of 2 works better. Hashing isn't much affected by a few bad cases and I think the fancy primes just hide the bad cases better. While technically equal, having KNF arg types internally and exposing POSIX ones via headers and docs as you suggest leading to harder eye catching, which makes things even worse when POSIX decide to change arg type. The kernel uses different (often wrong) types internally, and now misnames all entry points for syscalls with a sys_ prefix for bogus namespace reasons. clang -O uses single "idivl" calculating both quotient and reminder for current code, but for ldiv(3) case call/storage/additional calcs overhead will be added. ldiv(3) does not reduce anything, see stdlib/ldiv.c The extern functions give lots of overhead. Sometimes they get replaced automatically by builtins, so it is unclear when the extern functions are actually used. ldiv.c compiles to idivl for the division part but has lots of extras for the fixups. The fixups do nothing except waste time on most hardware/cstd combinations. IIRC, C99 requires direct division to have the same poor rounding rules as idiv(), so idiv() is not really needed in C99 and the fixups in it are negatively needed. The builtins know what is needed so they don't do the fixups in the usual case that the hardware follows the poor rounding rules. Bruce ___ svn-src-all@freebsd.org mailing list https://lists.freebsd.org/mailman/listinfo/svn-src-all To unsubscribe, send any mail to "svn-src-all-unsubscr...@freebsd.org"
Re: svn commit: r300956 - head/lib/libc/stdlib
(Interesting discussion)
On 05/29/16 21:17, Bruce Evans wrote:
On Sun, 29 May 2016, Andrey A. Chernov wrote:
Log:
1) Unifdef USE_WEAK_SEEDING since it is too obsolete to support and
makes
reading hard.
Good.
2) Instead of doing range transformation in each and every function
here,
do it single time directly in do_rand(). One "mod" operation overhead
is not
a big deal, but the code looks nicer and possible future functions
additions
or PRNG change do not miss range transformations neither have
unneeded ones.
The whole implementation is silly. It is manually optimized for 1980's
compilers. More below.
3) Use POSIX argument types for visible functions (cosmetic).
Not sure I like type changes.
Modified: head/lib/libc/stdlib/rand.c
==
--- head/lib/libc/stdlib/rand.cSun May 29 12:21:54 2016(r300955)
+++ head/lib/libc/stdlib/rand.cSun May 29 13:57:06 2016(r300956)
@@ -48,14 +48,6 @@ __FBSDID("$FreeBSD$");
static int
do_rand(unsigned long *ctx)
{
-#ifdef USE_WEAK_SEEDING
-/*
- * Historic implementation compatibility.
- * The random sequences do not vary much with the seed,
- * even with overflowing.
- */
-return ((*ctx = *ctx * 1103515245 + 12345) % ((u_long)RAND_MAX +
1));
This is a good implementation of a not very good LCG, made very bad by
botching RAND_MAX. The magic numbers except for RAND_MAX are copied
from the example in the C90 spec. I think they are good enough there.
The comment in at least the C99 spec says "// RAND_MAX assumed to be
32767". This means that these magic numbers were chosen to work with
this value of RAND_MAX. (unsigned) longs are used to give a period
much longer than RAND_MAX and for technical reasons. Taking the modulo
to many fewer bits than the minimum of 32 for an unsigned long then
disguises the linearity. The BSD version almost completly breaks this
on arches with 32 bit longs by taking the modulo to 31 bits (mod 32 bits
would give complete breakage). Arches with 64-bit longs accidentally
work a bit better, by the coefficients are poorly chosen -- they should
be 64 bits and the arithmetic 128 bits.
FWIW, there are coefficients for a 64 bit LCG here:
http://nuclear.llnl.gov/CNP/rng/rngman/node4.html
It would be interesting to have a LCG optimized for modern platforms,
even if we cannot produce more than 31 bits due to the standards.
Pedro.
___
svn-src-all@freebsd.org mailing list
https://lists.freebsd.org/mailman/listinfo/svn-src-all
To unsubscribe, send any mail to "svn-src-all-unsubscr...@freebsd.org"
Re: svn commit: r300956 - head/lib/libc/stdlib
On 30.05.2016 5:17, Bruce Evans wrote:
> On Sun, 29 May 2016, Andrey A. Chernov wrote:
>
>> Log:
>> 1) Unifdef USE_WEAK_SEEDING since it is too obsolete to support and
>> makes
>> reading hard.
>
> Good.
>
>> 2) Instead of doing range transformation in each and every function
>> here,
>> do it single time directly in do_rand(). One "mod" operation overhead
>> is not
>> a big deal, but the code looks nicer and possible future functions
>> additions
>> or PRNG change do not miss range transformations neither have
>> unneeded ones.
>
> The whole implementation is silly. It is manually optimized for 1980's
> compilers. More below.
>
>> 3) Use POSIX argument types for visible functions (cosmetic).
>
> Not sure I like type changes.
>
>> Modified: head/lib/libc/stdlib/rand.c
>> ==
>>
>> --- head/lib/libc/stdlib/rand.cSun May 29 12:21:54 2016(r300955)
>> +++ head/lib/libc/stdlib/rand.cSun May 29 13:57:06 2016(r300956)
>> @@ -48,14 +48,6 @@ __FBSDID("$FreeBSD$");
>> static int
>> do_rand(unsigned long *ctx)
>> {
>> -#ifdef USE_WEAK_SEEDING
>> -/*
>> - * Historic implementation compatibility.
>> - * The random sequences do not vary much with the seed,
>> - * even with overflowing.
>> - */
>> -return ((*ctx = *ctx * 1103515245 + 12345) % ((u_long)RAND_MAX +
>> 1));
>
> This is a good implementation of a not very good LCG, made very bad by
> botching RAND_MAX. The magic numbers except for RAND_MAX are copied
> from the example in the C90 spec. I think they are good enough there.
> The comment in at least the C99 spec says "// RAND_MAX assumed to be
> 32767". This means that these magic numbers were chosen to work with
> this value of RAND_MAX. (unsigned) longs are used to give a period
> much longer than RAND_MAX and for technical reasons. Taking the modulo
> to many fewer bits than the minimum of 32 for an unsigned long then
> disguises the linearity. The BSD version almost completly breaks this
> on arches with 32 bit longs by taking the modulo to 31 bits (mod 32 bits
> would give complete breakage). Arches with 64-bit longs accidentally
> work a bit better, by the coefficients are poorly chosen -- they should
> be 64 bits and the arithmetic 128 bits.
>
>> -#else /* !USE_WEAK_SEEDING */
>> /*
>> * Compute x = (7^5 * x) mod (2^31 - 1)
>> * without overflowing 31 bits:
>
> These coefficients are probably better, but they are still basically
> 32-bit ones and thus not very good for more than a 15-bit RAND_MAX,
> and the details of the calculation are excessively optimized for 1980's
> compilers and 32-bit uintmax_t.
>
> This can be written as x = (1687 * x) % 2147483647 (with some care about
> type sizes and signedness and overflow. It then looks like an even worse
> LCG than the botched C90 one, at least with the botch making its internals
> more visible. E.g., when x = 1, the first couple of iterations don't even
> involve the linear term in 31 bits.
>
> Even 1980's compiler technology was not far from reducing the division
> to a multiplication. The LCG expression would then reduce to
> (uintN_t)(A * x + B) where N is either 32 or 64. Perhaps N needs to
> be 64 even with the small coeefficients, due to the divisor being large
> and not a power of 2. But if we have 64-bit arithmetic, then we can
> choose much better coefficients than the C90 32-bit ones or the ACM
> barely 16-bit ones, and uses A * x + B directly, giving a 64-bit period,
> and have a chance of our 31-bit RAND_MAX finally working.
>
>> @@ -66,48 +58,34 @@ do_rand(unsigned long *ctx)
>> */
>> long hi, lo, x;
>>
>> -/* Must be in [1, 0x7ffe] range at this point. */
>> -hi = *ctx / 127773;
>> -lo = *ctx % 127773;
>> +/* Transform to [1, 0x7ffe] range. */
>> +x = (*ctx % 0x7ffe) + 1;
>> +hi = x / 127773;
>> +lo = x % 127773;
>> x = 16807 * lo - 2836 * hi;
>> if (x < 0)
>> x += 0x7fff;
>
> This does the division more magically but more slowly than the compiler
> would do. It uses one division and one remainder, and doesn't use
> the newfangled (late 1980's) ldiv() function to explicitly try to
> reduce these to one hardware divrem operation. But compilers can
> easily do this reduction. I think compilers can't easily (or perhaps
> at all) reduce to an A * x + B expression. It isn't clear if using
> signed long here makes things easier or harder for compilers. The
> algorithm is special to avoid overflow with signed longs, but it the
> compiler might not understand this. Then -fwrapv would inhibit it
> from doing much reduction, and -fno-wrapv is just complicated.
>
>> -*ctx = x;
>> /* Transform to [0, 0x7ffd] range. */
>> -return (x - 1);
>> -#endif /* !USE_WEAK_SEEDING */
>> +x--;
>> +*ctx = x;
>> +return (x);
>> }
>>
>>
>> int
>> -rand_r(unsigned int *ctx)
>> +rand_r(unsigned *ctx)
>
> You didn't change the type, but
Re: svn commit: r300956 - head/lib/libc/stdlib
On Sun, 29 May 2016, Andrey A. Chernov wrote:
Log:
1) Unifdef USE_WEAK_SEEDING since it is too obsolete to support and makes
reading hard.
Good.
2) Instead of doing range transformation in each and every function here,
do it single time directly in do_rand(). One "mod" operation overhead is not
a big deal, but the code looks nicer and possible future functions additions
or PRNG change do not miss range transformations neither have unneeded ones.
The whole implementation is silly. It is manually optimized for 1980's
compilers. More below.
3) Use POSIX argument types for visible functions (cosmetic).
Not sure I like type changes.
Modified: head/lib/libc/stdlib/rand.c
==
--- head/lib/libc/stdlib/rand.c Sun May 29 12:21:54 2016(r300955)
+++ head/lib/libc/stdlib/rand.c Sun May 29 13:57:06 2016(r300956)
@@ -48,14 +48,6 @@ __FBSDID("$FreeBSD$");
static int
do_rand(unsigned long *ctx)
{
-#ifdef USE_WEAK_SEEDING
-/*
- * Historic implementation compatibility.
- * The random sequences do not vary much with the seed,
- * even with overflowing.
- */
- return ((*ctx = *ctx * 1103515245 + 12345) % ((u_long)RAND_MAX + 1));
This is a good implementation of a not very good LCG, made very bad by
botching RAND_MAX. The magic numbers except for RAND_MAX are copied
from the example in the C90 spec. I think they are good enough there.
The comment in at least the C99 spec says "// RAND_MAX assumed to be
32767". This means that these magic numbers were chosen to work with
this value of RAND_MAX. (unsigned) longs are used to give a period
much longer than RAND_MAX and for technical reasons. Taking the modulo
to many fewer bits than the minimum of 32 for an unsigned long then
disguises the linearity. The BSD version almost completly breaks this
on arches with 32 bit longs by taking the modulo to 31 bits (mod 32 bits
would give complete breakage). Arches with 64-bit longs accidentally
work a bit better, by the coefficients are poorly chosen -- they should
be 64 bits and the arithmetic 128 bits.
-#else /* !USE_WEAK_SEEDING */
/*
* Compute x = (7^5 * x) mod (2^31 - 1)
* without overflowing 31 bits:
These coefficients are probably better, but they are still basically
32-bit ones and thus not very good for more than a 15-bit RAND_MAX,
and the details of the calculation are excessively optimized for 1980's
compilers and 32-bit uintmax_t.
This can be written as x = (1687 * x) % 2147483647 (with some care about
type sizes and signedness and overflow. It then looks like an even worse
LCG than the botched C90 one, at least with the botch making its internals
more visible. E.g., when x = 1, the first couple of iterations don't even
involve the linear term in 31 bits.
Even 1980's compiler technology was not far from reducing the division
to a multiplication. The LCG expression would then reduce to
(uintN_t)(A * x + B) where N is either 32 or 64. Perhaps N needs to
be 64 even with the small coeefficients, due to the divisor being large
and not a power of 2. But if we have 64-bit arithmetic, then we can
choose much better coefficients than the C90 32-bit ones or the ACM
barely 16-bit ones, and uses A * x + B directly, giving a 64-bit period,
and have a chance of our 31-bit RAND_MAX finally working.
@@ -66,48 +58,34 @@ do_rand(unsigned long *ctx)
*/
long hi, lo, x;
- /* Must be in [1, 0x7ffe] range at this point. */
- hi = *ctx / 127773;
- lo = *ctx % 127773;
+ /* Transform to [1, 0x7ffe] range. */
+ x = (*ctx % 0x7ffe) + 1;
+ hi = x / 127773;
+ lo = x % 127773;
x = 16807 * lo - 2836 * hi;
if (x < 0)
x += 0x7fff;
This does the division more magically but more slowly than the compiler
would do. It uses one division and one remainder, and doesn't use
the newfangled (late 1980's) ldiv() function to explicitly try to
reduce these to one hardware divrem operation. But compilers can
easily do this reduction. I think compilers can't easily (or perhaps
at all) reduce to an A * x + B expression. It isn't clear if using
signed long here makes things easier or harder for compilers. The
algorithm is special to avoid overflow with signed longs, but it the
compiler might not understand this. Then -fwrapv would inhibit it
from doing much reduction, and -fno-wrapv is just complicated.
- *ctx = x;
/* Transform to [0, 0x7ffd] range. */
- return (x - 1);
-#endif /* !USE_WEAK_SEEDING */
+ x--;
+ *ctx = x;
+ return (x);
}
int
-rand_r(unsigned int *ctx)
+rand_r(unsigned *ctx)
You didn't change the type, but fixed a style bug (the verbose spelling
of "unsigned") :-). It is interesting that the style bug is missing in
POSIX. I thought that standards mostly got this wrong. Actually, POSIX
almost never uses the verbose spelling. In the a 2001 draft, it has just
2
svn commit: r300956 - head/lib/libc/stdlib
Author: ache
Date: Sun May 29 13:57:06 2016
New Revision: 300956
URL: https://svnweb.freebsd.org/changeset/base/300956
Log:
1) Unifdef USE_WEAK_SEEDING since it is too obsolete to support and makes
reading hard.
2) Instead of doing range transformation in each and every function here,
do it single time directly in do_rand(). One "mod" operation overhead is not
a big deal, but the code looks nicer and possible future functions additions
or PRNG change do not miss range transformations neither have unneeded ones.
3) Use POSIX argument types for visible functions (cosmetic).
MFC after: 1 week
Modified:
head/lib/libc/stdlib/rand.c
Modified: head/lib/libc/stdlib/rand.c
==
--- head/lib/libc/stdlib/rand.c Sun May 29 12:21:54 2016(r300955)
+++ head/lib/libc/stdlib/rand.c Sun May 29 13:57:06 2016(r300956)
@@ -48,14 +48,6 @@ __FBSDID("$FreeBSD$");
static int
do_rand(unsigned long *ctx)
{
-#ifdef USE_WEAK_SEEDING
-/*
- * Historic implementation compatibility.
- * The random sequences do not vary much with the seed,
- * even with overflowing.
- */
- return ((*ctx = *ctx * 1103515245 + 12345) % ((u_long)RAND_MAX + 1));
-#else /* !USE_WEAK_SEEDING */
/*
* Compute x = (7^5 * x) mod (2^31 - 1)
* without overflowing 31 bits:
@@ -66,48 +58,34 @@ do_rand(unsigned long *ctx)
*/
long hi, lo, x;
- /* Must be in [1, 0x7ffe] range at this point. */
- hi = *ctx / 127773;
- lo = *ctx % 127773;
+ /* Transform to [1, 0x7ffe] range. */
+ x = (*ctx % 0x7ffe) + 1;
+ hi = x / 127773;
+ lo = x % 127773;
x = 16807 * lo - 2836 * hi;
if (x < 0)
x += 0x7fff;
- *ctx = x;
/* Transform to [0, 0x7ffd] range. */
- return (x - 1);
-#endif /* !USE_WEAK_SEEDING */
+ x--;
+ *ctx = x;
+ return (x);
}
int
-rand_r(unsigned int *ctx)
+rand_r(unsigned *ctx)
{
u_long val;
int r;
-#ifdef USE_WEAK_SEEDING
val = *ctx;
-#else
- /* Transform to [1, 0x7ffe] range. */
- val = (*ctx % 0x7ffe) + 1;
-#endif
r = do_rand();
-
-#ifdef USE_WEAK_SEEDING
- *ctx = (unsigned int)val;
-#else
- *ctx = (unsigned int)(val - 1);
-#endif
+ *ctx = (unsigned)val;
return (r);
}
-static u_long next =
-#ifdef USE_WEAK_SEEDING
-1;
-#else
-2;
-#endif
+static u_long next = 1;
int
rand(void)
@@ -116,13 +94,9 @@ rand(void)
}
void
-srand(u_int seed)
+srand(unsigned seed)
{
next = seed;
-#ifndef USE_WEAK_SEEDING
- /* Transform to [1, 0x7ffe] range. */
- next = (next % 0x7ffe) + 1;
-#endif
}
@@ -144,10 +118,6 @@ sranddev(void)
mib[0] = CTL_KERN;
mib[1] = KERN_ARND;
sysctl(mib, 2, (void *), , NULL, 0);
-#ifndef USE_WEAK_SEEDING
- /* Transform to [1, 0x7ffe] range. */
- next = (next % 0x7ffe) + 1;
-#endif
}
___
svn-src-all@freebsd.org mailing list
https://lists.freebsd.org/mailman/listinfo/svn-src-all
To unsubscribe, send any mail to "svn-src-all-unsubscr...@freebsd.org"
