crp/OpenMC
1
0
Fork 0
mirror of https://github.com/openmc-dev/openmc.git synced 2026-07-26 21:25:36 -04:00
Code Issues Projects Releases Packages Wiki Activity Actions

Simplify neighbor list implementation

Browse source
This commit is contained in:
Sterling Harper 2018-12-18 15:02:07 -05:00
parent 066f698da7
commit 35238259c0
6 changed files with 24 additions and 168 deletions
Download patch file Download diff file Expand all files Collapse all files

2
include/openmc/cell.h
View file

@ -105,7 +105,7 @@ public:
bool simple_; //!< Does the region contain only intersections?
//! \brief Neighboring cells in the same universe.
NeighborList neighbors;
NeighborList neighbors_;
Position translation_ {0, 0, 0}; //!< Translation vector for filled universe

9
include/openmc/geometry_aux.h
View file

@ -88,15 +88,6 @@ distribcell_path(int32_t target_cell, int32_t map, int32_t target_offset);
extern "C" int maximum_levels(int32_t univ);
//==============================================================================
//! Perform any geometry operations that should be done between generations.
//!
//! Currently, this just includes converting cell neighbor linked lists into
//! vectors, but similar optimizations could be added here in the future.
//==============================================================================
void geometry_finalize_generation();
//==============================================================================
//! Deallocates global vectors and maps for cells, universes, and lattices.
//==============================================================================

159
include/openmc/neighbor_list.h
View file

@ -10,174 +10,51 @@
namespace openmc {
// Forward declare the neighbor list iterator type.
template <typename T_value, typename T_prefix_iter, typename T_suffix_iter>
class NeighborListIter;
//==============================================================================
//! An efficient, threadsafe, dynamic container for listing neighboring cells.
//! A threadsafe, dynamic container for listing neighboring cells.
//
//! This container is a minor improvement upon a linked list. The linked list
//! allows for threadsafe dynamic growth; any number of threads can safely
//! read from the list without locks or reference counting. Write access must
//! be protected with a lock, but write events are fairly rare for neighbor
//! lists. This container also provides the option to flush the linked list
//! into a vector at some convenient time (like between generations where the
//! neighbor lists are not being read) for contiguous data. The resulting
//! performance benefit is small but consistent.
//! This container is a reduced interface for a linked list with an added OpenMP
//! lock for write operations. It allows for threadsafe dynamic growth; any
//! number of threads can safely read data without locks or reference counting.
//==============================================================================
class NeighborList
{
public:
using value_type = int32_t;
using prefix_iter = std::vector<value_type>::iterator;
using suffix_iter = std::forward_list<value_type>::iterator;
using iterator = NeighborListIter<value_type, prefix_iter, suffix_iter>;
using const_iterator = std::forward_list<value_type>::const_iterator;
//! Attempt to add an element.
// Attempt to add an element.
//
//! If the relevant OpenMP lock is currently owned by another thread, this
//! function will return without actually modifying the data. It is assumed
//! that returning the transport calculation and possibly re-adding the
//! element later is faster than waiting on the lock to be released.
// If the relevant OpenMP lock is currently owned by another thread, this
// function will return without actually modifying the data. It has been
// found that returning the transport calculation and possibly re-adding the
// element later is slightly faster than waiting on the lock to be released.
void push(int new_elem)
{
// Try to acquire the lock.
std::unique_lock<ThreadMutex> lock(mutex_, std::try_to_lock);
if (lock) {
// It is possible another thread already added this element to the suffix
// It is possible another thread already added this element to the list
// while this thread was searching for a cell so make sure the given
// element isn't a duplicate before adding it.
if (std::find(suffix_.cbegin(), suffix_.cend(), new_elem)
== suffix_.cend()) {
suffix_.push_front(new_elem);
if (std::find(list_.cbegin(), list_.cend(), new_elem) == list_.cend()) {
list_.push_front(new_elem);
}
}
}
//! Move data from the linked-list suffix to the consecutive vector prefix.
//
//! The consecutive data slightly improves runtime (likely due to cache
//! locality). Note that this function is not guaranteed threadsafe---the
//! caller is responsible for making sure only one thread at a time calls this
//! function.
void make_consecutive()
{
while (!suffix_.empty()) {
prefix_.push_back(suffix_.front());
suffix_.pop_front();
}
}
const_iterator cbegin() const
{return list_.cbegin();}
iterator begin();
const_iterator cend() const
{return list_.cend();}
iterator end();
private:
std::vector<value_type> prefix_;
std::forward_list<value_type> suffix_;
std::forward_list<value_type> list_;
ThreadMutex mutex_;
friend class NeighborListIter<value_type, prefix_iter, suffix_iter>;
};
//==============================================================================
template <typename T_value, typename T_prefix_iter, typename T_suffix_iter>
class NeighborListIter
{
public:
// Construct from a prefix iterator.
NeighborListIter(NeighborList* nl, T_prefix_iter it)
{
// If we were given an iterator to the end of the prefix, immediately switch
// over to suffix mode.
base_ = nl;
if (it != base_->prefix_.end()) {
in_prefix_ = true;
prefix_iter_ = it;
} else {
in_prefix_ = false;
suffix_iter_ = base_->suffix_.begin();
}
}
// Construct from a suffix iterator.
NeighborListIter(NeighborList* nl, T_suffix_iter it)
{
in_prefix_ = false;
suffix_iter_ = it;
}
T_value operator*()
{
if (in_prefix_) {
return *prefix_iter_;
} else {
return *suffix_iter_;
}
}
bool operator==(const NeighborListIter& other)
{
if (in_prefix_ != other.in_prefix_) return false;
if (in_prefix_) {
return prefix_iter_ == other.prefix_iter_;
} else {
return suffix_iter_ == other.suffix_iter_;
}
}
bool operator!=(const NeighborListIter& other)
{return !(*this == other);}
NeighborListIter& operator++()
{
if (in_prefix_) {
// We are in the prefix so increment the prefix iterator.
++prefix_iter_;
// If we've reached the end of the prefix, switch to the suffix iterator.
if (prefix_iter_ == base_->prefix_.end()) {
in_prefix_ = false;
suffix_iter_ = base_->suffix_.begin();
}
} else {
// We are in the suffix so increment the suffix iterator.
++suffix_iter_;
}
return *this;
}
private:
NeighborList* base_;
// This type essentially wraps the implementation for two different external
// iterators. A union is used to contain one iterator or the other, and the
// in_prefix_ flag indicates which version of that union is valid.
union {
T_prefix_iter prefix_iter_;
T_suffix_iter suffix_iter_;
};
bool in_prefix_;
};
//==============================================================================
inline NeighborList::iterator
NeighborList::begin()
{
return iterator(this, prefix_.begin());
}
inline NeighborList::iterator
NeighborList::end()
{
return iterator(this, suffix_.end());
}
} // namespace openmc
#endif // OPENMC_NEIGHBOR_LIST_H

10
src/geometry.cpp
View file

@ -63,14 +63,14 @@ check_cell_overlap(Particle* p)
//==============================================================================
bool
find_cell_inner(Particle* p, NeighborList* neighbor_list)
find_cell_inner(Particle* p, const NeighborList* neighbor_list)
{
// Find which cell of this universe the particle is in. Use the neighbor list
// to shorten the search if one was provided.
bool found = false;
int32_t i_cell;
if (neighbor_list) {
for (auto it = neighbor_list->begin(); it != neighbor_list->end(); ++it) {
for (auto it = neighbor_list->cbegin(); it != neighbor_list->cend(); ++it) {
i_cell = *it;
// Make sure the search cell is in the same universe.
@ -91,7 +91,7 @@ find_cell_inner(Particle* p, NeighborList* neighbor_list)
} else {
int i_universe = p->coord[p->n_coord-1].universe;
const auto& cells {model::universes[i_universe]->cells_};
for (auto it = cells.begin(); it != cells.end(); it++) {
for (auto it = cells.cbegin(); it != cells.cend(); it++) {
i_cell = *it;
// Make sure the search cell is in the same universe.
@ -290,14 +290,14 @@ find_cell(Particle* p, bool use_neighbor_lists)
// Search for the particle in that cell's neighbor list. Return if we
// found the particle.
bool found = find_cell_inner(p, &c.neighbors);
bool found = find_cell_inner(p, &c.neighbors_);
if (found) return found;
// The particle could not be found in the neighbor list. Try searching all
// cells in this universe, and update the neighbor list if we find a new
// neighboring cell.
found = find_cell_inner(p, nullptr);
if (found) c.neighbors.push(p->coord[coord_lvl].cell);
if (found) c.neighbors_.push(p->coord[coord_lvl].cell);
return found;
} else {

8
src/geometry_aux.cpp
View file

@ -432,14 +432,6 @@ maximum_levels(int32_t univ)
//==============================================================================
void
geometry_finalize_generation()
{
for (Cell* c : model::cells) c->neighbors.make_consecutive();
}
//==============================================================================
void
free_memory_geometry_c()
{

4
src/simulation.cpp
View file

@ -4,7 +4,6 @@
#include "openmc/container_util.h"
#include "openmc/eigenvalue.h"
#include "openmc/error.h"
#include "openmc/geometry_aux.h"
#include "openmc/message_passing.h"
#include "openmc/output.h"
#include "openmc/particle.h"
@ -450,9 +449,6 @@ void finalize_generation()
// For fixed-source mode, we need to sample the external source
fill_source_bank_fixedsource();
}
// Clean up neighbor lists
geometry_finalize_generation();
}
void initialize_history(Particle* p, int64_t index_source)