Skip to content
GitLab
Commit 87ba288b authored by BRAMAS Berenger's avatar BRAMAS Berenger
Browse files

Merge branch 'philox-random-gen' into 'master' 

Philox random gen

See merge request bramas/spetabaru!26
parents 123c779e 32402b9d
Hide whitespace changes
Inline Side-by-side
Examples/MC/mc.cpp
View file @ 87ba288b
......@@ -12,7 +12,7 @@
#include "Tasks/SpTask.hpp"
#include "Runtimes/SpRuntime.hpp"
#include "Random/SpMTGenerator.hpp"
#include "Random/SpPhiloxGenerator.hpp"
#include "Utils/small_vector.hpp"
#include "mcglobal.hpp"
......@@ -67,7 +67,7 @@ int main(){
SpTimer timerSpecNoCons[MaxidxConsecutiveSpec];
if(runSeq){
SpMTGenerator<double> randGen(0);
SpPhiloxGenerator<double> randGen(0);
small_vector<Domain<double>> domains = InitDomains<double>(NbDomains, NbParticlesPerDomain, BoxWidth, randGen);
always_assert(randGen.getNbValuesGenerated() == 3 * static_cast<size_t>(NbDomains) * static_cast<size_t>(NbParticlesPerDomain));
......@@ -132,7 +132,7 @@ int main(){
if(runTask){
SpRuntime runtime(NumThreads);
SpMTGenerator<double> randGen(0);
SpPhiloxGenerator<double> randGen(0);
timerTask.start();
small_vector<Domain<double>> domains = InitDomains<double>(NbDomains, NbParticlesPerDomain, BoxWidth, randGen);
......@@ -214,7 +214,7 @@ int main(){
return true;
});
SpMTGenerator<double> randGen(0);
SpPhiloxGenerator<double> randGen(0);
timerSpec.start();
small_vector<Domain<double>> domains = InitDomains<double>(NbDomains, NbParticlesPerDomain, BoxWidth, randGen);
......@@ -310,7 +310,7 @@ int main(){
return true;
});
SpMTGenerator<double> randGen(0);
SpPhiloxGenerator<double> randGen(0);
timerSpecNoCons[idxConsecutiveSpec].start();
small_vector<Domain<double>> domains = InitDomains<double>(NbDomains, NbParticlesPerDomain, BoxWidth, randGen);
......@@ -410,7 +410,7 @@ int main(){
return true;
});
SpMTGenerator<double> randGen(0);
SpPhiloxGenerator<double> randGen(0);
timerSpecAllReject.start();
......
Examples/MC/mcCollide.cpp
View file @ 87ba288b
......@@ -11,7 +11,7 @@
#include "Tasks/SpTask.hpp"
#include "Runtimes/SpRuntime.hpp"
#include "Random/SpMTGenerator.hpp"
#include "Random/SpPhiloxGenerator.hpp"
#include "Utils/small_vector.hpp"
#include "mcglobal.hpp"
......@@ -42,7 +42,7 @@ int main(){
const double collisionLimit = 0.00001;
if(runSeqMove){
SpMTGenerator<double> randGen(0);
SpPhiloxGenerator<double> randGen(0);
small_vector<Domain<double>> domains = InitDomains<double>(NbDomains, NbParticlesPerDomain, BoxWidth, randGen);
always_assert(randGen.getNbValuesGenerated() == 3 * NbDomains * NbParticlesPerDomain);
......@@ -118,7 +118,7 @@ int main(){
if(runTaskMove){
SpRuntime runtime(NumThreads);
SpMTGenerator<double> randGen(0);
SpPhiloxGenerator<double> randGen(0);
small_vector<Domain<double>> domains = InitDomains<double>(NbDomains, NbParticlesPerDomain, BoxWidth, randGen);
always_assert(randGen.getNbValuesGenerated() == 3 * NbDomains * NbParticlesPerDomain);
......@@ -229,7 +229,7 @@ int main(){
if(runSpecMove){
SpRuntime runtime(NumThreads);
SpMTGenerator<double> randGen(0);
SpPhiloxGenerator<double> randGen(0);
small_vector<Domain<double>> domains = InitDomains<double>(NbDomains, NbParticlesPerDomain, BoxWidth, randGen);
always_assert(randGen.getNbValuesGenerated() == 3 * NbDomains * NbParticlesPerDomain);
......
Examples/MC/remc.cpp
View file @ 87ba288b
......@@ -12,7 +12,7 @@
#include "Tasks/SpTask.hpp"
#include "Runtimes/SpRuntime.hpp"
#include "Random/SpMTGenerator.hpp"
#include "Random/SpPhiloxGenerator.hpp"
#include "Utils/small_vector.hpp"
#include "mcglobal.hpp"
......@@ -79,7 +79,7 @@ int main(){
SpTimer timerSpecNoCons[MaxidxConsecutiveSpec];
if(runSeq){
std::array<SpMTGenerator<double>,NbReplicas> replicaRandGen;
std::array<SpPhiloxGenerator<double>,NbReplicas> replicaRandGen;
std::array<small_vector<Domain<double>>,NbReplicas> replicaDomains;
std::array<Matrix<double>,NbReplicas> replicaEnergyAll;
std::array<size_t,NbReplicas> replicaCptGenerated;
......@@ -87,12 +87,12 @@ int main(){
timerSeq.start();
for(int idxReplica = 0 ; idxReplica < NbReplicas ; ++idxReplica){
SpMTGenerator<double>& randGen = replicaRandGen[idxReplica];
SpPhiloxGenerator<double>& randGen = replicaRandGen[idxReplica];
auto& domains = replicaDomains[idxReplica];
Matrix<double>& energyAll = replicaEnergyAll[idxReplica];
size_t& cptGenerated = replicaCptGenerated[idxReplica];
randGen = SpMTGenerator<double>(0/*idxReplica*/);
randGen = SpPhiloxGenerator<double>(0/*idxReplica*/);
domains = InitDomains<double>(NbDomains, NbParticlesPerDomain, BoxWidth, randGen);
always_assert(randGen.getNbValuesGenerated() == 3 * static_cast<size_t>(NbDomains) * static_cast<size_t>(NbParticlesPerDomain));
......@@ -107,7 +107,7 @@ int main(){
for(int idxLoop = 0 ; idxLoop < NbLoops ; idxLoop += NbInnerLoops){
const int NbInnerLoopsLimit = std::min(NbInnerLoops, NbLoops-idxLoop) + idxLoop;
for(int idxReplica = 0 ; idxReplica < NbReplicas ; ++idxReplica){
SpMTGenerator<double>& randGen = replicaRandGen[idxReplica];
SpPhiloxGenerator<double>& randGen = replicaRandGen[idxReplica];
auto& domains = replicaDomains[idxReplica];
Matrix<double>& energyAll = replicaEnergyAll[idxReplica];
size_t& cptGenerated = replicaCptGenerated[idxReplica];
......@@ -187,18 +187,18 @@ int main(){
if(runTask){
SpRuntime runtime(NumThreads);
std::array<SpMTGenerator<double>,NbReplicas> replicaRandGen;
std::array<SpPhiloxGenerator<double>,NbReplicas> replicaRandGen;
std::array<small_vector<Domain<double>>,NbReplicas> replicaDomains;
std::array<Matrix<double>,NbReplicas> replicaEnergyAll;
timerTask.start();
for(int idxReplica = 0 ; idxReplica < NbReplicas ; ++idxReplica){
SpMTGenerator<double>& randGen = replicaRandGen[idxReplica];
SpPhiloxGenerator<double>& randGen = replicaRandGen[idxReplica];
auto& domains = replicaDomains[idxReplica];
Matrix<double>& energyAll = replicaEnergyAll[idxReplica];
randGen = SpMTGenerator<double>(0/*idxReplica*/);
randGen = SpPhiloxGenerator<double>(0/*idxReplica*/);
domains = InitDomains<double>(NbDomains, NbParticlesPerDomain, BoxWidth, randGen);
......@@ -214,7 +214,7 @@ int main(){
for(int idxLoop = 0 ; idxLoop < NbLoops ; idxLoop += NbInnerLoops){
const int NbInnerLoopsLimit = std::min(NbInnerLoops, NbLoops-idxLoop) + idxLoop;
for(int idxReplica = 0 ; idxReplica < NbReplicas ; ++idxReplica){
SpMTGenerator<double>& randGen = replicaRandGen[idxReplica];
SpPhiloxGenerator<double>& randGen = replicaRandGen[idxReplica];
auto& domains = replicaDomains[idxReplica];
Matrix<double>& energyAll = replicaEnergyAll[idxReplica];
const double& Temperature = temperatures[idxReplica];
......@@ -272,7 +272,7 @@ int main(){
int* nbExchanges = new int(0);
const int startExchangeIdx = ((idxLoop/NbInnerLoops)&1);
for(int idxReplica = startExchangeIdx ; idxReplica+1 < NbReplicas ; idxReplica += 2){
SpMTGenerator<double>& randGen0 = replicaRandGen[idxReplica];
SpPhiloxGenerator<double>& randGen0 = replicaRandGen[idxReplica];
auto& domains0 = replicaDomains[idxReplica];
Matrix<double>& energyAll0 = replicaEnergyAll[idxReplica];
auto& domains1 = replicaDomains[idxReplica+1];
......@@ -332,18 +332,18 @@ int main(){
return true;
});
std::array<SpMTGenerator<double>,NbReplicas> replicaRandGen;
std::array<SpPhiloxGenerator<double>,NbReplicas> replicaRandGen;
std::array<small_vector<Domain<double>>,NbReplicas> replicaDomains;
std::array<Matrix<double>,NbReplicas> replicaEnergyAll;
timerSpec.start();
for(int idxReplica = 0 ; idxReplica < NbReplicas ; ++idxReplica){
SpMTGenerator<double>& randGen = replicaRandGen[idxReplica];
SpPhiloxGenerator<double>& randGen = replicaRandGen[idxReplica];
auto& domains = replicaDomains[idxReplica];
Matrix<double>& energyAll = replicaEnergyAll[idxReplica];
randGen = SpMTGenerator<double>(0/*idxReplica*/);
randGen = SpPhiloxGenerator<double>(0/*idxReplica*/);
domains = InitDomains<double>(NbDomains, NbParticlesPerDomain, BoxWidth, randGen);
......@@ -359,7 +359,7 @@ int main(){
for(int idxLoop = 0 ; idxLoop < NbLoops ; idxLoop += NbInnerLoops){
const int NbInnerLoopsLimit = std::min(NbInnerLoops, NbLoops-idxLoop) + idxLoop;
for(int idxReplica = 0 ; idxReplica < NbReplicas ; ++idxReplica){
SpMTGenerator<double>& randGen = replicaRandGen[idxReplica];
SpPhiloxGenerator<double>& randGen = replicaRandGen[idxReplica];
auto& domains = replicaDomains[idxReplica];
Matrix<double>& energyAll = replicaEnergyAll[idxReplica];
const double& Temperature = temperatures[idxReplica];
......@@ -426,7 +426,7 @@ int main(){
// TODO int* nbExchanges = new int(0);
const int startExchangeIdx = ((idxLoop/NbInnerLoops)&1);
for(int idxReplica = startExchangeIdx ; idxReplica+1 < NbReplicas ; idxReplica += 2){
SpMTGenerator<double>& randGen0 = replicaRandGen[idxReplica];
SpPhiloxGenerator<double>& randGen0 = replicaRandGen[idxReplica];
auto& domains0 = replicaDomains[idxReplica];
Matrix<double>& energyAll0 = replicaEnergyAll[idxReplica];
auto& domains1 = replicaDomains[idxReplica+1];
......@@ -495,18 +495,18 @@ int main(){
return true;
});
std::array<SpMTGenerator<double>,NbReplicas> replicaRandGen;
std::array<SpPhiloxGenerator<double>,NbReplicas> replicaRandGen;
std::array<small_vector<Domain<double>>,NbReplicas> replicaDomains;
std::array<Matrix<double>,NbReplicas> replicaEnergyAll;
timerSpecNoCons[idxConsecutiveSpec].start();
for(int idxReplica = 0 ; idxReplica < NbReplicas ; ++idxReplica){
SpMTGenerator<double>& randGen = replicaRandGen[idxReplica];
SpPhiloxGenerator<double>& randGen = replicaRandGen[idxReplica];
auto& domains = replicaDomains[idxReplica];
Matrix<double>& energyAll = replicaEnergyAll[idxReplica];
randGen = SpMTGenerator<double>(0/*idxReplica*/);
randGen = SpPhiloxGenerator<double>(0/*idxReplica*/);
domains = InitDomains<double>(NbDomains, NbParticlesPerDomain, BoxWidth, randGen);
......@@ -524,7 +524,7 @@ int main(){
for(int idxLoop = 0 ; idxLoop < NbLoops ; idxLoop += NbInnerLoops){
const int NbInnerLoopsLimit = std::min(NbInnerLoops, NbLoops-idxLoop) + idxLoop;
for(int idxReplica = 0 ; idxReplica < NbReplicas ; ++idxReplica){
SpMTGenerator<double>& randGen = replicaRandGen[idxReplica];
SpPhiloxGenerator<double>& randGen = replicaRandGen[idxReplica];
auto& domains = replicaDomains[idxReplica];
Matrix<double>& energyAll = replicaEnergyAll[idxReplica];
const double& Temperature = temperatures[idxReplica];
......@@ -594,7 +594,7 @@ int main(){
// TODO int* nbExchanges = new int(0);
const int startExchangeIdx = ((idxLoop/NbInnerLoops)&1);
for(int idxReplica = startExchangeIdx ; idxReplica+1 < NbReplicas ; idxReplica += 2){
SpMTGenerator<double>& randGen0 = replicaRandGen[idxReplica];
SpPhiloxGenerator<double>& randGen0 = replicaRandGen[idxReplica];
auto& domains0 = replicaDomains[idxReplica];
Matrix<double>& energyAll0 = replicaEnergyAll[idxReplica];
auto& domains1 = replicaDomains[idxReplica+1];
......@@ -664,18 +664,18 @@ int main(){
return true;
});
std::array<SpMTGenerator<double>,NbReplicas> replicaRandGen;
std::array<SpPhiloxGenerator<double>,NbReplicas> replicaRandGen;
std::array<small_vector<Domain<double>>,NbReplicas> replicaDomains;
std::array<Matrix<double>,NbReplicas> replicaEnergyAll;
timerSpecAllReject.start();
for(int idxReplica = 0 ; idxReplica < NbReplicas ; ++idxReplica){
SpMTGenerator<double>& randGen = replicaRandGen[idxReplica];
SpPhiloxGenerator<double>& randGen = replicaRandGen[idxReplica];
auto& domains = replicaDomains[idxReplica];
Matrix<double>& energyAll = replicaEnergyAll[idxReplica];
randGen = SpMTGenerator<double>(0/*idxReplica*/);
randGen = SpPhiloxGenerator<double>(0/*idxReplica*/);
domains = InitDomains<double>(NbDomains, NbParticlesPerDomain, BoxWidth, randGen);
......@@ -691,7 +691,7 @@ int main(){
for(int idxLoop = 0 ; idxLoop < NbLoops ; idxLoop += NbInnerLoops){
const int NbInnerLoopsLimit = std::min(NbInnerLoops, NbLoops-idxLoop) + idxLoop;
for(int idxReplica = 0 ; idxReplica < NbReplicas ; ++idxReplica){
SpMTGenerator<double>& randGen = replicaRandGen[idxReplica];
SpPhiloxGenerator<double>& randGen = replicaRandGen[idxReplica];
auto& domains = replicaDomains[idxReplica];
Matrix<double>& energyAll = replicaEnergyAll[idxReplica];
......@@ -767,7 +767,7 @@ int main(){
// TODO int* nbExchanges = new int(0);
const int startExchangeIdx = ((idxLoop/NbInnerLoops)&1);
for(int idxReplica = startExchangeIdx ; idxReplica+1 < NbReplicas ; idxReplica += 2){
SpMTGenerator<double>& randGen0 = replicaRandGen[idxReplica];
SpPhiloxGenerator<double>& randGen0 = replicaRandGen[idxReplica];
auto& domains0 = replicaDomains[idxReplica];
Matrix<double>& energyAll0 = replicaEnergyAll[idxReplica];
auto& domains1 = replicaDomains[idxReplica+1];
......
Examples/MC/remcCollide.cpp
View file @ 87ba288b
......@@ -11,7 +11,7 @@
#include "Tasks/SpTask.hpp"
#include "Runtimes/SpRuntime.hpp"
#include "Random/SpMTGenerator.hpp"
#include "Random/SpPhiloxGenerator.hpp"
#include "Utils/small_vector.hpp"
#include "mcglobal.hpp"
......@@ -56,18 +56,18 @@ int main(){
const double collisionLimit = 0.00001;
if(runSeqMove){
std::array<SpMTGenerator<double>,NbReplicas> replicaRandGen;
std::array<SpPhiloxGenerator<double>,NbReplicas> replicaRandGen;
std::array<small_vector<Domain<double>>,NbReplicas> replicaDomains;
std::array<Matrix<double>,NbReplicas> replicaEnergyAll;
std::array<size_t,NbReplicas> replicaCptGenerated;
for(int idxReplica = 0 ; idxReplica < NbReplicas ; ++idxReplica){
SpMTGenerator<double>& randGen = replicaRandGen[idxReplica];
SpPhiloxGenerator<double>& randGen = replicaRandGen[idxReplica];
auto& domains = replicaDomains[idxReplica];
Matrix<double>& energyAll = replicaEnergyAll[idxReplica];
size_t& cptGenerated = replicaCptGenerated[idxReplica];
randGen = SpMTGenerator<double>(0/*idxReplica*/);
randGen = SpPhiloxGenerator<double>(0/*idxReplica*/);
domains = InitDomains<double>(NbDomains, NbParticlesPerDomain, BoxWidth, randGen);
always_assert(randGen.getNbValuesGenerated() == 3 * NbDomains * NbParticlesPerDomain);
......@@ -82,7 +82,7 @@ int main(){
for(int idxLoop = 0 ; idxLoop < NbLoops ; idxLoop += NbInnerLoops){
const int NbInnerLoopsLimit = std::min(NbInnerLoops, NbLoops-idxLoop) + idxLoop;
for(int idxReplica = 0 ; idxReplica < NbReplicas ; ++idxReplica){
SpMTGenerator<double>& randGen = replicaRandGen[idxReplica];
SpPhiloxGenerator<double>& randGen = replicaRandGen[idxReplica];
auto& domains = replicaDomains[idxReplica];
Matrix<double>& energyAll = replicaEnergyAll[idxReplica];
size_t& cptGenerated = replicaCptGenerated[idxReplica];
......@@ -178,16 +178,16 @@ int main(){
if(runTaskMove){
SpRuntime runtime(NumThreads);
std::array<SpMTGenerator<double>,NbReplicas> replicaRandGen;
std::array<SpPhiloxGenerator<double>,NbReplicas> replicaRandGen;
std::array<small_vector<Domain<double>>,NbReplicas> replicaDomains;
std::array<Matrix<double>,NbReplicas> replicaEnergyAll;
for(int idxReplica = 0 ; idxReplica < NbReplicas ; ++idxReplica){
SpMTGenerator<double>& randGen = replicaRandGen[idxReplica];
SpPhiloxGenerator<double>& randGen = replicaRandGen[idxReplica];
auto& domains = replicaDomains[idxReplica];
Matrix<double>& energyAll = replicaEnergyAll[idxReplica];
randGen = SpMTGenerator<double>(0/*idxReplica*/);
randGen = SpPhiloxGenerator<double>(0/*idxReplica*/);
domains = InitDomains<double>(NbDomains, NbParticlesPerDomain, BoxWidth, randGen);
......@@ -202,7 +202,7 @@ int main(){
for(int idxLoop = 0 ; idxLoop < NbLoops ; idxLoop += NbInnerLoops){
const int NbInnerLoopsLimit = std::min(NbInnerLoops, NbLoops-idxLoop) + idxLoop;
for(int idxReplica = 0 ; idxReplica < NbReplicas ; ++idxReplica){
SpMTGenerator<double>& randGen = replicaRandGen[idxReplica];
SpPhiloxGenerator<double>& randGen = replicaRandGen[idxReplica];
auto& domains = replicaDomains[idxReplica];
Matrix<double>& energyAll = replicaEnergyAll[idxReplica];
const double& Temperature = temperatures[idxReplica];
......@@ -295,7 +295,7 @@ int main(){
int* nbExchanges = new int(0);
const int startExchangeIdx = ((idxLoop/NbInnerLoops)&1);
for(int idxReplica = startExchangeIdx ; idxReplica+1 < NbReplicas ; idxReplica += 2){
SpMTGenerator<double>& randGen0 = replicaRandGen[idxReplica];
SpPhiloxGenerator<double>& randGen0 = replicaRandGen[idxReplica];
auto& domains0 = replicaDomains[idxReplica];
Matrix<double>& energyAll0 = replicaEnergyAll[idxReplica];
auto& domains1 = replicaDomains[idxReplica+1];
......@@ -343,16 +343,16 @@ int main(){
if(runSpecMove){
SpRuntime runtime(NumThreads);
std::array<SpMTGenerator<double>,NbReplicas> replicaRandGen;
std::array<SpPhiloxGenerator<double>,NbReplicas> replicaRandGen;
std::array<small_vector<Domain<double>>,NbReplicas> replicaDomains;
std::array<Matrix<double>,NbReplicas> replicaEnergyAll;
for(int idxReplica = 0 ; idxReplica < NbReplicas ; ++idxReplica){
SpMTGenerator<double>& randGen = replicaRandGen[idxReplica];
SpPhiloxGenerator<double>& randGen = replicaRandGen[idxReplica];
auto& domains = replicaDomains[idxReplica];
Matrix<double>& energyAll = replicaEnergyAll[idxReplica];
randGen = SpMTGenerator<double>(0/*idxReplica*/);
randGen = SpPhiloxGenerator<double>(0/*idxReplica*/);
domains = InitDomains<double>(NbDomains, NbParticlesPerDomain, BoxWidth, randGen);
......@@ -367,7 +367,7 @@ int main(){
for(int idxLoop = 0 ; idxLoop < NbLoops ; idxLoop += NbInnerLoops){
const int NbInnerLoopsLimit = std::min(NbInnerLoops, NbLoops-idxLoop) + idxLoop;
for(int idxReplica = 0 ; idxReplica < NbReplicas ; ++idxReplica){
SpMTGenerator<double>& randGen = replicaRandGen[idxReplica];
SpPhiloxGenerator<double>& randGen = replicaRandGen[idxReplica];
auto& domains = replicaDomains[idxReplica];
Matrix<double>& energyAll = replicaEnergyAll[idxReplica];
const double& Temperature = temperatures[idxReplica];
......@@ -464,7 +464,7 @@ int main(){
int* nbExchanges = new int(0);
const int startExchangeIdx = ((idxLoop/NbInnerLoops)&1);
for(int idxReplica = startExchangeIdx ; idxReplica+1 < NbReplicas ; idxReplica += 2){
SpMTGenerator<double>& randGen0 = replicaRandGen[idxReplica];
SpPhiloxGenerator<double>& randGen0 = replicaRandGen[idxReplica];
auto& domains0 = replicaDomains[idxReplica];
Matrix<double>& energyAll0 = replicaEnergyAll[idxReplica];
auto& domains1 = replicaDomains[idxReplica+1];
......
Src/Random/SpPhiloxGenerator.hpp 0 → 100644
View file @ 87ba288b
///////////////////////////////////////////////////////////////////////////
// Thomas Millot (c), Unistra, 2020
// Under LGPL Licence, please you must read the LICENCE file.
///////////////////////////////////////////////////////////////////////////
#ifndef SPPHILOXGENERATOR_HPP
#define SPPHILOXGENERATOR_HPP
#include <iostream>
#include <random>
#include <array>
// Implementation of the Philox algorithm to generate random numbers in parallel.
// http://www.thesalmons.org/john/random123/papers/random123sc11.pdf
// Also based on an implementation of the algorithm by Tensorflow.
// https://github.com/tensorflow/tensorflow/blob/master/tensorflow/core/lib/random/philox_random.h
// It's the Philox-4×32 version, meaning 4 32bits random numbers each time.
// The number of cycles can be user defined, by default it's 10.
// The engine satisfies C++ named requirements RandomNumberDistribution,
// so it can be used with std::uniform_real_distribution for example.
// Exactly the same interface as SpMTGenerator, so it can be interchangeable
template<class RealType = double>
class SpPhiloxGenerator {
class philox4x32 {
typedef uint_fast32_t uint32;
typedef uint_fast64_t uint64;
static constexpr int DEFAULT_CYCLES = 10;
public:
// An array of four uint32, the results of the philox4 engine
using Result = std::array<uint32, 4>;
// 64-bit seed stored in two uint32
using Key = std::array<uint32, 2>;
philox4x32() = default;
explicit philox4x32(uint64 seed, int cycles = DEFAULT_CYCLES) {
// Splitting the seed in two
key_[0] = static_cast<uint32>(seed);
key_[1] = static_cast<uint32>(seed >> 32);
counter_.fill(0);
temp_results_.fill(0);
cycles_ = cycles;
operatorPPcounter = 0;
}
// Returns the minimum value productible by the engine
static constexpr uint32 min() { return _Min; }
// Returns the maximum value productible by the engine
static constexpr uint32 max() { return _Max; }
// Skip the specified number of steps
void Skip(uint64 count) {
if(count > 0) {
skipOnWrap();
const auto position = temp_counter_ + count;
if (position > 3) {
force_compute_ = true;
temp_counter_ = position % 4;
const auto nbOfCounterIncrements = position / 4;
const auto count_lo = static_cast<uint32>(nbOfCounterIncrements);
auto count_hi = static_cast<uint32>(nbOfCounterIncrements >> 32);
// 128 bit add
counter_[0] += count_lo;
if (counter_[0] < count_lo) {
++count_hi;
}
counter_[1] += count_hi;
if (counter_[1] < count_hi) {
if (++counter_[2] == 0) {
++counter_[3];
}
}
} else {
temp_counter_ = position;
}
}
}
// Returns an random number using the philox engine
uint32 operator()() {
operatorPPcounter += 1;
skipOnWrap();
if (force_compute_) {
force_compute_ = false;
temp_results_ = counter_;
ExecuteRounds();
}
uint32 value = temp_results_[temp_counter_];
temp_counter_++;
return value;
}
int getOperatorPPCounter() const{
return operatorPPcounter;
}
private:
// Using the same constants as recommended in the original paper.
static constexpr uint32 kPhiloxW32A = 0x9E3779B9;
static constexpr uint32 kPhiloxW32B = 0xBB67AE85;
static constexpr uint32 kPhiloxM4x32A = 0xD2511F53;
static constexpr uint32 kPhiloxM4x32B = 0xCD9E8D57;
// The minimum return value
static constexpr uint32 _Min = 0;
// The maximum return value
static constexpr uint32 _Max = UINT_FAST32_MAX;
// The counter for the current state of the engine
Result counter_;
// Keeping the last to results to improve performances during consecutive call
Result temp_results_;
// The split seed
Key key_;
// To iterate through the temp_results_
uint64 temp_counter_ = 0;
// The number of cycles used to generate randomness
int cycles_;
// To force the engine to compute the rounds to populates temp_results_
bool force_compute_ = true;
// The number of times operator () is called to ensure that the STL
// always call it once
int operatorPPcounter;
// Skip one step
void SkipOne() {
// 128 bit increment
if (++counter_[0] == 0) {
if (++counter_[1] == 0) {
if (++counter_[2] == 0) {
++counter_[3];
}
}
}
}
// Helper function to return the lower and higher 32-bits from two 32-bit integer multiplications.
static void MultiplyHighLow(uint32 a, uint32 b, uint32 *result_low, uint32 *result_high) {
const uint64 product = static_cast<uint64>(a) * b;
*result_low = static_cast<uint32>(product);
*result_high = static_cast<uint32>(product >> 32);
}
void ExecuteRounds() {
Key key = key_;
// Run the single rounds for ten times.
for (int i = 0; i < cycles_; ++i) {
temp_results_ = ComputeSingleRound(temp_results_, key);
RaiseKey(&key);
}
}
// Helper function for a single round of the underlying Philox algorithm.
static Result ComputeSingleRound(const Result &counter, const Key &key) {
uint32 lo0;
uint32 hi0;
MultiplyHighLow(kPhiloxM4x32A, counter[0], &lo0, &hi0);
uint32 lo1;
uint32 hi1;
MultiplyHighLow(kPhiloxM4x32B, counter[2], &lo1, &hi1);
Result result;
result[0] = hi1 ^ counter[1] ^ key[0];
result[1] = lo1;
result[2] = hi0 ^ counter[3] ^ key[1];
result[3] = lo0;
return result;
}
void RaiseKey(Key *key) {
(*key)[0] += kPhiloxW32A;
(*key)[1] += kPhiloxW32B;
}
void skipOnWrap() {
if(temp_counter_ > 3) {
temp_counter_ = 0;
SkipOne();