Added radioactive decay of 2 nuclei example, added simulation interface

ernestyalumni · Jan 31, 2019 · 6c2a2cb · 6c2a2cb
1 parent 2fa8544
commit 6c2a2cb
Show file tree

Hide file tree

Showing 7 changed files with 333 additions and 19 deletions.
diff --git a/MonteCarlo/Integrate1.h b/MonteCarlo/Integrate1.h
@@ -69,6 +69,8 @@ class Integrate1 : public Integrate<Domain, Codomain>
 
     Integrate1(const F& f, const long total_number_of_samples);
 
+    Integrate1(const F& f, const long total_number_of_samples, const long seed);
+
     void reset(const long total_number_of_samples)
     {
       assert(total_number_of_samples > 0);
@@ -138,11 +140,24 @@ Integrate1<Domain, Codomain, RandomNumberGenerator>::Integrate1(
   assert(N_ > 0);
 }
 
+template <typename Domain, typename Codomain, class RandomNumberGenerator>
+Integrate1<Domain, Codomain, RandomNumberGenerator>::Integrate1(
+  const F& f,
+  const long total_number_of_samples,
+  const long seed
+  ):
+  rng_{seed},
+  f_{f},
+  N_{total_number_of_samples},
+  running_total_{},
+  sum_sigma_{}
+{
+  assert(N_ > 0);
+}
+
 template <class Domain, typename Codomain, class RandomNumberGenerator>
 void Integrate1<Domain, Codomain, RandomNumberGenerator>::run_Monte_Carlo()
 {
-  long idum {-1};
-
   for (int i {0}; i < N_; ++i)
   {
     // Partially for debug purposes only (can be combined)
@@ -156,7 +171,7 @@ void Integrate1<Domain, Codomain, RandomNumberGenerator>::run_Monte_Carlo()
     //domain_runs_.push_back(x);
     //codomain_runs_.push_back(f_x);
 
-    const Codomain f_x {f_(rng_(&idum))};
+    const Codomain f_x {f_(rng_())};
     running_total_ += f_x;
     sum_sigma_ += f_x * f_x;
   }

diff --git a/MonteCarlo/Integrate1_main.cpp b/MonteCarlo/Integrate1_main.cpp
@@ -53,7 +53,7 @@ int main()
 
     for (int i {0}; i < 30; ++i)
     {
-      std::cout << minimal_park_miller(&idum) << ' ';
+      std::cout << minimal_park_miller() << ' ';
     }
   }
 

diff --git a/MonteCarlo/ParticlesInABox.h b/MonteCarlo/ParticlesInABox.h
@@ -6,6 +6,9 @@
 /// \url
 /// \ref
 /// \details Particles in a box.
+/// There is public inheritance of the interface (interface inheritance),
+/// separating the data collection and statistics from the actual Monte Carlo
+/// run.
 /// \copyright If you find this code useful, feel free to donate directly
 /// (username ernestyalumni or email address above), going directly to:
 ///
@@ -29,6 +32,7 @@
 #define _MONTE_CARLO_PARTICLES_IN_A_BOX_H_
 
 #include "RandomNumberGenerators/RandomNumberGenerators.h"
+#include "TimeEvolutionStateTransitions/Simulation.h"
 #include "TimeEvolutionStateTransitions/TimeEvolution.h"
 
 #include <cassert>
@@ -60,13 +64,16 @@ class ParticleTimeEvolution :
       N_{N}
     {}
 
+    ParticleTimeEvolution(const unsigned long N, const long seed):
+      rng_{seed},
+      N_{N}
+    {}
+
     unsigned long operator()(unsigned long& n_l)
     {
       assert(n_l <= N_);
 
-      long idum {-1};
-
-      const RR pr {rng_(&idum)}; // pr = probability
+      const RR pr {rng_()}; // pr = probability
 
       // boundary conditions
       if ((n_l != N_) && (n_l != 0))
@@ -95,7 +102,7 @@ class ParticleTimeEvolution :
 /// \class ParticlesInABox
 //------------------------------------------------------------------------------
 template <class RR, class RandomNumberGenerator>
-class ParticlesInABox
+class ParticlesInABox : public TimeEvolutionStateTransitions::Simulation
 {
   public:
 
@@ -112,6 +119,21 @@ class ParticlesInABox
       assert(T > 0);
     }
 
+    ParticlesInABox(
+      const unsigned long N,
+      const unsigned long T,
+      const long seed
+      ):
+      runs_{},
+      time_evolution_{N, seed},
+      n_l_{N},
+      N_{N},
+      T_{T}
+    {
+      assert(N > 0);
+      assert(T > 0);
+    }
+
     void run()
     {
       for (unsigned long t {0}; t < T_; ++t)