diff --git a/include/openmc/bank.h b/include/openmc/bank.h index 544f6a4da9..303c1e8be7 100644 --- a/include/openmc/bank.h +++ b/include/openmc/bank.h @@ -17,7 +17,7 @@ namespace simulation { extern std::vector source_bank; -extern Particle::Bank* fission_bank; +extern std::unique_ptr fission_bank; extern int64_t fission_bank_length; extern int64_t fission_bank_max; extern std::vector progeny_per_particle; diff --git a/include/openmc/particle.h b/include/openmc/particle.h index 555120e94b..2f5e96394c 100644 --- a/include/openmc/particle.h +++ b/include/openmc/particle.h @@ -172,7 +172,6 @@ public: Type particle; int64_t parent_id; int64_t progeny_id; - bool operator<(const Bank& bank) const { return (parent_id < bank.parent_id); } }; //! Saved ("banked") state of a particle, for nu-fission tallying diff --git a/openmc/lib/core.py b/openmc/lib/core.py index f26ddc5259..ffc8090e78 100644 --- a/openmc/lib/core.py +++ b/openmc/lib/core.py @@ -21,7 +21,8 @@ class _Bank(Structure): ('wgt', c_double), ('delayed_group', c_int), ('particle', c_int), - ('parent_id', c_int64)] + ('parent_id', c_int64), + ('progeny_id', c_int64)] # Define input type for numpy arrays that will be passed into C++ functions diff --git a/src/bank.cpp b/src/bank.cpp index 7a4ae244c4..b9bf9056d5 100644 --- a/src/bank.cpp +++ b/src/bank.cpp @@ -16,7 +16,13 @@ namespace simulation { std::vector source_bank; -Particle::Bank* fission_bank {nullptr}; +// The fission bank is allocated as a pointer, rather than a vector, as it will +// be shared by all threads in the simulation. It will be allocated to a fixed +// maximum capacity in the init_fission_bank() function. Then, Elements will be +// added to it by performing atomic incrementing captures on the fission_bank_length +// variable. A vector is not appropriate, as use of push_back() and size() +// functions would cause undefined or unintended behavior. +std::unique_ptr fission_bank; int64_t fission_bank_length {0}; int64_t fission_bank_max; @@ -34,24 +40,30 @@ std::vector progeny_per_particle; void free_memory_bank() { simulation::source_bank.clear(); - if (simulation::fission_bank) - delete[] simulation::fission_bank; - simulation::fission_bank = nullptr; + simulation::fission_bank.reset(); simulation::fission_bank_length = 0; } void init_fission_bank(int64_t max) { simulation::fission_bank_max = max; - simulation::fission_bank = new Particle::Bank[max]; + simulation::fission_bank = std::make_unique(max); simulation::fission_bank_length = 0; simulation::progeny_per_particle.resize(simulation::work_per_rank); } +// Performs an O(n) sort on the fission bank, by leveraging +// the parent_id and progeny_id fields of banked particles. See the following +// paper for more details: +// "Reproducibility and Monte Carlo Eigenvalue Calculations," F.B. Brown and +// T.M. Sutton, 1992 ANS Annual Meeting, Transactions of the American Nuclear +// Society, Volume 65, Page 235. void sort_fission_bank() { - if( simulation::progeny_per_particle.size() == 0 ) + // Ensure we don't read off the end of the array if we ran with 0 particles + if (simulation::progeny_per_particle.size() == 0) { return; + } // Perform exclusive scan summation to determine starting indices in fission // bank for each parent particle id @@ -62,8 +74,11 @@ void sort_fission_bank() tmp = simulation::progeny_per_particle[i]; simulation::progeny_per_particle[i] = value; } + + // TODO: C++17 introduces the exclusive_scan() function which could be + // used to replace everything above this point in this function. - // Create temporary scratch vector for permutation + // Create temporary scratch vector for permuting the fission bank std::vector sorted_bank(simulation::fission_bank_length); // Use parent and progeny indices to sort fission bank @@ -101,7 +116,7 @@ extern "C" int openmc_fission_bank(void** ptr, int64_t* n) set_errmsg("Fission bank has not been allocated."); return OPENMC_E_ALLOCATE; } else { - *ptr = simulation::fission_bank; + *ptr = simulation::fission_bank.get(); *n = simulation::fission_bank_length; return 0; } diff --git a/src/eigenvalue.cpp b/src/eigenvalue.cpp index b9f8a80435..c2008da937 100644 --- a/src/eigenvalue.cpp +++ b/src/eigenvalue.cpp @@ -506,7 +506,7 @@ void shannon_entropy() // Get source weight in each mesh bin bool sites_outside; xt::xtensor p = simulation::entropy_mesh->count_sites( - simulation::fission_bank, simulation::fission_bank_length, + simulation::fission_bank.get(), simulation::fission_bank_length, &sites_outside); // display warning message if there were sites outside entropy box diff --git a/src/physics.cpp b/src/physics.cpp index 10ec59d563..99f3b57249 100644 --- a/src/physics.cpp +++ b/src/physics.cpp @@ -192,7 +192,7 @@ create_fission_sites(Particle* p, int i_nuclide, const Reaction* rx) skipped++; break; } - site = simulation::fission_bank + idx; + site = &simulation::fission_bank[idx]; } else { // Create new bank site and get reference to last element auto& bank = p->secondary_bank_; diff --git a/src/physics_mg.cpp b/src/physics_mg.cpp index d3630ad19b..3a025a80eb 100644 --- a/src/physics_mg.cpp +++ b/src/physics_mg.cpp @@ -138,7 +138,7 @@ create_fission_sites(Particle* p) skipped++; break; } - site = simulation::fission_bank + idx; + site = &simulation::fission_bank[idx]; } else { // Create new bank site and get reference to last element auto& bank = p->secondary_bank_;