Date: Sat, 5 Oct 1996 21:52:27 +0400 (MSD) From: =?KOI8-R?Q?=E1=CE=C4=D2=C5=CA_=FE=C5=D2=CE=CF=D7?= (Andrey A. Chernov) <ache@nagual.ru> To: bde@zeta.org.au (Bruce Evans) Cc: joerg_wunsch@uriah.heep.sax.de, current@FreeBSD.org, freebsd-hackers@FreeBSD.org Subject: Re: I plan to change random() for -current (was Re: rand() and random()) Message-ID: <199610051752.VAA03726@nagual.ru> In-Reply-To: <199610051115.VAA21012@godzilla.zeta.org.au> from "Bruce Evans" at "Oct 5, 96 09:15:05 pm"
next in thread | previous in thread | raw e-mail | index | archive | help
> >IMHO we need to change our random() as suggested. > > How do you know that the suggested method is better? Well, I am not research random generators area personally, but given formulae is 1) well-known as good one, 2) give good results with proposed test, 3) already used in libkern/random.c, 4) better than previous one in all parameters. > ISO C example: next = 1103515245 * next + 12345; > return (unsigned int)(next / 65536) % 32768; > > I.e., it returns bits 16-31 of the current state (right shifted 16). This > is said to be better. Folklore says that someone broke rand() by not > discarding the low bits when ints became 32 bits. It is not good enough to live due to unrandom nature of original formulae. I.e. _all_ bits are unrandom, and lower ones are more unrandom. This formulae not good enough even for 16bits. At this moment I worry about random() only, lets consider rand() things after it. > >! /* > >! * Implements the following, without overflowing 31 bits: > >! * > >! * state[i] = (16807 * state[i - 1]) % 2147483647; > >! * > >! * 2^31-1 (prime) = 2147483647 = 127773*16807+2836 > >! */ > >! long int hi = state[i-1] / 127773; > >! long int lo = state[i-1] % 127773; > >! long int test = 16807*lo - 2836*hi; > >! state[i] = test + (test<0 ? 2147483647 : 0); > > This method is also found in the BSD4.4Lite libkern/rand.c. I guess it > can be trusted (as much as the BSD rand.c :-). So we can safely import it even not from this posting (which can be GNU diseased) but from our own libkern/rand.c Here proposed patch, initial table regenerated to conform default seed value (1): *** random.c.orig Sat Oct 5 20:41:57 1996 --- random.c Sat Oct 5 21:45:01 1996 *************** *** 122,128 **** /* * 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 --- 122,128 ---- /* * 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 *************** *** 135,146 **** static long randtbl[DEG_3 + 1] = { TYPE_3, ! 0x9a319039, 0x32d9c024, 0x9b663182, 0x5da1f342, 0xde3b81e0, 0xdf0a6fb5, ! 0xf103bc02, 0x48f340fb, 0x7449e56b, 0xbeb1dbb0, 0xab5c5918, 0x946554fd, ! 0x8c2e680f, 0xeb3d799f, 0xb11ee0b7, 0x2d436b86, 0xda672e2a, 0x1588ca88, ! 0xe369735d, 0x904f35f7, 0xd7158fd6, 0x6fa6f051, 0x616e6b96, 0xac94efdc, ! 0x36413f93, 0xc622c298, 0xf5a42ab8, 0x8a88d77b, 0xf5ad9d0e, 0x8999220b, ! 0x27fb47b9, }; /* --- 135,146 ---- static long randtbl[DEG_3 + 1] = { TYPE_3, ! 0x991539b1, 0x16a5bce3, 0x6774a4cd, 0x3e01511e, 0x4e508aaa, 0x61048c05, ! 0xf5500617, 0x846b7115, 0x6a19892c, 0x896a97af, 0xdb48f936, 0x14898454, ! 0x37ffd106, 0xb58bff9c, 0x59e17104, 0xcf918a49, 0x09378c83, 0x52c7a471, ! 0x8d293ea9, 0x1f4fc301, 0xc3db71be, 0x39b44e1c, 0xf8a44ef9, 0x4c8b80b1, ! 0x19edc328, 0x87bf4bdd, 0xc9b240e5, 0xe9ee4b1b, 0x4382aee7, 0x535b6b41, ! 0xf3bec5da }; /* *************** *** 192,206 **** srandom(x) unsigned int x; { ! register int i, j; if (rand_type == TYPE_0) state[0] = x; else { - j = 1; state[0] = x; ! for (i = 1; i < rand_deg; i++) ! state[i] = 1103515245 * state[i - 1] + 12345; fptr = &state[rand_sep]; rptr = &state[0]; for (i = 0; i < 10 * rand_deg; i++) --- 192,221 ---- srandom(x) unsigned int x; { ! register int i; if (rand_type == TYPE_0) state[0] = x; else { state[0] = x; ! for (i = 1; i < rand_deg; i++) { ! /* ! * Compute state[i] = (7^5 * state[i - 1]) mod (2^31 - 1) ! * wihout overflowing 31 bits: ! * (2^31 - 1) = 127773 * (7^5) + 2836 ! * From "Random number generators: good ones are hard to find", ! * Park and Miller, Communications of the ACM, vol. 31, no. 10, ! * October 1988, p. 1195. ! */ ! long hi, lo, t; ! ! hi = state[i - 1] / 127773; ! lo = state[i - 1] % 127773; ! t = 16807 * lo - 2836 * hi; ! if (t <= 0) ! t += 0x7fffffff; ! state[i] = t; ! } fptr = &state[rand_sep]; rptr = &state[0]; for (i = 0; i < 10 * rand_deg; i++) *************** *** 348,356 **** { long i; ! if (rand_type == TYPE_0) ! i = state[0] = (state[0] * 1103515245 + 12345) & 0x7fffffff; ! else { *fptr += *rptr; i = (*fptr >> 1) & 0x7fffffff; /* chucking least random bit */ if (++fptr >= end_ptr) { --- 363,386 ---- { long i; ! if (rand_type == TYPE_0) { ! /* ! * Compute i = (7^5 * state[0]) mod (2^31 - 1) ! * wihout overflowing 31 bits: ! * (2^31 - 1) = 127773 * (7^5) + 2836 ! * From "Random number generators: good ones are hard to find", ! * Park and Miller, Communications of the ACM, vol. 31, no. 10, ! * October 1988, p. 1195. ! */ ! long hi, lo; ! ! hi = state[0] / 127773; ! lo = state[0] % 127773; ! i = 16807 * lo - 2836 * hi; ! if (i <= 0) ! i += 0x7fffffff; ! state[0] = i; ! } else { *fptr += *rptr; i = (*fptr >> 1) & 0x7fffffff; /* chucking least random bit */ if (++fptr >= end_ptr) { -- Andrey A. Chernov <ache@nagual.ru> http://www.nagual.ru/~ache/
Want to link to this message? Use this URL: <https://mail-archive.FreeBSD.org/cgi/mid.cgi?199610051752.VAA03726>