smpcg32rng.hh

// This Source Code Form is licensed MPL-2.0: http://mozilla.org/MPL/2.0
// Author: 2014, Tim Janik, see http://testbit.org/keccak

#ifndef SPECTMORPH_PCG32_RNG_HH
#define SPECTMORPH_PCG32_RNG_HH

#include <glib.h>
#include <stdint.h>

namespace SpectMorph
{

class Pcg32Rng {
  uint64_t increment_;          // must be odd, allows for 2^63 distinct random sequences
  uint64_t accu_;               // can contain all 2^64 possible values
  static constexpr const uint64_t A = 6364136223846793005ULL; // from C. E. Hayness, see TAOCP by D. E. Knuth, 3.3.4, table 1, line 26.
  static inline constexpr uint32_t
  ror32 (const uint32_t bits, const uint32_t offset)
  {
    // bitwise rotate-right pattern recognized by gcc & clang iff 32==sizeof (bits)
    return (bits >> offset) | (bits << (32 - offset));
  }
  static inline constexpr uint32_t
  pcg_xsh_rr (const uint64_t input)
  {
    // Section 6.3.1. 32-bit Output, 64-bit State: PCG-XSH-RR
    // http://www.pcg-random.org/pdf/toms-oneill-pcg-family-v1.02.pdf
    return ror32 ((input ^ (input >> 18)) >> 27, input >> 59);
  }
public:
  template<class SeedSeq>
  explicit Pcg32Rng  (SeedSeq &seed_sequence) : increment_ (0), accu_ (0) { seed (seed_sequence); }
  explicit Pcg32Rng  (uint64_t offset, uint64_t sequence) :
    increment_ (0), accu_ (0)
  {
    seed (offset, sequence);
  }
  explicit Pcg32Rng  () :
    increment_ (0), accu_ (0)
  {
    auto_seed();
  }

  void auto_seed ()
  {
    //AutoSeeder seeder;
    //seed (seeder);
    seed (g_random_int(), g_random_int());
  }
  void seed (uint64_t offset, uint64_t sequence)
  {
    accu_ = sequence;
    increment_ = (sequence << 1) | 1;    // force increment_ to be odd
    accu_ += offset;
    accu_ = A * accu_ + increment_;
  }
  template<class SeedSeq> void
  seed (SeedSeq &seed_sequence)
  {
    uint64_t seeds[2];
    seed_sequence.generate (&seeds[0], &seeds[2]);
    seed (seeds[0], seeds[1]);
  }
  uint32_t
  random ()
  {
    const uint64_t lcgout = accu_;      // using the *last* state as ouput helps with CPU pipelining
    accu_ = A * accu_ + increment_;
    return pcg_xsh_rr (lcgout);         // PCG XOR-shift + random rotation
  }
};

}

#endif