diff --git a/include/openmc/event.h b/include/openmc/event.h index 44c08b27b..ff57121bf 100644 --- a/include/openmc/event.h +++ b/include/openmc/event.h @@ -5,6 +5,7 @@ //! \brief Event-based data structures and methods #include "openmc/particle.h" +#include "openmc/shared_array.h" namespace openmc { @@ -49,23 +50,14 @@ namespace simulation { // coordinated with atomics. This means that normal vector methods (e.g., // push_back(), size()) would cause undefined or unintended behavior. Rather, // adding particles to queues will be done via the enqueue_particle() function. -extern std::unique_ptr calculate_fuel_xs_queue; -extern std::unique_ptr calculate_nonfuel_xs_queue; -extern std::unique_ptr advance_particle_queue; -extern std::unique_ptr surface_crossing_queue; -extern std::unique_ptr collision_queue; +extern SharedArray calculate_fuel_xs_queue; +extern SharedArray calculate_nonfuel_xs_queue; +extern SharedArray advance_particle_queue; +extern SharedArray surface_crossing_queue; +extern SharedArray collision_queue; -// Event queue lengths -extern int64_t calculate_fuel_xs_queue_length; -extern int64_t calculate_nonfuel_xs_queue_length; -extern int64_t advance_particle_queue_length; -extern int64_t surface_crossing_queue_length; -extern int64_t collision_queue_length; - -// Particle buffer. This is an allocated pointer rather than a vector as it -// will be shared between threads, so most vector methods would result in -// undefined or unintended behavior. -extern std::unique_ptr particles; +// Particle buffer +extern std::vector particles; } // namespace simulation @@ -82,11 +74,9 @@ void free_event_queues(void); //! Atomically adds a particle to the specified queue //! \param queue The queue to append the particle to -//! \param length A reference to the length variable for the queue //! \param p A pointer to the particle //! \param buffer_idx The particle's actual index in the particle buffer -void enqueue_particle(QueueItem* queue, int64_t& length, const Particle* p, - int64_t buffer_idx); +void enqueue_particle(SharedArray& queue, const Particle* p, int64_t buffer_idx); //! Enqueues a particle based on if it is in fuel or a non-fuel material //! \param buffer_idx The particle's actual index in the particle buffer @@ -98,9 +88,8 @@ void dispatch_xs_event(int64_t buffer_idx); void process_init_events(int64_t n_particles, int64_t source_offset); //! Executes the calculate XS event for all particles in this event's buffer -//! \param queue The XS lookup queue to use -//! \param n_particles The number of particles in this queue -void process_calculate_xs_events(QueueItem* queue, int64_t n_particles); +//! \param queue A reference to the desired XS lookup queue +void process_calculate_xs_events(SharedArray& queue); //! Executes the advance particle event for all particles in this event's buffer void process_advance_particle_events(); diff --git a/include/openmc/shared_array.h b/include/openmc/shared_array.h index 872621ece..e45568ceb 100644 --- a/include/openmc/shared_array.h +++ b/include/openmc/shared_array.h @@ -29,7 +29,6 @@ public: SharedArray(int64_t capacity) : capacity_(capacity) { data_ = std::make_unique(capacity); - size_ = 0; } //========================================================================== @@ -45,7 +44,7 @@ public: //! Increases the size of the SharedArray by one and returns an index to the //! last element of the array. - int64_t append() + int64_t thread_safe_append() { // Atomically capture the index we want to write to int64_t idx; @@ -70,9 +69,13 @@ public: } int64_t size() {return size_;} + + int64_t resize(int64_t size) {size_ = size;} int64_t capacity() {return capacity_;} + T* data() {return data_.get();} + }; } // namespace openmc diff --git a/src/event.cpp b/src/event.cpp index d3dde9c71..ed61349a9 100644 --- a/src/event.cpp +++ b/src/event.cpp @@ -12,19 +12,13 @@ namespace openmc { namespace simulation { -std::unique_ptr calculate_fuel_xs_queue; -std::unique_ptr calculate_nonfuel_xs_queue; -std::unique_ptr advance_particle_queue; -std::unique_ptr surface_crossing_queue; -std::unique_ptr collision_queue; +SharedArray calculate_fuel_xs_queue; +SharedArray calculate_nonfuel_xs_queue; +SharedArray advance_particle_queue; +SharedArray surface_crossing_queue; +SharedArray collision_queue; -int64_t calculate_fuel_xs_queue_length {0}; -int64_t calculate_nonfuel_xs_queue_length {0}; -int64_t advance_particle_queue_length {0}; -int64_t surface_crossing_queue_length {0}; -int64_t collision_queue_length {0}; - -std::unique_ptr particles; +std::vector particles; } // namespace simulation @@ -34,30 +28,29 @@ std::unique_ptr particles; void init_event_queues(int64_t n_particles) { - simulation::calculate_fuel_xs_queue = std::make_unique(n_particles); - simulation::calculate_nonfuel_xs_queue = std::make_unique(n_particles); - simulation::advance_particle_queue = std::make_unique(n_particles); - simulation::surface_crossing_queue = std::make_unique(n_particles); - simulation::collision_queue = std::make_unique(n_particles); - simulation::particles = std::make_unique(n_particles); + simulation::calculate_fuel_xs_queue.reserve(n_particles); + simulation::calculate_nonfuel_xs_queue.reserve(n_particles); + simulation::advance_particle_queue.reserve(n_particles); + simulation::surface_crossing_queue.reserve(n_particles); + simulation::collision_queue.reserve(n_particles); + + simulation::particles.resize(n_particles); } void free_event_queues(void) { - simulation::calculate_fuel_xs_queue.reset(); - simulation::calculate_nonfuel_xs_queue.reset(); - simulation::advance_particle_queue.reset(); - simulation::surface_crossing_queue.reset(); - simulation::collision_queue.reset(); - simulation::particles.reset(); + simulation::calculate_fuel_xs_queue.clear(); + simulation::calculate_nonfuel_xs_queue.clear(); + simulation::advance_particle_queue.clear(); + simulation::surface_crossing_queue.clear(); + simulation::collision_queue.clear(); + + simulation::particles.clear(); } -void enqueue_particle(QueueItem* queue, int64_t& length, const Particle* p, - int64_t buffer_idx) +void enqueue_particle(SharedArray& queue, const Particle* p, int64_t buffer_idx) { - int64_t idx; - #pragma omp atomic capture - idx = length++; + int64_t idx = queue.thread_safe_append(); queue[idx].idx = buffer_idx; queue[idx].E = p->E_; @@ -68,13 +61,10 @@ void enqueue_particle(QueueItem* queue, int64_t& length, const Particle* p, void dispatch_xs_event(int64_t buffer_idx) { Particle* p = &simulation::particles[buffer_idx]; - if (p->material_ == MATERIAL_VOID || - !model::materials[p->material_]->fissionable_) { - enqueue_particle(simulation::calculate_nonfuel_xs_queue.get(), - simulation::calculate_nonfuel_xs_queue_length, p, buffer_idx); + if (p->material_ == MATERIAL_VOID || !model::materials[p->material_]->fissionable_) { + enqueue_particle(simulation::calculate_nonfuel_xs_queue, p, buffer_idx); } else { - enqueue_particle(simulation::calculate_fuel_xs_queue.get(), - simulation::calculate_fuel_xs_queue_length, p, buffer_idx); + enqueue_particle(simulation::calculate_fuel_xs_queue, p, buffer_idx); } } @@ -89,7 +79,7 @@ void process_init_events(int64_t n_particles, int64_t source_offset) simulation::time_event_init.stop(); } -void process_calculate_xs_events(QueueItem* queue, int64_t n_particles) +void process_calculate_xs_events(SharedArray& queue) { simulation::time_event_calculate_xs.start(); @@ -99,21 +89,24 @@ void process_calculate_xs_events(QueueItem* queue, int64_t n_particles) // to C++17, std::sort is a serial only operation, which in this case // makes it too slow to be practical for most test problems. // - // std::sort(std::execution::par_unseq, queue, queue+n); + // std::sort(std::execution::par_unseq, queue.data(), queue.data() + n); + + int64_t offset = simulation::advance_particle_queue.size();; #pragma omp parallel for schedule(runtime) - for (int64_t i = 0; i < n_particles; i++) { + for (int64_t i = 0; i < queue.size(); i++) { Particle* p = &simulation::particles[queue[i].idx]; p->event_calculate_xs(); // After executing a calculate_xs event, particles will // always require an advance event. Therefore, we don't need to use // the protected enqueuing function. - int64_t offset = simulation::advance_particle_queue_length + i; - simulation::advance_particle_queue[offset] = queue[i]; + simulation::advance_particle_queue[offset + i] = queue[i]; } - simulation::advance_particle_queue_length += n_particles; + simulation::advance_particle_queue.resize(offset + queue.size()); + + queue.resize(0); simulation::time_event_calculate_xs.stop(); } @@ -123,19 +116,19 @@ void process_advance_particle_events() simulation::time_event_advance_particle.start(); #pragma omp parallel for schedule(runtime) - for (int64_t i = 0; i < simulation::advance_particle_queue_length; i++) { + for (int64_t i = 0; i < simulation::advance_particle_queue.size(); i++) { int64_t buffer_idx = simulation::advance_particle_queue[i].idx; Particle* p = &simulation::particles[buffer_idx]; p->event_advance(); if (p->collision_distance_ > p->boundary_.distance) { - enqueue_particle(simulation::surface_crossing_queue.get(), - simulation::surface_crossing_queue_length, p, buffer_idx); + enqueue_particle(simulation::surface_crossing_queue, p, buffer_idx); } else { - enqueue_particle(simulation::collision_queue.get(), - simulation::collision_queue_length, p, buffer_idx); + enqueue_particle(simulation::collision_queue, p, buffer_idx); } } + simulation::advance_particle_queue.resize(0); + simulation::time_event_advance_particle.stop(); } @@ -144,7 +137,7 @@ void process_surface_crossing_events() simulation::time_event_surface_crossing.start(); #pragma omp parallel for schedule(runtime) - for (int64_t i = 0; i < simulation::surface_crossing_queue_length; i++) { + for (int64_t i = 0; i < simulation::surface_crossing_queue.size(); i++) { int64_t buffer_index = simulation::surface_crossing_queue[i].idx; Particle* p = &simulation::particles[buffer_index]; p->event_cross_surface(); @@ -152,6 +145,8 @@ void process_surface_crossing_events() if (p->alive_) dispatch_xs_event(buffer_index); } + + simulation::surface_crossing_queue.resize(0); simulation::time_event_surface_crossing.stop(); } @@ -161,7 +156,7 @@ void process_collision_events() simulation::time_event_collision.start(); #pragma omp parallel for schedule(runtime) - for (int64_t i = 0; i < simulation::collision_queue_length; i++) { + for (int64_t i = 0; i < simulation::collision_queue.size(); i++) { int64_t buffer_index = simulation::collision_queue[i].idx; Particle* p = &simulation::particles[buffer_index]; p->event_collide(); @@ -170,6 +165,8 @@ void process_collision_events() dispatch_xs_event(buffer_index); } + simulation::collision_queue.resize(0); + simulation::time_event_collision.stop(); } diff --git a/src/physics.cpp b/src/physics.cpp index 30925d6dc..ca7dd8701 100644 --- a/src/physics.cpp +++ b/src/physics.cpp @@ -16,6 +16,7 @@ #include "openmc/secondary_uncorrelated.h" #include "openmc/search.h" #include "openmc/settings.h" +#include "openmc/shared_array.h" #include "openmc/simulation.h" #include "openmc/string_utils.h" #include "openmc/thermal.h" diff --git a/src/simulation.cpp b/src/simulation.cpp index 5569358c3..b23f1d2e2 100644 --- a/src/simulation.cpp +++ b/src/simulation.cpp @@ -626,32 +626,25 @@ void transport_event_based() while (true) { // Determine which event kernel has the longest queue int64_t max = std::max({ - simulation::calculate_fuel_xs_queue_length, - simulation::calculate_nonfuel_xs_queue_length, - simulation::advance_particle_queue_length, - simulation::surface_crossing_queue_length, - simulation::collision_queue_length}); + simulation::calculate_fuel_xs_queue.size(), + simulation::calculate_nonfuel_xs_queue.size(), + simulation::advance_particle_queue.size(), + simulation::surface_crossing_queue.size(), + simulation::collision_queue.size()}); // Execute event with the longest queue if (max == 0) { break; - } else if (max == simulation::calculate_fuel_xs_queue_length) { - process_calculate_xs_events(simulation::calculate_fuel_xs_queue.get(), - simulation::calculate_fuel_xs_queue_length); - simulation::calculate_fuel_xs_queue_length = 0; - } else if (max == simulation::calculate_nonfuel_xs_queue_length) { - process_calculate_xs_events(simulation::calculate_nonfuel_xs_queue.get(), - simulation::calculate_nonfuel_xs_queue_length); - simulation::calculate_nonfuel_xs_queue_length = 0; - } else if (max == simulation::advance_particle_queue_length) { + } else if (max == simulation::calculate_fuel_xs_queue.size()) { + process_calculate_xs_events(simulation::calculate_fuel_xs_queue); + } else if (max == simulation::calculate_nonfuel_xs_queue.size()) { + process_calculate_xs_events(simulation::calculate_nonfuel_xs_queue); + } else if (max == simulation::advance_particle_queue.size()) { process_advance_particle_events(); - simulation::advance_particle_queue_length = 0; - } else if (max == simulation::surface_crossing_queue_length) { + } else if (max == simulation::surface_crossing_queue.size()) { process_surface_crossing_events(); - simulation::surface_crossing_queue_length = 0; - } else if (max == simulation::collision_queue_length) { + } else if (max == simulation::collision_queue.size()) { process_collision_events(); - simulation::collision_queue_length = 0; } }