[perl #115928] rand() on Windows only uses 15 bits of entropy

[Perl @ ResonatorSoft . org]

# New Ticket Created by Pe...@ResonatorSoft.org
# Please include the string: [perl #115928]
# in the subject line of all future correspondence about this issue.
# <URL: https://rt.perl.org:443/rt3/Ticket/Display.html?id=115928 >


This is a bug report for perl from Pe...@ResonatorSoft.org,
generated with the help of perlbug 1.39 running under perl 5.16.0.

-----------------------------------------------------------------
[Please describe your issue here]

The Win32 version of Perl doesn't seem to generate random numbers
very well at all, especially in comparison to the Linux version.

Here's an example script detailing the problem:

my %nums;
my $cnt = 0;
foreach my $i (1 .. 1_000_000) {
my $num = rand;
$cnt++ if ($nums{$num});
$nums{$num} = 1;
}
print "$cnt out of 1,000,000\n";

The purpose of the script is to reveal how many duplicate numbers
were generated out of a random sample of one million total numbers.
For Windows, this will output "967232 out of 1,000,000" every time.
In other words, there are only 32768 possible combinations of
"between 0 and 1" floats that rand() will generate.

If I test this same script on a Debian Linux box, I actually get no
dupes, meaning that it has a MUCH higher entropy.

[Please do not change anything below this line]
-----------------------------------------------------------------
---
Flags:
category=core
severity=medium
---
Site configuration information for perl 5.16.0:

Configured by strawberry-perl at Mon May 21 22:08:30 2012.

Summary of my perl5 (revision 5 version 16 subversion 0) configuration:

Platform:
osname=MSWin32, osvers=4.0, archname=MSWin32-x86-multi-thread
uname='Win32 strawberry-perl 5.16.0.1 #1 Mon May 21 22:07:30 2012 i386'
config_args='undef'
hint=recommended, useposix=true, d_sigaction=undef
useithreads=define, usemultiplicity=define
useperlio=define, d_sfio=undef, uselargefiles=define, usesocks=undef
use64bitint=undef, use64bitall=undef, uselongdouble=undef
usemymalloc=n, bincompat5005=undef
Compiler:
cc='gcc', ccflags =' -s -O2 -DWIN32 -DPERL_TEXTMODE_SCRIPTS
-DPERL_IMPLICIT_CONTEXT -DPERL_IMPLICIT_SYS -fno-strict-aliasing
-mms-bitfields',
optimize='-s -O2',
cppflags='-DWIN32'
ccversion='', gccversion='4.6.3', gccosandvers=''
intsize=4, longsize=4, ptrsize=4, doublesize=8, byteorder=1234
d_longlong=undef, longlongsize=8, d_longdbl=define, longdblsize=12
ivtype='long', ivsize=4, nvtype='double', nvsize=8, Off_t='long
long', lseeksize=8
alignbytes=8, prototype=define
Linker and Libraries:
ld='g++', ldflags ='-s -L"C:\strawberry\perl\lib\CORE"
-L"C:\strawberry\c\lib"'
libpth=C:\strawberry\c\lib C:\strawberry\c\i686-w64-mingw32\lib
libs=-lmoldname -lkernel32 -luser32 -lgdi32 -lwinspool -lcomdlg32
-ladvapi32 -lshell32 -lole32 -loleaut32 -lnetapi32 -luuid -lws2_32
-lmpr -lwinmm -lversion -lodbc32 -lodbccp32 -lcomctl32
perllibs=-lmoldname -lkernel32 -luser32 -lgdi32 -lwinspool
-lcomdlg32 -ladvapi32 -lshell32 -lole32 -loleaut32 -lnetapi32 -luuid
-lws2_32 -lmpr -lwinmm -lversion -lodbc32 -lodbccp32 -lcomctl32
libc=, so=dll, useshrplib=true, libperl=libperl516.a
gnulibc_version=''
Dynamic Linking:
dlsrc=dl_win32.xs, dlext=dll, d_dlsymun=undef, ccdlflags=' '
cccdlflags=' ', lddlflags='-mdll -s
-L"C:\strawberry\perl\lib\CORE" -L"C:\strawberry\c\lib"'

Locally applied patches:


---
@INC for perl 5.16.0:
C:/STRAWBERRY/perl/site/lib/MSWin32-x86-multi-thread
C:/STRAWBERRY/perl/site/lib
C:/STRAWBERRY/perl/vendor/lib
C:/STRAWBERRY/perl/lib
.

---
Environment for perl 5.16.0:
HOME (unset)
LANG (unset)
LANGUAGE (unset)
LD_LIBRARY_PATH (unset)
LOGDIR (unset)
PATH=C:\Windows\SYSTEM32;C:\Windows;C:\Windows\SYSTEM32\WBEM;C:\Windows\SYSTEM32\WINDOWSPOWERSHELL\V1.0\;C:\STRAWBERRY\C\BIN;C:\STRAWBERRY\PERL\SITE\BIN;C:\STRAWBERRY\PERL\BIN;C:\CYGWIN\BIN;C:\CYGWIN\USR\BIN;C:\Program
Files\TortoiseSVN\bin;C:\Program Files\Graphviz 2.28\bin;C:\Program
Files\TortoiseGit\bin;c:\Program Files\Microsoft SQL
Server\100\Tools\Binn\VSShell\Common7\IDE\;c:\Program Files\Microsoft
SQL Server\100\Tools\Binn\;c:\Program Files\Microsoft Visual Studio
9.0\Common7\IDE\PrivateAssemblies\;c:\Program Files\Microsoft SQL
Server\100\DTS\Binn\;C:\Program Files\Nmap
PERL_BADLANG (unset)
SHELL (unset)

[Christian Walde via RT]

I replicated the behavior with ActivePerl v5.12.4 (MSWin32-x86-multi-
thread).

---
via perlbug: queue: perl5 status: new
https://rt.perl.org:443/rt3/Ticket/Display.html?id=115928

[Leon Timmermans]

On Tue, Nov 27, 2012 at 3:59 PM, Pe...@ResonatorSoft.org
<perlbug-...@perl.org> wrote:
> The Win32 version of Perl doesn't seem to generate random numbers
> very well at all, especially in comparison to the Linux version.
>
> Here's an example script detailing the problem:
>
> my %nums;
> my $cnt = 0;
> foreach my $i (1 .. 1_000_000) {
> my $num = rand;
> $cnt++ if ($nums{$num});
> $nums{$num} = 1;
> }
> print "$cnt out of 1,000,000\n";
>
> The purpose of the script is to reveal how many duplicate numbers
> were generated out of a random sample of one million total numbers.
> For Windows, this will output "967232 out of 1,000,000" every time.
> In other words, there are only 32768 possible combinations of
> "between 0 and 1" floats that rand() will generate.
>
> If I test this same script on a Debian Linux box, I actually get no
> dupes, meaning that it has a MUCH higher entropy.

Perl already chooses between a bunch of random number sources. If
there's something obviously better available we could use that.

That said, I don't really think this is really severity=medium. rand()
doesn't make any promises on strength, and if you need something
stronger there are modules on CPAN that can provide just that.

Leon

[Christian Walde via RT]

On Tue Nov 27 07:19:53 2012, LeonT wrote:
> If
> there's something obviously better available we could use that.

On Windows XP and higher, there is RtlGenRandom, would that fit the bill?

http://msdn.microsoft.com/en-us/library/aa387694(v=vs.80).aspx
http://blogs.msdn.com/b/michael_howard/archive/2005/01/14/353379.aspx

---
via perlbug: queue: perl5 status: open
https://rt.perl.org:443/rt3/Ticket/Display.html?id=115928

[Leon Timmermans]

On Tue, Nov 27, 2012 at 4:40 PM, Christian Walde via RT
<perlbug-...@perl.org> wrote:
>> If there's something obviously better available we could use that.
>
> On Windows XP and higher, there is RtlGenRandom, would that fit the bill?
>
> http://msdn.microsoft.com/en-us/library/aa387694(v=vs.80).aspx
> http://blogs.msdn.com/b/michael_howard/archive/2005/01/14/353379.aspx

Looks workable, though it's linking is horrible even by Windows standards.

Leon

[Brendan Byrd via RT]

On Tue Nov 27 07:49:28 2012, LeonT wrote:
>
> Looks workable, though it's linking is horrible even by Windows
standards.

BTW, what does Perl currently use for random functionality on Windows?
Honestly, I'm surprised that Windows library routines have that few
combinations.

[Christian Walde via RT]

On Tue Nov 27 07:49:28 2012, LeonT wrote:

Well, as the doc says, if you don't mind a bit of extra memory use, you
can use CryptGenRandom [ http://msdn.microsoft.com/en-us/library/
aa379942(v=vs.85).aspx ] as a source. Am i remembering that right that
Windows XP is the oldest Windows supported by Perl? Because if so, then
that should be no issue.

[Christian Walde via RT]

On Tue Nov 27 07:53:01 2012, SineSwiper wrote:
> On Tue Nov 27 07:49:28 2012, LeonT wrote:
> >
> > Looks workable, though it's linking is horrible even by Windows
> standards.
>
> BTW, what does Perl currently use for random functionality on Windows?
> Honestly, I'm surprised that Windows library routines have that few
> combinations.

Windows currently supplies a rand() in the stdlib that uses 15 bits and
is only meant as a sort of placeholder, since you are meant to select
your own random generator if you actually care about it. Looks like Perl
simply uses that straight-up.

[Craig A. Berry]

On Tue, Nov 27, 2012 at 10:11 AM, Christian Walde via RT
<perlbug-...@perl.org> wrote:
> On Tue Nov 27 07:53:01 2012, SineSwiper wrote:
>> On Tue Nov 27 07:49:28 2012, LeonT wrote:
>> >
>> > Looks workable, though it's linking is horrible even by Windows
>> standards.
>>
>> BTW, what does Perl currently use for random functionality on Windows?
>> Honestly, I'm surprised that Windows library routines have that few
>> combinations.
>
> Windows currently supplies a rand() in the stdlib that uses 15 bits and
> is only meant as a sort of placeholder, since you are meant to select
> your own random generator if you actually care about it. Looks like Perl
> simply uses that straight-up.

It might be as simple as changing the configuration code (e.g., in
win32/config_H.vc and friends) that looks like:

#define Drand01()
(rand()/(double)((unsigned)1<<RANDBITS)) /**/
#define Rand_seed_t unsigned /**/
#define seedDrand01(x) srand((Rand_seed_t)x) /**/
#define RANDBITS 15 /**/

to point to something more capable.

[kmx]

On 27.11.2012 16:40, Christian Walde via RT wrote:
> On Tue Nov 27 07:19:53 2012, LeonT wrote:
>> If
>> there's something obviously better available we could use that.
> On Windows XP and higher, there is RtlGenRandom, would that fit the bill?
>
> http://msdn.microsoft.com/en-us/library/aa387694(v=vs.80).aspx
> http://blogs.msdn.com/b/michael_howard/archive/2005/01/14/353379.aspx

What about CryptGenRandom (available on WinXP+)

{
BYTE randomdata[256];
HCRYPTPROV hCryptProv = 0;
if (CryptAcquireContextW(&hCryptProv, NULL, NULL, PROV_RSA_FULL,
CRYPT_VERIFYCONTEXT)) {
if (CryptGenRandom(hCryptProv, sizeof(randomdata), randomdata)) {
/* success */
}
else {
/* failure */
}
CryptReleaseContext(hCryptProv, 0);
}
}

http://msdn.microsoft.com/en-us/library/windows/desktop/aa379942%28v=vs.85%29.aspx

--
kmx

[demerphq]

Just curious but wouldnt something like that be pretty slow?

I would think that ideally the Win32 version of Perl should "out of
the box" have a reasonably fast RNG.

Yves

--
perl -Mre=debug -e "/just|another|perl|hacker/"

[kmx]

[resending as I forget to Cc: the list]

CryptGenRandom itself is IMO pretty fast the trouble might be with the time spent in CryptAcquireContextW/CryptReleaseContext.

--
kmx

[bulk88 via RT]

On Tue Nov 27 08:11:01 2012, Mithaldu wrote:
> On Tue Nov 27 07:53:01 2012, SineSwiper wrote:
> > On Tue Nov 27 07:49:28 2012, LeonT wrote:
> > >
> > > Looks workable, though it's linking is horrible even by Windows
> > standards.
> >
> > BTW, what does Perl currently use for random functionality on Windows?
> > Honestly, I'm surprised that Windows library routines have that few
> > combinations.
>
> Windows currently supplies a rand() in the stdlib that uses 15 bits and
> is only meant as a sort of placeholder, since you are meant to select
> your own random generator if you actually care about it. Looks like Perl
> simply uses that straight-up.

Why not call the CRT's rand twice or thrice, once for low 15 bits, 2nd
for high 15 bits, and 3rd for remaining 2?

This is the implementation of rand in MS's CRT
http://code.google.com/searchframe#S3vzerue4i0/trunk/reactos/lib/sdk/crt/math/rand.c&type=cs&l=10
. ntdll's RtlRandom and RtlRandomEx are similar to CRT's rand. They make
no further calls. Going through RtlGenRandom or anything else in advapi
will be atleast a dozen if not a 100 times slower (or more, I didn't
benchmark it, I know it is 100s of calls tho and eventually an
ioctl/DeviceIoControl to ksecdd.sys driver for the randomness) by my
educated guess. A faster idea is to use QueryPerformanceCounter each
time, or save a QueryPerformanceCounter as a seed. MS itself uses tick
counts (GetTickCount makes no calls, just polls a global/static updated
by the kernel on every context switch) or QueryPerformanceCounter as
seeds for RtlRandom and RtlRandomEx. Also see
http://code.google.com/searchframe#S3vzerue4i0/trunk/reactos/lib/sdk/crt/startup/gs_support.c&q=__security_init_cookie%20package:reactos-mirror\.googlecode\.com&l=50

--
bulk88 ~ bulk88 at hotmail.com

[Christian Walde via RT]

On Wed Nov 28 00:15:00 2012, bulk88 wrote:
>
> Why not call the CRT's rand twice or thrice, once for low 15 bits, 2nd
> for high 15 bits, and 3rd for remaining 2?

That sounds like a very reasonable solution. Can you maybe implement a
very small c program to generate a million numbers like that and dump
them to STDOUT to check how many collisions result from that algorithm?

[demerphq]

On 28 November 2012 14:43, Christian Walde via RT

<perlbug-...@perl.org> wrote:
> On Wed Nov 28 00:15:00 2012, bulk88 wrote:
>>
>> Why not call the CRT's rand twice or thrice, once for low 15 bits, 2nd
>> for high 15 bits, and 3rd for remaining 2?
>
> That sounds like a very reasonable solution. Can you maybe implement a
> very small c program to generate a million numbers like that and dump
> them to STDOUT to check how many collisions result from that algorithm?

I am just curious why would that be preferable to just coding your own
longer period LCG?

Afaik they amount to a handful of lines of C.

[Christian Walde via RT]

On Wed Nov 28 00:15:00 2012, bulk88 wrote:

> Why not call the CRT's rand twice or thrice, once for low 15 bits, 2nd
> for high 15 bits, and 3rd for remaining 2?

Alright, i talked this through on #p5p and here's the summary as i
understood it:

First off, the offending property of the windows RNG is granularity in
this case. In this case the granularity is a maximum of 15 bits, which
provides a number of unique results that is far below even a modest sample
size.

Further, in Perl the case is that the random bits gotten from the system
are stuffed into a double float, bounded as [0,1[. This further reduces
the actual number of bits used from the RNG, as such a float can only
contain 53 bits.

On linux, drand48 is used, which provides 48 bits of entropy.

demerphq suggested that it might be possible to select an RNG that could
be implemented for Perl in C, and then be used for both Linux and Windows.
Since the granularity of the current linux generator is below the
granularity that a double can provide, this seems to be a viable path.
However i think this is something that should be carefully considered by
people with the appropiate amount of knowledge and experience, which will
take some time and thus might unnecessarily delay a resolution for
windows. As such i think this should definitely be considered in a
separate ticket.

To the issue at hand, granularity on windows, a solution that seems useful
and provides appropiate granularity would be this: https://
gist.github.com/956faeafb71b2ac735aa

I am unsure if it could be extended to provide 53 bits of entropy in a
double even on 32 bit systems and how it would be patched into Perl as it
is. bulk88, i've been told by demerphq that you have the chops and the
platform for this. Could you please look into this?

[Craig A. Berry]

The FreeBSD implementation of drand48() is indeed short and sweet even
if you follow the links a couple of levels deep:

<http://fxr.watson.org/fxr/source/gen/drand48.c?v=FREEBSD-LIBC;im=bigexcerpts>

It might just work verbatim on Windows.

[kmx]

Another candidate (long period, good speed, good randomness) that can be
quite easy to implement is MT19937
https://en.wikipedia.org/wiki/Mersenne_twister

The actual implementation (approx. 100 lines of code) can look like this:
https://metacpan.org/source/FANGLY/Math-Random-MT-1.16/_mt.c

--
kmx

[bulk88 via RT]

On Wed Nov 28 05:43:55 2012, Mithaldu wrote:
> On Wed Nov 28 00:15:00 2012, bulk88 wrote:
> >
> > Why not call the CRT's rand twice or thrice, once for low 15 bits, 2nd
> > for high 15 bits, and 3rd for remaining 2?
>
> That sounds like a very reasonable solution. Can you maybe implement a
> very small c program to generate a million numbers like that and dump
> them to STDOUT to check how many collisions result from that algorithm?

I am going to try to fix this by spreading the bits of whatever rand a
perl is using across the integer bits of a NV (usually 53 bits for a 64
bit double). I plan to test it as XSUBs first and I will post those here.

--
bulk88 ~ bulk88 at hotmail.com

[Green, Paul]

FWIW, the C-1990 and C-1999 standards define the rand() function as follows (relevant language extracted):

"The rand function computes a sequence of pseudo-random integers in the range 0 to RAND_MAX. ... The value of the RAND_MAX shall be at least 32767."

The Stratus OpenVOS implementation of rand() also limits itself to a RAND_MAX of 32767; I suspect this is due to a desire to maintain compatibility with legacy C programs. Once these limits are defined and widely used, it is impossible to change them without breaking working code.

PG

[h...@crypt.org]

"Christian Walde via RT" <perlbug-...@perl.org> wrote:

:To the issue at hand, granularity on windows, a solution that seems useful

:and provides appropiate granularity would be this: https://
:gist.github.com/956faeafb71b2ac735aa

If the RNG has a period of 2^15, I believe this solution will have
the same period. (If the RNG had a nearby period that happened to be
a multiple of 3, I believe this solution would have a period only
a third as long.)

Hugo

[Tony Cook]

That doesn't prevent the implementation storing more than 15 bits of
entropy.

It could simply return the lower 15 bits of the result from a 32 or
larger LCRNG implementation, as the simplest example.

Tony

[Brendan Byrd via RT]

Perl isn't C. We document that rand returns "a random fractional number
greater than or equal to 0 and less than the value of EXPR". While bits
of entropy isn't documented here, I think we can all agree that 15 bits
is too low.

I like bulk77's solution. It's fast and portably achieves the maximum
entropy based on IV length. Any crypto-security issues or special algorithms can be covered in other modules.

[Aristotle Pagaltzis]

* bulk88 via RT <perlbug-...@perl.org> [2012-11-28 09:20]:

> Why not call the CRT's rand twice or thrice, once for low 15 bits, 2nd
> for high 15 bits, and 3rd for remaining 2?

That doesn’t increase the entropy available, it just draws from the pool
faster. I don’t have my randomness down well enough to reason it through
with the necessary rigour but I’m pretty sure that at best that approach
does nothing to increase the quality of the random numbers and likely it
will decrease it.

Don’t.

Randomness is a treacherous sea – easy to get it wrong, hard to notice
when you do, and with far-reaching consequences. (Debian patched some
random number generator and caused half the world to have to mint new
SSH keys a few years down the line.)

Find a syscall that draws from a deeper pool or put in some other known
algorithm with a greater amount of hidden state or really just about any
other approach except this one.

Or just leave it as it is. If it’s bad, just let it look bad. Being bad
is not a problem (absent promises about quality – and perl makes none)
but looking better than it is would be.

The documentation already covers everything it needs to, so nothing to
be done there either.

Regards,
--
Aristotle Pagaltzis // <http://plasmasturm.org/>

[Aristotle Pagaltzis]

* Aristotle Pagaltzis <paga...@gmx.de> [2012-11-29 03:01]:

> Randomness is a treacherous sea – easy to get it wrong, hard to notice
> when you do, and with far-reaching consequences. (Debian patched some
> random number generator and caused half the world to have to mint new
> SSH keys a few years down the line.)

And just as I send this mail, I notice that last week’s LWN has come out
from under the paywall:

LCE: Don’t play dice with random numbers
https://lwn.net/Articles/525459/

[Jesse Luehrs]

As was pointed out on IRC, there are actually two orthogonal concerns
here - how long until the sequence of random numbers repeats, and how
many distinct numbers are able to be generated. The proposed solution
here will increase the second while potentially (although not
necessarily) decreasing the first. This isn't necessarily a bad
tradeoff, depending on what the period actually is (it could easily be
far greater than the number of distinct numbers that can be generated,
if the generator keeps hidden internal state). There's nothing
inherently problematic about this idea, as long as you're aware that
it's only going to increase the potential range of numbers generated,
not increase the randomness itself.

That said, I do tend to agree with Yves here - if we're actually
concerned about the quality of the random numbers that are avaialble,
good rngs are pretty trivial to implement (or steal), so we may as well
just do that if want to go down this road.

-doy

[bulk88 via RT]

On Wed Nov 28 23:13:32 2012, bulk88 wrote:
> I guess I have to somehow cache the Perl_pow
> results in global statics and calculate them once per process startup.
> The whole loop can be unrolled actually, but it wouldn't be portable if
> it is unrolled by hand.

On the collisions test posted by Brendan Byrd, the new rand returned 0.

[bulk88 via RT]

On Wed Nov 28 09:41:59 2012, bulk88 wrote:
> On Wed Nov 28 05:43:55 2012, Mithaldu wrote:
> > On Wed Nov 28 00:15:00 2012, bulk88 wrote:
> > >
> > > Why not call the CRT's rand twice or thrice, once for low 15 bits, 2nd
> > > for high 15 bits, and 3rd for remaining 2?
> >
> > That sounds like a very reasonable solution. Can you maybe implement a
> > very small c program to generate a million numbers like that and dump
> > them to STDOUT to check how many collisions result from that algorithm?
>
> I am going to try to fix this by spreading the bits of whatever rand a
> perl is using across the integer bits of a NV (usually 53 bits for a 64
> bit double). I plan to test it as XSUBs first and I will post those here.

Windows only right now.
______________________________________________________--

NV
RunRand(...)
PREINIT:
LARGE_INTEGER my_beg;
LARGE_INTEGER my_end;
LARGE_INTEGER my_freq;
NV value;
int io;
CODE:
if (!PL_srand_called) {
(void)seedDrand01((Rand_seed_t)seed());
PL_srand_called = TRUE;
}
RETVAL = 0.0;
QueryPerformanceFrequency(&my_freq);

QueryPerformanceCounter(&my_beg);
for(io=0; io < 10000000; io++)
{
int i;
NV randnv = 0.0;
const int nvofbits = NV_OVERFLOW_BITS;
const NV overflowat = NV_OVERFLOWS_INTEGERS_AT;
NV mulby;
value = 1.0;
for(i=0; i < nvofbits; i += RANDBITS){
NV randraw = rand();
if (nvofbits - i < RANDBITS) {
/* In the last loop iteration, we have
to cut off the high bits, so we don't overflow the
mantissa
*/
/* scalb is slower than pow */
NV mask = Perl_pow(2, nvofbits - i > RANDBITS ?
RANDBITS : nvofbits - i) - 1;
randraw = Perl_fmod(randraw, mask); /* & the bits */

}
randraw *= Perl_pow(2, i); /* << the bits */
randnv += randraw; /* | the bits */
}
mulby = (randnv / overflowat); /* >> so randnv is between
0.0 and 1.0 */
value *= mulby;
RETVAL += value;
}
QueryPerformanceCounter(&my_end);
printf("new time=%f, opt=" OPTIMIZE "\n",
((double)(my_end.QuadPart-my_beg.QuadPart))/(double)my_freq.QuadPart);

QueryPerformanceCounter(&my_beg);
for(io=0; io < 10000000; io++){
value = 1.0;
value *= Drand01();
RETVAL += value;
}
QueryPerformanceCounter(&my_end);
printf("old time=%f, opt=" OPTIMIZE "\n",
((double)(my_end.QuadPart-my_beg.QuadPart))/(double)my_freq.QuadPart);
OUTPUT:
RETVAL
_____________________________________________________

new time=6.786747, opt= -O1 -G7 -GL -Oi -Og
old time=0.167446, opt= -O1 -G7 -GL -Oi -Og
_____________________________________________________
That is 42 times slower than the previous implementation. I already
moved the bit masking into a conditional, that shaved off 4 seconds.
scalb made things .2 to .5 second worse. A hand written pow loop in a
function took 8 seconds instead of the 6. A header file is included that
RunRand requires. nvoverflowbits.h was created using a script, but
Visual C had a "parser stack overflow" (C1026) when upto 256 bits where
in the conditional tree. I guess I have to somehow cache the Perl_pow

results in global statics and calculate them once per process startup.
The whole loop can be unrolled actually, but it wouldn't be portable if
it is unrolled by hand.

[kmx]

On 29.11.2012 8:13, bulk88 via RT wrote:
> On Wed Nov 28 09:41:59 2012, bulk88 wrote:
>> On Wed Nov 28 05:43:55 2012, Mithaldu wrote:
>>> On Wed Nov 28 00:15:00 2012, bulk88 wrote:
>>>> Why not call the CRT's rand twice or thrice, once for low 15 bits, 2nd
>>>> for high 15 bits, and 3rd for remaining 2?
>>> That sounds like a very reasonable solution. Can you maybe implement a
>>> very small c program to generate a million numbers like that and dump
>>> them to STDOUT to check how many collisions result from that algorithm?
>> I am going to try to fix this by spreading the bits of whatever rand a
>> perl is using across the integer bits of a NV (usually 53 bits for a 64
>> bit double). I plan to test it as XSUBs first and I will post those here.
> Windows only right now.

> ...

bulk88, could you please check in your benchmark how slow will be the
following straightforward implementation?

double drand53()
{
unsigned int a = rand() & 0x1FFF; /* 13 bits */
unsigned int b = rand() & 0x1FFF; /* 13 bits */
unsigned int c = rand() & 0x1FFF; /* 13 bits */
unsigned int d = rand() & 0x3FFF; /* 14 bits */
/* (a * (2^40) + b * (2^27) + c * (2^14) + d) / (2^53) */
return (a*1099511627776.0+b*134217728.0+c*16384.0+d)/9007199254740992.0;
}

--
kmx

[Vincent Pit]

> As was pointed out on IRC, there are actually two orthogonal concerns
> here - how long until the sequence of random numbers repeats, and how
> many distinct numbers are able to be generated. The proposed solution
> here will increase the second while potentially (although not
> necessarily) decreasing the first.

If the period of the RNG is not coprime with the number of calls needed
to make a 32 bits random integer with that algorithm, then not only you
effectively divide the period of the RNG by that number, but the
resulting distribution is also not uniform. And of course you have no
reliable way to know what the period is, even though it's likely to be a
power of 2.

My opinion on this is : just leave it as is. The current implementation
is enough for basic uses, and implementing a better RNG for rand() is
not useful for people that want random numbers with good statistical
properties because that would be yet another black box (like nobody with
serious accuracy concerns will ever use the complex type from C99).


Vincent

[Christian Walde via RT]

On Thu Nov 29 03:59:27 2012, pe...@profvince.com wrote:

> My opinion on this is : just leave it as is. The current
> implementation is enough for basic uses

I strongly disagree with this statement and would like you to explain
exactly what you consider to be "basic uses".

Just to be clear about this: Anytime a script tries to select items
randomly from the entirety of a collection of 32769 items, it will not
work, because only 32768 of those can be accessed with rand(). Worse,
this failure will be silent.

[kmx]

You are probably right if talking about perl's rand function as rand() on
MS Windows is perhaps so poor that creating ugly workarounds is not worth
the effort (especially when implementing or "stealing" some of the well
established PRNG is IMO comparable amount of work).

The other thing (perhaps not related to this RT) is the fact that the same
15bit-low-entropy-short-period-not-so-excelent rand() is on MS Windows used
also in Perl_get_hash_seed(..)

--
kmx

[bulk88 via RT]

On Thu Nov 29 03:04:29 2012, k...@atlas.cz wrote:
>
> bulk88, could you please check in your benchmark how slow will be the
> following straightforward implementation?
>
> double drand53()
> {
> unsigned int a = rand() & 0x1FFF; /* 13 bits */
> unsigned int b = rand() & 0x1FFF; /* 13 bits */
> unsigned int c = rand() & 0x1FFF; /* 13 bits */
> unsigned int d = rand() & 0x3FFF; /* 14 bits */
> /* (a * (2^40) + b * (2^27) + c * (2^14) + d) / (2^53) */
> return
(a*1099511627776.0+b*134217728.0+c*16384.0+d)/9007199254740992.0;
> }
>
> --
> kmx
>

VC 2003 32 bit x86
___________________________________________________________________________

unsigned int a = rand() & 0x1FFF; /* 13 bits */
unsigned int b = rand() & 0x1FFF; /* 13 bits */
unsigned int c = rand() & 0x1FFF; /* 13 bits */
unsigned int d = rand() & 0x3FFF; /* 14 bits */
/* (a * (2^40) + b * (2^27) + c * (2^14) + d) / (2^53) */

value = 1.0;
value *=

(a*1099511627776.0+b*134217728.0+c*16384.0+d)/9007199254740992.0;
RETVAL += value;
}
QueryPerformanceCounter(&my_end);
printf("new kmx time=%f, opt=" OPTIMIZE "\n",

((double)(my_end.QuadPart-my_beg.QuadPart))/(double)my_freq.QuadPart);
QueryPerformanceCounter(&my_beg);
for(io=0; io < 10000000; io++){

/* (a * (2^40) + b * (2^27) + c * (2^14) + d) / (2^53) */

value = 1.0;
value *= ((rand() & 0x1FFF)*1099511627776.0
+(rand() & 0x1FFF)*134217728.0
+(rand() & 0x1FFF)*16384.0
+(rand() & 0x3FFF))/9007199254740992.0;
RETVAL += value;
}
QueryPerformanceCounter(&my_end);
printf("new kmx bulk88 1 time=%f, opt=" OPTIMIZE "\n",

((double)(my_end.QuadPart-my_beg.QuadPart))/(double)my_freq.QuadPart);
QueryPerformanceCounter(&my_beg);
for(io=0; io < 10000000; io++){

unsigned short a = rand() & 0x1FFF; /* 13 bits */
unsigned short b = rand() & 0x1FFF; /* 13 bits */
unsigned short c = rand() & 0x1FFF; /* 13 bits */
unsigned short d = rand() & 0x3FFF; /* 14 bits */

/* (a * (2^40) + b * (2^27) + c * (2^14) + d) / (2^53) */

value = 1.0;
value *=

(a*1099511627776.0+b*134217728.0+c*16384.0+d)/9007199254740992.0;
RETVAL += value;
}
QueryPerformanceCounter(&my_end);
printf("new kmx bulk88 2 time=%f, opt=" OPTIMIZE "\n",
((double)(my_end.QuadPart-my_beg.QuadPart))/(double)my_freq.QuadPart);
OUTPUT:
RETVAL
____________________________________________________________________
new time=6.196008, opt= -O1 -G7 -GL -Oi -Og -arch:SSE2
old time=0.159792, opt= -O1 -G7 -GL -Oi -Og -arch:SSE2
new kmx time=0.668545, opt= -O1 -G7 -GL -Oi -Og -arch:SSE2
new kmx bulk88 1 time=0.638748, opt= -O1 -G7 -GL -Oi -Og -arch:SSE2
new kmx bulk88 2 time=0.846463, opt= -O1 -G7 -GL -Oi -Og -arch:SSE2
____________________________________________________________________
new time=6.336704, opt= -O1 -G7 -GL -Oi -Og
old time=0.159920, opt= -O1 -G7 -GL -Oi -Og
new kmx time=0.665955, opt= -O1 -G7 -GL -Oi -Og
new kmx bulk88 1 time=0.680728, opt= -O1 -G7 -GL -Oi -Og
new kmx bulk88 2 time=0.820330, opt= -O1 -G7 -GL -Oi -Og
____________________________________________________________________

Explicit shorts were slower due to extra movzx instructions. The other 2
kmx/kmx derived implementations I didn't investigate. The kmx
implementation is windows specific. There is one more implementation
that is windows specific that I think will be faster. I will post it
later today or tomorrow.

[Brendan Byrd]

On Tue, Nov 27, 2012 at 10:19 AM, Leon Timmermans via RT
<perlbug-...@perl.org> wrote:
>
> Perl already chooses between a bunch of random number sources. If

> there's something obviously better available we could use that.
>

> That said, I don't really think this is really severity=medium. rand()
> doesn't make any promises on strength, and if you need something
> stronger there are modules on CPAN that can provide just that.

True, but the entropy seems to be abysmally low here, especially since
the -V settings indicate that the float bytes aren't the bottleneck
here. I would at least expect 24 to 32 bits of entropy. Not looking
for crypto security here.

I don't care too much about the severity, but only having 32K possible
random numbers still strikes me as a bug, not a "feature".

--
Brendan Byrd <Pe...@ResonatorSoft.org>
Brendan Byrd <BB...@CPAN.org>

[kmx]

On 30.11.2012 2:44, bulk88 via RT wrote:

> ... I never checked if they actually distribute the bits correctly.
> bulk88 3 I checked by watching hex digits but I didn't analyze it on a
> per bit level, line to line (I'll do that eventually).

Throwing random 8 bytes over raw double representation does not work as you
expect.

I have tested your code and resulting double values are from range:
0.06-0.24 (not even uniformly distributed)

If you want to construct really random double from random bytes without
involving FPU have a look for example at function
int_pair_to_real_inclusive() in
http://svn.ruby-lang.org/cgi-bin/viewvc.cgi/trunk/random.c?view=markup

(I apologise in advance if it is too offensive to link ruby source code in
p5p list)

--
kmx

[Aristotle Pagaltzis]

* kmx <k...@atlas.cz> [2012-11-30 12:25]:

> (I apologise in advance if it is too offensive to link ruby source
> code in p5p list)

Why?

How else do you steal?

[Nicholas Clark]

On Fri, Nov 30, 2012 at 10:17:24AM -0600, Craig A. Berry wrote:

> On Tue, Nov 27, 2012 at 1:52 PM, demerphq <deme...@gmail.com> wrote:
> > On 27 November 2012 20:43, kmx <k...@atlas.cz> wrote:
> >>
> >> On 27.11.2012 16:40, Christian Walde via RT wrote:
> >>>
> >>> On Tue Nov 27 07:19:53 2012, LeonT wrote:
>
> >>>
> >>> On Windows XP and higher, there is RtlGenRandom, would that fit the bill?
>

> And the CRT provides this as rand_s(). See

Which seems great. But

> <http://msdn.microsoft.com/en-us/library/sxtz2fa8%28v=vs.110%29.aspx>.


The rand_s function writes a pseudorandom integer in the range 0 to
UINT_MAX to the input pointer. The rand_s function uses the operating
system to generate cryptographically secure random numbers. It does not
use the seed generated by the srand function, nor does it affect the
random number sequence used by rand.


which means that this would bust the documented behaviour for srand:

Sets and returns the random number seed for the C<rand> operator.



So, what are we trying to achieve? And what guarantees do we want to keep?

> > I would think that ideally the Win32 version of Perl should "out of
> > the box" have a reasonably fast RNG.
>

> I agree, but more randomness is likely to cost *something* in run time
> and/or maintenance time.


I agree too.

I think that the answer to the first question is fairly unequivocably
"better random numbers from rand() on Win32", and, *specifically*, more than
15 bits.


I personally think that it's very useful to maintain the relationship
between srand() and rand() - ie that rand() is a *pseudo*-random number
generator, and that it's possible to fix and repeat the sequence by using
srand(). (I've used this in the past)


If that is true (and we don't want to re-write perlfunc.pod to remove
long standing functionality from srand()) then sadly that Win32 CRT function
is *not* useful.



I'm also very very wary of trying to roll *anything* ourselves. There will
be plenty of ways to really screw this up down the line. So, I think the
way to go is either

1) for Win32, adopt an existing BSD-licensed PRG which provides 32, if not
48 bits of result (and hopefully a lot more entropy internally)
What is in OpenBSD's libc?

2) for everything, switch to Mersenne Twister. If it's the default for PHP,
Python and Ruby, surely it's good enough for Perl?

Nicholas Clark

[Christian Walde via RT]

On Fri Nov 30 08:18:19 2012, craig....@gmail.com wrote:

> Yes, rand_s() takes about 6 times as long as rand(), using this as a
> benchmark:

I don't think that's the correct speed difference to test. Instead rand
should be compared between windows and linux on similar hardware. If a
better perl rand() on windows is 6 times slower than the current rand(),
but the current rand() on windows is 12 times faster than the linux rand
(), then the difference does not really matter. :)

[Craig A. Berry]

On Fri, Nov 30, 2012 at 10:34 AM, Nicholas Clark <ni...@ccl4.org> wrote:

> I personally think that it's very useful to maintain the relationship
> between srand() and rand() - ie that rand() is a *pseudo*-random number
> generator, and that it's possible to fix and repeat the sequence by using
> srand(). (I've used this in the past)

Microsoft has other APIs where the input value is used as part of the seed:

<http://msdn.microsoft.com/en-us/library/windows/desktop/aa379942(v=vs.85).aspx>

It's possible rand_s() does that too, though it's not documented.

> I'm also very very wary of trying to roll *anything* ourselves. There will
> be plenty of ways to really screw this up down the line.

Yes. I think this thread has been highly information on how not to
implement an RNG :-).

>So, I think the
> way to go is either
>
> 1) for Win32, adopt an existing BSD-licensed PRG which provides 32, if not
> 48 bits of result (and hopefully a lot more entropy internally)
> What is in OpenBSD's libc?

Earlier in the thread, I said:

The FreeBSD implementation of drand48() is indeed short and sweet even
if you follow the links a couple of levels deep:

<http://fxr.watson.org/fxr/source/gen/drand48.c?v=FREEBSD-LIBC;im=bigexcerpts>

It might just work verbatim on Windows.

> 2) for everything, switch to Mersenne Twister. If it's the default for PHP,
> Python and Ruby, surely it's good enough for Perl?

Perhaps TinyMT (<http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/TINYMT/index.html>)
would be adequate for a basic (and portable) random number generator.

[bulk88 via RT]

On Fri Nov 30 03:22:59 2012, k...@atlas.cz wrote:
>
> Throwing random 8 bytes over raw double representation does not work
as you
> expect.
>
> I have tested your code and resulting double values are from range:
> 0.06-0.24 (not even uniformly distributed)
>
> If you want to construct really random double from random bytes without
> involving FPU have a look for example at function
> int_pair_to_real_inclusive() in
> http://svn.ruby-lang.org/cgi-bin/viewvc.cgi/trunk/random.c?view=markup
>
> (I apologise in advance if it is too offensive to link ruby source
code in
> p5p list)
>
> --
> kmx
>

I had a feeling the first 2 digits were not as random as later digits
prior to you identifying. You proved me correct. The algorithm was
supposed to generate a number between .9999 and ~.0125 but not smaller
than .0125000000, but a bug probably narrowed that range. There was a
design issue whether to make the exponent random, but then
int(rand(100)) nearly always returns 0 rather than an equal distribution
between 0 to 100. I created a new version of it which uses FP
subtraction to cut off 1st unrandom digit, and writing up some digit
distribution tests, the old 15 bit rand has a very visible difference
with my tests. I will post it after some more testing and benchmarking
on my end.

[Brendan Byrd via RT]

On Wed Nov 28 18:01:44 2012, aristotle wrote:
>
> Don’t.
>
> Randomness is a treacherous sea – easy to get it wrong, hard to notice
> when you do, and with far-reaching consequences. (Debian patched some
> random number generator and caused half the world to have to mint new
> SSH keys a few years down the line.)

Math hard, so don't bother? Sorry, but I strongly disagree. A word of
caution is perfectly valid, as RNGs do need some careful thought for
implementation. But, it can be done.

Also, it's already broken right now.

> That doesn’t increase the entropy available, it just draws from the pool
> faster. I don’t have my randomness down well enough to reason it through
> with the necessary rigour but I’m pretty sure that at best that approach
> does nothing to increase the quality of the random numbers and likely it
> will decrease it.

Eh? What does it matter if it draws from the pool faster? We're not talking
about seed prediction (or crypto-level RNGs).

And we're talking about bitshifting on a binary boundary. (In other words, the
range is from 0 to 15-bits full.) As long as that is true, appending the RNG
results is perfectly valid, giving you a completely random number that uses all
bits in all places.

--

I am glad that this bug is getting quite a bit of attention. Keep up the hard
work!

[Aristotle Pagaltzis]

* Brendan Byrd via RT <perlbug-...@perl.org> [2012-12-01 05:50]:

> On Wed Nov 28 18:01:44 2012, aristotle wrote:
> > Don’t.
> >
> > Randomness is a treacherous sea – easy to get it wrong, hard to
> > notice when you do, and with far-reaching consequences. (Debian
> > patched some random number generator and caused half the world to
> > have to mint new SSH keys a few years down the line.)
>
> Math hard, so don't bother? Sorry, but I strongly disagree. A word of
> caution is perfectly valid, as RNGs do need some careful thought for
> implementation. But, it can be done.

Where did you see me saying that? What I said is, this math harder than
you think, so don’t just follow the first intuition that pops into your
head: if you want to do something about it, approach the topic with the
necessary humility. Of course it can be done, just not in 15 minutes.

(And you can always copy the work of someone who has already sweat the
details. I won’t raise any concerns over that. I would like it in fact.)

If you just meant the “Don’t”, it was only about one specific approach,
per preceding context.

> Also, it's already broken right now.

Yes, so don’t smear lipstick on it – fix it or leave it.

> > That doesn’t increase the entropy available, it just draws from the
> > pool faster. I don’t have my randomness down well enough to reason
> > it through with the necessary rigour but I’m pretty sure that at
> > best that approach does nothing to increase the quality of the
> > random numbers and likely it will decrease it.
>
> Eh? What does it matter if it draws from the pool faster? We're not
> talking about seed prediction (or crypto-level RNGs).
>
> And we're talking about bitshifting on a binary boundary. (In other
> words, the range is from 0 to 15-bits full.) As long as that is true,
> appending the RNG results is perfectly valid, giving you a completely
> random number that uses all bits in all places.

If the RNG cycle length shares prime factors with the number of times
you call the RNG to construct a single random number, then all you did
is shorten the effective cycle by those factors. And by definition the
length of the cycle is the upper bound on how many distinct numbers you
can generate. So the result would be a much wider spread of many fewer
distinct numbers.

If you’re throwing away any bits you got from the RNG, you can also end
up exposing a non-uniformity in the distribution that was not previously
apparent. (I think there is also some risk of this happening in general
but I’d have to bone up further on my math to substantiate this.)

[Brendan Byrd via RT]

On Sat, Dec 1, 2012 at 2:46 AM, A. Pagaltzis via RT <perlbug-...@perl.org> wrote:
>
> (And you can always copy the work of someone who has already sweat the
> details. I won’t raise any concerns over that. I would like it in fact.)

I'm open to that as long as there isn't a noticeable slowdown to the original
function, and the RNG code isn't specific to a type of need.

> Yes, so don’t smear lipstick on it – fix it or leave it.

Agreed.

> If the RNG cycle length shares prime factors with the number of times
> you call the RNG to construct a single random number, then all you did
> is shorten the effective cycle by those factors. And by definition the
> length of the cycle is the upper bound on how many distinct numbers you
> can generate. So the result would be a much wider spread of many fewer
> distinct numbers.

I'll admit that I only know enough about RNGs to be dangerous, but I'm still
not seeing the problem here.

Let's say that I have a RNG that picks a single digit between 0 and 9. It has
an equal chance (1:10) to pick each digit. I call that 50 times, append it to
the growing number, and then turn it into a float by putting a dot at the
beginning.

Bulk88's approach is essentially the same thing, except it's doing it in binary
with binary boundaries.

> If you’re throwing away any bits you got from the RNG, you can also end
> up exposing a non-uniformity in the distribution that was not previously
> apparent. (I think there is also some risk of this happening in general
> but I’d have to bone up further on my math to substantiate this.)

Yes, this is valid concern, which is why it's important to do this on a binary
boundary with the right information on the bits of entropy. If the number of
values is, say, 0x00 to 0xFF, and we try this on a 16-bit boundary, then
we would be adding extra binary zeros, polluting the RNG.

--
Brendan Byrd <Pe...@ResonatorSoft.org>
Brendan Byrd <BB...@CPAN.org>

[Jesse Luehrs]

How the values are distributed is another concern. A common issue with
bad RNGs is that they have much less entropy in the low bits than in the
high bits. This isn't typically noticeable when you use it in the way
perl does, since the low bits are only seen if you care about a lot of
precision in the random floating point numbers (the first couple decimal
places would still be reasonably distributed), but if you start
appending values in this way, it can move the problem up to a point
where it could be noticeable (if, say, you took the values from the
lowest two bits + a full 15 bits + a full 15 bits to get the full random
number). These are the kinds of issues people are talking about. Using
an existing, already written and tested implementation instead would
sidestep this whole issue.

-doy

[Peter Martini]

On Dec 1, 2012 8:31 AM, "Brendan Byrd via RT" <perlbug-...@perl.org> wrote:
>
> On Sat, Dec 1, 2012 at 2:46 AM, A. Pagaltzis via RT <perlbug-...@perl.org> wrote:
> > If the RNG cycle length shares prime factors with the number of times
> > you call the RNG to construct a single random number, then all you did
> > is shorten the effective cycle by those factors. And by definition the
> > length of the cycle is the upper bound on how many distinct numbers you
> > can generate. So the result would be a much wider spread of many fewer
> > distinct numbers.
>
> I'll admit that I only know enough about RNGs to be dangerous, but I'm still
> not seeing the problem here.
>
> Let's say that I have a RNG that picks a single digit between 0 and 9.  It has
> an equal chance (1:10) to pick each digit.  I call that 50 times, append it to
> the growing number, and then turn it into a float by putting a dot at the
> beginning.
>
> Bulk88's approach is essentially the same thing, except it's doing it in binary
> with binary boundaries.
>

The problem is that the kind of RNG perl can use isn't truly random, it just generates a sequence which *looks* random.  After a certain amount of numbers returned - the cycle length - the same sequence will start all over.

Let's say our RNG generates the sequence 0, 1, 2, 3, 4, 5, 6, 7, 8, 9.  The cycle length is 10.  If we then make a 4 digit number just by combining digits, we'll get 1234, 5678, 9012, 3456, 7890 and then repeat.  We've generated far larger numbers, but only half as many.

[demerphq]

FWIW as far as I know this is more or less defeated by letting the RNG
choose how many RNG "throws" to throw away each iteration.

Yves

--
perl -Mre=debug -e "/just|another|perl|hacker/"

[Aristotle Pagaltzis]

* Brendan Byrd via RT <perlbug-...@perl.org> [2012-12-01 14:35]:

> On Sat, Dec 1, 2012 at 2:46 AM, A. Pagaltzis via RT <perlbug-...@perl.org> wrote:
> > If the RNG cycle length shares prime factors with the number of
> > times you call the RNG to construct a single random number, then all
> > you did is shorten the effective cycle by those factors. And by
> > definition the length of the cycle is the upper bound on how many
> > distinct numbers you can generate. So the result would be a much
> > wider spread of many fewer distinct numbers.
>
> I'll admit that I only know enough about RNGs to be dangerous, but I'm
> still not seeing the problem here.

OK, I’ll try to illustrate.

First, we need to figure out cycle length here.

Well, necessarily for a PRNG you are using a deterministic equation, and
you are iterating on a seed with finite range, yes? Then inevitably the
seed must eventually reach a value it has had before, and since the PRNG
equation is deterministic, inevitably its output must from then on be
identical to the sequence it produced after the last time the seed had
that value.

So every PRNG must have a number of calls after which it its sequence of
outputs starts repeating.

Then it is also obvious that the number of different numbers the PRNG
can produce is at most as large as its cycle length, yes? After all, if
you repeat yourself after 17 outputs, then you can never produce more
than 17 different outputs, otherwise you wouldn’t be repeating yourself.

Now:

> Let's say that I have a RNG that picks a single digit between 0 and 9.
> It has an equal chance (1:10) to pick each digit. I call that 50
> times, append it to the growing number, and then turn it into a float
> by putting a dot at the beginning.
>
> Bulk88's approach is essentially the same thing, except it's doing it
> in binary with binary boundaries.

OK, let’s get to the illustration. To keep the numbers manageable let’s
say we have a PRNG that can generate a digit between 0 and 2. And let’s
say that because of its equation and size of seed, it will have a cycle
length of 18. Let’s say it produces this sequence, endlessly:

0 1 2 1 0 2 2 1 0 0 0 1 2 1 2 2 0 1

Now you say you want more range in its output. So you’ll just take two
random numbers at a time and you will return r1+3*r2 which will give you
a range of 0 to 8. Now the output looks like this:

3 5 6 5 0 3 5 8 3 …

At the “…” you have exactly exhausted your original PRNG’s cycle. So now
you must repeat yourself because your inputs will repeat themselves. The
cycle length is now 9, half as long as the original – because you were
drawing numbers from it twice as fast and its cycle length had one prime
factor of 2. Notice that even though the range of outputs is 0 to 8, the
actual number of different outputs is only 5, not 9.

It just so happens (badly chosen off the cuff example, sigh) that in
this case, the range of outputs and the cycle length are identical in
size so theoretically this PRNG could output its full range. But you can
see that if you tried to do the same on top of *this* PRNG in a way that
shortens the cycle, you would get much more range than 0-8, but actually
far fewer distinct outputs.

(Even if you didn’t exactly exhaust your PRNG, in fact even if the cycle
length and the number of outputs you draw from it at once are coprime,
you cannot truly increase the cycle length, you only hide it behind
a pattern.)

[demerphq]

I seem to recall the cycle length gets longer with an LCG if you do
something like this:

sub rand {
my $iters= int 1 + rand(100);
$rand= rand($_[0]) while $iters--;
return $rand;

[kmx]

Anyway the code trying to make things better with MS Windows' rand() is
getting quite complicated.

The unoptimized drand48 implementation looks like this:

static unsigned long long state;

void seed(unsigned int s)
{
state = ((unsigned long long)s << 16) | 0x330E;
}

double drand48()
{
state = (25214903917 * state + 11) % 281474976710656;
return (double)state / 281474976710656.0; /* state / 2^48 */
}

What about to simply implement drand48 - just for Win32 (as proposed by Craig)?

Could you please include drand48 implementation above to your benchmark and
show how slow it is?

--
kmx

[demerphq]

You need to support the reentrant version as well afaik.

[Aristotle Pagaltzis]

* demerphq <deme...@gmail.com> [2012-12-01 19:50]:

> I seem to recall the cycle length gets longer with an LCG if you do
> something like this:
>
> sub rand {
> my $iters= int 1 + rand(100);
> $rand= rand($_[0]) while $iters--;
> return $rand;
> }

That’ll work.

[kmx]

It'll be 50x slower

[Craig A. Berry]

On Sat, Dec 1, 2012 at 3:38 PM, demerphq <deme...@gmail.com> wrote:


>> Anyway the code trying to make things better with MS Windows' rand() is
>> getting quite complicated.
>>
>> The unoptimized drand48 implementation looks like this:
>>
>> static unsigned long long state;
>>
>> void seed(unsigned int s)
>> {
>> state = ((unsigned long long)s << 16) | 0x330E;
>> }
>>
>> double drand48()
>> {
>> state = (25214903917 * state + 11) % 281474976710656;
>> return (double)state / 281474976710656.0; /* state / 2^48 */
>> }
>>
>> What about to simply implement drand48 - just for Win32 (as proposed by
>> Craig)?
>>
>> Could you please include drand48 implementation above to your benchmark and
>> show how slow it is?

If "you" means me, I'll try to have a look but not sure when.

> You need to support the reentrant version as well afaik.

What does thread-safe mean in this context? Does each thread need to
have its own state variable such that threads cannot plant seeds in
each other? Or is it sufficient to have the drand48 call atomic,
i.e., just make sure the state doesn't change again while your'e in
the middle of computing something with it.

I suspect it's the former, which I suspect means a new interpreter variable.

[kmx]

On 1.12.2012 23:12, Craig A. Berry wrote:
> On Sat, Dec 1, 2012 at 3:38 PM, demerphq <deme...@gmail.com> wrote:
>
>
>>> Anyway the code trying to make things better with MS Windows' rand() is
>>> getting quite complicated.
>>>
>>> The unoptimized drand48 implementation looks like this:
>>>
>>> static unsigned long long state;
>>>
>>> void seed(unsigned int s)
>>> {
>>> state = ((unsigned long long)s << 16) | 0x330E;
>>> }
>>>
>>> double drand48()
>>> {
>>> state = (25214903917 * state + 11) % 281474976710656;
>>> return (double)state / 281474976710656.0; /* state / 2^48 */
>>> }
>>>
>>> What about to simply implement drand48 - just for Win32 (as proposed by
>>> Craig)?
>>>
>>> Could you please include drand48 implementation above to your benchmark and
>>> show how slow it is?
> If "you" means me, I'll try to have a look but not sure when.

I am sorry, I meant bulk88 and his benchmark with various rand()
workarounds for Win32

--
kmx

[Leon Timmermans]

On Sat, Dec 1, 2012 at 11:12 PM, Craig A. Berry <craig....@gmail.com> wrote:
> What does thread-safe mean in this context? Does each thread need to
> have its own state variable such that threads cannot plant seeds in
> each other? Or is it sufficient to have the drand48 call atomic,
> i.e., just make sure the state doesn't change again while your'e in
> the middle of computing something with it.
>
> I suspect it's the former, which I suspect means a new interpreter variable.

Yes, the former. The latter would be useless as sometimes one is
relaying on the deterministic nature of a PRNG. We already have
PL_reentrant_buffer for this kind of purpose, though I suspect we're
currently not using it on Windows.

Leon

[Paul Johnson]

Random number generation is hard. Let's build rockets.

If I understand your suggestion correctly, using it on this example
gives:

Step 1. Throw away 0 iterations:
Choose 1234
Step 2. Throw away 5 iterations: 6,7,8,9,0
Choose 1234
Step 3. Throw away 5 iterations: 6,7,8,9,0
Choose 1234

Oops. That's about the same quality as http://xkcd.com/221/ but
probably not as fast.

I think I'd recommend against trying to get better random numbers based
on an existing random number generator unless someone absolutely,
completely and fully understands the implications of so doing. Nicking
something from someone who has that understanding is probably the better
option.

I pretty much agree with Nicholas' conclusions and think we should use
Mersenne Twister by default on every platform.

--
Paul Johnson - pa...@pjcj.net
http://www.pjcj.net

[demerphq]

Afaik it already exists.

The api is to add an extra param which reference the state buffer.

see dword48_r() in linux docs for one implementations api.

[demerphq]

Well, on one hand point taken. On other hand your example doesn't
really match what I said, nor what I get if I actually implement it
even using such a bad source sequence.

> I think I'd recommend against trying to get better random numbers based
> on an existing random number generator unless someone absolutely,
> completely and fully understands the implications of so doing. Nicking
> something from someone who has that understanding is probably the better
> option.

Sure, this is definitely true. And from that point of view I shouldn't
have said anything. :-)

> I pretty much agree with Nicholas' conclusions and think we should use
> Mersenne Twister by default on every platform.

Yes this would be nice.

cheers,

[bulk88 via RT]

On Fri Nov 30 08:18:19 2012, craig....@gmail.com wrote:
>

> And the CRT provides this as rand_s(). See

> <http://msdn.microsoft.com/en-us/library/sxtz2fa8%28v=vs.110%29.aspx>.
> I've attached a proof-of-concept patch using rand_s() that out of
> laziness just pastes:
>
> static double win32_drand01(void) {
> unsigned int val;
> (void)rand_s(&val);
> return (double)val;
> }
>
> into perl.h. Obviously that would have to be moved to win32/win32.c
> and made visible, but this was enough to do a quick build and see what
> it gets us. The obvious advantage is that it's almost no code to
> maintain and immediately provides 32 bits of entropy instead of 15.
> The obvious disadvantage is ...
>
> > Just curious but wouldnt something like that be pretty slow?

>
> Yes, rand_s() takes about 6 times as long as rand(), using this as a
> benchmark:
>

> C:\perlgit>type randbench.pl
> use Benchmark;
>
> my $cnt = 0;
> my $t0 = new Benchmark;
> my $num = 0;
> foreach my $i (1 .. 10_000_000) {
> $num = rand;
> $cnt++;
> }
> my $t1 = new Benchmark;
> my $td = timediff($t1, $t0);
> print "$cnt iterations in ", timestr($td), "\n";
> print "The last random # generated was $num\n";

>
> > I would think that ideally the Win32 version of Perl should "out of
> > the box" have a reasonably fast RNG.
>
> I agree, but more randomness is likely to cost *something* in run time
> and/or maintenance time.

rand_s is VC 2005 (msvcr80.dll) or newer. Not VC 2003 or VC6 or
msvcrt.dll and no Windows 2000 (all of which perl supports). The MS
implementation in VC 2008 similar to
http://code.google.com/searchframe#S3vzerue4i0/trunk/reactos/lib/sdk/crt/math/rand.c&type=cs&l=33
but with many DecodePointers/EncodePointers and saving the func pointer
and not doing a loadlibrary each time. I will try to add rand_s to my
benchmark for the record if I can get the manifested newer CRT to load
into the process. Win32 C's rand_s also breaks perl's srand as mentioned
in this ticket.

--
bulk88 ~ bulk88 at hotmail.com

[Brendan Byrd via RT]

I think we should get a few randomization tests into the main Perl test
suite, similar to mine and bulk88's. Tests would include:

1. The number of numerical duplicates within a large batch (2**32)
2. The number of gaps in that batch, within numbers (2**16), if they
were all promoted to an integer.
3. Digit variation (should be uniform among all of the digits)
4. Random function is reasonably fast.

Tests 1&2 should define "large batch" as 2**32. I think bulk88's above
post covers what I was saying, but it should be an official test, ran
every time Perl is compiled. (Else we risk repeating the same Debian
OpenSSH mistake. They obviously didn't have a good randomization test
platform.) More test cases are welcome.

I am also starting to agree with rjbs' opinion of just putting MT as the
core rand() function. If one platform's random function isn't up to par
(like Windows), then we can't trust any of them. What's the difficulty
of taking the basic barebones functions within Math::Random::MT::Auto
and making it XS?

[kmx]

The original MT implementation is available at
http://www.math.sci.hiroshima-u.ac.jp/~m-MAT/MT/MT2002/emt19937ar.html

However I am not sure whether its BSD license is ok when talking about
including into perl core.

--
kmx

[Jesse Luehrs]

The BSD license (specifically, the 3-clause license used in that link)
is fine.

-doy

[Dr.Ruud]

On 2012-12-03 18:40, Brendan Byrd via RT wrote:

> 3. Digit variation (should be uniform among all of the digits)

Hmm, uniform, how would you define that,
"the same digit distribution as the same amount of consecutive numbers"?

--
Ruud

[Dana Jacobsen via RT]

Here is an issue I just had with Math-BigInt-Random-OO, which was
returning only even numbers as length(15) or length(80) random numbers.
It came down to an issue mentioned in the documentation -- drand48 on
many systems has an alternating last bit. This ends up interacting
quite badly with that module.

Given a system with randbits = 48:

perl -E 'my $rm = 281474976710656; $n=50000; my %freq; for (1..$n) { my
$num = int CORE::rand 1; $freq{ int(CORE::rand $rm) % 10 }++ }
printf("%4d %6.3f%%\n", $_, 100.0*$freq{$_}/$n) for sort {$a<=>$b} keys
%freq;'
0 19.992%
2 20.078%
4 20.132%
6 19.944%
8 19.854%

We only get 0 as the last bit from the second call. I've reproduced
this on Fedora x86_64 64-bit, Debian i386 32-bit, Debian ia64 (Itanium
2) 64-bit, Debian sparc (US T1) 64-bit, Fedora ppc64 (Power7) 64-bit,
NetBSD amd64 64-bit, Cygwin i686 64-bit. That's a pretty wide variety
of systems. I did write a little C program using drand48 that showed it
happening there, so this is an issue with the system RNG, not Perl code.


I think having Perl move to an internal standard RNG like MT or ISAAC
would be great. The O/S specific portion would then be figuring out how
to best get a default seed (and that's where we can use good random
sources like /dev/random or the Windows crypto stuff). It would still
allow seeds to be specified, though the state for these RNGs is no
longer a single integer.

Of the three options:

1) offer no default solution. Pick something explicitly!
2) offer a cruddy solution that may or may not work well
3) offer a decent default solution

I don't think (2) is a good answer, and that's what this ticket seems to
be about. (1) seems unlikely. Hence 3.

I'd also like to see a standard function like irand(n) that would return
a uniformly distributed UV between 0 and n inclusive, allowing n as
large as ~0, since that seems like a common case. But that's a separate
issue.

Here's an interesting report on some RNGs:

http://citeseerx.ist.psu.edu/viewdoc/summary?doi=10.1.1.141.2849

[demerphq]

So what should we do?

Is there a library that is portable that can solve this problem?

If there is point us at it and I will look into merging it one of these days.

Do we have license to use the mt code? Can we integrate the clan module into core?

Enough talk. Lets get some patches going. :-)

Yves

From an ipad, excuse the bad formatting and brevity and typos... Sigh

[bulk88 via RT]

On Wed Dec 12 13:58:25 2012, demerphq wrote:
> So what should we do?
>
> Is there a library that is portable that can solve this problem?
>
> If there is point us at it and I will look into merging it one of these
> days.
>
> Do we have license to use the mt code? Can we integrate the clan module
> into core?
>
> Enough talk. Lets get some patches going. :-)
>
> Yves

I am leaning towards drand48 on all platforms. Its fast and stuffs
doubles pretty good. A note should be placed in perlfunc not to use
beyond the 17th or 18th decimal digit (look at my earlier posts) of the
result of rand(). I assume a NV is atleast 8 bytes on all platforms
(correct me if I am wrong). A PRNG should be fast for quickly filling
out data with psuedo randomess, not crypto secure. A crypto secure RNG
is beyond the scope of the interp and is best left to platform specific
modules that use certified/audited entropy sources and algoritms. And
the srand replay must remain working. The quality of the seed is not the
issue in this ticket. The conversion of the <1 FP into an integer
through multiplication causing bad distribution is what this ticket is
about. IMO there are 3 ways to fix the problem, drand48 on Windows,
drand48 on all platforms, or call rand 4 times on Windows. I (and kmx I
think) already have a bunch of candidate code in this ticket.

--
bulk88 ~ bulk88 at hotmail.com

[bulk88 via RT]

On Wed Dec 12 01:55:08 2012, danaj wrote:
> Here is an issue I just had with Math-BigInt-Random-OO, which was
> returning only even numbers as length(15) or length(80) random numbers.
> It came down to an issue mentioned in the documentation -- drand48 on
> many systems has an alternating last bit. This ends up interacting
> quite badly with that module.
>
> Given a system with randbits = 48:

Random guess, that is the implied 1 bit in the FP number or a double can
store only 2^53 (or is that 2^52) before starting to round the low
digits, meanwhile a int 64 is 2^64.

>
> I'd also like to see a standard function like irand(n) that would return
> a uniformly distributed UV between 0 and n inclusive, allowing n as
> large as ~0, since that seems like a common case. But that's a separate
> issue.

I wrote something similar in
https://rt.perl.org/rt3/Ticket/Display.html?id=115928#txn-1175200 where
randomness was added to a NV in a portable way without looking inside
the NV, it was very slow due to non CPU native FP math. N/I/UVs are
different sizes on different platforms, so its not that simple.

--
bulk88 ~ bulk88 at hotmail.com

[Nicholas Clark]

OpenBSD source code is in src/lib/libc/stdlib:

http://www.openbsd.org/cgi-bin/cvsweb/src/lib/libc/stdlib/

drand48.c, erand48.c, _rand48.c (at least)

The BSD licence is "acceptable" :-)


Whilst that's a heck of a lot better than what we have generally, my gut
feeling is still that I'd prefer the Mersenne Twister, if it's not to hard
to get right. As I commented before, PHP Python and Ruby are all using it,
and "no news is good news" suggests that it doesn't cause any problems.

There is C source code linked from the Mersenne Twister Home page:

http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/emt.html

I've not looked at their code. Their licence is "acceptable" :-)

Nicholas Clark

[demerphq]

The one concern I have about this is what do we do with srand()?

I would guess there is code out there that expects to be able to use
srand() to get a predictable sequence, if we switch to MT we would
need to do /something/ about srand(), unfortunately given how MT is
seeded its not clear what.

cheers
Yves

[Craig A. Berry]

On Fri, Dec 14, 2012 at 8:55 AM, Nicholas Clark <ni...@ccl4.org> wrote:

> There is C source code linked from the Mersenne Twister Home page:
>
> http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/emt.html

Specifically the source for the original version is at:

<http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/MT2002/CODES/mt19937ar.c>

and looks pretty compact and portable. Probably best to start with
the readme at:

<http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/MT2002/CODES/readme-mt.txt>

> I've not looked at their code. Their licence is "acceptable" :-)

There are newer, purportedly faster, versions by the original authors.
They all look significantly more complicated and less portable than
the original, "classic" version I've pointed to above.

MT Classic stores state in an array of 624, 32-bit integers. I guess
this would have to be allocated per interpreter, probably the first
time the seed function is called. And then subsequent seed calls
would reinitialize (but not reallocate) it. So there would have to be
a new global PL_mt_state or similar.

Or should a single thread be allowed to have multiple RNGs going at
once? If so, the seed function would need to allocate a new state
vector every time and return it to the caller, who would then have to
pass it to the rand function. Which I think would be tricky to make
work for the current Perl operators without leaking state vectors.

Another design question is whether to boil MT down to a single seed
function and a single rand function in support of the existing Perl
operators, or to include the full MT API, which has the ability to
supply larger seeds up to a full state vector and return a random
number six different ways with various data types and randomness
characteristics. Of course there's nothing wrong with doing the
simplest possible thing first and (possibly) extending it later.

[kmx]

srand() is absolutely no problem with MT - the function you are looking for
is void init_genrand(unsigned long s) in MT's reference implementation.

In fact we have just opposite problem: How to initially seed PRNG when perl
starts, in other words where to collect some reasonable entropy we need for
initial MT seeding. It will be definitely some kind of platform dependent
kung-fu, e.g. use /dev/random if available, use CryptGenRandom on Win32,
call system's drandXY() if available, collect "something" in other cases.

My opinion (if it matters): let us implement something serious like MT or
perhaps consider ISAAC (proposed by Dana)

1/ MT: we probably want "Mersenne Twister with improved initialization"
website: http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/MT2002/emt19937ar.html
reference implementation:
http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/MT2002/CODES/mt19937ar.c
license of reference implementation: BSD 3-clause

2/ ISAAC: should be even a bit faster than MT with excellent
characteristics (some even better than MT)
website: http://burtleburtle.net/bob/rand/isaacafa.html
reference implementation (unoptimized version):
http://burtleburtle.net/bob/c/readable.c
license of reference implementation: public domain

The actual PRNG algorithm is just one part to be solved. The other part (in
my opinion much more complicated) is about how to glue it together with
perl core - how to initialize PRNG at startup, where to store state, how to
play nice with ithreads etc.

--
kmx

[demerphq]

Oh cool.

> In fact we have just opposite problem: How to initially seed PRNG when perl
> starts, in other words where to collect some reasonable entropy we need for
> initial MT seeding. It will be definitely some kind of platform dependent
> kung-fu, e.g. use /dev/random if available, use CryptGenRandom on Win32,
> call system's drandXY() if available, collect "something" in other cases.

We have most of this wrapped up in the routine we use for our current srand().

> My opinion (if it matters):

It seems to be useful. :-)

> let us implement something serious like MT or
> perhaps consider ISAAC (proposed by Dana)
>
> 1/ MT: we probably want "Mersenne Twister with improved initialization"
> website:
> http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/MT2002/emt19937ar.html
> reference implementation:
> http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/MT2002/CODES/mt19937ar.c
> license of reference implementation: BSD 3-clause
>
> 2/ ISAAC: should be even a bit faster than MT with excellent characteristics
> (some even better than MT)
> website: http://burtleburtle.net/bob/rand/isaacafa.html
> reference implementation (unoptimized version):
> http://burtleburtle.net/bob/c/readable.c
> license of reference implementation: public domain
>
> The actual PRNG algorithm is just one part to be solved. The other part (in
> my opinion much more complicated) is about how to glue it together with perl
> core - how to initialize PRNG at startup, where to store state, how to play
> nice with ithreads etc.

I have a bit of experience with the perl implementation issues, but
not with MT or ISAAC.

I am willing to look into implementing it if someone like you is
supporting me on it.

[Nicholas Clark]

On Fri, Dec 14, 2012 at 10:06:30AM -0600, Craig A. Berry wrote:

> Specifically the source for the original version is at:
>
> <http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/MT2002/CODES/mt19937ar.c>
>
> and looks pretty compact and portable. Probably best to start with
> the readme at:
>
> <http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/MT2002/CODES/readme-mt.txt>
>
> > I've not looked at their code. Their licence is "acceptable" :-)
>
> There are newer, purportedly faster, versions by the original authors.
> They all look significantly more complicated and less portable than
> the original, "classic" version I've pointed to above.

I think one codebase, portable, is the way to go.

> MT Classic stores state in an array of 624, 32-bit integers. I guess
> this would have to be allocated per interpreter, probably the first
> time the seed function is called. And then subsequent seed calls
> would reinitialize (but not reallocate) it. So there would have to be
> a new global PL_mt_state or similar.

Gulp. Something over 2k.

> Or should a single thread be allowed to have multiple RNGs going at
> once? If so, the seed function would need to allocate a new state
> vector every time and return it to the caller, who would then have to
> pass it to the rand function. Which I think would be tricky to make
> work for the current Perl operators without leaking state vectors.
>
> Another design question is whether to boil MT down to a single seed
> function and a single rand function in support of the existing Perl
> operators, or to include the full MT API, which has the ability to
> supply larger seeds up to a full state vector and return a random
> number six different ways with various data types and randomness
> characteristics. Of course there's nothing wrong with doing the
> simplest possible thing first and (possibly) extending it later.

The original bug report is that the win32 C library is minimally ANSI
conformant and only generates 15 bits of entropy. We're not trying to
solve anything more than that. So I think go for the simplest thing that
could possibly solve that (and not need re-working soon after)

In that code you link to, this looks pretty much like it provide a backend
for srand():

/* initializes mt[N] with a seed */
void init_genrand(unsigned long s)


(but change the default self-seed to call init_by_array(), using
/dev/urandom to fill the buffer, or whatever the current probing
code can pick up as a source of entropy)


Let CPAN authors worry about what nice Perl-esque interface to provide for
the Mersenne Twister.

Nicholas Clark

[Nicholas Clark]

On Fri, Dec 14, 2012 at 05:21:35PM +0100, kmx wrote:

> 2/ ISAAC: should be even a bit faster than MT with excellent
> characteristics (some even better than MT)
> website: http://burtleburtle.net/bob/rand/isaacafa.html
> reference implementation (unoptimized version):
> http://burtleburtle.net/bob/c/readable.c
> license of reference implementation: public domain

Also only 1K of state?

Right. I'm really claiming not to be an expert on this.


state (bytes) speed cycle length
Mersenne Twister: 2496 slow? 2**19937-1
ISAAC: 1024 faster at least 2**40, avg. 2**8295


Given we've been happy with drand48(), we don't "need" the massive cycle
length of the Mersenne Twister. The other characteristics of ISAAC seem to
be more suitable.

So, assuming it's not weak in some messy way (which it doesn't seem to be,
and I'm not going to be the one to bust it), and we have both

a) a good way to seed it "repeatably" with a parameter to srand()
b) a good way to seed it properly with platform specific Random data

it seems like a better fit to our needs.

[Mersenne Twister seems to have code to do (a), using a simpler RNG in a loop,
so if all else fails, steal that?]

Sorry if I managed to ignore it if it was mentioned earlier.

Nicholas Clark

[Craig A. Berry]

On Fri, Dec 14, 2012 at 10:28 AM, Nicholas Clark <ni...@ccl4.org> wrote:
> On Fri, Dec 14, 2012 at 10:06:30AM -0600, Craig A. Berry wrote:

>> MT Classic stores state in an array of 624, 32-bit integers.

> Gulp. Something over 2k.

My reaction as well. There is TinyMT, which I mentioned way back in
this thread, that uses only 127 bits of state:
<http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/TINYMT/index.html>.
Of course it doesn't provide as much randomness as the full MT. I
have no idea how much randomness is enough or how TinyMT compares to
drand48(). If it compares favorably, that might be the way to go.

[demerphq]

On 14 December 2012 18:24, Craig A. Berry <craig....@gmail.com> wrote:
> On Fri, Dec 14, 2012 at 10:28 AM, Nicholas Clark <ni...@ccl4.org> wrote:
>> On Fri, Dec 14, 2012 at 10:06:30AM -0600, Craig A. Berry wrote:
>
>>> MT Classic stores state in an array of 624, 32-bit integers.
>
>> Gulp. Something over 2k.
>
> My reaction as well.

Er. Why does it matter?

Yves

[demerphq]

On 14 December 2012 18:29, demerphq <deme...@gmail.com> wrote:
> On 14 December 2012 18:24, Craig A. Berry <craig....@gmail.com> wrote:
>> On Fri, Dec 14, 2012 at 10:28 AM, Nicholas Clark <ni...@ccl4.org> wrote:
>>> On Fri, Dec 14, 2012 at 10:06:30AM -0600, Craig A. Berry wrote:
>>
>>>> MT Classic stores state in an array of 624, 32-bit integers.
>>
>>> Gulp. Something over 2k.
>>
>> My reaction as well.
>
> Er. Why does it matter?

Seriously. A hash with a few hundred elements takes about 2k, and IMO
nobody would think twice about it. Why should we care at the C level?

[Leon Timmermans]

On Fri, Dec 14, 2012 at 6:24 PM, Craig A. Berry <craig....@gmail.com> wrote:
> I have no idea how much randomness is enough or how TinyMT compares
> to drand48().

An MT and an LCG are fundamentally different things; they have
different kinds of pseudo-random output. I'm not sure you can compare
them that easily. It's not just granularity and periodicy, but also
higher-scale patterns (e.g. the hyperplane issue). It's not that easy
to compare.

Leon

[demerphq]

It should be easier to compare now that I pushed smoke-me/tinymt32

:-)

Yves

[demerphq]

I just pushed smoke-me/tinymt32 which merges in the code for the 32
bit implementation of TinyMT.

Personally I think ripping out support for the system RNG would a big
improvement. It was way easier to sidestep the existing logic than it
was to integrate using tinymt32 separately.

There may be missing bits. I have a feeling I have to do something
during interpreter clone that Im not doing yet.

cheers,
Yves

commit e84331211cf28a7345869e4aef4520754a272683
Author: Yves Orton <deme...@gmail.com>
Date: Sat Dec 15 18:55:28 2012 +0100

Enable TINYMT32

Make perl use Tiny MT (32bit) RNG generator.

commit c5edafd44fe5a818fca631fc8a9c03c780b0ab4e
Author: Yves Orton <deme...@gmail.com>
Date: Sat Dec 15 18:40:33 2012 +0100

Add support for using Tiny MT, a 128 bit state RNG

Perl traditionally uses the systems random number generator
for the rand() function. This means that the behavior of code
using it is platform dependent. Worse, sometimes the random
number generator is of poor quality.

This patch adds support for using the TinyMT random number
generator (32 bit version) instead. This generator is of relatively high
quality, uses 128 bits of state and is relatively fast and should be
portable to all of our build targets.

Note that this patch does not *enable* using tinym32, to do that
you need to define TINYMT32 during configure.

Full details of the algorithm can be obtained at
http://www.math.sci.hiroshima-u.ac.jp/~m-...@math.sci.hiroshima-u.ac.jp/MT/TINYMT/

Derived from code by Mutsuo Saito and Makoto Matsumoto published under the
BSD 3-clause license.

Copyright (c) 2011 Mutsuo Saito, Makoto Matsumoto, Hiroshima
University and The University of Tokyo. All rights reserved.

[Craig A. Berry]

On Sat, Dec 15, 2012 at 12:02 PM, demerphq <deme...@gmail.com> wrote:

> I just pushed smoke-me/tinymt32 which merges in the code for the 32
> bit implementation of TinyMT.

Thanks very much for getting some code out there.

> Personally I think ripping out support for the system RNG would a big
> improvement. It was way easier to sidestep the existing logic than it
> was to integrate using tinymt32 separately.

I don't actually understand why what you did was easier than using the
existing configuration system so that the lines in config.h that
currently look something like:

#define Drand01() drand48() /**/
#define Rand_seed_t long int /**/
#define seedDrand01(x) srand48((Rand_seed_t)x) /**/
#define RANDBITS 48 /**/

would end up looking like:

#define Drand01() tinymt32_generate_double() /**/
#define Rand_seed_t U32 /**/
#define seedDrand01(x) tinymt32_init((Rand_seed_t)x) /**/
#define RANDBITS 32 /**/

I guess the question is whether people would want a way to go back to
using their system RNG for some reason.

Also, if you don't do anything to set RANDBITS, Windows will still get
15 bits, which will cause the wrong path to be taken in
util.c:Perl_seed(), where it says:

#if RANDBITS > 16
# define SEED_C1 1000003
#define SEED_C4 73819
#else
# define SEED_C1 25747
#define SEED_C4 20639
#endif

> There may be missing bits. I have a feeling I have to do something
> during interpreter clone that Im not doing yet.

Possibly something with PL_rand_state right after the line in sv.c
that looks like:

PL_srand_called = proto_perl->Isrand_called;

except that line seems wrong to me. Isn't that line saying that the
newly cloned interpreter will get the same initialization state for
rand() as the interpreter being cloned? I had thought a new
interpreter should have its own RNG and not share it with the
interpreter from which it was cloned.

I guess if the existing rand() implementation is not thread safe (such
as currently on Windows), that line does prevent us from
reinitializing it and potentially reducing randomness for the original
interpreter.

So I'd think you'd want something like (untested):

#if defined(TINYMT32) || defined(HAS_RANDOM_R) || defined(HAS_DRAND48_4)
PL_srand_called = FALSE;
Zero(&PL_rand_state, 1, tinymt32_t /* or appropriate type */);
#else
PL_srand_called = proto_perl->Isrand_called;
#endif

but I'm not really sure about that because I don't know if zeroing the
initial state is really the right thing to do for TinyMT.

Another way to do it would be to have PL_rand_state be a pointer that
gets set to null when cloning and then just malloc another tinymt_32
struct in the init function.

Other comments:

You don't need to include inttypes.h. You're using U32 rather than
uint32_t, so you don't need it, and I think it's a C99 thing rather
than C89 thing so you can't depend on its being there. Not sure if
there's any reason for stdint.h either.

I think these two functions:

+PERL_CALLCONV void Perl_tinymt32_next_state(pTHX);
+PERL_CALLCONV U32 Perl_tinymt32_temper(pTHX);

are for internal use only by TinyMT and probably shouldn't be part of
Perl's public API. They could probably be macros or static inline.

But none of this gets in the way of testing things out, so thanks
again for kick starting this.

[Craig A. Berry]

I was originally it would be part of the interpreter struct and thus
copied or initialized while cloning. If we keep adding things to the
interpreter struct 2k at a time, pretty soon we're talking about real
memory.

[demerphq]

On 15 December 2012 22:05, Craig A. Berry <craig....@gmail.com> wrote:
> On Sat, Dec 15, 2012 at 12:02 PM, demerphq <deme...@gmail.com> wrote:
>
>> I just pushed smoke-me/tinymt32 which merges in the code for the 32
>> bit implementation of TinyMT.
>
> Thanks very much for getting some code out there.

No problem. Glad to help. It is very educational too. :-)

>> Personally I think ripping out support for the system RNG would a big
>> improvement. It was way easier to sidestep the existing logic than it
>> was to integrate using tinymt32 separately.
>
> I don't actually understand why what you did was easier than using the
> existing configuration system so that the lines in config.h that
> currently look something like:
>
> #define Drand01() drand48() /**/
> #define Rand_seed_t long int /**/
> #define seedDrand01(x) srand48((Rand_seed_t)x) /**/
> #define RANDBITS 48 /**/
>
> would end up looking like:
>
> #define Drand01() tinymt32_generate_double() /**/
> #define Rand_seed_t U32 /**/
> #define seedDrand01(x) tinymt32_init((Rand_seed_t)x) /**/
> #define RANDBITS 32 /**/

Because I couldnt figure out how to make it actually end up like that.
And there are more macros defintions involved. Look for drand01, and
structs related to drand48_r and related macros. Do some greps for
things like Drand01 (case insensitive).

The support for existing reentrant versions (like we really want)
currently involves gunk all over the place, in embed.h, and the
rentrant struct etc. Not to mention every platform specific configure
file, etc.

So basically I couldnt figure out how to elegantly get all the
different builds to actually use the new tinymt code properly.

What I envisage is the ability to select one of several prebundled RNG
implementations. The state for these RNG's should live in
PL_random_state, and the logic to handling cloning it or using it
should be relatively generic.

Anyway, maybe its just my inexperience with this part of the
infrastructure talking but I was getting really frustrated doing it
the "right way" and then realized it was pretty simple to bypass the
frustrating stuff. :-)

> I guess the question is whether people would want a way to go back to
> using their system RNG for some reason.

Yeah exactly.

> Also, if you don't do anything to set RANDBITS, Windows will still get
> 15 bits, which will cause the wrong path to be taken in
> util.c:Perl_seed(), where it says:
>
> #if RANDBITS > 16
> # define SEED_C1 1000003
> #define SEED_C4 73819
> #else
> # define SEED_C1 25747
> #define SEED_C4 20639
> #endif

Oh, good catch. It reminds me that ideally TinyMT would be initialized
using the array interface, so we could supply 127 bits worth of seed
instead of just 32 bits.

>
>> There may be missing bits. I have a feeling I have to do something
>> during interpreter clone that Im not doing yet.
>
> Possibly something with PL_rand_state right after the line in sv.c
> that looks like:
>
> PL_srand_called = proto_perl->Isrand_called;

Yes right. Ill get on that pronto.

> except that line seems wrong to me. Isn't that line saying that the
> newly cloned interpreter will get the same initialization state for
> rand() as the interpreter being cloned? I had thought a new
> interpreter should have its own RNG and not share it with the
> interpreter from which it was cloned.

I think this is hysterical raisons.

We cant tell the difference between "they called srand() explicitly"
and "we called srand() implicitly" so any choice here is going to be
wrong some of the time. Currently I believe it is relatively well
known (although perhaps undocumented) that the flag is NOT reset at
fork, children get the same RNG sequence as the parent. This could be
construed as a feature. :-)

> I guess if the existing rand() implementation is not thread safe (such
> as currently on Windows), that line does prevent us from
> reinitializing it and potentially reducing randomness for the original
> interpreter.
>
> So I'd think you'd want something like (untested):
>
> #if defined(TINYMT32) || defined(HAS_RANDOM_R) || defined(HAS_DRAND48_4)
> PL_srand_called = FALSE;
> Zero(&PL_rand_state, 1, tinymt32_t /* or appropriate type */);
> #else
> PL_srand_called = proto_perl->Isrand_called;
> #endif
>
> but I'm not really sure about that because I don't know if zeroing the
> initial state is really the right thing to do for TinyMT.

No, it definitely isn't. MT style RNG's usually use LCG's to generate
their seed, so if we reinitialize we must call tinymt32_init().

IMO we should just copy the state as that will replicate existing
behaviour. If they call srand() (no args) then we will do the right
kind of initialization.

> Another way to do it would be to have PL_rand_state be a pointer that
> gets set to null when cloning and then just malloc another tinymt_32
> struct in the init function.

Nod.

> Other comments:
>
> You don't need to include inttypes.h. You're using U32 rather than
> uint32_t, so you don't need it, and I think it's a C99 thing rather
> than C89 thing so you can't depend on its being there. Not sure if
> there's any reason for stdint.h either.

Yeah, i ripped them out in the latest version of the smoke-me/tinymt32
branch. (I rebased so you may have to do a hard reset to the new
version).

> I think these two functions:
>
> +PERL_CALLCONV void Perl_tinymt32_next_state(pTHX);
> +PERL_CALLCONV U32 Perl_tinymt32_temper(pTHX);
>
> are for internal use only by TinyMT and probably shouldn't be part of
> Perl's public API. They could probably be macros or static inline.

Yes right. Another good catch. I meant to ask about what signature
these subs should have. I wasnt sure what would have to be public and
what not, or if they should be visible but private, or what.

> But none of this gets in the way of testing things out, so thanks
> again for kick starting this.

Thanks for pointing out the code and the patch review. If I have not
already acted on anything you mentioned that is important please
nudge me again.

cheers,
Yves

[demerphq]

On 15 December 2012 22:09, Craig A. Berry <craig....@gmail.com> wrote:
> On Fri, Dec 14, 2012 at 11:44 AM, demerphq <deme...@gmail.com> wrote:
>> On 14 December 2012 18:29, demerphq <deme...@gmail.com> wrote:
>>> On 14 December 2012 18:24, Craig A. Berry <craig....@gmail.com> wrote:
>>>> On Fri, Dec 14, 2012 at 10:28 AM, Nicholas Clark <ni...@ccl4.org> wrote:
>>>>> On Fri, Dec 14, 2012 at 10:06:30AM -0600, Craig A. Berry wrote:
>>>>
>>>>>> MT Classic stores state in an array of 624, 32-bit integers.
>>>>
>>>>> Gulp. Something over 2k.
>>>>
>>>> My reaction as well.
>>>
>>> Er. Why does it matter?
>>
>> Seriously. A hash with a few hundred elements takes about 2k, and IMO
>> nobody would think twice about it. Why should we care at the C level?
>
> I was originally it would be part of the interpreter struct and thus
> copied or initialized while cloning.

Yes I think it would have to be. So good point.

> If we keep adding things to the
> interpreter struct 2k at a time, pretty soon we're talking about real
> memory.

Ok, fair enough. I kinda think we should leave those kinds of
decisions to people building perl, but choosing "sensible" defaults
make sense too. Anyway, obviously I went with your counter proposal.
:-)

cheers,

[bulk88]

2KB for MT would double the interp in size. MT should not be in the
interp and should be left for being CPAN/XS module.

http://www.nntp.perl.org/group/perl.perl5.porters/2012/10/msg194823.html

[bulk88 via RT]

On Sat Dec 15 13:06:09 2012, craig....@gmail.com wrote:

> I guess if the existing rand() implementation is not thread safe (such
> as currently on Windows), that line does prevent us from
> reinitializing it and potentially reducing randomness for the original
> interpreter.
>

rand on windows uses TLS internally and is thread safe aslong as you
used a multithreaded clib/crt. Nobody uses anything but the
multithreaded crt on windows.

--
bulk88 ~ bulk88 at hotmail.com

[bulk88 via RT]

On Sat Dec 15 09:58:46 2012, demerphq wrote:
>
> It should be easier to compare now that I pushed smoke-me/tinymt32
>
> :-)
>
> Yves
>
>

I will benchmark that soon.

--
bulk88 ~ bulk88 at hotmail.com

[Craig A. Berry]

On Sat, Dec 15, 2012 at 8:25 PM, bulk88 via RT
<perlbug-...@perl.org> wrote:
> On Sat Dec 15 13:06:09 2012, craig....@gmail.com wrote:
>
>> I guess if the existing rand() implementation is not thread safe (such
>> as currently on Windows), that line does prevent us from
>> reinitializing it and potentially reducing randomness for the original
>> interpreter.
>>
>
> rand on windows uses TLS internally and is thread safe aslong as you
> used a multithreaded clib/crt. Nobody uses anything but the
> multithreaded crt on windows.

OK, thanks. I can't remember where I got that impression but I'm glad
to be corrected.

[demerphq]

On 15 December 2012 22:33, demerphq <deme...@gmail.com> wrote:
> Thanks for pointing out the code and the patch review. If I have not
> already acted on anything you mentioned that is important please
> nudge me again.

I have rebased and repushed all of these changes as
smoke-me/builtin_rng, along with the addition of another RNG, WellRNG
512.

If you feel up to redoing your code review id be grateful. :-)

cheers

[Craig A. Berry]

On Sun, Dec 16, 2012 at 4:01 PM, demerphq <deme...@gmail.com> wrote:
> On 15 December 2012 22:33, demerphq <deme...@gmail.com> wrote:
>> Thanks for pointing out the code and the patch review. If I have not
>> already acted on anything you mentioned that is important please
>> nudge me again.
>
> I have rebased and repushed all of these changes as
> smoke-me/builtin_rng, along with the addition of another RNG, WellRNG
> 512.
>
> If you feel up to redoing your code review id be grateful. :-)

It might be preferable to have slightly more descriptive macro names
saying which RNG we're using, such as PERL_RAND_IS_WELLRNG512A,
PERL_RAND_IS_TINYMT32, etc. Or USE_xxx_RANDFUNC or something. Other
than that, no problems with the code jumped out at me.

I had already done a basic build with TinyMT32 on OS X (Intel) and on
VMS (Itanium) and encountered no troubles.

There are a number of TODOs. And no, I'm not saying Yves should be
the one to do them (nor volunteering to do them myself!).

How do we test these things? I think ideally we should be able to
prove that the RNGs display the same characteristics described by the
authors in their original papers. And we should be able to prove it
for big endian as well as little endian, quadword and octaword
alignment as well as longword alignment, etc. I don't know offhand
how to do that. The comments at the top of t/op/rand.t indicate it
was written in 1997. I would guess it's not completely suitable for
the new regime.

Sorting out what to do in Configure will have to be done before we
ship so that people don't check $Config{randfunc} and get 'drand48'
and then wonder why the numbers are so different from what drand48()
provides on their system.

There are a couple of intersecting concerns related to Perl
interpreter state and RNG state. If we ever want to change the
hysterical raisins that make a newly cloned interpreter start off with
but then diverge from the RNG state of the parent thread (such as by
initializing its own state), now would be the time to do it. I have
a feeling that this is yet another aspect of "threads happened" and
weren't exactly planned for from the start.

But even if this current scheme is desired for some reason, adding a
large amount of data that gets cloned may be undesirable. TinyMT
nominally has 7 longwords in its state struct and WellRNG512 has 17,
but with padding for alignment, actual size will usually be somewhat
more than that. Maybe that's not large, especially compared to the
624 longwords of the original MT. But maybe in general the cost of
allocating and initializing the state of a RNG should happen when
asking for a random number rather than when asking for a new thread.

[Nicholas Clark]

Yes, I was assuming this too. 2k per interpreter struct feels like a lot
of "overhead". Also (and I've not checked the code for how much it does
this), having to inspect 2k of state to generate each random number
potentially means quite a bit of L1 cache churn, for a lot more quality
than most code is going to need.

Nicholas Clark

[Nicholas Clark]

On Sat, Dec 15, 2012 at 10:33:09PM +0100, demerphq wrote:
> On 15 December 2012 22:05, Craig A. Berry <craig....@gmail.com> wrote:

> > except that line seems wrong to me. Isn't that line saying that the
> > newly cloned interpreter will get the same initialization state for
> > rand() as the interpreter being cloned? I had thought a new
> > interpreter should have its own RNG and not share it with the
> > interpreter from which it was cloned.
>
> I think this is hysterical raisons.

Right now, sv.c does this in perl_clone():

PL_srand_called = proto_perl->Isrand_called;

but threads.xs fights it with this:

SAVEBOOL(PL_srand_called); /* Save this so it becomes the correct value */
PL_srand_called = FALSE; /* Set it to false so we can detect if it gets
set during the clone */


which seems to be a somewhat obfuscated way of not writing this post-clone:

thread->interp->Isrand_called = FALSE;


Upshot is that srand is marked as "not yet called" in any child ithreads.
(But the state is preserved in any pseudo-forked threads)


This is consistent with how C code would behave it written to use re-entrant
variants of rand()-like functions and a per-thread state. C code written
like this would have a different pseudo-random sequence in each thread,
but the state would be the same in any child process forked.

> We cant tell the difference between "they called srand() explicitly"
> and "we called srand() implicitly" so any choice here is going to be
> wrong some of the time. Currently I believe it is relatively well
> known (although perhaps undocumented) that the flag is NOT reset at
> fork, children get the same RNG sequence as the parent. This could be
> construed as a feature. :-)

Well, we can't unless we change the code :-)

> > I guess if the existing rand() implementation is not thread safe (such
> > as currently on Windows), that line does prevent us from
> > reinitializing it and potentially reducing randomness for the original
> > interpreter.
> >
> > So I'd think you'd want something like (untested):
> >
> > #if defined(TINYMT32) || defined(HAS_RANDOM_R) || defined(HAS_DRAND48_4)
> > PL_srand_called = FALSE;
> > Zero(&PL_rand_state, 1, tinymt32_t /* or appropriate type */);
> > #else
> > PL_srand_called = proto_perl->Isrand_called;
> > #endif
> >
> > but I'm not really sure about that because I don't know if zeroing the
> > initial state is really the right thing to do for TinyMT.
>
> No, it definitely isn't. MT style RNG's usually use LCG's to generate
> their seed, so if we reinitialize we must call tinymt32_init().
>
> IMO we should just copy the state as that will replicate existing
> behaviour. If they call srand() (no args) then we will do the right
> kind of initialization.

To be consistent with what we currently do, we need to copy the state for
pseudofork, but reset it for ithread creation.

> > Another way to do it would be to have PL_rand_state be a pointer that
> > gets set to null when cloning and then just malloc another tinymt_32
> > struct in the init function.
>
> Nod.

This seems more space efficient, if most threads don't use rand(). But
without trying it an measuring it, it's not clear how much it complicates
the code.


Thanks for starting to condense the ideas into working code.

Nicholas Clark

[demerphq]

Indeed. I think that there is some kind of difference between someone
calling srand() explicitly, and us calling it implicitly.

I think that an explicit call to srand() should cause all forks/child
threads to use the SAME srand() value.

I think an implicit call to srand() should cause forks or child
threads to re-seed on the first use of rand.

>
>> > I guess if the existing rand() implementation is not thread safe (such
>> > as currently on Windows), that line does prevent us from
>> > reinitializing it and potentially reducing randomness for the original
>> > interpreter.
>> >
>> > So I'd think you'd want something like (untested):
>> >
>> > #if defined(TINYMT32) || defined(HAS_RANDOM_R) || defined(HAS_DRAND48_4)
>> > PL_srand_called = FALSE;
>> > Zero(&PL_rand_state, 1, tinymt32_t /* or appropriate type */);
>> > #else
>> > PL_srand_called = proto_perl->Isrand_called;
>> > #endif
>> >
>> > but I'm not really sure about that because I don't know if zeroing the
>> > initial state is really the right thing to do for TinyMT.
>>
>> No, it definitely isn't. MT style RNG's usually use LCG's to generate
>> their seed, so if we reinitialize we must call tinymt32_init().
>>
>> IMO we should just copy the state as that will replicate existing
>> behaviour. If they call srand() (no args) then we will do the right
>> kind of initialization.
>
> To be consistent with what we currently do, we need to copy the state for
> pseudofork, but reset it for ithread creation.

I am not quite getting you here. Can you check what I did in the
smoke-me/builtin_rng and let me know if I got it right?

We copy the state. If the flag saying that srand was called is false
we will reinitialize it. Im assuming that I do not need to do anything
beyond ensure that the state is copied properly, and that the way the
flag is managed takes care of the rest.

>> > Another way to do it would be to have PL_rand_state be a pointer that
>> > gets set to null when cloning and then just malloc another tinymt_32
>> > struct in the init function.
>>
>> Nod.
>
> This seems more space efficient, if most threads don't use rand(). But
> without trying it an measuring it, it's not clear how much it complicates
> the code.

The biggest objection is that a naive implementation would make
non-threaded perls pay a price every startup to avoid a cost that only
affects threaded builds at thread spawn time. It seems to me like
worrying about this is silly until there is evidence that it actually
matters. Its seems pretty clear that 2000 bytes worth of reference
data is too much, but it does not seem so clear to me that we should
worry about 20 to 64 bytes. (the size of the state for TintMT or
WellRNG512a).

> Thanks for starting to condense the ideas into working code.

You are welcome.

cheers

[demerphq]

Oh, one point worth keeping in mind is that the FreeBSD version of
dword48 uses 7 * U16, four are statically declared data structures
holding constants, but the srand() implementation updates all of them.
Thus the current dword48 implementation is using 112 bits of storage.
Which isnt that much less than the 128 bits of storage used by TinyMT.

Anyway, I am nearly ready with code to put the FreeBSD LCG based
dword48 directly into Perl itself.

[demerphq]

On 17 December 2012 04:11, Craig A. Berry <craig....@gmail.com> wrote:
> On Sun, Dec 16, 2012 at 4:01 PM, demerphq <deme...@gmail.com> wrote:
>> On 15 December 2012 22:33, demerphq <deme...@gmail.com> wrote:
>>> Thanks for pointing out the code and the patch review. If I have not
>>> already acted on anything you mentioned that is important please
>>> nudge me again.
>>
>> I have rebased and repushed all of these changes as
>> smoke-me/builtin_rng, along with the addition of another RNG, WellRNG
>> 512.
>>
>> If you feel up to redoing your code review id be grateful. :-)
>
> It might be preferable to have slightly more descriptive macro names
> saying which RNG we're using, such as PERL_RAND_IS_WELLRNG512A,
> PERL_RAND_IS_TINYMT32, etc. Or USE_xxx_RANDFUNC or something. Other
> than that, no problems with the code jumped out at me.

Good idea, I will sort that out. Thanks. I was also thinking that the
stuff with PERL_RANDOM_STATE_TYPE could be simplified, there is an
extra layer of indirection there that could be avoiding, by just
having a random_state_t which we typedef appropriately. I guess I was
keeping the door open to selecting the random number engine, possibly
via use, or something like that. (Think how we support multiple sort
algorithms).

> I had already done a basic build with TinyMT32 on OS X (Intel) and on
> VMS (Itanium) and encountered no troubles.

Cool.

> There are a number of TODOs. And no, I'm not saying Yves should be
> the one to do them (nor volunteering to do them myself!).

:-)

> How do we test these things?

Once we have a correct implementation, we should, for any given seed,
always get the same results. So we can do some form of test that
things are working as expected. I will look for test input/output for
the published algorithms.

> I think ideally we should be able to
> prove that the RNGs display the same characteristics described by the
> authors in their original papers.
> And we should be able to prove it
> for big endian as well as little endian, quadword and octaword
> alignment as well as longword alignment, etc. I don't know offhand
> how to do that.

Here are some RNG testers:

http://www.phy.duke.edu/~rgb/General/dieharder.php
http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/ADAPTIVE/
http://www.iro.umontreal.ca/~simardr/testu01/tu01.html

But these kind of things will be unsuitable for general testing
purposes. They would have to be run by hand, and might be hard to
build on all of our target platforms.

> The comments at the top of t/op/rand.t indicate it
> was written in 1997. I would guess it's not completely suitable for
> the new regime.

It looks to me like its not so bad. The new regime will misrepresent
RANDBITS but beyond that the tests seem reasonable. Its pretty similar
to what I was doing by hand for quick testing. (Like rolling 2d6 100k
times and checking the distribution).

> Sorting out what to do in Configure will have to be done before we
> ship so that people don't check $Config{randfunc} and get 'drand48'
> and then wonder why the numbers are so different from what drand48()
> provides on their system.

I think we should be careful here. I think the current Configure logic
should mostly stay the same. It exposes what RNG is offered by the
/system/. We haven't changed that, we have merely added the option to
use something other than the system RNG (whatever it might be). So
instead of changing $Config{randfunc} to include say tinymt32, we
should add a new var, whose current legal value would be one of
"$Config{randfunc}", "tinymt32" or "wellrng512a". (Add some salt and
handwaving here)

cheers,

[Craig A. Berry]

I sympathize with wanting to simply step over the puddle of goo, but
I'm not sure we can.

Porting/Glossary says:

4277 randfunc (randfunc.U):
4278 Indicates the name of the random number function to use.
4279 Values include drand48, random, and rand. In C programs,
4280 the 'Drand01' macro is defined to generate uniformly distributed
4281 random numbers over the range [0., 1.[ (see drand01 and nrand).

nrand appears not to exist anymore, so that looks like a broken cross
reference. drand01 is documented to be a synonym for Drand01 and to
give you randfunc. XS code could be using Drand01 or drand01 and
expecting to get the RNG that was selected at Configure-time.

But as far as HAS_RANDOM_R and so forth, I agree with you -- those are
indications of system capability rather than what Perl is using and
can be left alone.

[demerphq]

Maybe we could clean stuff up and still keep Drand01 and drand01
defines for back compat?

I really think the way that the flag about whether the seed is set or
not, and the lack of explicit state arguments in the interfaces for
Drand01 and srand() functions makes life more difficult for everybody
than it should be. If we clean it up and restructure it a bit we
could probably still maintain back compat but have something saner to
use in the future.

> But as far as HAS_RANDOM_R and so forth, I agree with you -- those are
> indications of system capability rather than what Perl is using and
> can be left alone.

Nod.

[demerphq]

On 18 December 2012 00:25, bulk88 via RT <perlbug-...@perl.org> wrote:

> On Sat Dec 15 18:28:47 2012, bulk88 wrote:
>> On Sat Dec 15 09:58:46 2012, demerphq wrote:
>> >
>> > It should be easier to compare now that I pushed smoke-me/tinymt32
>> >
>> > :-)
>> >
>> > Yves
>> >
>> >
>>
>> I will benchmark that soon.
>>
>

> VC2008, activeperl 5.12.3 32 bit, opteron 875 Server 2003 x64, nothing
> inlined, tinymt32, wellrng, and drand48 were all func calls in asm in
> RunRand, Perl_tinymt32_generate_U32 was NOT inlined away (unlike VC 2003
> later on)
> ______________________________________________________________________
> new time=14.107506, opt=-MD -Zi -DNDEBUG -O1 4970233.211236
> old time=0.639507, opt=-MD -Zi -DNDEBUG -O1 9997076.043274
> new kmx time=2.266837, opt=-MD -Zi -DNDEBUG -O1 10000148.315117
> new kmx bulk88 1 time=2.649125, opt=-MD -Zi -DNDEBUG -O1 10001037.278379
> new kmx bulk88 2 time=2.376385, opt=-MD -Zi -DNDEBUG -O1 9999965.239603
> new kmx bulk88 3 time=3.024081, opt=-MD -Zi -DNDEBUG -O1 9999970.013225
> just rand time=2.130050, opt=-MD -Zi -DNDEBUG -O1 1310540572828.000000
> rand_s time=16.052782, opt=-MD -Zi -DNDEBUG -O1 4999226.450462
> drand48 time=0.510524, opt=-MD -Zi -DNDEBUG -O1 10000680.628264
> TINYMT32 time=0.747539, opt=-MD -Zi -DNDEBUG -O1 9998574.703049
> WELLRNG512A time=0.465683, opt=-MD -Zi -DNDEBUG -O1 9998797.404216
> ______________________________________________________________________
> VC 2008, x64 Perl 5.17 (93a641ae382638ffd1980378be4810244d04f4b0),
> opteron 875 Server 2003 x64, asm shows everything in the DLL was inlined
> automatically by the compiler (drand48, tinymt32 and wellrng and all
> children calls were
> inlined). rand and rand_s can't be inlined for obvious reasons), this is
> apples to oranges, so I will have to repeat this with disabling the inlining
> ______________________________________________________________________
> new time=10.479223, opt=-MD -Zi -DNDEBUG -O1 -GL -GS- -favor:AMD64
> -fp:precise 4
> 970314.444337
> old time=0.427892, opt=-MD -Zi -DNDEBUG -O1 -GL -GS- -favor:AMD64
> -fp:precise 99
> 98732.212921
> new kmx time=1.629771, opt=-MD -Zi -DNDEBUG -O1 -GL -GS- -favor:AMD64
> -fp:precis
> e 9999588.845403
> new kmx bulk88 1 time=1.711441, opt=-MD -Zi -DNDEBUG -O1 -GL -GS-
> -favor:AMD64 -
> fp:precise 10001166.772142
> new kmx bulk88 2 time=1.639529, opt=-MD -Zi -DNDEBUG -O1 -GL -GS-
> -favor:AMD64 -
> fp:precise 9999228.837134
> new kmx bulk88 3 time=1.484076, opt=-MD -Zi -DNDEBUG -O1 -GL -GS-
> -favor:AMD64 -
> fp:precise 10000031.193598
> just rand time=1.429209, opt=-MD -Zi -DNDEBUG -O1 -GL -GS- -favor:AMD64
> -fp:prec
> ise 1310703320700.000000
> rand_s time=10.888868, opt=-MD -Zi -DNDEBUG -O1 -GL -GS- -favor:AMD64
> -fp:precis
> e 5000939.396292
> drand48 time=0.068243, opt=-MD -Zi -DNDEBUG -O1 -GL -GS- -favor:AMD64
> -fp:precis
> e 10001741.884865
> TINYMT32 time=0.427861, opt=-MD -Zi -DNDEBUG -O1 -GL -GS- -favor:AMD64
> -fp:preci
> se 9997611.127257
> WELLRNG512A time=0.254936, opt=-MD -Zi -DNDEBUG -O1 -GL -GS-
> -favor:AMD64 -fp:pr
> ecise 10001033.492300
> ______________________________________________________________________
> VC 2003, 32 bit Perl 5.17 93a641ae382638ffd1980378be4810244d04f4b0,
> Windows XP 32 bit, core 2 duo t7200, nothing was inlined in asm in
> RunRand, Perl_tinymt32_generate_U32 was inlined away.
> ______________________________________________________________________
> new time=6.801278, opt=-O1 -GL -arch:SSE2 -GF 4970279.404719
> old time=0.342173, opt=-O1 -GL -arch:SSE2 -GF 9998786.194519
> new kmx time=1.348873, opt=-O1 -GL -arch:SSE2 -GF 10001408.239714
> new kmx bulk88 1 time=1.284122, opt=-O1 -GL -arch:SSE2 -GF 9998050.657198
> new kmx bulk88 2 time=1.723125, opt=-O1 -GL -arch:SSE2 -GF 10001843.021335
> new kmx bulk88 3 time=1.273438, opt=-O1 -GL -arch:SSE2 -GF 10000463.515876
> just rand time=1.175337, opt=-O1 -GL -arch:SSE2 -GF 1310685646068.000000
> rand_s time=16.356017, opt=-O1 -GL -arch:SSE2 -GF 9997707.487858
> drand48 time=0.498862, opt=-O1 -GL -arch:SSE2 -GF 9999970.685309
> TINYMT32 time=0.601499, opt=-O1 -GL -arch:SSE2 -GF 10001623.470601
> WELLRNG512A time=0.488246, opt=-O1 -GL -arch:SSE2 -GF 9999946.450606
> _______________________________________________________________________
>
> I'm not sure how to interpret the benchmark results. Since Wellrng512a
> and drand48 are the same on 32 bit. drand48 on 32 bit windows makes a
> call to __allmul, Perl_wellrng512a_generate_double makes no calls. On
> x64 drand48 is many times faster than wellrng512a (both were inlined, I
> will try to stop that in the near future). wellrng512a takes 68 bytes
> (according to the asm in x64 dll), tinymt32 takes 28 bytes (according to
> the asm in x64 dll), drand48 takes 8 bytes. x64 Windows isn't aligning
> the U32s to U64. IDK if any other 64 bit platforms will try align all
> those U32s to U64s. drand48 is in a OS's clib. Do any OSes use
> wellrng512a and tinymt32 in their clibs?
>
> If someone is interested, I will write a patch to use drand48 on Windows
> only.

I will push a patch so everyone can use the same "drand48".

> One thing not addressed in this thread so far, what happens on platforms
> where an NV is a float (32), or a long double (80), or a quad
> float/__float128? None of the solutions in this ticket address that
> other my uselessly slow first try (with the Perl_pow and Perl_fmod calls).

I picked 32 bit based RNG's. There are other RNG's that are more
suitable for 64bit. I was thinking for instance of supporting TinyMT64
as well as TinyMT32. Also there are other RNG's to try. :-)

[Nicholas Clark]

On Mon, Dec 17, 2012 at 03:25:15PM -0800, bulk88 via RT wrote:

> One thing not addressed in this thread so far, what happens on platforms
> where an NV is a float (32), or a long double (80), or a quad
> float/__float128? None of the solutions in this ticket address that
> other my uselessly slow first try (with the Perl_pow and Perl_fmod calls).

NV is never a (32 bit) float.

As to the others, I don't know, but

a) they all have more precision than a double, so problems are going to be
more subtle than losing information.
b) all the current system random number functions we've used haven't had
problems.


We seem to be loosing track of the original bug - that 15 bits of random
number is not enough. 32 seemed generally to be considered just fine.

Nicholas Clark

[Nicholas Clark]

On Mon, Dec 17, 2012 at 08:54:19PM +0100, demerphq wrote:

> Good idea, I will sort that out. Thanks. I was also thinking that the
> stuff with PERL_RANDOM_STATE_TYPE could be simplified, there is an
> extra layer of indirection there that could be avoiding, by just
> having a random_state_t which we typedef appropriately. I guess I was
> keeping the door open to selecting the random number engine, possibly
> via use, or something like that. (Think how we support multiple sort
> algorithms).

You seem to have caught second system effect. I don't think that antibiotics
work on that. :-(

On Mon, Dec 17, 2012 at 10:26:10PM +0100, demerphq wrote:

> Maybe we could clean stuff up and still keep Drand01 and drand01
> defines for back compat?
>
> I really think the way that the flag about whether the seed is set or
> not, and the lack of explicit state arguments in the interfaces for
> Drand01 and srand() functions makes life more difficult for everybody
> than it should be. If we clean it up and restructure it a bit we
> could probably still maintain back compat but have something saner to
> use in the future.

I only counted 12 distributions on CPAN using Drand01() outside of the core.
IIRC two of them are dual-life.*

I agree that the lack of state makes it hard to do anything flexible, but
I'm not convinced that that's our target market. The XS APIs tend to assume
an interpreter struct is kicking around to provide state, and that it's
one god**-object of state per interpreter thread, take it or leave it.

I don't think that it's worth spending effort creating a better random
number framework, either for Perl code or for XS code. If Perl code wants
a more flexible framework than the core provides, it can use something from
CPAN. If XS code wants a more flexible C interface, then it has plenty of
choice of code to (re)use.

Nicholas Clark

* csearch 'use\s+sort.*(stable|qsort|mergesort)'
suggests that 22 distributions outside the core use sort.
Of those, only 3 use anything other than 'stable'.
It's really popular.
** our god object is neither very god-like nor very object-like.
If it even had feet of clay that would be progress.