libcruft-util/rand/pcg.hpp

87 lines
2.8 KiB
C++

/*
* 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/.
*
* Copyright 2019 Danny Robson <danny@nerdcruft.net>
*/
#pragma once
#include "../bitwise.hpp"
#include "../std.hpp"
#include <limits>
#include <random>
// Implementations of algorithms from the paper: "PCG: A Family of Simple
// Fast Space-Efficient Statistically Good Algorithms for Random Number
// Generation"
namespace cruft::rand {
// Implements the PCG-XSH-RR algorithm; xorshift-high, random rotation.
template <
typename StateT = u64,
typename OutputT = u32
>
class pcg_xsh_rr {
public:
using result_type = OutputT;
explicit pcg_xsh_rr (std::seed_seq seed)
: state (increment)
{
// Overlap a series of unsigned words with state words. We might
// over-provision but that shouldn't be an issue with the sizes we
// support.
union {
unsigned u[sizeof (StateT) / sizeof (unsigned) + 1];
StateT s;
} words;
seed.generate (std::begin (words.u), std::end (words.u));
state += words.s;
}
explicit pcg_xsh_rr (StateT seed) noexcept
: state (increment + seed)
{ (*this)(); }
static constexpr auto ShiftV = 5u;
static constexpr auto state_bits = sizeof (StateT ) * 8u;
static constexpr auto state_shift = sizeof (StateT ) * 8u - ShiftV;
static constexpr auto output_shift = sizeof (OutputT) * 8u - ShiftV;
static constexpr auto xshift = (state_bits - output_shift) / 2u;
result_type operator() (void) noexcept
{
StateT x = state;
StateT const rrotate = x >> state_shift;
state = x * multiplier + increment;
x ^= x >> xshift;
// Be very careful to cast to result_type before doing the
// rotation otherwise we'll rotate the zeroed upper bits into the
// returned value.
return rotater (result_type (x >> output_shift), rrotate);
}
static constexpr auto max () noexcept { return std::numeric_limits<result_type>::max (); }
static constexpr auto min () noexcept { return std::numeric_limits<result_type>::min (); }
private:
// Initialiser from the PCG example code.
StateT state = 0x4d595df4d0f33173;
// Using multiplier and increment from the PCG paper, which appear to
// also correspond to constants chosen by Knuth.
static constexpr StateT const multiplier = 6364136223846793005u;
static constexpr StateT const increment = 1442695040888963407u;
};
}