diff --git a/include/openmc/cell.h b/include/openmc/cell.h index 032475ce982..d01020f8e7a 100644 --- a/include/openmc/cell.h +++ b/include/openmc/cell.h @@ -29,7 +29,6 @@ namespace openmc { enum class Fill { MATERIAL, UNIVERSE, LATTICE }; -// TODO: Convert to enum constexpr int32_t OP_LEFT_PAREN {std::numeric_limits::max()}; constexpr int32_t OP_RIGHT_PAREN {std::numeric_limits::max() - 1}; constexpr int32_t OP_COMPLEMENT {std::numeric_limits::max() - 2}; diff --git a/include/openmc/constants.h b/include/openmc/constants.h index 605ae1839d8..a83a4a07d4b 100644 --- a/include/openmc/constants.h +++ b/include/openmc/constants.h @@ -350,6 +350,10 @@ enum class RandomRaySourceShape { FLAT, LINEAR, LINEAR_XY }; enum class GeometryType { CSG, DAG }; +// a surface token cannot be zero due to the unsigned nature of zero for integer +// representations. This value represents no surface. +constexpr int32_t SURFACE_NONE {0}; + } // namespace openmc #endif // OPENMC_CONSTANTS_H diff --git a/include/openmc/particle_data.h b/include/openmc/particle_data.h index 164148cce10..66137c5f688 100644 --- a/include/openmc/particle_data.h +++ b/include/openmc/particle_data.h @@ -185,10 +185,14 @@ struct CacheDataMG { struct BoundaryInfo { double distance {INFINITY}; //!< distance to nearest boundary - int surface_index {0}; //!< if boundary is surface, index in surfaces vector - int coord_level; //!< coordinate level after crossing boundary + int surface { + SURFACE_NONE}; //!< surface token, non-zero if boundary is surface + int coord_level; //!< coordinate level after crossing boundary array lattice_translation {}; //!< which way lattice indices will change + + // TODO: off-by-one + int surface_index() const { return std::abs(surface) - 1; } }; /* @@ -226,7 +230,7 @@ class GeometryState { void init_from_r_u(Position r_a, Direction u_a) { clear(); - surface() = 0; + surface() = SURFACE_NONE; material() = C_NONE; r() = r_a; u() = u_a; @@ -296,10 +300,17 @@ class GeometryState { Direction& u_local() { return coord_[n_coord_ - 1].u; } const Direction& u_local() const { return coord_[n_coord_ - 1].u; } - // Surface that the particle is on + // Surface token for the surface that the particle is currently on int& surface() { return surface_; } const int& surface() const { return surface_; } + // Surface index based on the current value of the surface_ attribute + int surface_index() const + { + // TODO: off-by-one + return std::abs(surface_) - 1; + } + // Boundary information BoundaryInfo& boundary() { return boundary_; } @@ -337,7 +348,8 @@ class GeometryState { Position r_last_; //!< previous coordinates Direction u_last_; //!< previous direction coordinates - int surface_ {0}; //!< index for surface particle is on + int surface_ { + SURFACE_NONE}; //!< surface token for surface the particle is currently on BoundaryInfo boundary_; //!< Info about the next intersection diff --git a/src/boundary_condition.cpp b/src/boundary_condition.cpp index b58054dce8c..7216ac89649 100644 --- a/src/boundary_condition.cpp +++ b/src/boundary_condition.cpp @@ -129,8 +129,7 @@ TranslationalPeriodicBC::TranslationalPeriodicBC(int i_surf, int j_surf) void TranslationalPeriodicBC::handle_particle( Particle& p, const Surface& surf) const { - // TODO: off-by-one on surface indices throughout this function. - int i_particle_surf = std::abs(p.surface()) - 1; + int i_particle_surf = p.surface_index(); // Figure out which of the two BC surfaces were struck then find the // particle's new location and surface. @@ -255,8 +254,7 @@ RotationalPeriodicBC::RotationalPeriodicBC(int i_surf, int j_surf) void RotationalPeriodicBC::handle_particle( Particle& p, const Surface& surf) const { - // TODO: off-by-one on surface indices throughout this function. - int i_particle_surf = std::abs(p.surface()) - 1; + int i_particle_surf = p.surface_index(); // Figure out which of the two BC surfaces were struck to figure out if a // forward or backward rotation is required. Specify the other surface as diff --git a/src/dagmc.cpp b/src/dagmc.cpp index 134e31ecf3c..94c220850f3 100644 --- a/src/dagmc.cpp +++ b/src/dagmc.cpp @@ -847,7 +847,7 @@ void check_dagmc_root_univ() int32_t next_cell(int32_t surf, int32_t curr_cell, int32_t univ) { - auto surfp = dynamic_cast(model::surfaces[surf - 1].get()); + auto surfp = dynamic_cast(model::surfaces[surf].get()); auto cellp = dynamic_cast(model::cells[curr_cell].get()); auto univp = static_cast(model::universes[univ].get()); diff --git a/src/geometry.cpp b/src/geometry.cpp index 445b19faac1..9e975c00dee 100644 --- a/src/geometry.cpp +++ b/src/geometry.cpp @@ -421,15 +421,15 @@ BoundaryInfo distance_to_boundary(GeometryState& p) // have to explicitly check which half-space the particle would be // traveling into if the surface is crossed if (c.is_simple() || d == INFTY) { - info.surface_index = level_surf_cross; + info.surface = level_surf_cross; } else { Position r_hit = r + d_surf * u; Surface& surf {*model::surfaces[std::abs(level_surf_cross) - 1]}; Direction norm = surf.normal(r_hit); if (u.dot(norm) > 0) { - info.surface_index = std::abs(level_surf_cross); + info.surface = std::abs(level_surf_cross); } else { - info.surface_index = -std::abs(level_surf_cross); + info.surface = -std::abs(level_surf_cross); } } @@ -441,7 +441,7 @@ BoundaryInfo distance_to_boundary(GeometryState& p) } else { if (d == INFINITY || (d - d_lat) / d >= FP_REL_PRECISION) { d = d_lat; - info.surface_index = 0; + info.surface = SURFACE_NONE; info.lattice_translation = level_lat_trans; info.coord_level = i + 1; } diff --git a/src/output.cpp b/src/output.cpp index a430fe9a6c6..8494450c409 100644 --- a/src/output.cpp +++ b/src/output.cpp @@ -200,9 +200,9 @@ void print_particle(Particle& p) } // Display miscellaneous info. - if (p.surface() != 0) { + if (p.surface() != SURFACE_NONE) { // Surfaces identifiers are >= 1, but indices are >= 0 so we need -1 - const Surface& surf {*model::surfaces[std::abs(p.surface()) - 1]}; + const Surface& surf {*model::surfaces[p.surface_index()]}; fmt::print(" Surface = {}\n", (p.surface() > 0) ? surf.id_ : -surf.id_); } fmt::print(" Weight = {}\n", p.wgt()); diff --git a/src/particle.cpp b/src/particle.cpp index 0ea8650143d..65e8065fded 100644 --- a/src/particle.cpp +++ b/src/particle.cpp @@ -108,7 +108,7 @@ void Particle::from_source(const SourceSite* src) { // Reset some attributes clear(); - surface() = 0; + surface() = SURFACE_NONE; cell_born() = C_NONE; material() = C_NONE; n_collision() = 0; @@ -277,7 +277,7 @@ void Particle::event_cross_surface() n_coord_last() = n_coord(); // Set surface that particle is on and adjust coordinate levels - surface() = boundary().surface_index; + surface() = boundary().surface; n_coord() = boundary().coord_level; if (boundary().lattice_translation[0] != 0 || @@ -291,7 +291,7 @@ void Particle::event_cross_surface() } else { // Particle crosses surface // TODO: off-by-one - const auto& surf {model::surfaces[std::abs(surface()) - 1].get()}; + const auto& surf {model::surfaces[surface_index()].get()}; // If BC, add particle to surface source before crossing surface if (surf->surf_source_ && surf->bc_) { add_surf_source_to_bank(*this, *surf); @@ -328,7 +328,7 @@ void Particle::event_collide() score_surface_tally(*this, model::active_meshsurf_tallies); // Clear surface component - surface() = 0; + surface() = SURFACE_NONE; if (settings::run_CE) { collision(*this); @@ -549,7 +549,7 @@ void Particle::cross_surface(const Surface& surf) #ifdef DAGMC // in DAGMC, we know what the next cell should be if (surf.geom_type() == GeometryType::DAG) { - int32_t i_cell = next_cell(std::abs(surface()), cell_last(n_coord() - 1), + int32_t i_cell = next_cell(surface_index(), cell_last(n_coord() - 1), lowest_coord().universe) - 1; // save material and temp @@ -587,7 +587,7 @@ void Particle::cross_surface(const Surface& surf) // the particle is really traveling tangent to a surface, if we move it // forward a tiny bit it should fix the problem. - surface() = 0; + surface() = SURFACE_NONE; n_coord() = 1; r() += TINY_BIT * u(); diff --git a/src/plot.cpp b/src/plot.cpp index 85131d67a01..b36ed6f5d93 100644 --- a/src/plot.cpp +++ b/src/plot.cpp @@ -1301,7 +1301,7 @@ void ProjectionPlot::create_output() const while (intersection_found) { bool inside_cell = false; - int32_t i_surface = std::abs(p.surface()) - 1; + int32_t i_surface = p.surface_index(); if (i_surface > 0 && model::surfaces[i_surface]->geom_type() == GeometryType::DAG) { #ifdef DAGMC @@ -1334,13 +1334,13 @@ void ProjectionPlot::create_output() const this_line_segments[tid][horiz].emplace_back( color_by_ == PlotColorBy::mats ? p.material() : p.lowest_coord().cell, - dist.distance, std::abs(dist.surface_index)); + dist.distance, std::abs(dist.surface)); // Advance particle for (int lev = 0; lev < p.n_coord(); ++lev) { p.coord(lev).r += dist.distance * p.coord(lev).u; } - p.surface() = dist.surface_index; + p.surface() = dist.surface; p.n_coord_last() = p.n_coord(); p.n_coord() = dist.coord_level; if (dist.lattice_translation[0] != 0 || diff --git a/src/tallies/filter_musurface.cpp b/src/tallies/filter_musurface.cpp index 58fe39b9113..340149d4cff 100644 --- a/src/tallies/filter_musurface.cpp +++ b/src/tallies/filter_musurface.cpp @@ -12,7 +12,7 @@ void MuSurfaceFilter::get_all_bins( const Particle& p, TallyEstimator estimator, FilterMatch& match) const { // Get surface normal (and make sure it is a unit vector) - const auto surf {model::surfaces[std::abs(p.surface()) - 1].get()}; + const auto surf {model::surfaces[p.surface_index()].get()}; auto n = surf->normal(p.r()); n /= n.norm(); diff --git a/src/tallies/filter_surface.cpp b/src/tallies/filter_surface.cpp index a84d1772822..1fbf8d44e38 100644 --- a/src/tallies/filter_surface.cpp +++ b/src/tallies/filter_surface.cpp @@ -47,7 +47,7 @@ void SurfaceFilter::set_surfaces(gsl::span surfaces) void SurfaceFilter::get_all_bins( const Particle& p, TallyEstimator estimator, FilterMatch& match) const { - auto search = map_.find(std::abs(p.surface()) - 1); + auto search = map_.find(p.surface_index()); if (search != map_.end()) { match.bins_.push_back(search->second); if (p.surface() < 0) {