Bug 1206356: Add mfbt/Random.h, implementing the xorshift128+ random number generator. r=waldo

mfbt/XorShift128PlusRNG.h

@@ -0,0 +1,95 @@
+/* -*- Mode: C++; tab-width: 8; indent-tabs-mode: nil; c-basic-offset: 2 -*- */
+/* vim: set ts=8 sts=2 et sw=2 tw=80: */
+/* This Source Code Form is subject to the terms of the Mozilla Public
+ * License, v. 2.0. If a copy of the MPL was not distributed with this
+ * file, You can obtain one at http://mozilla.org/MPL/2.0/. */
+
+/* The xorshift128+ pseudo-random number generator. */
+
+#ifndef mozilla_XorShift128Plus_h
+#define mozilla_XorShift128Plus_h
+
+#include "mozilla/Assertions.h"
+#include "mozilla/FloatingPoint.h"
+
+#include <math.h>
+#include <memory.h>
+#include <stdint.h>
+
+namespace mozilla {
+
+/*
+ * A stream of pseudo-random numbers generated using the xorshift+ technique
+ * described here:
+ *
+ * Vigna, Sebastiano (2014). "Further scramblings of Marsaglia's xorshift
+ * generators". arXiv:1404.0390 (http://arxiv.org/abs/1404.0390)
+ *
+ * That paper says:
+ *
+ *     In particular, we propose a tightly coded xorshift128+ generator that
+ *     does not fail systematically any test from the BigCrush suite of TestU01
+ *     (even reversed) and generates 64 pseudorandom bits in 1.10 ns on an
+ *     Intel(R) Core(TM) i7-4770 CPU @3.40GHz (Haswell). It is the fastest
+ *     generator we are aware of with such empirical statistical properties.
+ *
+ * This generator is not suitable as a cryptographically secure random number
+ * generator.
+ */
+class XorShift128PlusRNG {
+  uint64_t mState[2];
+
+ public:
+  /*
+   * Construct a xorshift128+ pseudo-random number stream using |aInitial0| and
+   * |aInitial1| as the initial state. These may not both be zero; ideally, they
+   * should have an almost even mix of zero and one bits.
+   */
+  XorShift128PlusRNG(uint64_t aInitial0, uint64_t aInitial1) {
+    setState(aInitial0, aInitial1);
+  }
+
+  /* Return a pseudo-random 64-bit number. */
+  uint64_t next() {
+    uint64_t s1 = mState[0];
+    const uint64_t s0 = mState[1];
+    mState[0] = s0;
+    s1 ^= s1 << 23;
+    mState[1] = s1 ^ s0 ^ (s1 >> 17) ^ (s0 >> 26);
+    return mState[1] + s0;
+  }
+
+  /*
+   * Return a pseudo-random floating-point value in the range [0, 1).
+   * More precisely, choose an integer in the range [0, 2**53) and
+   * divide it by 2**53.
+   */
+  double nextDouble() {
+    /*
+     * Because the IEEE 64-bit floating point format stores the leading '1' bit
+     * of the mantissa implicitly, it effectively represents a mantissa in the
+     * range [0, 2**53) in only 52 bits. FloatingPoint<double>::kExponentShift
+     * is the width of the bitfield in the in-memory format, so we must add one
+     * to get the mantissa's range.
+     */
+    static const int kMantissaBits =
+      mozilla::FloatingPoint<double>::kExponentShift + 1;
+    uint64_t mantissa = next() & ((1ULL << kMantissaBits) - 1);
+    return ldexp(mantissa, -kMantissaBits);
+  }
+
+  /*
+   * Set the stream's current state to |aState0| and |aState1|. These must not
+   * both be zero; ideally, they should have an almost even mix of zero and one
+   * bits.
+   */
+  void setState(uint64_t aState0, uint64_t aState1) {
+    MOZ_ASSERT(aState0 || aState1);
+    mState[0] = aState0;
+    mState[1] = aState1;
+  }
+};
+
+} // namespace mozilla
+
+#endif // mozilla_XorShift128Plus_h

mfbt/moz.build

@@ -91,6 +91,7 @@ EXPORTS.mozilla = [
     'Variant.h',
     'Vector.h',
     'WeakPtr.h',
+    'XorShift128PlusRNG.h',
     'unused.h',
 ]
 

mfbt/tests/TestXorShift128PlusRNG.cpp

@@ -0,0 +1,110 @@
+/* -*- Mode: C++; tab-width: 8; indent-tabs-mode: nil; c-basic-offset: 2 -*- */
+/* vim: set ts=8 sts=2 et sw=2 tw=80: */
+/* This Source Code Form is subject to the terms of the Mozilla Public
+ * License, v. 2.0. If a copy of the MPL was not distributed with this file,
+ * You can obtain one at http://mozilla.org/MPL/2.0/. */
+
+#include <math.h>
+
+#include "mozilla/Assertions.h"
+#include "mozilla/PodOperations.h"
+#include "mozilla/XorShift128PlusRNG.h"
+
+static void
+TestDumbSequence()
+{
+  mozilla::XorShift128PlusRNG rng(1, 4);
+
+  // Calculated by hand following the algorithm given in the paper. The upper
+  // bits are mostly zero because we started with a poor seed; once it has run
+  // for a while, we'll get an even mix of ones and zeros in all 64 bits.
+  MOZ_RELEASE_ASSERT(rng.next() == 0x800049);
+  MOZ_RELEASE_ASSERT(rng.next() == 0x3000186);
+  MOZ_RELEASE_ASSERT(rng.next() == 0x400003001145);
+
+  // Using ldexp here lets us write out the mantissa in hex, so we can compare
+  // them with the results generated by hand.
+  MOZ_RELEASE_ASSERT(rng.nextDouble() == ldexp(0x1400003105049, -53));
+  MOZ_RELEASE_ASSERT(rng.nextDouble() == ldexp(0x2000802e49146, -53));
+  MOZ_RELEASE_ASSERT(rng.nextDouble() == ldexp(0x248300468544d, -53));
+}
+
+static size_t
+Population(uint64_t n)
+{
+  size_t pop = 0;
+
+  while (n > 0) {
+    n &= n-1; // Clear the rightmost 1-bit in n.
+    pop++;
+  }
+
+  return pop;
+}
+
+static void
+TestPopulation()
+{
+  mozilla::XorShift128PlusRNG rng(698079309544035222ULL, 6012389156611637584ULL);
+
+  // Give it some time to warm up; it should tend towards more
+  // even distributions of zeros and ones.
+  for (size_t i = 0; i < 40; i++)
+    rng.next();
+
+  for (size_t i = 0; i < 40; i++) {
+    size_t pop = Population(rng.next());
+    MOZ_RELEASE_ASSERT(24 <= pop && pop <= 40);
+  }
+}
+
+static void
+TestSetState()
+{
+  static const uint64_t seed[2] = { 1795644156779822404ULL, 14162896116325912595ULL };
+  mozilla::XorShift128PlusRNG rng(seed[0], seed[1]);
+
+  const size_t n = 10;
+  uint64_t log[n];
+
+  for (size_t i = 0; i < n; i++)
+    log[i] = rng.next();
+
+  rng.setState(seed[0], seed[1]);
+
+  for (size_t i = 0; i < n; i++)
+    MOZ_RELEASE_ASSERT(log[i] == rng.next());
+}
+
+static void
+TestDoubleDistribution()
+{
+  mozilla::XorShift128PlusRNG rng(0xa207aaede6859736, 0xaca6ca5060804791);
+
+  const size_t n = 100;
+  size_t bins[n];
+  mozilla::PodArrayZero(bins);
+
+  // This entire file runs in 0.006s on my laptop. Generating
+  // more numbers lets us put tighter bounds on the bins.
+  for (size_t i = 0; i < 100000; i++) {
+    double d = rng.nextDouble();
+    MOZ_RELEASE_ASSERT(0.0 <= d && d < 1.0);
+    bins[(int) (d * n)]++;
+  }
+
+  for (size_t i = 0; i < n; i++) {
+    MOZ_RELEASE_ASSERT(900 <= bins[i] && bins[i] <= 1100);
+  }
+}
+
+int
+main()
+{
+  TestDumbSequence();
+  TestPopulation();
+  TestSetState();
+  TestDoubleDistribution();
+
+  return 0;
+}

mfbt/tests/moz.build

@@ -40,6 +40,7 @@ CppUnitTests([
     'TestVariant',
     'TestVector',
     'TestWeakPtr',
+    'TestXorShift128PlusRNG',
 ])
 
 if not CONFIG['MOZ_ASAN']: