randomGenerator: Global flag

A random generator can now be constructed with a global flag. If the
flag is false then the provided seed will be randomised across the
different processes. If the flag is true then the synchronicity of the
generators state will be checked when performing certain operations in
debug mode.

applications/test/parallel/Test-parallel.C

@@ -56,7 +56,7 @@ int main(int argc, char *argv[])
 
     if (true)
     {
-        randomGenerator rndGen(43544*Pstream::myProcNo());
+        randomGenerator rndGen(43544);
 
         // Generate random data.
         List<Tuple2<label, List<scalar>>> complexData(100);

applications/test/syncTools/Test-syncTools.C

@@ -593,7 +593,7 @@ int main(int argc, char *argv[])
     #include "createPolyMesh.H"
 
 
-    randomGenerator rndGen(5341*(Pstream::myProcNo()+1));
+    randomGenerator rndGen(5341);
 
 
     // Face sync

src/OpenFOAM/primitives/randomGenerator/randomGenerator.C

@@ -28,17 +28,75 @@ License
 
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
-Foam::randomGenerator::randomGenerator(Istream& is)
+Foam::randomGenerator::randomGenerator(Istream& is, const bool global)
 :
+    global_(global),
     x_(pTraits<uint64_t>(is))
+{
+    checkSync();
+}
+
+
+Foam::randomGenerator::randomGenerator
+(
+    const word& name,
+    const dictionary& dict,
+    randomGenerator&& defaultRndGen
+)
+:
+    global_(defaultRndGen.global_),
+    x_
+    (
+        dict.found(name)
+      ? dict.lookup<uint64_t>(name)
+      : dict.found(name + "Seed")
+      ? seed(dict.lookup<label>(name + "Seed")).x(global_)
+      : defaultRndGen.x_
+    )
 {}
 
 
+Foam::randomGenerator::randomGenerator
+(
+    const word& name,
+    const dictionary& dict,
+    const seed defaultS,
+    const bool global
+)
+:
+    randomGenerator(name, dict, randomGenerator(defaultS, global))
+{}
+
+
+// * * * * * * * * * * * * * * Member Operators  * * * * * * * * * * * * * * //
+
+void Foam::randomGenerator::operator=(const randomGenerator& rndGen)
+{
+    if (global_ != rndGen.global_)
+    {
+        FatalErrorInFunction
+            << "Attempted assignment of a " << (global_ ? "" : "non-")
+            << "global random generator to a " << (rndGen.global_ ? "" : "non-")
+            << "global random generator"
+            << exit(FatalError);
+    }
+    x_ = rndGen.x_;
+    checkSync();
+}
+
+
+void Foam::randomGenerator::operator=(randomGenerator&& rndGen)
+{
+    *this = static_cast<const randomGenerator&>(rndGen);
+}
+
+
 // * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * * //
 
 Foam::Istream& Foam::operator>>(Istream& is, randomGenerator& rndGen)
 {
     is >> rndGen.x_;
+    rndGen.checkSync();
     is.check("operator>>(Istream& is, randomGenerator& rndGen)");
     return is;
 }
@@ -46,10 +104,19 @@ Foam::Istream& Foam::operator>>(Istream& is, randomGenerator& rndGen)
 
 Foam::Ostream& Foam::operator<<(Ostream& os, const randomGenerator& rndGen)
 {
+    rndGen.checkSync();
     os << rndGen.x_;
     os.check("operator<<(Ostream& os, const randomGenerator& rndGen)");
     return os;
 }
 
 
+// * * * * * * * * * * * * * * * Global Functions  * * * * * * * * * * * * * //
+
+void Foam::writeEntry(Ostream& os, const randomGenerator& rndGen)
+{
+    os << rndGen;
+}
+
+
 // ************************************************************************* //

src/OpenFOAM/primitives/randomGenerator/randomGenerator.H

@@ -67,6 +67,43 @@ Ostream& operator<<(Ostream&, const randomGenerator&);
 
 class randomGenerator
 {
+public:
+
+    // Public Classes
+
+        //- Seed class
+        class seed
+        {
+            // Private Data
+
+                //- The seed value
+                const uint64_t s_;
+
+
+            // Private Member Functions
+
+                //- Return the initial integer
+                inline uint64_t x(const bool global) const;
+
+
+        public:
+
+            //- Allow randomGenerator to access the seed value
+            friend class randomGenerator;
+
+
+            // Constructors
+
+                //- Construct from a label
+                inline seed(const label s);
+
+                //- Construct from a word
+                inline seed(const word& s);
+        };
+
+
+private:
+
     // Private Static Data
 
         //- The parameters of the linear congruential iteration
@@ -75,42 +112,78 @@ class randomGenerator
 
     // Private Data
 
+        //- Is this generator global?
+        const bool global_;
+
         //- The stored integer
         uint64_t x_;
 
 
     // Private Member Functions
 
+        //- Check the state is synchronised
+        void checkSync() const;
+
         //- Advance the state and return an integer sample
         inline uint64_t sample();
 
 
-protected:
+        // Scalars
 
-    // Protected Constructors
+            //- Return a scalar uniformly distributed between zero and one.
+            //  Don't check synchronisation.
+            inline scalar scalar01NoCheckSync();
 
-        //- Construct from a 64-bit unsigned integer seed
+            //- Return a scalar uniformly distributed between zero and one.
-        inline randomGenerator(const pTraits<uint64_t> s);
+            //  Don't check synchronisation.
+            inline scalar scalarABNoCheckSync(const scalar a, const scalar b);
 
 
-    // Protected Member Functions
+        // Other Types
 
-        //- Access the stored integer
+            //- Return a type with components uniformly distributed between
-        uint64_t x() const;
+            //  zero and one. Don't check synchronisation.
+            template<class Type>
+            inline Type sample01NoCheckSync();
+
+            //- Return a type with components uniformly distributed between two
+            //  limits. Don't check synchronisation.
+            template<class Type>
+            inline Type sampleABNoCheckSync(const Type& a, const Type& b);
 
 
 public:
 
     // Constructors
 
-        //- Construct from a label seed
+        //- Construct from a seed
-        inline randomGenerator(const label s);
+        inline randomGenerator(const seed s, const bool global = false);
 
-        //- Construct from a word seed
+        //- Copy construct
-        inline randomGenerator(const word& s);
+        inline randomGenerator(const randomGenerator&);
+
+        //- Move construct
+        inline randomGenerator(randomGenerator&&);
 
         //- Construct from a stream
-        randomGenerator(Istream& is);
+        randomGenerator(Istream& is, const bool global = false);
+
+        //- Construct from a dictionary
+        randomGenerator
+        (
+            const word& name,
+            const dictionary& dict,
+            randomGenerator&& defaultRndGen
+        );
+
+        //- Construct from a dictionary
+        randomGenerator
+        (
+            const word& name,
+            const dictionary& dict,
+            const seed defaultS,
+            const bool global = false
+        );
 
 
     //- Destructor
@@ -171,6 +244,18 @@ public:
         template<class Container>
         inline void permute(Container& l);
 
+        //- Create a randomly seeded generator
+        inline randomGenerator generator();
+
+
+    // Member Operators
+
+        //- Copy-assignment
+        void operator=(const randomGenerator&);
+
+        //- Move-assignment
+        void operator=(randomGenerator&&);
+
 
     // IOstream Operators
 
@@ -183,16 +268,27 @@ public:
 
 
 template<>
-inline scalar randomGenerator::sample01();
+inline scalar randomGenerator::sample01NoCheckSync();
 
 template<>
-inline label randomGenerator::sample01();
+inline label randomGenerator::sample01NoCheckSync();
 
 template<>
-inline scalar randomGenerator::sampleAB(const scalar& a, const scalar& b);
+inline scalar randomGenerator::sampleABNoCheckSync
+(
+    const scalar& a,
+    const scalar& b
+);
 
 template<>
-inline label randomGenerator::sampleAB(const label& a, const label& b);
+inline label randomGenerator::sampleABNoCheckSync
+(
+    const label& a,
+    const label& b
+);
+
+
+void writeEntry(Ostream& os, const randomGenerator& rndGen);
 
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

src/OpenFOAM/primitives/randomGenerator/randomGeneratorI.H

@@ -28,6 +28,31 @@ License
 
 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
 
+inline uint64_t Foam::randomGenerator::seed::x(const bool global) const
+{
+    const uint64_t localS =
+        global ? s_ : s_ + (UINT64_MAX/Pstream::nProcs())*Pstream::myProcNo();
+
+    return (localS << 16) + 0x330E;
+}
+
+
+inline void Foam::randomGenerator::checkSync() const
+{
+    if (global_)
+    {
+        uint64_t xMaster = x_;
+        Pstream::scatter(xMaster);
+        if (xMaster != x_)
+        {
+            FatalErrorInFunction
+                << "Global random number generator is not synchronised"
+                << exit(FatalError);
+        }
+    }
+}
+
+
 inline uint64_t Foam::randomGenerator::sample()
 {
     x_ = (A*x_ + C) % M;
@@ -35,33 +60,116 @@ inline uint64_t Foam::randomGenerator::sample()
 }
 
 
-// * * * * * * * * * * * * * Protected Constructors  * * * * * * * * * * * * //
+inline Foam::scalar Foam::randomGenerator::scalar01NoCheckSync()
-
-inline Foam::randomGenerator::randomGenerator(const pTraits<uint64_t> s)
-:
-    x_((s << 16) + 0x330E)
-{}
-
-
-// * * * * * * * * * * * * Protected Member Functions  * * * * * * * * * * * //
-
-inline uint64_t Foam::randomGenerator::x() const
 {
-    return x_;
+    return scalar(sample())/(M >> 17);
+}
+
+
+inline Foam::scalar Foam::randomGenerator::scalarABNoCheckSync
+(
+    const scalar a,
+    const scalar b
+)
+{
+    return a + scalar01NoCheckSync()*(b - a);
+}
+
+
+template<class Type>
+inline Type Foam::randomGenerator::sample01NoCheckSync()
+{
+    Type value;
+    for (direction i = 0; i < pTraits<Type>::nComponents; ++ i)
+    {
+        value.component(i) = scalar01NoCheckSync();
+    }
+    return value;
+}
+
+
+template<>
+inline Foam::scalar Foam::randomGenerator::sample01NoCheckSync()
+{
+    return scalar01NoCheckSync();
+}
+
+
+template<>
+inline Foam::label Foam::randomGenerator::sample01NoCheckSync()
+{
+    return sample() % 2;
+}
+
+
+template<class Type>
+inline Type Foam::randomGenerator::sampleABNoCheckSync
+(
+    const Type& a,
+    const Type& b
+)
+{
+    return a + cmptMultiply(sample01NoCheckSync<Type>(), b - a);
+}
+
+
+template<>
+inline Foam::scalar Foam::randomGenerator::sampleABNoCheckSync
+(
+    const scalar& a,
+    const scalar& b
+)
+{
+    return scalarABNoCheckSync(a, b);
+}
+
+
+template<>
+inline Foam::label Foam::randomGenerator::sampleABNoCheckSync
+(
+    const label& a,
+    const label& b
+)
+{
+    return a + sample() % (b - a);
 }
 
 
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
-inline Foam::randomGenerator::randomGenerator(const label s)
+inline Foam::randomGenerator::seed::seed(const label s)
 :
-    randomGenerator(pTraits<uint64_t>(s))
+    s_(s)
 {}
 
 
-inline Foam::randomGenerator::randomGenerator(const word& s)
+inline Foam::randomGenerator::seed::seed(const word& s)
 :
-    randomGenerator(pTraits<uint64_t>(Hash<word>()(s)))
+    s_(Hash<word>()(s))
+{}
+
+
+inline Foam::randomGenerator::randomGenerator(const seed s, const bool global)
+:
+    global_(global),
+    x_(s.x(global))
+{
+    checkSync();
+}
+
+
+inline Foam::randomGenerator::randomGenerator(const randomGenerator& rndGen)
+:
+    global_(rndGen.global_),
+    x_(rndGen.x_)
+{
+    checkSync();
+}
+
+
+inline Foam::randomGenerator::randomGenerator(randomGenerator&& rndGen)
+:
+    randomGenerator(static_cast<const randomGenerator&>(rndGen))
 {}
 
 
@@ -75,7 +183,11 @@ inline Foam::randomGenerator::~randomGenerator()
 
 inline Foam::scalar Foam::randomGenerator::scalar01()
 {
-    return scalar(sample())/(M >> 17);
+    #ifdef FULLDEBUG
+    checkSync();
+    #endif
+
+    return scalar01NoCheckSync();
 }
 
 
@@ -84,10 +196,14 @@ inline Foam::tmp<Foam::scalarField> Foam::randomGenerator::scalar01
     const label n
 )
 {
+    #ifdef FULLDEBUG
+    checkSync();
+    #endif
+
     tmp<scalarField> tValues(new scalarField(n));
     for (label i = 0; i < n; ++ i)
     {
-        tValues.ref()[i] = scalar01();
+        tValues.ref()[i] = scalar01NoCheckSync();
     }
     return tValues;
 }
@@ -99,7 +215,11 @@ inline Foam::scalar Foam::randomGenerator::scalarAB
     const scalar b
 )
 {
-    return a + scalar01()*(b - a);
+    #ifdef FULLDEBUG
+    checkSync();
+    #endif
+
+    return scalarABNoCheckSync(a, b);
 }
 
 
@@ -110,10 +230,14 @@ inline Foam::tmp<Foam::scalarField> Foam::randomGenerator::scalarAB
     const scalar b
 )
 {
+    #ifdef FULLDEBUG
+    checkSync();
+    #endif
+
     tmp<scalarField> tValues(new scalarField(n));
     for (label i = 0; i < n; ++ i)
     {
-        tValues.ref()[i] = scalarAB(a, b);
+        tValues.ref()[i] = scalarABNoCheckSync(a, b);
     }
     return tValues;
 }
@@ -122,26 +246,11 @@ inline Foam::tmp<Foam::scalarField> Foam::randomGenerator::scalarAB
 template<class Type>
 inline Type Foam::randomGenerator::sample01()
 {
-    Type value;
+    #ifdef FULLDEBUG
-    for (direction i = 0; i < pTraits<Type>::nComponents; ++ i)
+    checkSync();
-    {
+    #endif
-        value.component(i) = scalar01();
-    }
-    return value;
-}
 
-
+    return sample01NoCheckSync<Type>();
-template<>
-inline Foam::scalar Foam::randomGenerator::sample01()
-{
-    return scalar01();
-}
-
-
-template<>
-inline Foam::label Foam::randomGenerator::sample01()
-{
-    return sample() % 2;
 }
 
 
@@ -151,10 +260,14 @@ inline Foam::tmp<Foam::Field<Type>> Foam::randomGenerator::sample01
     const label n
 )
 {
+    #ifdef FULLDEBUG
+    checkSync();
+    #endif
+
     tmp<Field<Type>> tValues(new Field<Type>(n));
     for (label i = 0; i < n; ++ i)
     {
-        tValues.ref()[i] = sample01<Type>();
+        tValues.ref()[i] = sample01NoCheckSync<Type>();
     }
     return tValues;
 }
@@ -163,29 +276,11 @@ inline Foam::tmp<Foam::Field<Type>> Foam::randomGenerator::sample01
 template<class Type>
 inline Type Foam::randomGenerator::sampleAB(const Type& a, const Type& b)
 {
-    return a + cmptMultiply(sample01<Type>(), b - a);
+    #ifdef FULLDEBUG
-}
+    checkSync();
+    #endif
 
-
+    return sampleABNoCheckSync(a, b);
-template<>
-inline Foam::scalar Foam::randomGenerator::sampleAB
-(
-    const scalar& a,
-    const scalar& b
-)
-{
-    return scalarAB(a, b);
-}
-
-
-template<>
-inline Foam::label Foam::randomGenerator::sampleAB
-(
-    const label& a,
-    const label& b
-)
-{
-    return a + sample() % (b - a);
 }
 
 
@@ -197,10 +292,14 @@ inline Foam::tmp<Foam::Field<Type>> Foam::randomGenerator::sampleAB
     const Type& b
 )
 {
+    #ifdef FULLDEBUG
+    checkSync();
+    #endif
+
     tmp<Field<Type>> tValues(new Field<Type>(n));
     for (label i = 0; i < n; ++ i)
     {
-        tValues.ref()[i] = sampleAB(a, b);
+        tValues.ref()[i] = sampleABNoCheckSync(a, b);
     }
     return tValues;
 }
@@ -209,11 +308,21 @@ inline Foam::tmp<Foam::Field<Type>> Foam::randomGenerator::sampleAB
 template<class Container>
 inline void Foam::randomGenerator::permute(Container& l)
 {
+    #ifdef FULLDEBUG
+    checkSync();
+    #endif
+
     for (label i = 0; i < l.size(); ++ i)
     {
-        Swap(l[i], l[sampleAB<label>(i, l.size())]);
+        Swap(l[i], l[sampleABNoCheckSync<label>(i, l.size())]);
     }
 }
 
 
+inline Foam::randomGenerator Foam::randomGenerator::generator()
+{
+    return randomGenerator(sample(), global_);
+}
+
+
 // ************************************************************************* //

src/lagrangian/DSMC/clouds/Templates/DSMCCloud/DSMCCloud.C

@@ -650,7 +650,7 @@ Foam::DSMCCloud<ParcelType>::DSMCCloud
         mesh_
     ),
     constProps_(),
-    rndGen_(label(149382906) + 7183*Pstream::myProcNo()),
+    rndGen_(label(149382906)),
     boundaryT_
     (
         volScalarField
@@ -904,7 +904,7 @@ Foam::DSMCCloud<ParcelType>::DSMCCloud
         )
     ),
     constProps_(),
-    rndGen_(label(971501) + 1526*Pstream::myProcNo()),
+    rndGen_(label(971501)),
     boundaryT_
     (
         volScalarField

src/meshTools/patchIntersection/star/starTemplates.C

@@ -23,8 +23,7 @@ License
 
 \*---------------------------------------------------------------------------*/
 
-#include "OFstream.H"
+#include "star.H"
-
 
 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
 

src/parallel/decompose/decompositionMethods/random/random.C

@@ -79,7 +79,7 @@ Foam::labelList Foam::decompositionMethods::random::decompose
 {
     checkWeights(points, pointWeights);
 
-    randomGenerator rndGen(seed_);
+    randomGenerator rndGen(seed_, true);
 
     // Simple random integer. No guarantees about balance.
     /*

src/sampling/sampledSet/boundaryRandom/boundaryRandom.C

@@ -118,7 +118,7 @@ void Foam::sampledSets::boundaryRandom::calcSamples
     }
 
     // Generate the samples
-    randomGenerator rndGen(261782);
+    randomGenerator rndGen(261782, true);
     const label proci = Pstream::myProcNo();
     for (label i = 0; i < nPoints_; ++ i)
     {

src/sampling/sampledSet/circleRandom/circleRandom.C

@@ -55,7 +55,7 @@ void Foam::sampledSets::circleRandom::calcSamples
     DynamicList<label>& samplingFaces
 ) const
 {
-    randomGenerator rndGen(261782);
+    randomGenerator rndGen(261782, true);
 
     const vector radial1 = normalised(perpendicular(normal_));
     const vector radial2 = normalised(normal_ ^ radial1);

src/sampling/sampledSet/sphereRandom/sphereRandom.C

@@ -55,7 +55,7 @@ void Foam::sampledSets::sphereRandom::calcSamples
     DynamicList<label>& samplingFaces
 ) const
 {
-    randomGenerator rndGen(261782);
+    randomGenerator rndGen(261782, true);
 
     for (label i = 0; i < nPoints_; ++ i)
     {

src/waves/waveModels/irregular/irregular.C

@@ -135,7 +135,7 @@ Foam::waveModels::irregular::irregular
     }
 
     // Set random phases
-    randomGenerator rndGen(0);
+    randomGenerator rndGen(0, true);
     forAll(phases_, i)
     {
         phases_[i] = rndGen.scalarAB(0, constant::mathematical::twoPi);