From 230e054ec39134369fd908a1c516293eb449b8e4 Mon Sep 17 00:00:00 2001 From: Patrick Shriwise Date: Thu, 1 Apr 2021 23:24:06 -0500 Subject: [PATCH] Move DAGCell implementation into the dagmc files. --- include/openmc/cell.h | 23 ----------- include/openmc/dagmc.h | 22 ++++++++++- src/cell.cpp | 90 ------------------------------------------ src/dagmc.cpp | 90 ++++++++++++++++++++++++++++++++++++++++++ 4 files changed, 111 insertions(+), 114 deletions(-) diff --git a/include/openmc/cell.h b/include/openmc/cell.h index 54e06f782e..ad17b6f78b 100644 --- a/include/openmc/cell.h +++ b/include/openmc/cell.h @@ -260,28 +260,6 @@ protected: vector::iterator start, const vector& rpn); }; -//============================================================================== - -#ifdef DAGMC -class DAGCell : public Cell -{ -public: - DAGCell(); - - bool contains(Position r, Direction u, int32_t on_surface) const; - - std::pair - distance(Position r, Direction u, int32_t on_surface, Particle* p) const; - - BoundingBox bounding_box() const; - - void to_hdf5_inner(hid_t group_id) const; - - std::shared_ptr dagmc_ptr_; //!< Pointer to DagMC instance - int32_t dag_index_; //!< DagMC index of cell -}; -#endif - //============================================================================== //! Speeds up geometry searches by grouping cells in a search tree. // @@ -351,7 +329,6 @@ void read_cells(pugi::xml_node node); #ifdef DAGMC class DAGUniverse; -int32_t next_cell(DAGUniverse* dag_univ, DAGCell* cur_cell, DAGSurface* surf_xed); #endif } // namespace openmc diff --git a/include/openmc/dagmc.h b/include/openmc/dagmc.h index 76bd0a4ec3..1dc1c2bc37 100644 --- a/include/openmc/dagmc.h +++ b/include/openmc/dagmc.h @@ -10,13 +10,32 @@ extern "C" const bool DAGMC_ENABLED; #include "DagMC.hpp" #include "openmc/cell.h" +#include "openmc/particle.h" #include "openmc/position.h" +#include "openmc/surface.h" #include "openmc/xml_interface.h" class UWUW; namespace openmc { +class DAGCell : public Cell { +public: + DAGCell(); + + bool contains(Position r, Direction u, int32_t on_surface) const; + + std::pair + distance(Position r, Direction u, int32_t on_surface, Particle* p) const; + + BoundingBox bounding_box() const; + + void to_hdf5_inner(hid_t group_id) const; + + std::shared_ptr dagmc_ptr_; //!< Pointer to DagMC instance + int32_t dag_index_; //!< DagMC index of cell +}; + class DAGUniverse : public Universe { public: @@ -26,7 +45,6 @@ public: void initialize(); //!< Sets up the DAGMC instance and OpenMC internals void read_uwuw_materials(); //!< Reads UWUW materials and returns an ID map - bool uses_uwuw() const; int32_t implicit_complement_idx() const; @@ -52,7 +70,9 @@ public: bool adjust_material_ids_; }; +// DAGMC Functions void read_dagmc_universes(pugi::xml_node node); +int32_t next_cell(DAGUniverse* dag_univ, DAGCell* cur_cell, DAGSurface* surf_xed); } // namespace openmc diff --git a/src/cell.cpp b/src/cell.cpp index ad711573b2..23786727f1 100644 --- a/src/cell.cpp +++ b/src/cell.cpp @@ -825,81 +825,6 @@ CSGCell::contains_complex(Position r, Direction u, int32_t on_surface) const } } -//============================================================================== -// DAGMC Cell implementation -//============================================================================== -#ifdef DAGMC -DAGCell::DAGCell() : Cell{} { - geom_type_ = GeometryType::DAG; - simple_ = true; -}; - -std::pair -DAGCell::distance(Position r, Direction u, int32_t on_surface, Particle* p) const -{ - Expects(p); - // if we've changed direction or we're not on a surface, - // reset the history and update last direction - if (u != p->last_dir()) { p->last_dir() = u; p->history().reset(); } - if (on_surface == 0) { p->history().reset(); } - - const auto& univ = model::universes[p->coord_[p->n_coord_ - 1].universe]; - - DAGUniverse* dag_univ = static_cast(univ.get()); - if (!dag_univ) fatal_error("DAGMC call made for particle in a non-DAGMC universe"); - - moab::ErrorCode rval; - moab::EntityHandle vol = dagmc_ptr_->entity_by_index(3, dag_index_); - moab::EntityHandle hit_surf; - double dist; - double pnt[3] = {r.x, r.y, r.z}; - double dir[3] = {u.x, u.y, u.z}; - rval = dagmc_ptr_->ray_fire(vol, pnt, dir, hit_surf, dist, &p->history()); - MB_CHK_ERR_CONT(rval); - int surf_idx; - if (hit_surf != 0) { - surf_idx = dag_univ->surf_idx_offset_ + dagmc_ptr_->index_by_handle(hit_surf); - } else { - // indicate that particle is lost - surf_idx = -1; - dist = INFINITY; - if (!dagmc_ptr_->is_implicit_complement(vol) || model::universe_map[dag_univ->id_] == model::root_universe) { - p->mark_as_lost(fmt::format("No intersection found with DAGMC cell {}", id_)); - } - } - - return {dist, surf_idx}; -} - -bool DAGCell::contains(Position r, Direction u, int32_t on_surface) const -{ - moab::ErrorCode rval; - moab::EntityHandle vol = dagmc_ptr_->entity_by_index(3, dag_index_); - - int result = 0; - double pnt[3] = {r.x, r.y, r.z}; - double dir[3] = {u.x, u.y, u.z}; - rval = dagmc_ptr_->point_in_volume(vol, pnt, result, dir); - MB_CHK_ERR_CONT(rval); - return result; -} - -void DAGCell::to_hdf5_inner(hid_t group_id) const { - write_string(group_id, "geom_type", "dagmc", false); -} - -BoundingBox DAGCell::bounding_box() const -{ - moab::ErrorCode rval; - moab::EntityHandle vol = dagmc_ptr_->entity_by_index(3, dag_index_); - double min[3], max[3]; - rval = dagmc_ptr_->getobb(vol, min, max); - MB_CHK_ERR_CONT(rval); - return {min[0], max[0], min[1], max[1], min[2], max[2]}; -} - -#endif - //============================================================================== // UniversePartitioner implementation //============================================================================== @@ -1356,21 +1281,6 @@ openmc_extend_cells(int32_t n, int32_t* index_start, int32_t* index_end) return 0; } -#ifdef DAGMC -int32_t next_cell(DAGUniverse* dag_univ, DAGCell* cur_cell, DAGSurface* surf_xed) -{ - moab::EntityHandle surf = - surf_xed->dagmc_ptr_->entity_by_index(2, surf_xed->dag_index_); - moab::EntityHandle vol = - cur_cell->dagmc_ptr_->entity_by_index(3, cur_cell->dag_index_); - - moab::EntityHandle new_vol; - cur_cell->dagmc_ptr_->next_vol(surf, vol, new_vol); - - return cur_cell->dagmc_ptr_->index_by_handle(new_vol) + dag_univ->cell_idx_offset_; -} -#endif - extern "C" int cells_size() { return model::cells.size(); } } // namespace openmc diff --git a/src/dagmc.cpp b/src/dagmc.cpp index 635f05bd88..10d620d7d0 100644 --- a/src/dagmc.cpp +++ b/src/dagmc.cpp @@ -7,6 +7,7 @@ #include "openmc/file_utils.h" #include "openmc/geometry.h" #include "openmc/geometry_aux.h" +#include "openmc/hdf5_interface.h" #include "openmc/material.h" #include "openmc/string_utils.h" #include "openmc/settings.h" @@ -489,5 +490,94 @@ DAGUniverse::read_uwuw_materials() { } } + +//============================================================================== +// DAGMC Cell implementation +//============================================================================== + +DAGCell::DAGCell() : Cell{} { + geom_type_ = GeometryType::DAG; + simple_ = true; +}; + +std::pair +DAGCell::distance(Position r, Direction u, int32_t on_surface, Particle* p) const +{ + Expects(p); + // if we've changed direction or we're not on a surface, + // reset the history and update last direction + if (u != p->last_dir_) { p->last_dir_ = u; p->history_.reset(); } + if (on_surface == 0) { p->history_.reset(); } + + const auto& univ = model::universes[p->coord_[p->n_coord_ - 1].universe]; + + DAGUniverse* dag_univ = static_cast(univ.get()); + if (!dag_univ) fatal_error("DAGMC call made for particle in a non-DAGMC universe"); + + moab::ErrorCode rval; + moab::EntityHandle vol = dagmc_ptr_->entity_by_index(3, dag_index_); + moab::EntityHandle hit_surf; + double dist; + double pnt[3] = {r.x, r.y, r.z}; + double dir[3] = {u.x, u.y, u.z}; + rval = dagmc_ptr_->ray_fire(vol, pnt, dir, hit_surf, dist, &p->history_); + MB_CHK_ERR_CONT(rval); + int surf_idx; + if (hit_surf != 0) { + surf_idx = dag_univ->surf_idx_offset_ + dagmc_ptr_->index_by_handle(hit_surf); + } else { + // indicate that particle is lost + surf_idx = -1; + dist = INFINITY; + if (!dagmc_ptr_->is_implicit_complement(vol) || model::universe_map[dag_univ->id_] == model::root_universe) { + p->mark_as_lost(fmt::format("No intersection found with DAGMC cell {}", id_)); + } + } + + return {dist, surf_idx}; +} + +bool DAGCell::contains(Position r, Direction u, int32_t on_surface) const +{ + moab::ErrorCode rval; + moab::EntityHandle vol = dagmc_ptr_->entity_by_index(3, dag_index_); + + int result = 0; + double pnt[3] = {r.x, r.y, r.z}; + double dir[3] = {u.x, u.y, u.z}; + rval = dagmc_ptr_->point_in_volume(vol, pnt, result, dir); + MB_CHK_ERR_CONT(rval); + return result; +} + +void DAGCell::to_hdf5_inner(hid_t group_id) const { + write_string(group_id, "geom_type", "dagmc", false); +} + +BoundingBox DAGCell::bounding_box() const +{ + moab::ErrorCode rval; + moab::EntityHandle vol = dagmc_ptr_->entity_by_index(3, dag_index_); + double min[3], max[3]; + rval = dagmc_ptr_->getobb(vol, min, max); + MB_CHK_ERR_CONT(rval); + return {min[0], max[0], min[1], max[1], min[2], max[2]}; +} + +// DAGMC Functions +int32_t next_cell(DAGUniverse* dag_univ, DAGCell* cur_cell, DAGSurface* surf_xed) +{ + moab::EntityHandle surf = + surf_xed->dagmc_ptr_->entity_by_index(2, surf_xed->dag_index_); + moab::EntityHandle vol = + cur_cell->dagmc_ptr_->entity_by_index(3, cur_cell->dag_index_); + + moab::EntityHandle new_vol; + cur_cell->dagmc_ptr_->next_vol(surf, vol, new_vol); + + return cur_cell->dagmc_ptr_->index_by_handle(new_vol) + dag_univ->cell_idx_offset_; +} + + } #endif