parent
edb58bdc1d
commit
f7754a142e
Binary file not shown.
Binary file not shown.
File diff suppressed because one or more lines are too long
File diff suppressed because it is too large
Load Diff
File diff suppressed because one or more lines are too long
File diff suppressed because one or more lines are too long
File diff suppressed because it is too large
Load Diff
File diff suppressed because one or more lines are too long
File diff suppressed because one or more lines are too long
File diff suppressed because one or more lines are too long
@ -0,0 +1,597 @@
|
||||
<!DOCTYPE HTML PUBLIC "-//W3C//DTD HTML 4.01//EN" "http://www.w3.org/TR/html4/strict.dtd">
|
||||
<html>
|
||||
<head>
|
||||
<link href="https://fonts.googleapis.com/css?family=Open+Sans:400,400i,700,700i" rel="stylesheet">
|
||||
<meta http-equiv="Content-Type" content="text/html;charset=utf-8">
|
||||
<meta name="keywords" content="rng, prng, xoshiro, xoroshiro, xorshift, pseudorandom number generator, random number generator">
|
||||
<style type="text/css">
|
||||
@import "css/content.php";
|
||||
@import "css/layout.php";
|
||||
@import "css/tablesorter.css";
|
||||
|
||||
</style>
|
||||
<title>xoshiro/xoroshiro generators and the PRNG shootout</title>
|
||||
<script type="text/javascript" src="js/jquery.js"></script>
|
||||
<script type="text/javascript" src="js/tablesorter.js"></script>
|
||||
<script type="text/javascript" src="js/metadata.js"></script>
|
||||
<script type="text/javascript">
|
||||
$.tablesorter.defaults.widgets = ['zebra'];
|
||||
$(document).ready( function() {
|
||||
$("#prng").tablesorter({
|
||||
headers: {
|
||||
3: {
|
||||
sorter: false
|
||||
},
|
||||
4: {
|
||||
sorter: false
|
||||
}
|
||||
}
|
||||
});
|
||||
$("#vect").tablesorter({
|
||||
});
|
||||
$("#prngf").tablesorter({
|
||||
headers: {
|
||||
3: {
|
||||
sorter: false
|
||||
},
|
||||
4: {
|
||||
sorter: false
|
||||
}
|
||||
}
|
||||
});
|
||||
} );
|
||||
</script>
|
||||
|
||||
<script type="text/javascript" src="https://polyfill.io/v3/polyfill.min.js?features=es6"></script>
|
||||
<script type="text/javascript" src="https://cdn.jsdelivr.net/npm/mathjax@3/es5/tex-mml-chtml.js"></script>
|
||||
|
||||
</head>
|
||||
<body>
|
||||
<div id=header><code>xoshiro</code> / <code>xoroshiro</code> generators and the PRNG shootout</div>
|
||||
|
||||
<div id=left>
|
||||
<ul id="left-nav">
|
||||
<li><a href="http://vigna.di.unimi.it/">Home</a></li>
|
||||
<li><a href="http://vigna.di.unimi.it/papers.php">Papers<span class=arrow>➔</span></a></li>
|
||||
<li><a href="http://vigna.di.unimi.it/software.php">Software<span class=arrow>➔</span></a>
|
||||
<li><strong><a href="http://prng.di.unimi.it/">PRNG shootout</a></strong>
|
||||
|
||||
<ul>
|
||||
<li><a href="#intro">Introduction</A></li>
|
||||
<li><a href="#shootout">A PRNG Shootout</A></li>
|
||||
<li><a href="#speed">Speed</A></li>
|
||||
<li><a href="#quality">Quality</A></li>
|
||||
<li><a href="#remarks">Remarks</A></li>
|
||||
</ul>
|
||||
|
||||
|
||||
<li><a href="http://pcg.di.unimi.it/pcg.php">The wrap-up on PCG generators<span class=arrow>➔</span></a>
|
||||
<li><a href="http://fastutil.di.unimi.it/"><code>fastutil</code><span class=arrow>➔</span></a>
|
||||
<li><a href="http://dsiutils.di.unimi.it/">DSI utilities<span class=arrow>➔</span></a>
|
||||
<li><a href="http://webgraph.di.unimi.it/">WebGraph<span class=arrow>➔</span></a>
|
||||
<li><a href="http://sux.di.unimi.it/">Sux<span class=arrow>➔</span></a>
|
||||
<li><a href="http://vigna.di.unimi.it/music.php">Music<span class=arrow>➔</span></a>
|
||||
<li><a href="https://shrinkai.di.unimi.it/">Shrink AI<span class=arrow>➔</span></a>
|
||||
</ul>
|
||||
</div>
|
||||
|
||||
<div id="main">
|
||||
|
||||
|
||||
<div id="content">
|
||||
|
||||
<h1 class=first>Introduction</h1>
|
||||
|
||||
<p>This page describes some new pseudorandom number generators (PRNGs) we (David Blackman and I) have been working on recently, and
|
||||
a shootout comparing them with other generators. Details about the generators can
|
||||
be found in our <a
|
||||
href="http://vigna.di.unimi.it/papers.php#BlVSLPNG">paper</a>. Information about my previous <code>xorshift</code>-based
|
||||
generators can be found <a href="xorshift.php">here</a>, but they have been entirely superseded by the new ones, which
|
||||
are faster <em>and</em> better. As part of our study, we developed a very strong <a href="hwd.php">test for Hamming-weight dependencies</a>
|
||||
that gave a number of surprising results.
|
||||
|
||||
<h1>64-bit Generators</h1>
|
||||
|
||||
<P><a href="xoshiro256plusplus.c"><code>xoshiro256++</code></a>/<a href="xoshiro256starstar.c"><code>xoshiro256**</code></a>
|
||||
(XOR/shift/rotate) are our <strong>all-purpose</strong>
|
||||
generators (not <em>cryptographically secure</em> generators, though,
|
||||
like all PRNGs in these pages). They have excellent (sub-ns) speed, a state
|
||||
space (256 bits) that is large enough for any parallel application, and
|
||||
they pass all tests we are aware of. See the <a href="http://vigna.di.unimi.it/papers.php#BlVSLPNG">paper</a>
|
||||
for a discussion of their differences.
|
||||
|
||||
<p>If, however, one has to generate only 64-bit <strong>floating-point</strong> numbers
|
||||
(by extracting the upper 53 bits) <a
|
||||
href="xoshiro256plus.c"><code>xoshiro256+</code></a> is a slightly (≈15%)
|
||||
faster generator with analogous statistical properties. For general
|
||||
usage, one has to consider that its lowest bits have low linear
|
||||
complexity and will <a href="lowcomp.php">fail linearity tests</a>; however, low linear
|
||||
complexity of the lowest bits can have hardly any impact in practice, and certainly has no
|
||||
impact at all if you generate floating-point numbers using the upper bits (we computed a <a
|
||||
href="http://vigna.di.unimi.it/papers.php#BlVSLPNG">precise
|
||||
estimate</a> of the linear complexity of the lowest bits).
|
||||
|
||||
<p>If you are <strong>tight on space</strong>, <a
|
||||
href="xoroshiro128plusplus.c"><code>xoroshiro128++</code></a>/<a
|
||||
href="xoroshiro128starstar.c"><code>xoroshiro128**</code></a>
|
||||
(XOR/rotate/shift/rotate) and <a
|
||||
href="xoroshiro128plus.c"><code>xoroshiro128+</code></a> have the same
|
||||
speed and use half of the space; the same comments apply. They are suitable only for
|
||||
low-scale parallel applications; moreover, <code>xoroshiro128+</code>
|
||||
exhibits a mild dependency in Hamming weights that generates a failure
|
||||
after 5 TB of output in <a href="hwd.php">our test</a>. We believe
|
||||
this slight bias cannot affect any application.
|
||||
|
||||
<p>Finally, if for any reason (which reason?) you need <strong>more
|
||||
state</strong>, we provide in the same
|
||||
vein <a href="xoshiro512plusplus.c"><code>xoshiro512++</code></a> / <a href="xoshiro512starstar.c"><code>xoshiro512**</code></a> / <a href="xoshiro512plus.c"><code>xoshiro512+</code></a> and
|
||||
<a href="xoroshiro1024plusplus.c"><code>xoroshiro1024++</code></a> / <a href="xoroshiro1024starstar.c"><code>xoroshiro1024**</code></a> / <a href="xoroshiro1024star.c"><code>xoroshiro1024*</code></a> (see the <a
|
||||
href="http://vigna.di.unimi.it/papers.php#BlVSLPNG">paper</a>).
|
||||
|
||||
<p>All generators, being based on linear recurrences, provide <em>jump
|
||||
functions</em> that make it possible to simulate any number of calls to
|
||||
the next-state function in constant time, once a suitable <em>jump
|
||||
polynomial</em> has been computed. We provide ready-made jump functions for
|
||||
a number of calls equal to the square root of the period, to make it easy
|
||||
generating non-overlapping sequences for parallel computations, and equal
|
||||
to the cube of the fourth root of the period, to make it possible to
|
||||
generate independent sequences on different parallel processors.
|
||||
|
||||
<p>We suggest to use <a href="splitmix64.c"><span
|
||||
style="font-variant: small-caps">SplitMix64</span></a> to initialize
|
||||
the state of our generators starting from a 64-bit seed, as <a href="https://dl.acm.org/citation.cfm?doid=1276927.1276928">research
|
||||
has shown</a> that initialization must be performed with a generator
|
||||
radically different in nature from the one initialized to avoid
|
||||
correlation on similar seeds.
|
||||
|
||||
|
||||
<h1>32-bit Generators</h1>
|
||||
|
||||
<P><a href="xoshiro128plusplus.c"><code>xoshiro128++</code></a>/<a href="xoshiro128starstar.c"><code>xoshiro128**</code></a> are our
|
||||
<strong>32-bit</strong> all-purpose generators, whereas <a
|
||||
href="xoshiro128plus.c"><code>xoshiro128+</code></a> is
|
||||
for floating-point generation. They are the 32-bit counterpart of
|
||||
<code>xoshiro256++</code>, <code>xoshiro256**</code> and <code>xoshiro256+</code>, so similar comments apply.
|
||||
Their state is too small for
|
||||
large-scale parallelism: their intended usage is inside embedded
|
||||
hardware or GPUs. For an even smaller scale, you can use <a
|
||||
href="xoroshiro64starstar.c"><code>xoroshiro64**</code></a> and <a
|
||||
href="xoroshiro64star.c"><code>xoroshiro64*</code></a>. We not believe
|
||||
at this point in time 32-bit generator with a larger state can be of
|
||||
any use (but there are 32-bit <code>xoroshiro</code> generators of much larger size).
|
||||
|
||||
<p>All 32-bit generators pass all tests we are aware of, with the
|
||||
exception of linearity tests (binary rank and linear complexity) for
|
||||
<code>xoshiro128+</code> and <code>xoroshiro64*</code>: in this case,
|
||||
due to the smaller number of output bits the low linear complexity of the
|
||||
lowest bits is sufficient to trigger BigCrush tests when the output is bit-reversed. Analogously to
|
||||
the 64-bit case, generating 32-bit floating-point number using the
|
||||
upper bits will not use any of the bits with <a href="lowcomp.php">low linear complexity</a>.
|
||||
|
||||
<h1>16-bit Generators</h1>
|
||||
|
||||
<p>We do not suggest any particular 16-bit generator, but it is possible
|
||||
to design relatively good ones using our techniques. For example,
|
||||
Parallax has embedded in their <a href="https://www.parallax.com/propeller-2/">Propeller 2 microcontroller</a> multiple 16-bit
|
||||
<code>xoroshiro32++</code> generators.
|
||||
|
||||
<h1>Congruential Generators</h1>
|
||||
|
||||
<p>In case you are interested in 64-bit PRNGs based on congruential arithmetic, I provide
|
||||
three instances of a
|
||||
<a href="https://groups.google.com/forum/#!searchin/sci.stat.math/Yet$20another$20rng%7Csort:date/sci.stat.math/p7aLW3TsJys/QGb1kti6kN0J">Marsaglia's Multiply-With-Carry generators</a>,
|
||||
<a href="MWC128.c"><code>MWC128</code></a>, <a href="MWC192.c"><code>MWC192</code></a>, and <a href="MWC256.c"><code>MWC256</code></a>, for which I computed good constants. They are some
|
||||
of the fastest generator available, but they need 128-bit operations.
|
||||
|
||||
<p>Stronger theoretical guarantees are provided by the
|
||||
<a href="https://www.math.ias.edu/~goresky/pdf/p1-goresky.pdf">generalized multiply-with-carry generators defined by Goresky and Klapper</a>:
|
||||
also in this case I provide two instances, <a href="GMWC128.c"><code>GMWC128</code></a> and <a href="GMWC256.c"><code>GMWC256</code></a>, for which I computed good constants.
|
||||
This generators, however, are about twice slower than MWC generators.
|
||||
|
||||
<h1>JavaScript</h1>
|
||||
|
||||
<p><code>xorshift128+</code> is presently used in the JavaScript engines of
|
||||
<a href="http://v8project.blogspot.com/2015/12/theres-mathrandom-and-then-theres.html">Chrome</a>,
|
||||
<a href="https://nodejs.org/">Node.js</a>,
|
||||
<a href="https://bugzilla.mozilla.org/show_bug.cgi?id=322529#c99">Firefox</a>,
|
||||
<a href="https://bugs.webkit.org/show_bug.cgi?id=151641">Safari</a> and
|
||||
<a href="https://github.com/Microsoft/ChakraCore/commit/dbda0182dc0a983dfb37d90c05000e79b6fc75b0">Microsoft Edge</a>.
|
||||
|
||||
<h1>Rust</h1>
|
||||
<p>The <a href="https://docs.rs/rand/latest/rand/rngs/struct.SmallRng.html">SmallRng</a> from the <a href="https://docs.rs/rand/latest/rand/">rand</a>
|
||||
crate is <a HREF="xoshiro256plusplus.c"><code>xoshiro256++</code></a> or <a HREF="xoshiro128plusplus.c"><code>xoshiro128++</code></a>, depending
|
||||
on the platform.
|
||||
|
||||
<h1><code>java.util.random</code></h1>
|
||||
|
||||
<p>I worked with Guy Steele at the <a
|
||||
href="https://docs.oracle.com/en/java/javase/17/docs/api/java.base/java/util/random/package-summary.html">new
|
||||
family of PRNGs available in Java 17</a>. The family, called <a
|
||||
href="http://vigna.di.unimi.it/papers.php#StVLXM">LXM</a>, uses <a
|
||||
href="http://vigna.di.unimi.it/papers.php#StVCESGMCPNG">new, better
|
||||
tables of multipliers for LCGs with power-of-two moduli</a>. Moreover,
|
||||
<code>java.util.random</code> contains ready-to-use implementations of
|
||||
<a HREF="xoroshiro128plusplus.c"><code>xoroshiro128++</code></a> and <a
|
||||
HREF="xoshiro256plusplus.c"><code>xoshiro256++</code></a>.
|
||||
|
||||
<h1>.NET</h1>
|
||||
|
||||
<p>In version 6, Microsoft's .NET framework <a href="https://devblogs.microsoft.com/dotnet/performance-improvements-in-net-6/">has adopted</a>
|
||||
<a HREF="xoshiro256starstar.c"><code>xoshiro256**</code></a> and <a
|
||||
HREF="xoshiro128starstar.c"><code>xoshiro128**</code></a> as default PRNGs.
|
||||
|
||||
<h1>Erlang</h1>
|
||||
|
||||
<p>The parallel functional language <a href="https://www.erlang.org/">Erlang</a> implements <a href="https://www.erlang.org/doc/man/rand.html">several
|
||||
variants of <code>xorshift</code>/<code>xoroshiro</code>-based generators</a> adapted in collaboration with Raimo Niskanen for Erlang's
|
||||
58/59-bit arithmetic.
|
||||
|
||||
<h1>GNU FORTRAN</h1>
|
||||
<p>GNU's <a href="https://gcc.gnu.org/fortran/">implementation of the FORTRAN language</a> <a href="https://gcc.gnu.org/onlinedocs/gfortran/RANDOM_005fNUMBER.html">uses</a>
|
||||
<a HREF="xoshiro256starstar.c"><code>xoshiro256**</code></a> as default PRNG.
|
||||
|
||||
<h1>Julia</h1>
|
||||
<p>The <a href="https://julialang.org/">Julia programming language</a> <a href="https://docs.julialang.org/en/v1/stdlib/Random/">uses</a>
|
||||
<a HREF="xoshiro256plusplus.c"><code>xoshiro256++</code></a> as default PRNG.
|
||||
|
||||
<h1>Lua</h1>
|
||||
<p>The scripting language <a href="http://www.lua.org/">Lua</a> <a href="https://www.lua.org/manual/5.4/manual.html#pdf-math.random">uses</a> <a HREF="xoshiro256starstar.c"><code>xoshiro256**</code></a> as default PRNG.
|
||||
|
||||
<h1>IoT</h1>
|
||||
|
||||
<p>The IoT operating systems <a href="https://os.mbed.com/">Mbed</a> and <a href="https://www.zephyrproject.org/">Zephyr</a> use
|
||||
<a HREF="xoroshiro128plus.c"><code>xoroshiro128+</code></a> as default PRNG.
|
||||
|
||||
<h1><a name=shootout></a>A PRNG Shootout</h1>
|
||||
|
||||
<p>I provide here a shootout of a few recent 64-bit PRNGs that are quite widely used.
|
||||
The purpose is that of providing a consistent, reproducible assessment of two properties of the generators: speed and quality.
|
||||
The code used to perform the tests and all the output from statistical test suites is available for download.
|
||||
|
||||
<h2><a name=speed></a>Speed</h2>
|
||||
|
||||
<p>The speed reported in this page is the time required to emit 64
|
||||
random bits, and the number of clock cycles required to generate a byte (thanks to the <a href="http://icl.utk.edu/papi/">PAPI</a> library). If a generator is 32-bit in nature, I glue two
|
||||
consecutive outputs. Note that
|
||||
I do not report results using GPUs or SSE instructions, with an exception for the very common SFMT: for that to be
|
||||
meaningful, I should have implementations for all generators.
|
||||
Otherwise, with suitable hardware support I could just use AES in
|
||||
counter mode and get 64 secure bits in 0.56 ns (or just use <a href="https://github.com/google/randen">Randen</a>). The tests were performed on a
|
||||
12th Gen Intel® Core™ i7-12700KF @3.60GHz using <code>gcc</code> 12.2.1.
|
||||
|
||||
<p>A few <i>caveats</i>:
|
||||
<ul>
|
||||
<li>There is some looping overhead, but subtracting it from the timings is not going to
|
||||
be particularly meaningful due to instruction rescheduling, etc.
|
||||
<li>Relative speed might be different on different CPUs and on different scenarios.
|
||||
<li>I do not use <code>-march=native</code>, which can improve the timing of some generators
|
||||
by vectorization or special instructions, because those improvements might not be possible
|
||||
when the generator is embedded in user code.
|
||||
<li>Code has been compiled using <code>gcc</code>'s <code>-fno-unroll-loops</code>
|
||||
option. This options is essential to get a sensible result: without it, the compiler
|
||||
may perform different loop unrolling depending on the generator. Previosuly I was using also
|
||||
<code>-fno-move-loop-invariants</code>, which was essential in not giving generators using several
|
||||
large constants an advantage by preventing the compiler from loading them into registers. However,
|
||||
as of <code>gcc</code> 12.2.1 the compiler loads the constants into registers anyway, so the
|
||||
option is no longer used. Timings
|
||||
with <a href="http://clang.llvm.org/"><code>clang</code></a> at the time of this writing
|
||||
are very close to those obtained with <code>gcc</code>.
|
||||
If you find timings that are significantly better than those shown here on
|
||||
comparable hardware, they are likely to be due to compiler artifacts (e.g., vectorization).
|
||||
<li>Timings are taken running a generator for billions of times in a loop; but this is not the way you use generators. Register
|
||||
allocation might be very different when the generator is embedded in an application, leading to constants being reloaded
|
||||
or part of the state space being written to main memory at each iteration. These costs do not appear in the benchmarks below.
|
||||
</ul>
|
||||
|
||||
<p>To ease replicability, I distribute a <a href="harness.c"><em>harness</em></a> performing the measurement. You just
|
||||
have to define a <a href="xoroshiro128plus-speed.c"><code>next()</code></a> function and include the harness. But the only realistic
|
||||
suggestion is to try different generators in your application and see what happens.
|
||||
|
||||
<h2><a name=quality></a>Quality</h2>
|
||||
|
||||
<p>This is probably the more <a
|
||||
href="http://dilbert.com/strips/comic/2001-10-25/">elusive</a> property
|
||||
of a PRNG. Here quality is measured using the powerful
|
||||
BigCrush suite of tests. BigCrush is part of <a
|
||||
href="http://simul.iro.umontreal.ca/testu01/tu01.html">TestU01</a>,
|
||||
a monumental framework for testing PRNGs developed by Pierre L'Ecuyer
|
||||
and Richard Simard (“TestU01: A C library for empirical testing
|
||||
of random number generators”, <i>ACM Trans. Math. Softw.</i>
|
||||
33(4), Article 22, 2007).
|
||||
|
||||
<p>I run BigCrush starting from 100 equispaced points of the state space
|
||||
of the generator and collect <em>failures</em>—tests in which the
|
||||
<i>p</i>-value statistics is outside the interval [0.001..0.999]. A failure
|
||||
is <em>systematic</em> if it happens at all points.
|
||||
|
||||
<p>Note that TestU01 is a 32-bit test suite. Thus, two 32-bit integer values
|
||||
are passed to the test suite for each generated 64-bit value. Floating point numbers
|
||||
are generated instead by dividing the unsigned output of the generator by 2<sup>64</sup>.
|
||||
Since this implies a bias towards the high bits (which is anyway a known characteristic
|
||||
of TestU01), I run the test suite also on the <em>reverse</em>
|
||||
generator. More detail about the whole process can be found in this <a
|
||||
href="http://vigna.di.unimi.it/papers.php#VigEEMXGS">paper</a>.
|
||||
|
||||
<p>Beside BigCrush, I analyzed generators using a test for <a href="hwd.php">Hamming-weight dependencies</a>
|
||||
described in our <a href="http://vigna.di.unimi.it/papers.php#BlVNTHWD">paper</a>. As I already remarked, our only
|
||||
generator failing the test (but only after 5 TB of output) is <code>xoroshiro128+</code>.
|
||||
|
||||
<p>I report the period of each generator and its footprint in bits: a generator gives “bang-for-the-buck”
|
||||
if the base-2 logarithm of the period is close to the footprint. Note
|
||||
that the footprint has been always padded to a multiple of 64, and it can
|
||||
be significantly larger than expected because of padding and
|
||||
cyclic access indices.
|
||||
|
||||
<div style="align: center"><table id='prng' style='margin: 2em 0' class='tablesorter'>
|
||||
<thead><tr>
|
||||
<th>PRNG
|
||||
<th>Footprint (bits)
|
||||
<th class="{ sorter: 'metadata' }">Period
|
||||
<th><a href="http://simul.iro.umontreal.ca/testu01/tu01.html">BigCrush</a> Systematic Failures
|
||||
<th><a href="http://prng.di.unimi.it/hwd.php">HWD failure</a>
|
||||
<th>ns/64 bits
|
||||
<th>cycles/B
|
||||
<tbody>
|
||||
<tr><td><a href="xoroshiro128plus.c"><code>xoroshiro128+</code></a><td align=right>128 <td align=right class='{sortValue: 128}'>2<sup>128</sup> − 1<td align=right>—<td align=right>5 TB<td align=right>0.80<td align=right>0.36
|
||||
<tr><td><a href="xoroshiro128plusplus.c"><code>xoroshiro128++</code></a><td align=right>128 <td align=right class='{sortValue: 128}'>2<sup>128</sup> − 1<td align=right>—<td align=right>—<td align=right>0.90<td align=right>0.40
|
||||
<tr><td><a href="xoroshiro128starstar.c"><code>xoroshiro128**</code></a><td align=right>128 <td align=right class='{sortValue: 128}'>2<sup>128</sup> − 1<td align=right>—<td align=right>—<td align=right>0.78<td align=right>0.36
|
||||
<tr><td><a href="xoshiro256plus.c"><code>xoshiro256+</code></a><td align=right>256 <td align=right class='{sortValue: 256}'>2<sup>256</sup> − 1<td align=right>—<td align=right>—<td align=right>0.61<td align=right>0.27
|
||||
<tr><td><a href="xoshiro256plusplus.c"><code>xoshiro256++</code></a><td align=right>256 <td align=right class='{sortValue: 256}'>2<sup>256</sup> − 1<td align=right>—<td align=right>—<td align=right>0.75<td align=right>0.34
|
||||
<tr><td><a href="xoshiro256starstar.c"><code>xoshiro256**</code></a><td align=right>256 <td align=right class='{sortValue: 256}'>2<sup>256</sup> − 1<td align=right>—<td align=right>—<td align=right>0.75<td align=right>0.34
|
||||
<tr><td><a href="xoshiro512plus.c"><code>xoshiro512+</code></a><td align=right>512<td align=right class='{sortValue: 512}'>2<sup>512</sup> − 1<td align=right>—<td align=right>—<td align=right>0.68<td align=right>0.30
|
||||
<tr><td><a href="xoshiro512plusplus.c"><code>xoshiro512++</code></a><td align=right>512<td align=right class='{sortValue: 512}'>2<sup>512</sup> − 1<td align=right>—<td align=right>—<td align=right>0.79<td align=right>0.36
|
||||
<tr><td><a href="xoshiro512starstar.c"><code>xoshiro512**</code></a><td align=right>512<td align=right class='{sortValue: 512}'>2<sup>512</sup> − 1<td align=right>—<td align=right>—<td align=right>0.81<td align=right>0.37
|
||||
<tr><td><a href="xoroshiro1024star.c"><code>xoroshiro1024*</code></a><td align=right>1068<td align=right class='{sortValue: 1024}'>2<sup>1024</sup> − 1<td align=right>—<td align=right>—<td align=right>0.82<td align=right>0.37
|
||||
<tr><td><a href="xoroshiro1024plusplus.c"><code>xoroshiro1024++</code></a><td align=right>1068<td align=right class='{sortValue: 1024}'>2<sup>1024</sup> − 1<td align=right>—<td align=right>—<td align=right>1.01<td align=right>0.46
|
||||
<tr><td><a href="xoroshiro1024starstar.c"><code>xoroshiro1024**</code></a><td align=right>1068<td align=right class='{sortValue: 1024}'>2<sup>1024</sup> − 1<td align=right>—<td align=right>—<td align=right>0.98<td align=right>0.44
|
||||
<tr><td><a href="MWC128.c"><span style="font-variant: small-caps">MWC128</span></a><td align=right>128 <td align=right class='{sortValue: 127}'>≈2<sup>127</sup><td align=right>—<td align=right>—<td align=right>0.83<td align=right>0.37
|
||||
<tr><td><a href="MWC192.c"><span style="font-variant: small-caps">MWC192</span></a><td align=right>192 <td align=right class='{sortValue: 127}'>≈2<sup>191</sup><td align=right>—<td align=right>—<td align=right>1.42<td align=right>0.19
|
||||
<tr><td><a href="MWC256.c"><span style="font-variant: small-caps">MWC256</span></a><td align=right>256 <td align=right class='{sortValue: 255}'>≈2<sup>255</sup><td align=right>—<td align=right>—<td align=right>0.45<td align=right>0.20
|
||||
<tr><td><a href="GMWC128.c"><span style="font-variant: small-caps">GMWC128</span></a><td align=right>128 <td align=right class='{sortValue: 127}'>≈2<sup>127</sup><td align=right>—<td align=right>—<td align=right>1.84<td align=right>0.83
|
||||
<tr><td><a href="GMWC256.c"><span style="font-variant: small-caps">GMWC256</span></a><td align=right>256 <td align=right class='{sortValue: 255}'>≈2<sup>255</sup><td align=right>—<td align=right>—<td align=right>1.85<td align=right>0.83
|
||||
<tr><td><a href="http://pracrand.sourceforge.net/"><span style="font-variant: small-caps">SFC64</span></a><td align=right>256 <td align=right class='{sortValue: 64}'>≥2<sup>64</sup><td align=right>—<td align=right>—<td align=right>0.66<td align=right>0.30
|
||||
<tr><td><a href="splitmix64.c"><span style="font-variant: small-caps">SplitMix64</span></a><td align=right>64 <td align=right class='{sortValue: 64}'>2<sup>64</sup><td align=right>—<td align=right>—<td align=right>0.63<td align=right>0.29
|
||||
<tr><td><a href="http://pcg-random.org/">PCG 128 XSH RS 64 (LCG)</a> <td align=right>128 <td align=right class='{sortValue: 128}'>2<sup>128</sup><td align=right>—<td align=right>—<td align=right>1.70<td align=right>0.77
|
||||
<tr><td><a href="https://github.com/numpy/numpy">PCG64-DXSM (NumPy)</a> <td align=right>128 <td align=right class='{sortValue: 128}'>2<sup>128</sup><td align=right>—<td align=right>—<td align=right>1.11<td align=right>0.50
|
||||
<tr><td><a href="http://numerical.recipes/"><code>Ran</code></a><td align=right>192 <td align=right class='{sortValue: 191}'>≈2<sup>191</sup><td align=right>—<td align=right>—<td align=right>1.37<td align=right>0.62
|
||||
<tr><td><a href="http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/emt.html"><code>MT19937-64</code> (Mersenne Twister)</a><td align=right>20032 <td align=right class='{sortValue: 19937}'>2<sup>19937</sup> − 1<td align=right>LinearComp<td align=right>—<td align=right>1.36<td align=right>0.62
|
||||
<tr><td><a href="http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/SFMT/"><code>SFMT19937 (uses SSE2 instructions)</code></a><td align=right>20032 <td align=right class='{sortValue: 19937}'>2<sup>19937</sup> − 1<td align=right>LinearComp<td align=right>—<td align=right>0.93<td align=right>0.42
|
||||
<tr><td><a href="http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/SFMT/"><code>SFMT607 (uses SSE2 instructions)</code></a><td align=right>672 <td align=right class='{sortValue: 607}'>2<sup>607</sup> − 1<td align=right>MatrixRank, LinearComp<td align=right>400 MB<td align=right>0.78<td align=right>0.34
|
||||
<tr><td><a href="http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/TINYMT/index.html">Tiny Mersenne Twister</a> (64 bits)<td align=right>256<td align=right class='{sortValue: 127}'>2<sup>127</sup> − 1<td align=right>—<td align=right>90 TB→<td align=right>2.76<td align=right>1.25
|
||||
<tr><td><a href="http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/TINYMT/index.html">Tiny Mersenne Twister</a> (32 bits)<td align=right>224<td align=right class='{sortValue: 127}'>2<sup>127</sup> − 1<td align=right>CollisionOver, Run, SimPoker, AppearanceSpacings, MatrixRank, LinearComp, LongestHeadRun, Run of Bits (reversed)<td align=right>40 TB→<td align=right>4.27<td align=right>1.92
|
||||
<tr><td><a href="http://www.iro.umontreal.ca/~panneton/WELLRNG.html"><code>WELL512a</code></a><td align=right>544 <td align=right class='{sortValue: 512}'>2<sup>512</sup> − 1 <td align=right>MatrixRank, LinearComp<td align=right>3.5 PB<td align=right>5.42<td align=right>2.44
|
||||
<tr><td><a href="http://www.iro.umontreal.ca/~panneton/WELLRNG.html"><code>WELL1024a</code></a><td align=right>1056 <td align=right class='{sortValue: 1024}'>2<sup>1024</sup> − 1 <td align=right>MatrixRank, LinearComp<td align=right>—<td align=right>5.30<td align=right>2.38
|
||||
</table></div>
|
||||
|
||||
<p>The following table compares instead two ways of generating floating-point numbers, namely the 521-bit <a href="http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/SFMT/">dSFMT</a>, which
|
||||
generates directly floating-point numbers with 52 significant bits, and
|
||||
<a href="xoshiro256plus.c"><code>xoshiro256+</code></a> followed by a standard conversion of its upper bits to a floating-point number with 53 significant bits (see below).
|
||||
|
||||
<div style="align: center"><table id='prngf' style='margin: 2em 0' class='tablesorter'>
|
||||
<thead><tr>
|
||||
<th>PRNG
|
||||
<th>Footprint (bits)
|
||||
<th class="{ sorter: 'metadata' }">Period
|
||||
<th> <a href="http://simul.iro.umontreal.ca/testu01/tu01.html">BigCrush</a> Systematic Failures
|
||||
<th><a href="http://prng.di.unimi.it/hwd.php">HWD failure</a>
|
||||
<th>ns/double
|
||||
<th>cycles/B
|
||||
<tbody>
|
||||
<tr><td><a href="xoshiro256plus.c"><code>xoshiro256+</code></a> (returns 53 significant bits) <td align=right>256<td align=right class='{sortValue: 256}'>2<sup>256</sup> − 1<td align=right>—<td align=right>—<td align=right>0.92<td align=right>3.40
|
||||
<tr><td><a href="http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/SFMT/"><code>dSFMT</code></a> (uses SSE2 instructions, returns only 52 significant bits)<td align=right>704<td align=right class='{sortValue: 521}'>2<sup>521</sup> − 1<td align=right>MatrixRank, LinearComp<td align=right>6 TB<td align=right>0.85<td align=right>3.07
|
||||
</table></div>
|
||||
|
||||
<p><code>xoshiro256+</code> is ≈8% slower than the dSFMT, but it has a doubled range of output values, does not need any extra SSE instruction (can be programmed in Java, etc.),
|
||||
has a much smaller footprint, and its upper bits do not fail any test.
|
||||
|
||||
<h1><a name=remarks></a>Remarks</h1>
|
||||
|
||||
<h2>Vectorization</h2>
|
||||
|
||||
<p>Some of the generators can be very easily vectorized, so that multiple instances can be run in parallel to provide
|
||||
fast bulk generation. Thanks to an interesting <a href="https://github.com/JuliaLang/julia/issues/27614">discussion with the Julia developers</a>,
|
||||
I've become aware that AVX2 vectorizations of multiple instances of generators using the <code>+</code>/<code>++</code> scrambler are impressively fast (links
|
||||
below point at a speed test to be used with the <a href="harness.c">harness</a>, and the result will be multiplied by 1000):
|
||||
|
||||
<div style="align: center"><table id='vec' style='margin: 2em 0' class='tablesorter'>
|
||||
<thead><tr>
|
||||
<th>PRNG
|
||||
<th>ns/64 bits
|
||||
<th>cycles/B
|
||||
<tbody>
|
||||
<tr><td><a href="xoroshiro128+-vect-speed.c"><code>xoroshiro128+</code></a> (4 parallel instances)<td align=right>0.36<td align=right>0.14
|
||||
<tr><td><a href="xoroshiro128++-vect-speed.c"><code>xoroshiro128++</code></a> (4 parallel instances)<td align=right>0.45<td align=right>0.18
|
||||
<tr><td><a href="xoshiro256+-vect-speed.c"><code>xoshiro256+</code></a> (8 parallel instances)<td align=right>0.19<td align=right>0.08
|
||||
<tr><td><a href="xoshiro256++-vect-speed.c"><code>xoshiro256++</code></a> (8 parallel instances)<td align=right>0.26<td align=right>0.09
|
||||
</table></div>
|
||||
|
||||
<p>Note that sometimes convincing the compiler to vectorize is a
|
||||
slightly quirky process: for example, on <code>gcc</code> 12.2.1 I have to use <code>-O3 -fdisable-tree-cunrolli -march=native</code>
|
||||
to vectorize <code>xoshiro256</code>-based generators
|
||||
(<code>-O3</code> alone will not vectorize; thanks to to Chris Elrod for pointing me at <code>-fdisable-tree-cunrolli</code>).
|
||||
|
||||
<h2>A long period does not imply high quality</h2>
|
||||
|
||||
<p>This is a common misconception. The generator <code>x++</code> has
|
||||
period \(2^k\), for any \(k\geq0\), provided that <code>x</code> is
|
||||
represented using \(k\) bits: nonetheless, it is a horrible generator.
|
||||
The generator returning \(k-1\) zeroes followed by a one has period
|
||||
\(k\).
|
||||
|
||||
<p>It is however important that the period is long enough. A first heuristic rule of thumb
|
||||
is that if you need to use \(t\) values, you need a generator with period at least \(t^2\).
|
||||
Moreover, if you run \(n\) independent computations starting at random seeds,
|
||||
the sequences used by each computation should not overlap.
|
||||
|
||||
<p>Now, given a generator with period \(P\), the probability that \(n\) subsequences of length \(L\) starting at random points in the state space
|
||||
overlap <a href="http://vigna.di.unimi.it/papers.php#VigPORSPNG">is bounded by \(n^2L/P\)</a>. If your generator has period \(2^{256}\) and you run
|
||||
on \(2^{64}\) cores (you will never have them) a computation using \(2^{64}\) pseudorandom numbers (you will never have the time)
|
||||
the probability of overlap would be less than \(2^{-64}\).
|
||||
|
||||
<p>In other words: any generator with a period beyond
|
||||
\(2^{256}\) has a period that is
|
||||
sufficient for every imaginable application. Unless there are other motivations (e.g., provably
|
||||
increased quality), a generator with a larger period is only a waste of
|
||||
memory (as it needs a larger state), of cache lines, and of
|
||||
precious high-entropy random bits for seeding (unless you're using
|
||||
small seeds, but then it's not clear why you would want a very long
|
||||
period in the first place—the computation above is valid only if you seed all bits of the state
|
||||
with independent, uniformly distributed random bits).
|
||||
|
||||
<p>In case the generator provides a <em>jump function</em> that lets you skip through chunks of the output in constant
|
||||
time, even a period of \(2^{128}\) can be sufficient, as it provides \(2^{64}\) non-overlapping sequences of length \(2^{64}\).
|
||||
|
||||
<h2>Equidistribution</h2>
|
||||
|
||||
<p>Every 64-bit generator of ours with <var>n</var> bits of state scrambled
|
||||
with <code>*</code> or <code>**</code> is <var>n</var>/64-dimensionally
|
||||
equidistributed: every <var>n</var>/64-tuple of consecutive 64-bit
|
||||
values appears exactly once in the output, except for the zero tuple
|
||||
(and this is the largest possible dimension). Generators based on the
|
||||
<code>+</code> or <code>++</code> scramblers are however only (<var>n</var>/64 −
|
||||
1)-dimensionally equidistributed: every (<var>n</var>/64 −
|
||||
1)-tuple of consecutive 64-bit values appears exactly 2<sup>64</sup>
|
||||
times in the output, except for a missing zero tuple. The same considerations
|
||||
apply to 32-bit generators.
|
||||
|
||||
<h2>Generating uniform doubles in the unit interval</h2>
|
||||
|
||||
<p>A standard double (64-bit) floating-point number in
|
||||
<a href="https://en.wikipedia.org/wiki/IEEE_floating_point">IEEE floating point format</a> has 52 bits of
|
||||
significand, plus an implicit bit at the left of the significand. Thus,
|
||||
the representation can actually store numbers with <em>53</em> significant binary digits.
|
||||
|
||||
<p>Because of this fact, in C99 a 64-bit unsigned integer <code>x</code> should be converted to a 64-bit double
|
||||
using the expression
|
||||
<pre>
|
||||
#include <stdint.h>
|
||||
|
||||
(x >> 11) * 0x1.0p-53
|
||||
</pre>
|
||||
<p>In Java you can use almost the same expression for a (signed) 64-bit integer:
|
||||
<pre>
|
||||
(x >>> 11) * 0x1.0p-53
|
||||
</pre>
|
||||
|
||||
|
||||
<p>This conversion guarantees that all dyadic rationals of the form <var>k</var> / 2<sup>−53</sup>
|
||||
will be equally likely. Note that this conversion prefers the high bits of <code>x</code> (usually, a good idea), but you can alternatively
|
||||
use the lowest bits.
|
||||
|
||||
<p>An alternative, multiplication-free conversion is
|
||||
<pre>
|
||||
#include <stdint.h>
|
||||
|
||||
static inline double to_double(uint64_t x) {
|
||||
const union { uint64_t i; double d; } u = { .i = UINT64_C(0x3FF) << 52 | x >> 12 };
|
||||
return u.d - 1.0;
|
||||
}
|
||||
</pre>
|
||||
<p>The code above cooks up by bit manipulation
|
||||
a real number in the interval [1..2), and then subtracts
|
||||
one to obtain a real number in the interval [0..1). If <code>x</code> is chosen uniformly among 64-bit integers,
|
||||
<code>d</code> is chosen uniformly among dyadic rationals of the form <var>k</var> / 2<sup>−52</sup>. This
|
||||
is the same technique used by generators providing directly doubles, such as the
|
||||
<a href="http://dx.doi.org/10.1007/978-3-540-85912-3_26">dSFMT</a>.
|
||||
|
||||
<p>This technique is supposed to be fast, but on recent hardare it does not seem to give a significant advantage.
|
||||
More importantly, <em>you will be generating half the values you could actually generate</em>.
|
||||
The same problem plagues the dSFMT. All doubles generated will have the lowest significand bit set to zero (I must
|
||||
thank Raimo Niskanen from the Erlang team for making me notice this—a previous version of this site
|
||||
did not mention this issue).
|
||||
|
||||
<p>In Java you can obtain an analogous result using suitable static methods:
|
||||
<pre>
|
||||
Double.longBitsToDouble(0x3FFL << 52 | x >>> 12) - 1.0
|
||||
</pre>
|
||||
|
||||
<p>To adhere to the principle of least surprise, my implementations now use the multiplicative version, everywhere.
|
||||
|
||||
<p>Interestingly, these are not the only notions of “uniformity” you can come up with. Another possibility
|
||||
is that of generating 1074-bit integers, normalize and return the nearest value representable as a
|
||||
64-bit double (this is the theory—in practice, you will almost never
|
||||
use more than two integers per double as the remaining bits would not be representable). This approach guarantees that all
|
||||
representable doubles could be in principle generated, albeit not every
|
||||
returned double will appear with the same probability. A reference
|
||||
implementation can be found <a href="random_real.c">here</a>. Note that unless your generator has
|
||||
at least 1074 bits of state and suitable equidistribution properties, the code above will not do what you expect
|
||||
(e.g., it might <em>never</em> return zero).
|
||||
|
||||
|
||||
</div>
|
||||
|
||||
|
||||
</div>
|
||||
|
||||
<div id="right">
|
||||
|
||||
<!-- <h1>Download</h1>
|
||||
<p><ul>
|
||||
<li><a HREF="prng-1.2.tgz">source tarball</A>
|
||||
<li><a HREF="prng-data.tar.bz2">data tarball (large!)</a>
|
||||
</ul>
|
||||
-->
|
||||
<h1>C code (64 bits)</h1>
|
||||
<p><ul>
|
||||
<li><a HREF="xoshiro256plusplus.c"><code>xoshiro256++</code></a>
|
||||
<li><a HREF="xoshiro256starstar.c"><code>xoshiro256**</code></a>
|
||||
<li><a HREF="xoshiro256plus.c"><code>xoshiro256+</code></a>
|
||||
<li><a HREF="xoroshiro128plusplus.c"><code>xoroshiro128++</code></a>
|
||||
<li><a HREF="xoroshiro128starstar.c"><code>xoroshiro128**</code></a>
|
||||
<li><a HREF="xoroshiro128plus.c"><code>xoroshiro128+</code></a>
|
||||
<li><a HREF="https://github.com/vigna/MRG32k3a">Testless <code>MRG32k3a</code></a>
|
||||
<li><a HREF="MWC128.c"><code>MWC128</code></a> + <a HREF="mp.c"><code>mp.c</code></a>
|
||||
<li><a HREF="MWC192.c"><code>MWC192</code></a> + <a HREF="mp.c"><code>mp.c</code></a>
|
||||
<li><a HREF="MWC256.c"><code>MWC256</code></a> + <a HREF="mp.c"><code>mp.c</code></a>
|
||||
<li><a HREF="GMWC128.c"><code>GMWC128</code></a> + <a HREF="mp.c"><code>mp.c</code></a>
|
||||
<li><a HREF="GMWC256.c"><code>GMWC256</code></a> + <a HREF="mp.c"><code>mp.c</code></a>
|
||||
<!-- <li><a HREF="xoshiro512starstar.c"><code>xoshiro512**</code></a>
|
||||
<li><a HREF="xoshiro512plus.c"><code>xoshiro512+</code></a>
|
||||
<li><a HREF="xoroshiro1024starstar.c"><code>xoroshiro1024**</code></a>
|
||||
<li><a HREF="xoroshiro1024plus.c"><code>xoroshiro1024*</code></a>-->
|
||||
</ul>
|
||||
|
||||
<h1>C code (32 bits)</h1>
|
||||
<p><ul>
|
||||
<li><a HREF="xoshiro128plusplus.c"><code>xoshiro128++</code></a>
|
||||
<li><a HREF="xoshiro128starstar.c"><code>xoshiro128**</code></a>
|
||||
<li><a HREF="xoshiro128plus.c"><code>xoshiro128+</code></a>
|
||||
<li><a HREF="xoroshiro64starstar.c"><code>xoroshiro64**</code></a>
|
||||
<li><a HREF="xoroshiro64star.c"><code>xoroshiro64*</code></a>
|
||||
</ul>
|
||||
|
||||
<h1>Java code (<a HREF="https://github.com/openjdk/jdk17/tree/master/src/jdk.random/share/classes/jdk/random"><code>java.util.random</code></a>)</h1>
|
||||
|
||||
<h1>Java code (<a href="http://dsiutils.di.unimi.it">DSI utilities</a>)</h1>
|
||||
<p><ul>
|
||||
<li><a HREF="http://dsiutils.di.unimi.it/docs/it/unimi/dsi/util/package-summary.html">Overview</a>
|
||||
<li><a HREF="http://dsiutils.di.unimi.it/docs/it/unimi/dsi/util/XoShiRo256PlusPlusRandom.html"><code>xoshiro256++</code></a>
|
||||
<li><a HREF="http://dsiutils.di.unimi.it/docs/it/unimi/dsi/util/XoShiRo256StarStarRandom.html"><code>xoshiro256**</code></a>
|
||||
<li><a HREF="http://dsiutils.di.unimi.it/docs/it/unimi/dsi/util/XoShiRo256PlusRandom.html"><code>xoshiro256+</code></a>
|
||||
<li><a HREF="http://dsiutils.di.unimi.it/docs/it/unimi/dsi/util/XoRoShiRo128PlusPlusRandom.html"><code>xoroshiro128++</code></a>
|
||||
<li><a HREF="http://dsiutils.di.unimi.it/docs/it/unimi/dsi/util/XoRoShiRo128StarStarRandom.html"><code>xoroshiro128**</code></a>
|
||||
<li><a HREF="http://dsiutils.di.unimi.it/docs/it/unimi/dsi/util/XoRoShiRo128PlusRandom.html"><code>xoroshiro128+</code></a>
|
||||
<li><a HREF="https://github.com/vigna/MRG32k3a">Testless <code>MRG32k3a</code></a>
|
||||
</ul>
|
||||
|
||||
<h1>Java code (<a HREF="https://gitbox.apache.org/repos/asf?p=commons-rng.git">Apache Commons RNG implementations</a>)</h1>
|
||||
|
||||
<h1>Documentation</h1>
|
||||
<p><ul>
|
||||
<li>The <a href="http://vigna.di.unimi.it/papers.php#BlVSLPNG">paper</a> introducing <code>xoshiro</code>/<code>xoroshiro</code>.
|
||||
<li>The <a href="http://vigna.di.unimi.it/papers.php#BlVNTHWD">paper</a> describing our <a href="hwd.php">test for Hamming-weight dependencies</a>.
|
||||
<li>A <a href="http://vigna.di.unimi.it/papers.php#VigHTLGMT">paper</a> discussing the defects of the Mersenne Twister family of PRNGs.
|
||||
<li>A <a href="http://vigna.di.unimi.it/papers.php#VigPORSPNG">paper</a> discussing the probability of overlap of random subsequences.
|
||||
<li>A <a href="http://vigna.di.unimi.it/papers.php#StVCESGMCPNG">paper</a> with new tables of multipliers for LCGs with power-of-two moduli.
|
||||
<li>A <a href="http://vigna.di.unimi.it/papers.php#StVLXM">paper</a> presenting the family LXM of PRNGs.
|
||||
</ul>
|
||||
|
||||
<h1>Discussion</h1>
|
||||
|
||||
<p>There is a <a href="http://groups.google.com/group/prng">discussion group</a>
|
||||
about this page. You can join or <a href="mailto:prng@googlegroups.com">send a message</a>.
|
||||
<h1><a href="https://validator.w3.org/check/referer">This is valid HTML 4.01</a></h1>
|
||||
|
||||
</div>
|
||||
</body>
|
||||
</html>
|
||||
File diff suppressed because one or more lines are too long
@ -0,0 +1,577 @@
|
||||
<!DOCTYPE html>
|
||||
<html prefix="og: http://ogp.me/ns# article: http://ogp.me/ns/article# " lang="en">
|
||||
<head>
|
||||
<meta charset="utf-8">
|
||||
<meta name="viewport" content="width=device-width, initial-scale=1">
|
||||
<title>On Vigna's PCG Critique | PCG, A Better Random Number Generator</title>
|
||||
<link href="../assets/css/bootstrap.min.css" rel="stylesheet" type="text/css">
|
||||
<link href="../assets/css/html4css1.css" rel="stylesheet" type="text/css">
|
||||
<link href="../assets/css/nikola_rst.css" rel="stylesheet" type="text/css">
|
||||
<link href="../assets/css/code.css" rel="stylesheet" type="text/css">
|
||||
<link href="../assets/css/colorbox.css" rel="stylesheet" type="text/css">
|
||||
<link href="../assets/css/theme.css" rel="stylesheet" type="text/css">
|
||||
<meta name="theme-color" content="#5670d4">
|
||||
<meta name="generator" content="Nikola (getnikola.com)">
|
||||
<link rel="alternate" type="application/rss+xml" title="RSS" href="../rss.xml">
|
||||
<link rel="canonical" href="http://www.pcg-random.org/posts/on-vignas-pcg-critique.html">
|
||||
<!--[if lt IE 9]><script src="../assets/js/html5.js"></script><![endif]--><meta name="author" content="M.E. O'Neill">
|
||||
<link rel="prev" href="implausible-output-from-xoshiro256.html" title="Implausible Output from Xoshiro256**" type="text/html">
|
||||
<link rel="next" href="bob-jenkins-small-prng-passes-practrand.html" title="Bob Jenkins's Small PRNG Passes PractRand (And More!)" type="text/html">
|
||||
<meta property="og:site_name" content="PCG, A Better Random Number Generator">
|
||||
<meta property="og:title" content="On Vigna's PCG Critique">
|
||||
<meta property="og:url" content="http://www.pcg-random.org/posts/on-vignas-pcg-critique.html">
|
||||
<meta property="og:description" content="On 14 May 2018, Sebastiano Vigna added a page to his website (archived here) entitled “The wrap-up on PCG generators” that attempts to persuade readers to avoid various PCG generators.
|
||||
That day, he al">
|
||||
<meta property="og:type" content="article">
|
||||
<meta property="article:published_time" content="2018-05-25T16:49:25-07:00">
|
||||
<meta property="article:tag" content="pcg">
|
||||
<meta property="article:tag" content="practrand">
|
||||
<meta property="article:tag" content="splitmix">
|
||||
<meta property="article:tag" content="testing">
|
||||
<meta property="article:tag" content="xoroshiro">
|
||||
</head>
|
||||
<body>
|
||||
<a href="#content" class="sr-only sr-only-focusable">Skip to main content</a>
|
||||
|
||||
<!-- Menubar -->
|
||||
|
||||
<nav class="navbar navbar-inverse navbar-static-top"><div class="container">
|
||||
<!-- This keeps the margins nice -->
|
||||
<div class="navbar-header">
|
||||
<button type="button" class="navbar-toggle collapsed" data-toggle="collapse" data-target="#bs-navbar" aria-controls="bs-navbar" aria-expanded="false">
|
||||
<span class="sr-only">Toggle navigation</span>
|
||||
<span class="icon-bar"></span>
|
||||
<span class="icon-bar"></span>
|
||||
<span class="icon-bar"></span>
|
||||
</button>
|
||||
<a class="navbar-brand" href="http://www.pcg-random.org/">
|
||||
|
||||
<span id="blog-title">PCG, A Better Random Number Generator</span>
|
||||
</a>
|
||||
</div>
|
||||
<!-- /.navbar-header -->
|
||||
<div class="collapse navbar-collapse" id="bs-navbar" aria-expanded="false">
|
||||
<ul class="nav navbar-nav">
|
||||
<li>
|
||||
<a href="../download.html">Download</a>
|
||||
</li>
|
||||
<li>
|
||||
<a href="../using-pcg.html">Docs</a>
|
||||
</li>
|
||||
<li>
|
||||
<a href="../paper.html">Paper</a>
|
||||
</li>
|
||||
<li>
|
||||
<a href="stanford-colloquium-talk.html">Video</a>
|
||||
</li>
|
||||
<li>
|
||||
<a href="../blog/">Blog</a>
|
||||
|
||||
|
||||
</li>
|
||||
</ul>
|
||||
<ul class="nav navbar-nav navbar-right"></ul>
|
||||
</div>
|
||||
<!-- /.navbar-collapse -->
|
||||
</div>
|
||||
<!-- /.container -->
|
||||
</nav><!-- End of Menubar --><div class="container" id="content" role="main">
|
||||
<div class="body-content">
|
||||
<!--Body content-->
|
||||
<div class="row">
|
||||
|
||||
|
||||
<article class="post-text h-entry hentry postpage col-md-8" itemscope="itemscope" itemtype="http://schema.org/Article"><header><h1 class="p-name entry-title" itemprop="headline name"><a href="#" class="u-url">On Vigna's PCG Critique</a></h1>
|
||||
|
||||
<div class="metadata">
|
||||
<p class="byline author vcard"><span class="byline-name fn" itemprop="author">
|
||||
M.E. O'Neill
|
||||
</span></p>
|
||||
<p class="dateline"><a href="#" rel="bookmark"><time class="published dt-published" datetime="2018-05-25T16:49:25-07:00" itemprop="datePublished" title="2018-05-25 16:49">2018-05-25 16:49</time></a></p>
|
||||
|
||||
|
||||
</div>
|
||||
|
||||
|
||||
</header><div class="e-content entry-content" itemprop="articleBody text">
|
||||
<div>
|
||||
<p>On 14 May 2018, Sebastiano Vigna added <a href="http://pcg.di.unimi.it/pcg.php">a page to his website</a> (archived <a href="http://archive.is/VE0sX">here</a>) entitled “<em>The wrap-up on PCG generators</em>” that attempts to persuade readers to avoid various PCG generators.</p>
|
||||
<p>That day, he also submitted a link to his critique to <a href="https://www.reddit.com/r/programming/comments/8jbkgy/the_wrapup_on_pcg_generators/">Reddit</a> (archived <a href="http://archive.is/b76dd">here</a>). I think it is fair to say that his remarks did not get quite the reception he might have hoped for. Readers mostly seemed to infer a certain animosity in his tone and his criticisms gained little traction with that audience.</p>
|
||||
<p>Although I'm pleased to see readers of Reddit thinking critically about these things, it is worth taking the time to dive in and see what what lessons we can learn from all of this.</p>
|
||||
<!-- TEASER_END -->
|
||||
|
||||
<h3 id="background">Background</h3>
|
||||
<p>We have to feel a little sympathy for Vigna. On May 4, he updated <a href="http://xoshiro.di.unimi.it/">his website</a> to announce a new generation scheme, <em>Xoshiro</em> and accompanying <a href="http://vigna.di.unimi.it/ftp/papers/ScrambledLinear.pdf">paper</a>, the product of two years of work. He posted a link to his work <a href="https://www.reddit.com/r/programming/comments/8gx2d3/new_line_of_fast_prngs_released_by_the_author_of/">on Reddit</a> (archived <a href="http://archive.is/lv3js">here</a> and <a href="http://archive.is/1UULl">here</a>), and although he got some praise and thanks for his work, he ended up spending quite a lot of time talking not about his new work, but about flaws in his old work and about my work.</p>
|
||||
<p>Here is an example of the kind of remarks he had to contend with; Reddit user “TomatoCo” wrote:</p>
|
||||
<blockquote>
|
||||
<p>I liked xoroshiro a lot until I read all of the dire condemnations of it, so I switched to PCG. I'm not a mathematician, I can't understand your papers and PCG's write ups are a lot easier to understand. I'm sure that you've analyzed the shit out of your previous generator and I can see on your site you've come up with new techniques to measure if xoshiro suffers the same flaws. But once bitten, twice shy. Xoroshiro was defended as great with the sole exception of the lowest bit. But then it was "the lowest bit is just a LSFR, so don't use that. Well, actually, the other low bits are also just really long period LSFRs, well, actually," and new flaws were constantly appearing.
|
||||
Respectfully, I think you need to explain more and in simpler terms to earn everyone's trust back.</p>
|
||||
<p>The reason I picked PCG was because its author could, in plain language, describe its behavior and why some authors witnessed patterns in your RNG.</p>
|
||||
</blockquote>
|
||||
<p>I think it's quite understandable that Vigna would want to look for ways to take PCG (and me) down a peg or two, and in various comment replies he endeavored to express things he didn't like about PCG (and the PCG website).</p>
|
||||
<p>Most of the issues he raised were, I thought, adequately addressed and refuted in the Reddit discussion, but having gone to the effort already to try to articulate the things he did not like, even writing code to do so, it makes sense that he would want circulate these thoughts more broadly.</p>
|
||||
<h3 id="reddit-reaction-2">Reddit Reaction #2</h3>
|
||||
<p>Reddit's reaction to Vigna's new PCG-critique page was perhaps not what he hoped for. From what I can tell, pretty much none of the commenters were persuaded by his claims, and much was made of his tone.</p>
|
||||
<p>Regarding tone, user “notfancy” said:</p>
|
||||
<blockquote>
|
||||
<p>Take your feud somewhere else. […] theory and practice definitely belong here. The petty squabbling and the name calling definitely don't. Seeing that Vigna himself is posting links to his own site, this is to me self-promoting spam.</p>
|
||||
</blockquote>
|
||||
<p>and user “foofel” added:</p>
|
||||
<blockquote>
|
||||
<p>the style in which he presents his stuff is always full of hate and despise, that's not a good way to represent it and probably why people are fed up.</p>
|
||||
</blockquote>
|
||||
<p>and user “evand” added:</p>
|
||||
<blockquote>
|
||||
<p>I would describe a lot of it as written very... condescendingly. There's also a lot that is written to attack her and not PCG</p>
|
||||
</blockquote>
|
||||
<p>and user “AntiauthoritarianNow” chimed in, saying;</p>
|
||||
<blockquote>
|
||||
<p>Yeah, it's one thing to tease other researchers a little bit, but this guy has a real problem sticking to arguments on the merits rather than derailing into reddit-esque ad-hom.</p>
|
||||
</blockquote>
|
||||
<p>But the thread also had plenty of rebuttals. For just about every claim Vigna had made in his critique, there was a comment explaining why the claim was flawed.</p>
|
||||
<h3 id="my-reaction">My Reaction</h3>
|
||||
<p>I could settle back into my chair here, and say, “Thank you, Reddit, for keeping your wits about you!”, but since (at the time of writing) Vigna's page remains live with the same claims, it seems sensible for me to create my own writeup (this one) to address his claims directly.</p>
|
||||
<p>Moreover, I believe firmly that although it's never much fun to be on the receiving end for invective or personal attacks, in academia peer critique makes everything stronger. While much of what Vigna says about PCG doesn't hold up to closer scrutiny, it is worth trying to find value of some kind in every criticism. I believe in the approach taken in the world of improvisational comedy, known as <a href="https://en.wikipedia.org/wiki/Yes,_and...">“Yes, and…”</a>, which suggests that a participant should accept what another participant has stated (“yes”) and then expand on that line of thinking (“and”).</p>
|
||||
<p>Thus, in the subsequent sections, I'll look at each of Vigna's critiques, first give a defensive response, and then endeavor to find a way to say “Yes, and…” to each one.</p>
|
||||
<h3 id="correlations-due-to-contrived-seeding">Correlations Due to Contrived Seeding</h3>
|
||||
<p>Vigna's first two claims relate to creating two PCG generators whose outputs are correlated because he has specifically set them up to have internal states that would cause them to be correlated.</p>
|
||||
<h4 id="pcg-ext-variants-single-bit-change-to-the-extension-array">PCG <code>ext</code> Variants: Single Bit Change to the Extension Array</h4>
|
||||
<p>In the first claim, he modifies the code for the PCG of the extended generation scheme to so that he can flip a single bit in the extension array that adds <em>k</em>-dimensional equidistribution to a base generator.</p>
|
||||
<p>Vigna creates two <code>pcg64_k32</code> generators that are the same in all respects except for a single bit difference in one element of the 32-element extension array, and then observes that 31 of every 32 outputs will remain identical between the generators for some considerable time. Vigna clearly considers this behavior to be problematic and notes multiple LFSR-based PRNGs where such behavior would not occur.</p>
|
||||
<p>Vigna states</p>
|
||||
<blockquote>
|
||||
<p>Said otherwise, the whole sequence of the generator is made by an enormous number of strongly correlated, very short sequences. And this makes the correlation tests fail.</p>
|
||||
</blockquote>
|
||||
<p>Vigna concludes that no one should use generators like <code>pcg64_k32</code> as a result.</p>
|
||||
<h5 id="defensive-response">Defensive Response</h5>
|
||||
<p>Vigna actually created a <em>custom version</em> of the PCG code to effect his single bit change. The <code>pcg64_k32</code> generator has 2303 bits of state, 127 bits of LCG increment (which stays constant), 128 bits of LCG current state, and 32 64-bit words in the extension array. The odds of seeding two <code>pcg64_k32</code> generators each with 2303 bits of seed entropy and finding that they only differ by a single bit in the extension array is 1 in 2<sup>2292</sup>, an order of magnitude so vast that it cannot be represented as a floating point double.</p>
|
||||
<p>If the PRNG were properly initialized (e.g., using <code>std::seed_seq</code> or <code>pcg_extras::sed_seq_from<std::random_device></code>), Vigna's observed correlation would not have occurred. Likewise, had the single bit change been in the LCG side of the PRNG, it would also not have occurred.</p>
|
||||
<p>But what of Vigna's other claim, that PRNGs that are slow to diffuse single-bit changes to their internal state are necessarily bad? Vigna is right that for LFSR-based designs, the rate of bit diffusion (a.k.a. “<em>avalanche</em>”) matters a lot.</p>
|
||||
<p>However, numerous perfectly good designs for PRNGs would fail Vigna's criteria. All counter-based designs (e.g., SplitMix, Random123, Chacha) will preserve the single bit difference indefinitely if we examine their internal state. In fact, Vigna's collaborator, David Blackman, is author of <code>gjrand</code>, which also includes a counter whose internal state won't diverge significantly over time. But of these designs, only SplitMix would fail a test that looks for output correlations rather than similar internal states.</p>
|
||||
<p>The closest design to PCG's extension array is found in George Marsaglia's venerable <a href="https://en.wikipedia.org/wiki/Xorshift#xorwow">XorWow PRNG</a>, shown below (code taken from the Wikipedia page):</p>
|
||||
<pre><code>/* The state array must be initialized to not be all zero in the first four
|
||||
words */
|
||||
uint32_t xorwow(uint32_t state[static 5])
|
||||
{
|
||||
/* Algorithm "xorwow" from p. 5 of Marsaglia, "Xorshift RNGs" */
|
||||
uint32_t s, t = state[3];
|
||||
t ^= t >> 2;
|
||||
t ^= t << 1;
|
||||
state[3] = state[2]; state[2] = state[1]; state[1] = s = state[0];
|
||||
t ^= s;
|
||||
t ^= s << 4;
|
||||
state[0] = t;
|
||||
return t + (state[4] += 362437);
|
||||
}
|
||||
</code></pre>
|
||||
|
||||
<p>In Marsaglia's design, <code>state[4]</code> is a counter in much the same way that PCG's extension array is a “funky counter”. Marsaglia calls this counter a <em>Weyl sequence</em> after Hermann Weyl, who proved the equidistribution theorem in 1916.</p>
|
||||
<p>We can exactly reproduce Vigna's claim's about <code>pcg64_k32</code> producing similar output with XorWow. The program <a href="../downloads/snippets/uncxorwow.c"><code>uncxorwow.c</code></a> is a port of his demonstration program to XorWow. It fails if tested with PractRand, and, if we uncomment the <code>printf</code> statements, after 1 billion iterations we see that the outputs have not become uncorrelated. They continue to differ only in their high bit. And they will continue that way forever:</p>
|
||||
<pre><code>61b0be0f
|
||||
e1b0be0f
|
||||
c5a003d8
|
||||
45a003d8
|
||||
20e14479
|
||||
a0e14479
|
||||
5a5ebe42
|
||||
da5ebe42
|
||||
99ce85af
|
||||
19ce85af
|
||||
d2a1aabb
|
||||
52a1aabb
|
||||
6bf29670
|
||||
ebf29670
|
||||
948587d6
|
||||
148587d6
|
||||
e2c0f91c
|
||||
62c0f91c
|
||||
536fe7eb
|
||||
d36fe7eb
|
||||
</code></pre>
|
||||
|
||||
<p>Similarly, Vigna's complaint about “strongly correlated very short sequences” could likewise be applied to XorWow. It consists of 2<sup>64</sup> very similar sequences (differing only by a constant). It might seem bad at a glance to concatenate a number of very similar sequences but it is worth realizing that the nearest similar sequence is 2<sup>128</sup>-1 steps away. If Vigna would characterize 2<sup>128</sup>-1 as “very short”, he must be using a mathematician's sense of scale.</p>
|
||||
<p>Marsaglia's design of Xorwow quite deliberately uses a very simple and weak generator (a Weyl sequence) for a specific purpose. We could say “a counter isn't a very good random number generator”, but the key idea is that <em>it doesn't need to be</em>. It's not the whole story. It's a piece with a specific role to play, and it doesn't need to be any better than it is.</p>
|
||||
<p>PCG's extended generation scheme is a similar story. The extension array is a funky counter akin to a Weyl sequence (each array element is like a digit of a counter). It's slightly better than a Weyl sequence (a single bit change will quickly affect all the bits in the in that array element), but it is essentially the same idea.</p>
|
||||
<p>The <code>pcg64_k32_oneseq</code> and <code>pcg64_k32_fast</code> generators follow XorWow's scheme of just joining together the similar sequences. <code>pcg64_k32</code> swaps around chunks of size 2<sup>16</sup> from each similar sequence. In all cases, from any starting point you would need 2<sup>128</sup> outputs before the base linear congruential generator lined up to the same place again, and vastly more for the extension array to line up similarly. In short, for <code>pcg64_k32</code> the correlated states are quite literally unimaginably far away from each other.</p>
|
||||
<p>Talking about his contrived seedings, Vigna notes that, “This is all the decorrelation we get after a billion iterations, and it will not improve (not significantly before the thermodynamical death of the universe).” What he seems to have missed is the corollary to his statements—correlation and decorrelation are sides of the same coin. Two currently <em>uncorrelated</em> <code>pcg64_k32</code> states will not correlate before the heat death of the universe either.</p>
|
||||
<p>In short, Vigna contrived a seed to show correlation that would never arise in practice with normal seeding, nor could arise by advancing one generator. His critique is not unique to PCG, and should not be a concern for users of PCG.</p>
|
||||
<h5 id="yes-and-response">“Yes, and…” Response</h5>
|
||||
<p>A rather flippant “Yes, and…” response is that I'm perfectly happy for people to avoid <code>pcg64_k32</code>, as I'm not at all sure it is buying you anything meaningful over and above <code>pcg64</code>— it's a fair amount of added code complexity for something of dubious value. In fact, I didn't even bother to implement it in the C version and only a small number of people who have ported PCG have implemented it. As I see it, <em>k</em>-dimensional equidistribution sounds like a cool property, but the only use case I've found for such a property is <a href="http://www.pcg-random.org/party-tricks.html">performing party tricks</a>. But some people do like <em>k</em>-dimensional equidistribution, so let's press on…</p>
|
||||
<p>First, Vigna went to far too much trouble to create correlated states. He copied the entire C++ source for PCG and hacked it to make a private data member public so he could set a single bit. Had he been more familiar with the features the extended generators provide, he could instead have written.</p>
|
||||
<pre><code>pcg64_k32 rng0;
|
||||
pcg64_k32 rng1 = rng0;
|
||||
rng1.set(rng0() ^ 1);
|
||||
</code></pre>
|
||||
<p>This code uses <code>pcg64_k32</code>'s party-trick functionality to leap unimaginably huge distances across the state space to find exactly the correlated generator you want, one that is the same in every respect except for one differing output.</p>
|
||||
<p>In other words, what he sees as a deficiency, I've already highlighted as a feature.</p>
|
||||
<p>But whether it is achieved by the simple method above, or the more convoluted method Vigna used, we have the question of what to do if people are allowed to create very correlated generator states that would not normally arise in practice. One option is to just say “don't do that”, but a more “Yes, and…” perspective would be to allow people to create such states if they choose but provide a means to detect them. More on that in the next section.</p>
|
||||
<p>It's also worth asking whether the slowness with which a single bit change diffuses across the extension array is something inherent in the design of PCG's extended generation scheme, or mere happenstance. In fact, it is the latter.</p>
|
||||
<p>The only cleverness in the extended generation scheme isn't the idea of combining two generators, a strong one and a weaker-but-<em>k</em>-dimensionally-equidistributed one, it's the fact that we can do so without any <em>extra</em> state to keep track of what we're doing.</p>
|
||||
<p>I'm thus not wedded to the particular Weyl-sequence inspired method I used. If it's important that unimaginably distant similar generators do not stay correlated for long, that's a very easy feature to provide.</p>
|
||||
<p>When I designed how the extension array advances, I made a choice to make it “no better than it needs to be”. It doesn't need good avalanche properties, so that wasn't a design concern. But that doesn't mean it couldn't be tweaked to have good avalanche properties, so that a single bit change affects <em>all the bits</em> the next time the extension array advances. In fact, having designed <code>seed_seq_fe</code> for <code>randutils</code>, I'm aware of elegant and amply efficient ways to have better avalanche, so why not?</p>
|
||||
<p>It may not really be <em>necessary</em>, but I actually like this idea. So thanks, Sebastiano, I'll address this issue in a future update to PCG that provides some alternative schemes for updating the extension array!</p>
|
||||
<h4 id="pcg-regular-variants-contrived-seeds-for-inter-stream-correlations">PCG Regular Variants: Contrived seeds for Inter-Stream Correlations</h4>
|
||||
<p>In his next concern, Vigna uses makes correlated generators from two “random looking” seeds. He presents a program, <a href="../downloads/snippets/corrpcg.c"><code>corrpcg.c</code></a> that mixes together the two correlated generators and can then be fed into statistical tests (which will fail because of the correlation).</p>
|
||||
<h5 id="defensive-response_1">Defensive Response</h5>
|
||||
<p>We can devise bad seed pairs for just about any PRNG. Here are three example programs, <a href="../downloads/snippets/corrxoshiro.c"><code>corrxoshiro.c</code></a>, <a href="../downloads/snippets/corrsplitmix.c"><code>corrsplitmix.c</code></a>, and <a href="../downloads/snippets/corrxorwow.c"><code>corrxorwow.c</code></a>, which initialize generators with two “random looking” seeds but create correlated streams that will fail statistical tests if mixed.</p>
|
||||
<p>In all cases, despite being “random looking”, the seeds are carefully contrived. Seeds such as these would be vanishingly unlikely with proper seeding practice.</p>
|
||||
<p>As before, the concerns Vigna expresses apply to many prior generators. We can view XorWow's <code>state[4]</code> value as being a stream selection constant, but this time let's focus in on SplitMix. For SplitMix, different <code>gamma_</code> values constitute different streams.</p>
|
||||
<p>In <code>corrsplitmix.c</code> the implementation is hard-wired to use a single stream (<code>0x9e3779b97f4a7c15</code>), but in <a href="../downloads/snippets/corrsplitmix2.c"><code>corrsplitmix2.c</code></a> we mix two streams, (<code>0x9e3779b97f4a7c15</code> and <code>0xdf67d33dd518d407</code>) and observe correlations. Although these gamma values look random, they are not, they are carefully contrived. In particular, here <code>0xdf67d33dd518d407 * 3</code> = <code>0x9e3779b97f4a7c15</code> (in 64-bit arithmetic), which means that every third output from the second stream will exactly match an output from the first.</p>
|
||||
<p>Vigna's critique thus applies at least as strongly to SplitMix's streams as it does to PCG's.</p>
|
||||
<p>I have <a href="critiquing-pcg-streams.html">written at length about PCG's streams</a> (and discussed SplitMix's, too). I freely acknowledge that these streams exist in a space of trade-offs where we are choosing to do the cheap thing, leveraging the properties of the underlying LCG (or Weyl sequence for SplitMix). In that article, I say:</p>
|
||||
<blockquote>
|
||||
<p>Changing the increment parameter is just barely enough for streams that are actually useful. They aren't statistically independent, far from it, but they are distinct and they do help.</p>
|
||||
</blockquote>
|
||||
<p>No one should worry that PCG's streams makes anything worse.</p>
|
||||
<h5 id="yes-and-response_1">“Yes, and…” Response</h5>
|
||||
<p>Although it is vanishingly unlikely that two randomly seeded <code>pcg64</code> generators would be correlated (it would only happen with poor/adversarial seeding), it is reasonable to ask if this kind of correlation due to bad seeding can be detected.</p>
|
||||
<p>We can even argue that another <em>checklist feature</em> for a general-purpose PRNG is the ability to tell how independent the sequences from two seeds are likely to be. PCG goes some way towards this goal with its <code>-</code> operator that calculates the distance between two generators, but the functionality was originally designed for generators on the same stream. I've now updated that functionality so that for generators on different streams, it will calculate the distance to their point of closest approach (i.e., where the differences between successive values of the underlying LCG align).</p>
|
||||
<p>So it's now possible with PCG to compare two generators to see whether they have been badly seeded so that they correlate.</p>
|
||||
<p>Here's a short <a href="../downloads/snippets/strmdist.c">test program</a>:</p>
|
||||
<pre><code>#include "pcg_random.hpp"
|
||||
#include "pcg_extras.hpp"
|
||||
|
||||
#include <iostream>
|
||||
#include <iomanip>
|
||||
#include <random>
|
||||
|
||||
int main() {
|
||||
using namespace pcg_extras;
|
||||
|
||||
#if USE_VIGNA_CONTRIVED_SEEDS
|
||||
pcg64 x(PCG_128BIT_CONSTANT(0x83EED115C9CBCC30, 0x4C55E45838B75647),
|
||||
PCG_128BIT_CONSTANT(0x3E0897751B1A19E7, 0xD9D50DD3E3A454DC));
|
||||
pcg64 y(PCG_128BIT_CONSTANT(0x7C112EEA363433CF, 0xB3AA1BA7C748A9B9),
|
||||
PCG_128BIT_CONSTANT(0x41F7688AE4E5E618, 0x262AF22C1C5BAB23));
|
||||
#elif USE_PCG_UNIQUE
|
||||
pcg64_unique x,y;
|
||||
#elif USE_SMALL_SEEDS1
|
||||
pcg64 x(0), y(1);
|
||||
#elif USE_SMALL_SEEDS2
|
||||
pcg64 x(0,0), y(0,1);
|
||||
#elif USE_SMALL_SEEDS3
|
||||
pcg64 x(0,0), y(1,1);
|
||||
#elif USE_RANDOM_DEVICE
|
||||
pcg64 x(seed_seq_from<std::random_device>{}),
|
||||
y(seed_seq_from<std::random_device>{});
|
||||
#endif
|
||||
|
||||
std::cout << std::hex;
|
||||
for (int i = 0; i < 10; ++i) {
|
||||
std::cout << (x - y) << ": ";
|
||||
std::cout << x() << ", " << y() << "\n";
|
||||
}
|
||||
}
|
||||
</code></pre>
|
||||
|
||||
<p>And here are the results of running it (in each case, each line shows the distance between the streams and a value from each PRNG; the distance stays the same because the PRNGs are advancing together):</p>
|
||||
<pre><code>unix% c++ -Wall -std=c++11 -o strmdist strmdist.cpp -Iinclude -DUSE_RANDOM_DEVICE && ./strmdist
|
||||
a571d615b08fea47c84f39f0811f04f: c021049beac5efd0, ceaa3596f168e8b6
|
||||
a571d615b08fea47c84f39f0811f04f: 573371998db59a67, e5d84a00b37c3556
|
||||
a571d615b08fea47c84f39f0811f04f: bc4246c671ef9a1f, 1b13ad2f224707c7
|
||||
a571d615b08fea47c84f39f0811f04f: b1f3e4ffcfef569, 11b50b226a67cdbe
|
||||
a571d615b08fea47c84f39f0811f04f: 8a378ec693dc1e4, 903ccfd4dc769389
|
||||
a571d615b08fea47c84f39f0811f04f: 4799de5c580be6ab, 22d13ce52d83c9cb
|
||||
a571d615b08fea47c84f39f0811f04f: e8fdf041a93626e8, f24c8f49866b7b4e
|
||||
a571d615b08fea47c84f39f0811f04f: f29e3d08104d7630, b37e5b58ae91d45c
|
||||
a571d615b08fea47c84f39f0811f04f: 28f524ad8f57bedb, 52d41d39b1186616
|
||||
a571d615b08fea47c84f39f0811f04f: 9be8cb37ea8952b5, e6812ed8f0613d3
|
||||
|
||||
unix% c++ -Wall -std=c++11 -o strmdist strmdist.cpp -Iinclude -DUSE_RANDOM_DEVICE && ./strmdist
|
||||
25c3990ef6e7766ab543435aa25f4326: 2f76ab68249fd7f5, 4fbfc0ce19119391
|
||||
25c3990ef6e7766ab543435aa25f4326: 933845d6c7ad9396, 7572dae64b2cc5a
|
||||
25c3990ef6e7766ab543435aa25f4326: d7d1dc18bae0604a, 5b1f8310e1f0dc8a
|
||||
25c3990ef6e7766ab543435aa25f4326: 85cd1dcff8830ad5, a1cfea3c01314c8d
|
||||
25c3990ef6e7766ab543435aa25f4326: 543ba46266a0b6ba, 7217b15c05cba254
|
||||
25c3990ef6e7766ab543435aa25f4326: 5a3bd5d4d6c49a55, a243af7df5cfe287
|
||||
25c3990ef6e7766ab543435aa25f4326: 9f2dc30afc3dcead, deaa9d03f7ca1117
|
||||
25c3990ef6e7766ab543435aa25f4326: 5856b884c1298dc9, 67502e4490b77bae
|
||||
25c3990ef6e7766ab543435aa25f4326: 9b94ebb084cc6fdd, 2e07957697add77c
|
||||
25c3990ef6e7766ab543435aa25f4326: efe6b451c262a3fb, 2e94d782daae964d
|
||||
|
||||
unix% c++ -Wall -std=c++11 -o strmdist strmdist.cpp -Iinclude -DUSE_RANDOM_DEVICE && ./strmdist
|
||||
32982840d1ddcb5e7f1ed57a6d496525: 96ed26957ef938db, 568fe0aa7e9e8a26
|
||||
32982840d1ddcb5e7f1ed57a6d496525: 33270d80d24b0965, 44e42e1afc4db710
|
||||
32982840d1ddcb5e7f1ed57a6d496525: 6de9ac5272dd1193, 90696d1c4f52e71d
|
||||
32982840d1ddcb5e7f1ed57a6d496525: 43c5c899c7123e57, 337b9d25e00fb0de
|
||||
32982840d1ddcb5e7f1ed57a6d496525: 753954b73076704d, f4fce4c33756df7e
|
||||
32982840d1ddcb5e7f1ed57a6d496525: 3b5dc9402b56584d, fd7ae3c708355dc0
|
||||
32982840d1ddcb5e7f1ed57a6d496525: 15a9227305a442d8, 78fa04eb7f881590
|
||||
32982840d1ddcb5e7f1ed57a6d496525: b9e58872c3a299, 381a8f851acbc5f4
|
||||
32982840d1ddcb5e7f1ed57a6d496525: 1b624879e6cf5128, aa908d3a4f2d8f02
|
||||
32982840d1ddcb5e7f1ed57a6d496525: 79d4836bb5a56a77, 1650f74b3ef617f9
|
||||
|
||||
unix% c++ -Wall -std=c++11 -o strmdist strmdist.cpp -Iinclude -DUSE_SMALL_SEEDS1 && ./strmdist
|
||||
1c31b969dc65d7b0df636de659042bb1: 1070196e695f8f1, e175e32ed3507bfa
|
||||
1c31b969dc65d7b0df636de659042bb1: 703ec840c59f4493, c0bf922a0b283109
|
||||
1c31b969dc65d7b0df636de659042bb1: e54954914b3a44fa, 140bfa21e68785bb
|
||||
1c31b969dc65d7b0df636de659042bb1: 96130ff204b9285e, c5ec8bcc4fe35830
|
||||
1c31b969dc65d7b0df636de659042bb1: 7d9fdef535ceb21a, 4dd8ed1ca22869c5
|
||||
1c31b969dc65d7b0df636de659042bb1: 666feed42e1219a0, c9bffa29c802ef4c
|
||||
1c31b969dc65d7b0df636de659042bb1: 981f685721c8326f, 3aa09aa4e147478b
|
||||
1c31b969dc65d7b0df636de659042bb1: ad80710d6eab4dda, 1dfdf6222d06378c
|
||||
1c31b969dc65d7b0df636de659042bb1: e202c480b037a029, 5a05dacf4df61d4e
|
||||
1c31b969dc65d7b0df636de659042bb1: 5d3390eaedd907e2, 489650b1eb840a26
|
||||
|
||||
unix% c++ -Wall -std=c++11 -o strmdist strmdist.cpp -Iinclude -DUSE_SMALL_SEEDS2 && ./strmdist
|
||||
151361a7e7368c239a3988178df4d76d: d4feb4e5a4bcfe09, acdbf879b3c73375
|
||||
151361a7e7368c239a3988178df4d76d: e85a7fe071b026e6, 7ea754d074e8d88f
|
||||
151361a7e7368c239a3988178df4d76d: 3a5b9037fe928c11, f8fc7aec8ae6245a
|
||||
151361a7e7368c239a3988178df4d76d: 7b044380d100f216, 7d2ebc3c0b5bedb4
|
||||
151361a7e7368c239a3988178df4d76d: 1c7850a6b6d83e6a, cbaf666f55051666
|
||||
151361a7e7368c239a3988178df4d76d: 240b82fcc04f0926, 4eba9f04dfb9903b
|
||||
151361a7e7368c239a3988178df4d76d: 7e43df85bf9fba26, 4fab6bcf361bd63d
|
||||
151361a7e7368c239a3988178df4d76d: 43adf3380b1fe129, 257fcac1ed3817df
|
||||
151361a7e7368c239a3988178df4d76d: 3f0fb307287219c, bf6f5515988a494
|
||||
151361a7e7368c239a3988178df4d76d: 781f4b84f42a2df, 1081ed38c84c1c9d
|
||||
|
||||
unix% c++ -Wall -std=c++11 -o strmdist strmdist.cpp -Iinclude -DUSE_SMALL_SEEDS3 && ./strmdist
|
||||
edfe668df810de6e58b8e92e878fefa: d4feb4e5a4bcfe09, d4692f845d3a3706
|
||||
edfe668df810de6e58b8e92e878fefa: e85a7fe071b026e6, bb0f09b0eebab6ff
|
||||
edfe668df810de6e58b8e92e878fefa: 3a5b9037fe928c11, e26ac904ad283c09
|
||||
edfe668df810de6e58b8e92e878fefa: 7b044380d100f216, 83860212b5d92197
|
||||
edfe668df810de6e58b8e92e878fefa: 1c7850a6b6d83e6a, 1c3601ed5afd3f49
|
||||
edfe668df810de6e58b8e92e878fefa: 240b82fcc04f0926, 5e4fa027be29b47e
|
||||
edfe668df810de6e58b8e92e878fefa: 7e43df85bf9fba26, b930e28d59383019
|
||||
edfe668df810de6e58b8e92e878fefa: 43adf3380b1fe129, e0d61e1b074df835
|
||||
edfe668df810de6e58b8e92e878fefa: 3f0fb307287219c, f42c38b1aca3ac9d
|
||||
edfe668df810de6e58b8e92e878fefa: 781f4b84f42a2df, 19e9cc4fa58fd0ad
|
||||
|
||||
unix% c++ -Wall -std=c++11 -o strmdist strmdist.cpp -Iinclude -DUSE_PCG_UNIQUE && ./strmdist
|
||||
534a7c98f86b50b72fad6990038ba18: af8a07de4c8d67d1, d649257470c0180d
|
||||
534a7c98f86b50b72fad6990038ba18: 3789d12fe8e452b1, 1017152e85f732fc
|
||||
534a7c98f86b50b72fad6990038ba18: c3c4e780fd60901b, 91a9d78551f0c776
|
||||
534a7c98f86b50b72fad6990038ba18: e7257e02f7fa5b40, 46fb62417ebf2f13
|
||||
534a7c98f86b50b72fad6990038ba18: 3697948fa9aa8378, 60e44721c6fbc9d0
|
||||
534a7c98f86b50b72fad6990038ba18: 7bdbcc91de7efbcf, 21de9d1dc03e2ca6
|
||||
534a7c98f86b50b72fad6990038ba18: 9cf598a61c9ad958, 62e8c3dc421f4e58
|
||||
534a7c98f86b50b72fad6990038ba18: 5c8a6da6c91b7d35, 3cb08b7e59fd655a
|
||||
534a7c98f86b50b72fad6990038ba18: f55a8b190a85c9c0, 5a71766fac52ec8a
|
||||
534a7c98f86b50b72fad6990038ba18: 906b1a30904fe59, f71525dc1d91a06e
|
||||
|
||||
unix% c++ -Wall -std=c++11 -o strmdist strmdist.cpp -Iinclude -DUSE_PCG_UNIQUE && ./strmdist
|
||||
1b7a9a85b5ed2b6a2a92da9e093eba18: a11d6aa92efc9a79, e646943445e368a
|
||||
1b7a9a85b5ed2b6a2a92da9e093eba18: 35026a6e1a195a29, 906b9bed756e1667
|
||||
1b7a9a85b5ed2b6a2a92da9e093eba18: af1f1193515d9e7b, fe51967d5d532f70
|
||||
1b7a9a85b5ed2b6a2a92da9e093eba18: 61baa5620ceeff38, 644345c453ee3b11
|
||||
1b7a9a85b5ed2b6a2a92da9e093eba18: 71e88c9c27a7abbf, 1b6a254f565f6c70
|
||||
1b7a9a85b5ed2b6a2a92da9e093eba18: 1125753cd420e3c1, 8be4065858e93c57
|
||||
1b7a9a85b5ed2b6a2a92da9e093eba18: a53ce57ffaa57eb3, 7f1c546ae9bf7b61
|
||||
1b7a9a85b5ed2b6a2a92da9e093eba18: 4cf2c7c152326c4, ada2d31650f07ef8
|
||||
1b7a9a85b5ed2b6a2a92da9e093eba18: b731cbec3bfba773, 92ce80f0c8dc855f
|
||||
1b7a9a85b5ed2b6a2a92da9e093eba18: b8c449d4872f7971, 44ed4207442550da
|
||||
|
||||
unix% c++ -Wall -std=c++11 -o strmdist strmdist.cpp -Iinclude -DUSE_PCG_UNIQUE && ./strmdist
|
||||
360981a27aee6d34271feaa80270ba18: 5da8c0afa4330059, 67af26ab1d05ed52
|
||||
360981a27aee6d34271feaa80270ba18: ef0ef074871cc9a0, cda2688372cb72b7
|
||||
360981a27aee6d34271feaa80270ba18: 6a15c49d4ae8d89d, 3708ddd964f616fe
|
||||
360981a27aee6d34271feaa80270ba18: dd8f24112bcbf580, 69309c3ffa6cea2e
|
||||
360981a27aee6d34271feaa80270ba18: e8f252a4132fd0e3, e3ff9751773f6db
|
||||
360981a27aee6d34271feaa80270ba18: e23a1246ea5980be, 1161fd499cbecafa
|
||||
360981a27aee6d34271feaa80270ba18: 1d19a64904134065, a9e31a01b4c51a43
|
||||
360981a27aee6d34271feaa80270ba18: 2c3166d304f9dedf, fdd3f540a6859c19
|
||||
360981a27aee6d34271feaa80270ba18: 8f73778d1f6133ea, 13a54957b3c65205
|
||||
360981a27aee6d34271feaa80270ba18: c8d362ba3d62239, 66db0b2ae6908dc8
|
||||
|
||||
unix% c++ -Wall -std=c++11 -o strmdist strmdist.cpp -Iinclude -DUSE_PCG_UNIQUE && ./strmdist
|
||||
1266069359d4404d4fe77f291da43a18: 9994872b3cc3104c, 5582722b3f354f4b
|
||||
1266069359d4404d4fe77f291da43a18: cec9ae92f2f0a929, 7a2d534e7c3a7281
|
||||
1266069359d4404d4fe77f291da43a18: ce777879518e6169, c384bb65c1d4364b
|
||||
1266069359d4404d4fe77f291da43a18: 2cb082454d09aa19, 703c5ad7747a9b42
|
||||
1266069359d4404d4fe77f291da43a18: a581d3154c60654, b4b9369d997cda6e
|
||||
1266069359d4404d4fe77f291da43a18: 5ba66e3d99cd33c9, 80aa887fbb5fdef3
|
||||
1266069359d4404d4fe77f291da43a18: 1038e3281dcae11d, 54c304cf2a66182c
|
||||
1266069359d4404d4fe77f291da43a18: 9df3df9d27af7148, 7ddd385e114299b9
|
||||
1266069359d4404d4fe77f291da43a18: bf1656198867bd08, 7aeae9ba84a17dbe
|
||||
1266069359d4404d4fe77f291da43a18: 60aef1418aa1c6f1, 8a7196feda932f06
|
||||
|
||||
unix% c++ -Wall -std=c++11 -o strmdist strmdist.cpp -Iinclude -DUSE_VIGNA_CONTRIVED_SEEDS && ./strmdist
|
||||
0: e1e4e4b44cca9ade, 43dc3c9c96899953
|
||||
0: a3ef563648055140, 2b8a051f7ab1b24
|
||||
0: 7aa3dc341221459a, 1a0960a2cd3d51ee
|
||||
0: cfa0d055fbe9f476, a0abf5d3e8ed9f41
|
||||
0: b69403f2c93f3fce, 807e58a7e7f9d6d2
|
||||
0: a2550ed76e8d9ae, 144aa1daedd1b35e
|
||||
0: a1f898a64347533b, c532263a99dd0fc4
|
||||
0: d483377a20c295f0, bbd10614af86a019
|
||||
0: 5c6469b1053d2ce1, 9c2b8c8d2e20a7a5
|
||||
0: 5f91b4bd64d5eeb1, 58afc8da4eb26af7
|
||||
</code></pre>
|
||||
|
||||
<p>As we can see in the last example, Vigna contrived seeds that had the streams exactly aligned. The values from each stream are distinct in this case, but a statistical test will see that they are correlated.</p>
|
||||
<p>Interpreting the distance value is easy in this case, but not every user will be able to do so, and some distances (e.g., a just a single high bit set) would also be bad, so better detection of contrived seeds probably demands a new function, <code>independence_score()</code>, based on this distance metric.</p>
|
||||
<p>Beyond these functions, there is also the question of whether it is wise to allow users to seed generators where they can specify the entire internal state. Vigna's generators (and all basically LFSR generators) must avoid the all-zeros state and do not like states with low hamming weight (so <code>{ seed, 0, 0, 0 }</code> is also a poor choice). With these issues in mind, perhaps we should deny users the ability to seed the entire state. That might prevent some contrived seedings like the one Vigna used. I'm not fully sold on this idea, but it is a widely-used approach use by other generators (e.g., Blackman's gjrand) and worth considering.</p>
|
||||
<p>Although Vigna's contrived seeding was a bit silly, his example has helped me improve the PCG distance metric, given us another checklist feature that some people might want (detecting bad seed pairs), got me thinking about future features, and returned me to the topic of good seeding. All in all, we can call this a positive contribution. Thanks, Sebastiano!</p>
|
||||
<h3 id="prediction-difficulty">Prediction Difficulty</h3>
|
||||
<p>The next two sections relate to predicting PCG.</p>
|
||||
<h4 id="predicting-pcg32_oneseq">Predicting <code>pcg32_oneseq</code>
|
||||
</h4>
|
||||
<p>Vigna writes:</p>
|
||||
<blockquote>
|
||||
<p>To me, it has been always evident that PCG generators are very easy to predict. Of course, no security expert ever tried to to that: it would be like beating 5-year-old kid on a race. It would be embarrassing.</p>
|
||||
<p>So we had this weird chicken-and-egg situation: nobody who could easily predict a PCG generator would write the code, because they are too easy to predict; but since nobody was predicting a PCG generator, Melissa O'Neill kept on the absurd claim that they were challenging to predict.</p>
|
||||
</blockquote>
|
||||
<p>Vigna then goes on to show code to predict <code>pcg32_oneseq</code>, a 64-bit PRNG with 32-bit output.</p>
|
||||
<h5 id="defensive-response_2">Defensive Response</h5>
|
||||
<p>As one reddit observer wrote:</p>
|
||||
<blockquote>
|
||||
<p>[Vigna's] program needs to totally brute force half of the state, and then some additional overhead to brute force bits of the rest of the state, so runtime is 2n/2, exponential, not polynomial.</p>
|
||||
</blockquote>
|
||||
<p>Vigna has written an exponential algorithm to brute force 32 bits of state. I hope it was obvious to almost everyone that I never claimed that brute-forcing 32-bits of state was hard. In fact, I have already <a href="http://www.pcg-random.org/predictability.html#predictability-of-the-pcg-family">outlined how</a> to predict <code>pcg32</code> (more bits to figure out given the unknown stream). I observed that <code>pcg32</code> is predictable using established techniques (specifically the LLL algorithm), and I have even linked to <a href="https://github.com/mariuslp/PCG_attack">an implementation</a> of those ideas by Marius Lombard-Platet.</p>
|
||||
<p>I characterize <code>pcg32_oneseq</code> as easy to brute force, and <code>pcg32</code> as annoying (as Marius Lombard-Platet discovered). Only when we get to <code>pcg64</code> do we have something where there is a meaningful challenge.</p>
|
||||
<p>If Vigna really believes that <em>all</em> members of the PCG family are easy to predict, he should have predicted <code>pcg64</code> or <code>pcg64_c32</code>.</p>
|
||||
<h5 id="yes-and-response_2">“Yes, and…” Response</h5>
|
||||
<p>The best part of Vigna's critique are these lines:</p>
|
||||
<blockquote>
|
||||
<p>Writing the function that performs the prediction, <code>recover()</code>, took maybe half an hour of effort. It's a couple of loops, a couple of if's and a few logical operations. Less than 10 lines of code (of course this can be improved, made faster, etc.).</p>
|
||||
</blockquote>
|
||||
<p>and the source code comment that reads:</p>
|
||||
<blockquote>
|
||||
<p>Pass an initial state (in decimal or hexadecimal), see it recovered from the output in a few seconds.</p>
|
||||
</blockquote>
|
||||
<p>So, here Vigna is essentially endorsing all the <em>practical</em> aspects I've previously noted regarding trivial predictability. Specifically, he's noting that with little <em>time</em> or <em>effort</em>, he can write a <em>simple</em> program that <em>quickly</em> predicts a PRNG and has actually done so. This is very different from taking a purely theoretical perspective (e.g., noting that techniques exist to solve a problem in polynomial time without ever implementing them).</p>
|
||||
<p>In other words, clearly <em>ease of prediction</em> matters to Vigna. So we both <em>agree</em>—<code>pcg32_oneseq</code> is easy to predict.</p>
|
||||
<p>Now let's keep that characterization of easiness and move on to some of the other generators.</p>
|
||||
<p>Vigna and I would agree, I think, that <em>I</em> lack the necessary insight to develop fast prediction methods for <code>pcg64</code> or <code>pcg64_c32</code> (it's an instance of <a href="https://www.schneier.com/blog/archives/2011/04/schneiers_law.html">Schneier's Law</a>). Vigna is also right that, if it is tractable to predict, those who might have the necessary skill lack much incentive to try. For some years I have been intending to have a prediction contest with real prizes and I remain hopeful that I'll find the time to run such a contest this summer. When the contest finally launches, I hope he'll have a go—I'd be delighted to send him a prize.</p>
|
||||
<h4 id="predicting-pcg64_once_insecure">Predicting <code>pcg64_once_insecure</code>
|
||||
</h4>
|
||||
<p>Vigna also notes that he can invert the bijection that serves as the output function for <code>pcg64_once_insecure</code>, which reveals the underlying LCG with all its statistical flaws.</p>
|
||||
<h5 id="defensive-response_3">Defensive Response</h5>
|
||||
<p>I noted this exact issue in 2014 in the PCG paper. It's why <code>pcg64_once_insecure</code> has the name it does. I discourage its use as a general-purpose PRNG precisely because of its invertible output function.</p>
|
||||
<h5 id="yes-and-response_3">“Yes, and…” Response</h5>
|
||||
<p>Vigna is at least acknowledging that some people might care about this property.</p>
|
||||
<h3 id="speed-and-comparison-against-lcgs">Speed and Comparison against LCGs</h3>
|
||||
<p>Finally, Vigna develops a PCG variant using a traditional integer hash function based on MurmurHash (I would call it PCG XS M XS M XS). He claims it is faster than the PCG variants I recommend and notes that he doesn't consider PCG especially fast.</p>
|
||||
<h4 id="defensive-response_4">Defensive Response</h4>
|
||||
<p>I considered this exact idea in the 2014 PCG paper. In my tests, I found that a variant using a very similar general integer hash function was not as fast as the PCG permutations I used.</p>
|
||||
<p>Testing is a finicky business.</p>
|
||||
<h4 id="yes-and-response_4">“Yes, and…” Response</h4>
|
||||
<p>I absolutely agree with Vigna's claim that people should run their own speed tests.</p>
|
||||
<p>I also realized long ago that PCG probably won't have the speed crown, because it can't. A simple truncated 128-bit LCG passes all standard statistical tests once we get up to 128 bits, and beats everything, including Vigna's generators. Because <code>pcg64</code> is built from a 128-bit LCG, it can never beat it in speed.</p>
|
||||
<p>I should write a blog post on speed testing. But here's a taste. We'll use Vigna's hamming-weight test as our benchmark, because it is a real program that puts randomness to actual use but is coded with execution speed in mind.</p>
|
||||
<p>First, let's test the Mersenne Twister. Compiling with Clang, we get</p>
|
||||
<pre><code>processed 1.75e+11 bytes in 130 seconds (1.346 GB/s, 4.847 TB/h). Fri May 25 14:03:25 2018
|
||||
</code></pre>
|
||||
|
||||
<p>whereas compiling with GCC, we get</p>
|
||||
<pre><code>processed 1.75e+11 bytes in 73 seconds (2.397 GB/s, 8.631 TB/h). Fri May 25 14:05:44 2018
|
||||
</code></pre>
|
||||
|
||||
<p>With GCC, it runs almost twice as fast.</p>
|
||||
<p>Now let's contrast that result with this 128-bit MCG:</p>
|
||||
<pre><code>static uint128_t state = 1; // can be seeded to any odd number
|
||||
|
||||
static inline uint64_t next()
|
||||
{
|
||||
constexpr uint128_t MULTIPLIER =
|
||||
(uint128_t(0x0fc94e3bf4e9ab32ULL) << 64) | 0x866458cd56f5e605ULL;
|
||||
// Spectral test: M8 = 0.71005, M16 = 0.66094, M24 = 0.61455
|
||||
state *= MULTIPLIER;
|
||||
return state >> 64;
|
||||
}
|
||||
</code></pre>
|
||||
|
||||
<p>Compiling with Clang, we get</p>
|
||||
<pre><code>processed 1.75e+11 bytes in 39 seconds (4.488 GB/s, 16.16 TB/h). Fri May 25 14:16:25 2018
|
||||
</code></pre>
|
||||
|
||||
<p>whereas with GCC we get</p>
|
||||
<pre><code>processed 1.75e+11 bytes in 58 seconds (3.017 GB/s, 10.86 TB/h). Fri May 25 14:18:14 2018
|
||||
</code></pre>
|
||||
|
||||
<p>The GCC code is no slouch, but Clang's code here is <em>much</em> faster. Clang is apparently better at 128-bit math.</p>
|
||||
<p>If we really care about speed though, <em>this</em> 128-bit MCG (which uses a carefully chosen 64-bit multiplier instead of a more typical 128-bit multiplier) is even faster and still passes statistical tests:</p>
|
||||
<pre><code>static uint128_t state = 1; // can be seeded to any odd number
|
||||
|
||||
static inline uint64_t next()
|
||||
{
|
||||
return (state *= 0xda942042e4dd58b5ULL) >> 64;
|
||||
}
|
||||
</code></pre>
|
||||
|
||||
<p>Compiling with Clang, we get</p>
|
||||
<pre><code>processed 1.75e+11 bytes in 37 seconds (4.73 GB/s, 17.03 TB/h). Fri May 25 14:09:26 2018
|
||||
</code></pre>
|
||||
|
||||
<p>whereas with GCC we get</p>
|
||||
<pre><code>processed 1.75e+11 bytes in 44 seconds (3.978 GB/s, 14.32 TB/h). Fri May 25 14:11:40 2018
|
||||
</code></pre>
|
||||
|
||||
<p>Again, Clang takes the speed crown; its executable generates and checks 1 TB of randomness about every 3.5 minutes.</p>
|
||||
<p>If we test Vigna's latest generator, xoshiro256**, and compile with Clang, it gives us</p>
|
||||
<pre><code>processed 1.75e+11 bytes in 50 seconds (3.5 GB/s, 12.6 TB/h). Fri May 25 14:30:05 2018
|
||||
</code></pre>
|
||||
|
||||
<p>whereas with GCC we get</p>
|
||||
<pre><code>processed 1.75e+11 bytes in 43 seconds (4.07 GB/s, 14.65 TB/h). Fri May 25 14:31:52 2018
|
||||
</code></pre>
|
||||
|
||||
<p>This result is very fast, but not faster than either 128-bit MCG.</p>
|
||||
<p>Finally, let's look at PCG-style generators. First let's look at Vigna's proposed variant. Compiling with Clang, we get</p>
|
||||
<pre><code>processed 1.75e+11 bytes in 59 seconds (2.966 GB/s, 10.68 TB/h). Fri May 25 14:44:37 2018
|
||||
</code></pre>
|
||||
|
||||
<p>and with GCC we get</p>
|
||||
<pre><code>processed 1.75e+11 bytes in 62 seconds (2.823 GB/s, 10.16 TB/h). Fri May 25 14:46:42 2018
|
||||
</code></pre>
|
||||
|
||||
<p>This is one of the rare occasions where GCC and Clang actually turn in almost equivalent times.</p>
|
||||
<p>In contrast, with the general-purpose <code>pcg64</code> generator, compiling with Clang I see:</p>
|
||||
<pre><code>processed 1.75e+11 bytes in 57 seconds (3.07 GB/s, 11.05 TB/h). Fri May 25 14:57:02 2018
|
||||
</code></pre>
|
||||
|
||||
<p>whereas with GCC, I see</p>
|
||||
<pre><code>processed 1.75e+11 bytes in 64 seconds (2.735 GB/s, 9.844 TB/h). Fri May 25 14:59:07 2018
|
||||
</code></pre>
|
||||
|
||||
<p>Thus, depending on which compiler we choose, Vigna's variant is either slightly faster or slightly slower.</p>
|
||||
<p>Finally, if we look at <code>pcg64_fast</code>, compiling with Clang gives us</p>
|
||||
<pre><code>processed 1.75e+11 bytes in 49 seconds (3.572 GB/s, 12.86 TB/h). Fri May 25 15:00:45 2018
|
||||
</code></pre>
|
||||
|
||||
<p>and with GCC we get</p>
|
||||
<pre><code>processed 1.75e+11 bytes in 65 seconds (2.693 GB/s, 9.693 TB/h). Fri May 25 15:02:15 2018
|
||||
</code></pre>
|
||||
|
||||
<p>Again the performance of GCC is a bit disappointing; this MCG-based generator is actually running slower than the LCG-based one.</p>
|
||||
<p>From this small amount of testing, we can see that <code>pcg64</code> is not as fast as <code>xoshiro256**</code>, but a lot depends on the compiler you're using—if you're using Clang (which is the default compiler on OS X), <code>pcg64_fast</code> will beat xoshiro256**.</p>
|
||||
<p>There's plenty of room for speed improvement in PCG. My original goal was to be faster than the Mersenne Twister, which it is, but knowing that it'll always be beaten by the speed of its underlying LCG I haven't put a lot of effort into optimizing the code. I could have used the faster multiplier that I used above, and there is actually a completely different way of handling LCG increment that reduces dependences and enhances speed but implementing LCGs that way makes the code more opaque. If PCG's speed is an issue, these are design decisions are worth revisiting.</p>
|
||||
<p>But the speed winner is clearly a 128-bit MCG. It's actually what I use when speed is the primary criterion.</p>
|
||||
<h3 id="conclusion">Conclusion</h3>
|
||||
<p>None of Vigna's concerns raise any serious worries about PCG. But critique is useful, and helps spur us to do better.</p>
|
||||
<p>I'm sure Vigna has spent far longer thinking about PCG than he would like, so it is best to say a big thank you to him for all the thought and energy he has expended in these efforts. I'm pleased that I've mostly been able to put the critique to good use—it may be mostly specious for users, but it is certainly helpful for me. Reddit mostly saw vitriol and condescension, but I prefer to see it as a gift of his time and thought.</p>
|
||||
<p>Thanks, Sebastiano!</p>
|
||||
</div>
|
||||
</div>
|
||||
<aside class="postpromonav"><nav><ul itemprop="keywords" class="tags">
|
||||
<li><a class="tag p-category" href="../categories/pcg.html" rel="tag">pcg</a></li>
|
||||
<li><a class="tag p-category" href="../categories/practrand.html" rel="tag">practrand</a></li>
|
||||
<li><a class="tag p-category" href="../categories/splitmix.html" rel="tag">splitmix</a></li>
|
||||
<li><a class="tag p-category" href="../categories/testing.html" rel="tag">testing</a></li>
|
||||
<li><a class="tag p-category" href="../categories/xoroshiro.html" rel="tag">xoroshiro</a></li>
|
||||
</ul>
|
||||
<ul class="pager hidden-print">
|
||||
<li class="previous">
|
||||
<a href="implausible-output-from-xoshiro256.html" rel="prev" title="Implausible Output from Xoshiro256**">Previous post</a>
|
||||
</li>
|
||||
<li class="next">
|
||||
<a href="bob-jenkins-small-prng-passes-practrand.html" rel="next" title="Bob Jenkins's Small PRNG Passes PractRand (And More!)">Next post</a>
|
||||
</li>
|
||||
</ul></nav></aside></article>
|
||||
</div>
|
||||
<!--End of body content-->
|
||||
|
||||
<footer id="footer">
|
||||
Contents © 2018 <a href="mailto:oneill@pcg-random.org">M.E. O'Neill</a> - Powered by <a href="https://getnikola.com" rel="nofollow">Nikola</a>
|
||||
|
||||
</footer>
|
||||
</div>
|
||||
</div>
|
||||
|
||||
|
||||
<script src="../assets/js/jquery.min.js"></script><script src="../assets/js/bootstrap.min.js"></script><script src="../assets/js/moment-with-locales.min.js"></script><script src="../assets/js/fancydates.js"></script><script src="../assets/js/jquery.colorbox-min.js"></script><script>$('a.image-reference:not(.islink) img:not(.islink)').parent().colorbox({rel:"gal",maxWidth:"100%",maxHeight:"100%",scalePhotos:true});</script><!-- fancy dates --><script>
|
||||
moment.locale("en");
|
||||
fancydates(0, "YYYY-MM-DD HH:mm");
|
||||
</script><!-- end fancy dates -->
|
||||
</body>
|
||||
</html>
|
||||
@ -0,0 +1,19 @@
|
||||
$ cat squiggle.c | grep http | sed 's|.*http|http|g'
|
||||
|
||||
https://en.wikipedia.org/wiki/Xorshift
|
||||
https://stackoverflow.com/questions/53886131/how-does-xorshift32-works
|
||||
https://www.pcg-random.org/posts/on-vignas-pcg-critique.html
|
||||
https://prng.di.unimi.it/
|
||||
https://en.wikipedia.org/wiki/Box%E2%80%93Muller_transform
|
||||
https://stackoverflow.com/questions/20626994/how-to-calculate-the-inverse-of-the-normal-cumulative-distribution-function-in-p
|
||||
https://www.wolframalpha.com/input?i=N%5BInverseCDF%28normal%280%2C1%29%2C+0.05%29%2C%7B%E2%88%9E%2C100%7D%5D
|
||||
https://en.wikipedia.org/wiki/Normal_distribution?lang=en#Operations_on_a_single_normal_variable
|
||||
https://dl.acm.org/doi/pdf/10.1145/358407.358414
|
||||
https://en.wikipedia.org/wiki/Gamma_distribution
|
||||
https://dl.acm.org/doi/pdf/10.1145/358407.358414
|
||||
https://en.wikipedia.org/wiki/Gamma_distribution#Related_distributions
|
||||
https://en.wikipedia.org/wiki/Beta_distribution?lang=en#Rule_of_succession
|
||||
|
||||
$ cat squiggle_more.c | grep http | sed 's|.*http|http|g'
|
||||
https://en.wikipedia.org/wiki/Quickselect
|
||||
|
||||
Loading…
Reference in new issue