Tweaks to Universe::find_cell (#2662)

This commit is contained in:
Patrick Shriwise 2023-08-24 17:08:15 -05:00 committed by GitHub
parent ffbf1b9917
commit 06d64370c8
No known key found for this signature in database
GPG key ID: 4AEE18F83AFDEB23
6 changed files with 44 additions and 42 deletions

View file

@ -182,10 +182,11 @@ void DAGUniverse::init_geometry()
model::cell_map[c->id_] = model::cells.size();
} else {
warning(fmt::format("DAGMC Cell IDs: {}", dagmc_ids_for_dim(3)));
fatal_error(fmt::format("DAGMC Universe {} contains a cell with ID {}, which "
"already exists elsewhere in the geometry. Setting auto_geom_ids "
"to True when initiating the DAGMC Universe may "
"resolve this issue",
fatal_error(fmt::format(
"DAGMC Universe {} contains a cell with ID {}, which "
"already exists elsewhere in the geometry. Setting auto_geom_ids "
"to True when initiating the DAGMC Universe may "
"resolve this issue",
this->id_, c->id_));
}
@ -395,7 +396,7 @@ bool DAGUniverse::find_cell(Particle& p) const
// cells, place it in the implicit complement
bool found = Universe::find_cell(p);
if (!found && model::universe_map[this->id_] != model::root_universe) {
p.coord(p.n_coord() - 1).cell = implicit_complement_idx();
p.lowest_coord().cell = implicit_complement_idx();
found = true;
}
return found;
@ -581,7 +582,7 @@ std::pair<double, int32_t> DAGCell::distance(
p->history().reset();
}
const auto& univ = model::universes[p->coord(p->n_coord() - 1).universe];
const auto& univ = model::universes[p->lowest_coord().universe];
DAGUniverse* dag_univ = static_cast<DAGUniverse*>(univ.get());
if (!dag_univ)
@ -609,7 +610,8 @@ std::pair<double, int32_t> DAGCell::distance(
// into the implicit complement on the other side where no intersection will
// be found. Treating this as a lost particle is problematic when plotting.
// Instead, the infinite distance and invalid surface index are returned.
if (settings::run_mode == RunMode::PLOTTING) return {INFTY, -1};
if (settings::run_mode == RunMode::PLOTTING)
return {INFTY, -1};
// the particle should be marked as lost immediately if an intersection
// isn't found in a volume that is not the implicit complement. In the case
@ -739,7 +741,8 @@ void check_dagmc_root_univ()
}
}
int32_t next_cell(int32_t surf, int32_t curr_cell, int32_t univ) {
int32_t next_cell(int32_t surf, int32_t curr_cell, int32_t univ)
{
auto surfp = dynamic_cast<DAGSurface*>(model::surfaces[surf - 1].get());
auto cellp = dynamic_cast<DAGCell*>(model::cells[curr_cell].get());
auto univp = static_cast<DAGUniverse*>(model::universes[univ].get());
@ -750,11 +753,12 @@ int32_t next_cell(int32_t surf, int32_t curr_cell, int32_t univ) {
cellp->dagmc_ptr()->entity_by_index(3, cellp->dag_index());
moab::EntityHandle new_vol;
moab::ErrorCode rval = cellp->dagmc_ptr()->next_vol(surf_handle, curr_vol, new_vol);
if (rval != moab::MB_SUCCESS) return -1;
moab::ErrorCode rval =
cellp->dagmc_ptr()->next_vol(surf_handle, curr_vol, new_vol);
if (rval != moab::MB_SUCCESS)
return -1;
return cellp->dagmc_ptr()->index_by_handle(new_vol) +
univp->cell_idx_offset_;
return cellp->dagmc_ptr()->index_by_handle(new_vol) + univp->cell_idx_offset_;
}
} // namespace openmc

View file

@ -110,7 +110,7 @@ bool find_cell_inner(Particle& p, const NeighborList* neighbor_list)
i_cell = *it;
// Make sure the search cell is in the same universe.
int i_universe = p.coord(p.n_coord() - 1).universe;
int i_universe = p.lowest_coord().universe;
if (model::cells[i_cell]->universe_ != i_universe)
continue;
@ -119,7 +119,7 @@ bool find_cell_inner(Particle& p, const NeighborList* neighbor_list)
Direction u {p.u_local()};
auto surf = p.surface();
if (model::cells[i_cell]->contains(r, u, surf)) {
p.coord(p.n_coord() - 1).cell = i_cell;
p.lowest_coord().cell = i_cell;
found = true;
break;
}
@ -145,7 +145,7 @@ bool find_cell_inner(Particle& p, const NeighborList* neighbor_list)
// code below this conditional, we set i_cell back to C_NONE to indicate
// that.
if (i_cell == C_NONE) {
int i_universe = p.coord(p.n_coord() - 1).universe;
int i_universe = p.lowest_coord().universe;
const auto& univ {model::universes[i_universe]};
found = univ->find_cell(p);
}
@ -153,7 +153,7 @@ bool find_cell_inner(Particle& p, const NeighborList* neighbor_list)
if (!found) {
return found;
}
i_cell = p.coord(p.n_coord() - 1).cell;
i_cell = p.lowest_coord().cell;
// Announce the cell that the particle is entering.
if (found && (settings::verbosity >= 10 || p.trace())) {
@ -289,7 +289,7 @@ bool neighbor_list_find_cell(Particle& p)
bool exhaustive_find_cell(Particle& p)
{
int i_universe = p.coord(p.n_coord() - 1).universe;
int i_universe = p.lowest_coord().universe;
if (i_universe == C_NONE) {
p.coord(0).universe = model::root_universe;
p.n_coord() = 1;
@ -306,7 +306,7 @@ bool exhaustive_find_cell(Particle& p)
void cross_lattice(Particle& p, const BoundaryInfo& boundary)
{
auto& coord {p.coord(p.n_coord() - 1)};
auto& coord {p.lowest_coord()};
auto& lat {*model::lattices[coord.lattice]};
if (settings::verbosity >= 10 || p.trace()) {
@ -344,7 +344,7 @@ void cross_lattice(Particle& p, const BoundaryInfo& boundary)
} else {
// Find cell in next lattice element.
p.coord(p.n_coord() - 1).universe = lat[coord.lattice_i];
p.lowest_coord().universe = lat[coord.lattice_i];
bool found = exhaustive_find_cell(p);
if (!found) {
@ -472,7 +472,7 @@ extern "C" int openmc_find_cell(
return OPENMC_E_GEOMETRY;
}
*index = p.coord(p.n_coord() - 1).cell;
*index = p.lowest_coord().cell;
*instance = p.cell_instance();
return 0;
}

View file

@ -147,7 +147,7 @@ void Particle::event_calculate_xs()
// If the cell hasn't been determined based on the particle's location,
// initiate a search for the current cell. This generally happens at the
// beginning of the history and again for any secondary particles
if (coord(n_coord() - 1).cell == C_NONE) {
if (lowest_coord().cell == C_NONE) {
if (!exhaustive_find_cell(*this)) {
mark_as_lost(
"Could not find the cell containing particle " + std::to_string(id()));
@ -156,7 +156,7 @@ void Particle::event_calculate_xs()
// Set birth cell attribute
if (cell_born() == C_NONE)
cell_born() = coord(n_coord() - 1).cell;
cell_born() = lowest_coord().cell;
}
// Write particle track.
@ -392,7 +392,7 @@ void Particle::event_revive_from_secondary()
// Since the birth cell of the particle has not been set we
// have to determine it before the energy of the secondary particle can be
// removed from the pulse-height of this cell.
if (coord(n_coord() - 1).cell == C_NONE) {
if (lowest_coord().cell == C_NONE) {
if (!exhaustive_find_cell(*this)) {
mark_as_lost("Could not find the cell containing particle " +
std::to_string(id()));
@ -400,7 +400,7 @@ void Particle::event_revive_from_secondary()
}
// Set birth cell attribute
if (cell_born() == C_NONE)
cell_born() = coord(n_coord() - 1).cell;
cell_born() = lowest_coord().cell;
}
pht_secondary_particles();
}
@ -456,7 +456,7 @@ void Particle::pht_collision_energy()
// determine index of cell in pulse_height_cells
auto it = std::find(model::pulse_height_cells.begin(),
model::pulse_height_cells.end(), coord(n_coord() - 1).cell);
model::pulse_height_cells.end(), lowest_coord().cell);
if (it != model::pulse_height_cells.end()) {
int index = std::distance(model::pulse_height_cells.begin(), it);
@ -535,7 +535,7 @@ void Particle::cross_surface()
material_last() = material();
sqrtkT_last() = sqrtkT();
// set new cell value
coord(n_coord() - 1).cell = i_cell;
lowest_coord().cell = i_cell;
cell_instance() = 0;
material() = model::cells[i_cell]->material_[0];
sqrtkT() = model::cells[i_cell]->sqrtkT_[0];

View file

@ -58,7 +58,7 @@ void IdData::set_value(size_t y, size_t x, const Particle& p, int level)
}
// set material data
Cell* c = model::cells.at(p.coord(p.n_coord() - 1).cell).get();
Cell* c = model::cells.at(p.lowest_coord().cell).get();
if (p.material() == MATERIAL_VOID) {
data_(y, x, 2) = MATERIAL_VOID;
return;
@ -79,7 +79,7 @@ PropertyData::PropertyData(size_t h_res, size_t v_res)
void PropertyData::set_value(size_t y, size_t x, const Particle& p, int level)
{
Cell* c = model::cells.at(p.coord(p.n_coord() - 1).cell).get();
Cell* c = model::cells.at(p.lowest_coord().cell).get();
data_(y, x, 0) = (p.sqrtkT() * p.sqrtkT()) / K_BOLTZMANN;
if (c->type_ != Fill::UNIVERSE && p.material() != MATERIAL_VOID) {
Material* m = model::materials.at(p.material()).get();
@ -1327,9 +1327,8 @@ void ProjectionPlot::create_output() const
// edges on the model boundary for the same cell.
if (first_inside_model) {
this_line_segments[tid][horiz].emplace_back(
color_by_ == PlotColorBy::mats
? p.material()
: p.coord(p.n_coord() - 1).cell,
color_by_ == PlotColorBy::mats ? p.material()
: p.lowest_coord().cell,
0.0, first_surface);
first_inside_model = false;
}
@ -1339,7 +1338,7 @@ void ProjectionPlot::create_output() const
auto dist = distance_to_boundary(p);
this_line_segments[tid][horiz].emplace_back(
color_by_ == PlotColorBy::mats ? p.material()
: p.coord(p.n_coord() - 1).cell,
: p.lowest_coord().cell,
dist.distance, std::abs(dist.surface_index));
// Advance particle

View file

@ -72,7 +72,7 @@ void CellInstanceFilter::set_cell_instances(gsl::span<CellInstance> instances)
void CellInstanceFilter::get_all_bins(
const Particle& p, TallyEstimator estimator, FilterMatch& match) const
{
gsl::index index_cell = p.coord(p.n_coord() - 1).cell;
gsl::index index_cell = p.lowest_coord().cell;
gsl::index instance = p.cell_instance();
if (cells_.count(index_cell) > 0) {

View file

@ -41,18 +41,17 @@ bool Universe::find_cell(Particle& p) const
const auto& cells {
!partitioner_ ? cells_ : partitioner_->get_cells(p.r_local(), p.u_local())};
for (auto it = cells.begin(); it != cells.end(); it++) {
int32_t i_cell = *it;
int32_t i_univ = p.coord(p.n_coord() - 1).universe;
Position r {p.r_local()};
Position u {p.u_local()};
auto surf = p.surface();
int32_t i_univ = p.lowest_coord().universe;
for (auto i_cell : cells) {
if (model::cells[i_cell]->universe_ != i_univ)
continue;
// Check if this cell contains the particle;
Position r {p.r_local()};
Direction u {p.u_local()};
auto surf = p.surface();
// Check if this cell contains the particle
if (model::cells[i_cell]->contains(r, u, surf)) {
p.coord(p.n_coord() - 1).cell = i_cell;
p.lowest_coord().cell = i_cell;
return true;
}
}