Merge pull request #1666 from paulromano/tally-fixes

Several fixes/improvements for large tally cases

docs/source/io_formats/statepoint.rst

@@ -161,3 +161,5 @@ All values are given in seconds and are measured on the master process.
              source sites between processes for load balancing.
            - **accumulating tallies** (*double*) -- Time spent communicating
              tally results and evaluating their statistics.
+          - **writing statepoints** (*double*) -- Time spent writing statepoint
+            files

include/openmc/timer.h

@@ -22,6 +22,7 @@ extern Timer time_inactive;
 extern Timer time_initialize;
 extern Timer time_read_xs;
 extern Timer time_sample_source;
+extern Timer time_statepoint;
 extern Timer time_tallies;
 extern Timer time_total;
 extern Timer time_transport;

src/output.cpp

@@ -420,6 +420,7 @@ void print_runtime()
     show_time("SEND/RECV source sites", time_bank_sendrecv.elapsed(), 2);
   }
   show_time("Time accumulating tallies", time_tallies.elapsed(), 1);
+  show_time("Time writing statepoints", time_statepoint.elapsed(), 1);
   show_time("Total time for finalization", time_finalize.elapsed());
   show_time("Total time elapsed", time_total.elapsed());
 

src/settings.cpp

@@ -666,7 +666,7 @@ void read_settings_xml()
   // Check if the user has specified to not reduce tallies at the end of every
   // batch
   if (check_for_node(root, "no_reduce")) {
-    reduce_tallies = get_node_value_bool(root, "no_reduce");
+    reduce_tallies = !get_node_value_bool(root, "no_reduce");
   }
 
   // Check if the user has specified to use confidence intervals for

src/state_point.cpp

@@ -33,6 +33,8 @@ namespace openmc {
 extern "C" int
 openmc_statepoint_write(const char* filename, bool* write_source)
 {
+  simulation::time_statepoint.start();
+
   // Set the filename
   std::string filename_;
   if (filename) {
@@ -272,7 +274,9 @@ openmc_statepoint_write(const char* filename, bool* write_source)
   } else if (mpi::master) {
     // Write number of global realizations
     write_dataset(file_id, "n_realizations", simulation::n_realizations);
+  }
 
+  if (mpi::master) {
     // Write out the runtime metrics.
     using namespace simulation;
     hid_t runtime_group = create_group(file_id, "runtime");
@@ -294,6 +298,7 @@ openmc_statepoint_write(const char* filename, bool* write_source)
     }
     write_dataset(runtime_group, "accumulating tallies", time_tallies.elapsed());
     write_dataset(runtime_group, "total", time_total.elapsed());
+    write_dataset(runtime_group, "writing statepoints", time_statepoint.elapsed());
     close_group(runtime_group);
 
     file_close(file_id);
@@ -312,6 +317,8 @@ openmc_statepoint_write(const char* filename, bool* write_source)
     if (mpi::master || parallel) file_close(file_id);
   }
 
+  simulation::time_statepoint.stop();
+
   return 0;
 }
 
@@ -774,9 +781,6 @@ void write_tally_results_nr(hid_t file_id)
     // Write number of realizations
     write_dataset(file_id, "n_realizations", simulation::n_realizations);
 
-    // Write number of global tallies
-    write_dataset(file_id, "n_global_tallies", N_GLOBAL_TALLIES);
-
     tallies_group = open_group(file_id, "tallies");
   }
 
@@ -808,7 +812,7 @@ void write_tally_results_nr(hid_t file_id)
     if (!t->active_) continue;
     if (!t->writable_) continue;
 
-    if (mpi::master && !object_exists(file_id, "tallies_present")) {
+    if (mpi::master && !attribute_exists(file_id, "tallies_present")) {
       write_attribute(file_id, "tallies_present", 1);
     }
 
@@ -817,7 +821,7 @@ void write_tally_results_nr(hid_t file_id)
       xt::range(static_cast<int>(TallyResult::SUM), static_cast<int>(TallyResult::SUM_SQ) + 1));
 
     // Make copy of tally values in contiguous array
-    xt::xtensor<double, 2> values = values_view;
+    xt::xtensor<double, 3> values = values_view;
 
     if (mpi::master) {
       // Open group for tally
@@ -863,6 +867,8 @@ void write_tally_results_nr(hid_t file_id)
       // Indicate that tallies are off
       write_dataset(file_id, "tallies_present", 0);
     }
+
+    close_group(tallies_group);
   }
 }
 

src/tallies/tally.cpp

@@ -642,6 +642,7 @@ void Tally::accumulate()
     double norm = total_source / (settings::n_particles * settings::gen_per_batch);
 
     // Accumulate each result
+    #pragma omp parallel for
     for (int i = 0; i < results_.shape()[0]; ++i) {
       for (int j = 0; j < results_.shape()[1]; ++j) {
         double val = results_(i, j, TallyResult::VALUE) * norm;
@@ -718,29 +719,35 @@ void read_tallies_xml()
 #ifdef OPENMC_MPI
 void reduce_tally_results()
 {
-  for (int i_tally : model::active_tallies) {
-    // Skip any tallies that are not active
-    auto& tally {model::tallies[i_tally]};
+  // Don't reduce tally is no_reduce option is on
+  if (settings::reduce_tallies) {
+    for (int i_tally : model::active_tallies) {
+      // Skip any tallies that are not active
+      auto& tally {model::tallies[i_tally]};
 
-    // Get view of accumulated tally values
-    auto values_view = xt::view(tally->results_, xt::all(), xt::all(), static_cast<int>(TallyResult::VALUE));
+      // Get view of accumulated tally values
+      auto values_view = xt::view(tally->results_, xt::all(), xt::all(), static_cast<int>(TallyResult::VALUE));
 
-    // Make copy of tally values in contiguous array
-    xt::xtensor<double, 2> values = values_view;
-    xt::xtensor<double, 2> values_reduced = xt::empty_like(values);
+      // Make copy of tally values in contiguous array
+      xt::xtensor<double, 2> values = values_view;
+      xt::xtensor<double, 2> values_reduced = xt::empty_like(values);
 
-    // Reduce contiguous set of tally results
-    MPI_Reduce(values.data(), values_reduced.data(), values.size(),
-      MPI_DOUBLE, MPI_SUM, 0, mpi::intracomm);
+      // Reduce contiguous set of tally results
+      MPI_Reduce(values.data(), values_reduced.data(), values.size(),
+        MPI_DOUBLE, MPI_SUM, 0, mpi::intracomm);
 
-    // Transfer values on master and reset on other ranks
-    if (mpi::master) {
-      values_view = values_reduced;
-    } else {
-      values_view = 0.0;
+      // Transfer values on master and reset on other ranks
+      if (mpi::master) {
+        values_view = values_reduced;
+      } else {
+        values_view = 0.0;
+      }
     }
   }
 
+  // Note that global tallies are *always* reduced even when no_reduce option is
+  // on.
+
   // Get view of global tally values
   auto& gt = simulation::global_tallies;
   auto gt_values_view = xt::view(gt, xt::all(), static_cast<int>(TallyResult::VALUE));
@@ -774,11 +781,11 @@ accumulate_tallies()
 {
 #ifdef OPENMC_MPI
   // Combine tally results onto master process
-  if (settings::reduce_tallies) reduce_tally_results();
+  if (mpi::n_procs > 1) reduce_tally_results();
 #endif
 
   // Increase number of realizations (only used for global tallies)
-  simulation::n_realizations += settings::reduce_tallies ? 1 : mpi::n_procs;
+  simulation::n_realizations += 1;
 
   // Accumulate on master only unless run is not reduced then do it on all
   if (mpi::master || !settings::reduce_tallies) {

src/timer.cpp

@@ -17,6 +17,7 @@ Timer time_inactive;
 Timer time_initialize;
 Timer time_read_xs;
 Timer time_sample_source;
+Timer time_statepoint;
 Timer time_tallies;
 Timer time_total;
 Timer time_transport;
@@ -76,6 +77,7 @@ void reset_timers()
   simulation::time_initialize.reset();
   simulation::time_read_xs.reset();
   simulation::time_sample_source.reset();
+  simulation::time_statepoint.reset();
   simulation::time_tallies.reset();
   simulation::time_total.reset();
   simulation::time_transport.reset();