From d8868cd3b1bab598f1678bf41fb51dfd9102c7a7 Mon Sep 17 00:00:00 2001
From: Ye Luo <yeluo@anl.gov>
Date: Mon, 26 Nov 2018 18:56:32 -0600
Subject: [PATCH] Change RNG seeds generation scheme.

---
 src/OhmmsApp/RandomNumberControl.cpp | 84 ++++++++++------------------
 src/OhmmsApp/RandomNumberControl.h   |  6 +-
 src/Utilities/BoostRandom.h          |  8 +++
 3 files changed, 40 insertions(+), 58 deletions(-)

diff --git a/src/OhmmsApp/RandomNumberControl.cpp b/src/OhmmsApp/RandomNumberControl.cpp
index 8adc20c507..9c74cf2584 100755
--- a/src/OhmmsApp/RandomNumberControl.cpp
+++ b/src/OhmmsApp/RandomNumberControl.cpp
@@ -77,22 +77,30 @@ void RandomNumberControl::reset()
 }
 
 /// reset the generator
-void RandomNumberControl::make_seeds()
+void RandomNumberControl::make_seeds(bool init_from_time, uint_type time_seed)
 {
+  // step 1, get initial seed by time
+  if(init_from_time)
+  {
+    time_seed = std::chrono::system_clock::now().time_since_epoch().count();
+    app_summary() << "  Offset for the random number seeds based on time: " << time_seed << std::endl;
+  }
+  mpi::bcast(*OHMMS::Controller,time_seed);
+  // step 2, generate seeds for RNGs used by MC. Each rank needs omp_get_max_threads()+1 seeds
+  std::minstd_rand seed_generator(time_seed);
+  // jump ahead
   int pid = OHMMS::Controller->rank();
-  int nprocs = OHMMS::Controller->size();
-  uint_type iseed=static_cast<uint_type>(std::time(0))%1024;
-  mpi::bcast(*OHMMS::Controller,iseed);
-  //OHMMS::Controller->bcast(iseed);//broadcast the seed
-  Offset=iseed;
+  seed_generator.discard(pid*(omp_get_max_threads()+1));
+  // step 3, seed the per-rank RNG
+  Random.init(seed_generator());
+  // step 4, seed the per-thread RNG
   std::vector<uint_type> mySeeds;
-  RandomNumberControl::PrimeNumbers.get(Offset,nprocs*(omp_get_max_threads()+2), mySeeds);
-  Random.init(pid,nprocs,mySeeds[pid],Offset+pid);
-  //change children as well
-  make_children();
+  for(int tid=0; tid<omp_get_max_threads(); tid++)
+    mySeeds.push_back(seed_generator());
+  make_children(mySeeds);
 }
 
-void RandomNumberControl::make_children()
+void RandomNumberControl::make_children(std::vector<uint_type> &mySeeds)
 {
   int nthreads=omp_get_max_threads();
   int n=nthreads-Children.size();
@@ -101,16 +109,8 @@ void RandomNumberControl::make_children()
     Children.push_back(new RandomGenerator_t);
     n--;
   }
-  int rank=OHMMS::Controller->rank();
-  int nprocs=OHMMS::Controller->size();
-  int baseoffset=Offset+nprocs+nthreads*rank;
-  std::vector<uint_type> myprimes;
-  PrimeNumbers.get(baseoffset,nthreads,myprimes);
   for(int ip=0; ip<nthreads; ip++)
-  {
-    int offset=baseoffset+ip;
-    Children[ip]->init(rank,nprocs,myprimes[ip],offset);
-  }
+    Children[ip]->init(mySeeds[ip]);
 }
 
 xmlNodePtr
@@ -173,48 +173,20 @@ bool RandomNumberControl::put(xmlNodePtr cur)
 {
   if(NeverBeenInitialized)
   {
-    bool init_mpi = true;
-    int offset_in = -1; // default is to generate by Wall-clock
+    app_summary() << " Random Number" << std::endl;
+    app_summary() << " -------------" << std::endl;
+
+    bool seed_from_time = true;
+    size_t offset_in = 0; // default is to generate by Wall-clock
     if(cur != NULL)
     {
-      std::string pname("yes");
       OhmmsAttributeSet oAttrib;
-      oAttrib.add(pname,"parallel");
       oAttrib.add(offset_in,"seed");
       oAttrib.put(cur);
-      if(pname == "0" || pname == "false" || pname == "no")
-        init_mpi=false;
-    }
-    int nprocs = 1;
-    int pid = 0;
-    if(init_mpi)
-    {
-      pid = OHMMS::Controller->rank();
-      nprocs = OHMMS::Controller->size();
+      app_summary() << "  Offset for the random number seeds from input file: " << offset_in << std::endl;
+      seed_from_time = false;
     }
-    app_summary() << " Random Number" << std::endl;
-    app_summary() << " -------------" << std::endl;
-    if(offset_in<0)
-    {
-      offset_in=static_cast<int>(static_cast<uint_type>(std::time(0))%1024);
-      app_summary() << "  Offset for the random number seeds based on time: " << offset_in << std::endl;
-      mpi::bcast(*OHMMS::Controller,offset_in);
-    }
-    else
-    {
-      offset_in%=1024;
-      app_summary() << "  Offset for the random number seeds from input file (mod 1024): " << offset_in << std::endl;
-    }
-    app_summary() << std::endl;
-    Offset=offset_in;
-    std::vector<uint_type> mySeeds;
-    //allocate twice of what is required
-    PrimeNumbers.get(Offset,nprocs*(omp_get_max_threads()+2), mySeeds);
-    Random.init(pid,nprocs,mySeeds[pid],Offset+pid);
-    app_log() <<  "  Range of prime numbers to use as seeds over processors and threads = " << mySeeds[0] <<"-" << mySeeds[nprocs*omp_get_max_threads()] <<  std::endl;
-    app_log() << std::endl;
-
-    make_children();
+    make_seeds(seed_from_time, offset_in);
     NeverBeenInitialized = false;
     app_log() << std::endl;
   }
diff --git a/src/OhmmsApp/RandomNumberControl.h b/src/OhmmsApp/RandomNumberControl.h
index b82ffb8b7d..9fc85cc0b8 100755
--- a/src/OhmmsApp/RandomNumberControl.h
+++ b/src/OhmmsApp/RandomNumberControl.h
@@ -18,6 +18,8 @@
 #include "Utilities/RandomGenerator.h"
 #include "Utilities/PrimeNumberSet.h"
 #include <io/hdf_archive.h>
+#include <chrono>
+#include <random>
 
 class Communicate;
 
@@ -50,8 +52,8 @@ class RandomNumberControl : public OhmmsElementBase
   void reset();
   static void test();
 
-  static void make_seeds();
-  static void make_children();
+  static void make_seeds(bool init_from_time = true, uint_type time_seed = 0);
+  static void make_children(std::vector<uint_type> &mySeeds);
 
   xmlNodePtr initialize(xmlXPathContextPtr);
 
diff --git a/src/Utilities/BoostRandom.h b/src/Utilities/BoostRandom.h
index 036d2a7acc..1390eeacd1 100644
--- a/src/Utilities/BoostRandom.h
+++ b/src/Utilities/BoostRandom.h
@@ -92,6 +92,14 @@ class BoostRandom
     uni.engine().seed(baseSeed);
   }
 
+  /** initialize the generator with a seed
+   * @param iseed_in input seed
+   */
+  void init(uint_type iseed_in)
+  {
+    uni.engine().seed(iseed_in);
+  }
+
   ///get baseOffset
   inline int offset() const
   {