From 2f9dba6a5813505a5de60dbd0b735628632815ef Mon Sep 17 00:00:00 2001
From: John Tramm <john.tramm@gmail.com>
Date: Mon, 18 Nov 2019 22:28:19 +0000
Subject: [PATCH] removed print statements. Works in serial as expected.

---
 examples/xml/pincell/settings.xml |  4 ++--
 src/physics.cpp                   | 30 +++++++++++++-------------
 src/simulation.cpp                | 35 ++++++++++++++++++-------------
 3 files changed, 37 insertions(+), 32 deletions(-)

diff --git a/examples/xml/pincell/settings.xml b/examples/xml/pincell/settings.xml
index 5845898c8..7b0a985f1 100644
--- a/examples/xml/pincell/settings.xml
+++ b/examples/xml/pincell/settings.xml
@@ -3,8 +3,8 @@
 
   <!-- Define how many particles to run and for how many batches -->
   <run_mode>eigenvalue</run_mode>
-  <batches>100</batches>
-  <inactive>10</inactive>
+  <batches>10</batches>
+  <inactive>5</inactive>
   <particles>1000</particles>
 
   <!-- The starting source is a uniform distribution over the entire pin
diff --git a/src/physics.cpp b/src/physics.cpp
index 7dae2431c..ef15e866b 100644
--- a/src/physics.cpp
+++ b/src/physics.cpp
@@ -79,7 +79,7 @@ void collision(Particle* p)
 
 void sample_neutron_reaction(Particle* p)
 {
-	std::cout << "particle energy e = " << p->E_ << std::endl;
+	//std::cout << "particle energy e = " << p->E_ << std::endl;
   // Sample a nuclide within the material
   int i_nuclide = sample_nuclide(p);
 
@@ -128,12 +128,12 @@ void sample_neutron_reaction(Particle* p)
   }
   if (!p->alive_) return;
 
-	std::cout << "before scatter particle energy e = " << p->E_ << " with nuclide = " << i_nuclide <<  std::endl;
+	//std::cout << "before scatter particle energy e = " << p->E_ << " with nuclide = " << i_nuclide <<  std::endl;
   // Sample a scattering reaction and determine the secondary energy of the
   // exiting neutron
   scatter(p, i_nuclide);
 	
-  std::cout << "after scatter particle energy e = " << p->E_ << std::endl;
+  //std::cout << "after scatter particle energy e = " << p->E_ << std::endl;
 
   // Advance URR seed stream 'N' times after energy changes
   if (p->E_ != p->E_last_) {
@@ -602,7 +602,7 @@ void scatter(Particle* p, int i_nuclide)
     // Determine temperature
     double kT = nuc->multipole_ ? p->sqrtkT_*p->sqrtkT_ : nuc->kTs_[i_temp];
 
-	std::cout << "we are doing an elastic scatter" << std::endl;
+	//std::cout << "we are doing an elastic scatter" << std::endl;
     // Perform collision physics for elastic scattering
     elastic_scatter(i_nuclide, *nuc->reactions_[0], kT, p);
 
@@ -615,7 +615,7 @@ void scatter(Particle* p, int i_nuclide)
     // =======================================================================
     // S(A,B) SCATTERING
 
-	std::cout << "we are doing an SAB scatter" << std::endl;
+	//std::cout << "we are doing an SAB scatter" << std::endl;
     sab_scatter(i_nuclide, micro.index_sab, p);
 
     p->event_mt_ = ELASTIC;
@@ -623,7 +623,7 @@ void scatter(Particle* p, int i_nuclide)
   }
 
   if (!sampled) {
-	std::cout << "we are doing an Inelastic scatter" << std::endl;
+	//std::cout << "we are doing an Inelastic scatter" << std::endl;
     // =======================================================================
     // INELASTIC SCATTERING
 
@@ -650,9 +650,9 @@ void scatter(Particle* p, int i_nuclide)
 
     // Perform collision physics for inelastic scattering
     const auto& rx {nuc->reactions_[i]};
-	std::cout << "Energy before Inelastic scatter = " << p->E_ << std::endl;
+	//std::cout << "Energy before Inelastic scatter = " << p->E_ << std::endl;
     inelastic_scatter(nuc.get(), rx.get(), p);
-	std::cout << "Energy after Inelastic scatter = " << p->E_ << std::endl;
+	//std::cout << "Energy after Inelastic scatter = " << p->E_ << std::endl;
     p->event_mt_ = rx->mt_;
   }
 
@@ -1061,35 +1061,35 @@ void inelastic_scatter(const Nuclide* nuc, const Reaction* rx, Particle* p)
   double mu;
   // This gives a negative sometimes for some unknown reason
   rx->products_[0].sample(E_in, E, mu);
-std::cout << "E_in = " << E_in << "E = " << E << " mu = " << mu << std::endl;
+//std::cout << "E_in = " << E_in << "E = " << E << " mu = " << mu << std::endl;
 
   // if scattering system is in center-of-mass, transfer cosine of scattering
   // angle and outgoing energy from CM to LAB
   if (rx->scatter_in_cm_) {
     double E_cm = E;
-	std::cout << "E_cm = " << E_cm << std::endl;
+	//std::cout << "E_cm = " << E_cm << std::endl;
 
     // determine outgoing energy in lab
     double A = nuc->awr_;
-	std::cout << "A = " << nuc->awr_ << std::endl;
+//	std::cout << "A = " << nuc->awr_ << std::endl;
     E = E_cm + (E_in + 2.0*mu*(A + 1.0) * std::sqrt(E_in*E_cm))
           / ((A + 1.0)*(A + 1.0));
 	// here's where it goes wrong
-	std::cout << "E = " << E << std::endl;
+	//std::cout << "E = " << E << std::endl;
 
     // determine outgoing angle in lab
     mu = mu*std::sqrt(E_cm/E) + 1.0/(A+1.0) * std::sqrt(E_in/E);
-	std::cout << "mu = " << mu << std::endl;
+//	std::cout << "mu = " << mu << std::endl;
   }
 
   // Because of floating-point roundoff, it may be possible for mu to be
   // outside of the range [-1,1). In these cases, we just set mu to exactly -1
   // or 1
   if (std::abs(mu) > 1.0) mu = std::copysign(1.0, mu);
-	std::cout << "mu = " << mu << std::endl;
+//	std::cout << "mu = " << mu << std::endl;
 
   // Set outgoing energy and scattering angle
-	std::cout << "FINAL E = " << E << " mu = " << mu << std::endl;
+//	std::cout << "FINAL E = " << E << " mu = " << mu << std::endl;
   p->E_ = E;
   p->mu_ = mu;
 
diff --git a/src/simulation.cpp b/src/simulation.cpp
index f56dfb2ef..90772bd23 100644
--- a/src/simulation.cpp
+++ b/src/simulation.cpp
@@ -371,7 +371,7 @@ void process_collision_events()
   //for (auto& p : collision_queue) {
   for (int i = 0; i < collision_queue_length; i++) {
 	  Particle * p = particles + collision_queue[i];
-	  std::cout << "Beginning collision of particle id " << collision_queue[i] << " with energy E = " << p->E_ << std::endl;
+	  //std::cout << "Beginning collision of particle id " << collision_queue[i] << " with energy E = " << p->E_ << std::endl;
     // Score collision estimate of keff
     if (settings::run_mode == RUN_MODE_EIGENVALUE &&
         p->type_ == Particle::Type::neutron) {
@@ -386,7 +386,7 @@ void process_collision_events()
     if (!model::active_meshsurf_tallies.empty())
       score_surface_tally(p, model::active_meshsurf_tallies);
 	    
-	std::cout << "After surface tally of particle id " << collision_queue[i] << " with energy E = " << p->E_ << std::endl;
+	//std::cout << "After surface tally of particle id " << collision_queue[i] << " with energy E = " << p->E_ << std::endl;
 
     // Clear surface component
     p->surface_ = 0;
@@ -397,7 +397,7 @@ void process_collision_events()
       collision_mg(p);
     }
 
-	std::cout << "After collision() of particle id " << collision_queue[i] << " with energy E = " << p->E_ << std::endl;
+	//std::cout << "After collision() of particle id " << collision_queue[i] << " with energy E = " << p->E_ << std::endl;
 
     // Score collision estimator tallies -- this is done after a collision
     // has occurred rather than before because we need information on the
@@ -410,7 +410,7 @@ void process_collision_events()
         score_analog_tally_mg(p);
       }
     }
-	std::cout << "After analog tally of particle id " << collision_queue[i] << " with energy E = " << p->E_ << std::endl;
+	//std::cout << "After analog tally of particle id " << collision_queue[i] << " with energy E = " << p->E_ << std::endl;
 
     // Reset banked weight during collision
     p->n_bank_ = 0;
@@ -456,7 +456,7 @@ void process_collision_events()
     if (p->alive_)
 	{
 		dispatch_xs_event(collision_queue[i]);
-	    std::cout << "Ended collision of particle id " << collision_queue[i] << " with energy E = " << p->E_ << std::endl;
+	    //std::cout << "Ended collision of particle id " << collision_queue[i] << " with energy E = " << p->E_ << std::endl;
 		assert(std::isfinite(p->E_) );
 	}
   }
@@ -561,42 +561,46 @@ void transport()
 			dispatch_xs_event(i);
 		}
 
+		int event_kernel_executions = 0;
 		while (true) {
+			event_kernel_executions++;
+			/*
 			std::cout << "Fuel XS Lookups = " << calculate_fuel_xs_queue_length << std::endl;
 			std::cout << "Non Fuel XS Lookups = " << calculate_nonfuel_xs_queue_length << std::endl;
 			std::cout << "Advance Particles = " << advance_particle_queue_length << std::endl;
 			std::cout << "Surface Crossings = " << surface_crossing_queue_length << std::endl;
 			std::cout << "Collisions = " << collision_queue_length << std::endl;
+			*/
 			int max = std::max({calculate_fuel_xs_queue_length, calculate_nonfuel_xs_queue_length, advance_particle_queue_length, surface_crossing_queue_length, collision_queue_length});
 			check_energies();
 			if (max == 0) {
 				break;
 			} else if (max == calculate_fuel_xs_queue_length) {
-				std::cout << "pre fuel XS check..." << std::endl;
+				//std::cout << "pre fuel XS check..." << std::endl;
 				//check_energies(calculate_fuel_xs_queue, calculate_fuel_xs_queue_length);
-				std::cout << "Performing Fuel XS Lookups..." << std::endl;
+				//std::cout << "Performing Fuel XS Lookups..." << std::endl;
 				process_calculate_xs_events(calculate_fuel_xs_queue, calculate_fuel_xs_queue_length);
 				calculate_fuel_xs_queue_length = 0;
 			} else if (max == calculate_nonfuel_xs_queue_length) {
-				std::cout << "pre non fuel XS check..." << std::endl;
+				//std::cout << "pre non fuel XS check..." << std::endl;
 				//check_energies(calculate_nonfuel_xs_queue, calculate_nonfuel_xs_queue_length);
-				std::cout << "Performing Non Fuel XS Lookups..." << std::endl;
+			//	std::cout << "Performing Non Fuel XS Lookups..." << std::endl;
 				process_calculate_xs_events(calculate_nonfuel_xs_queue, calculate_nonfuel_xs_queue_length);
 				calculate_nonfuel_xs_queue_length = 0;
 			} else if (max == advance_particle_queue_length) {
-				std::cout << "pre advancing check..." << std::endl;
+				//std::cout << "pre advancing check..." << std::endl;
 				//check_energies(advance_particle_queue, advance_particle_queue_length);
-				std::cout << "Advancing Particles..." << std::endl;
+				//std::cout << "Advancing Particles..." << std::endl;
 				process_advance_particle_events();
 			} else if (max == surface_crossing_queue_length) {
-				std::cout << "pre surface crossing check..." << std::endl;
+				//std::cout << "pre surface crossing check..." << std::endl;
 				//check_energies(surface_crossing_queue, surface_crossing_queue_length);
-				std::cout << "Surface Crossings..." << std::endl;
+				//std::cout << "Surface Crossings..." << std::endl;
 				process_surface_crossing_events();
 			} else if (max == collision_queue_length) {
-				std::cout << "pre Colliding check..." << std::endl;
+				//std::cout << "pre Colliding check..." << std::endl;
 				//check_energies(collision_queue, collision_queue_length);
-				std::cout << "Colliding..." << std::endl;
+				//std::cout << "Colliding..." << std::endl;
 				process_collision_events();
 			}
 		}
@@ -615,6 +619,7 @@ void transport()
 		
 		// TODO: Do this only once per PI
 		free_event_queues();
+		std::cout << "Event kernels retired: " << event_kernel_executions << std::endl;
 	}
 }