Using the DAGUniverse constructor to setup the DAGMC instance.

include/openmc/cell.h

@@ -80,11 +80,11 @@ public:
 class DAGUniverse : public Universe {
 
 public:
-  explicit DAGUniverse(std::string filename);
+  explicit DAGUniverse(const std::string& filename);
 
 
   // Data Members
-  std::shared_ptr<moab::DagMC> dag_instance_; //! DAGMC Instance for this universe
+  std::shared_ptr<moab::DagMC> dagmc_instance_; //! DAGMC Instance for this universe
 };
 
 #endif

include/openmc/dagmc.h

@@ -9,8 +9,9 @@ extern "C" const bool DAGMC_ENABLED;
 #ifdef DAGMC
 
 #include "DagMC.hpp"
-#include "openmc/xml_interface.h"
+
 #include "openmc/position.h"
+#include "openmc/xml_interface.h"
 
 namespace openmc {
 

src/dagmc.cpp

@@ -151,7 +151,215 @@ void legacy_assign_material(std::string mat_string, DAGCell* c)
 
 void load_dagmc_geometry()
 {
-  int32_t univ_idx = create_dagmc_universe(dagmc_file());
+  model::universes.push_back(std::make_unique<DAGUniverse>(dagmc_file()));
+  model::universe_map[model::universes.back()->id_] = model::universes.size() - 1;
+}
+
+DAGUniverse::DAGUniverse(const std::string& filename) {
+
+  // determine the next cell id
+  int32_t next_cell_id = 0;
+  for (const auto& c : model::cells) {
+    if (c->id_ > next_cell_id) next_cell_id = c->id_;
+  }
+  next_cell_id++;
+
+  // determine the next surface id
+  int32_t next_surf_id = 0;
+  for (const auto& s : model::surfaces) {
+    if (s->id_ > next_surf_id) next_surf_id = s->id_;
+  }
+  next_surf_id++;
+
+  // determine the next universe id
+  int32_t next_univ_id = 0;
+  for (const auto& u : model::universes) {
+    if (u->id_ > next_univ_id) next_univ_id = u->id_;
+  }
+  next_univ_id++;
+
+  // set the universe id
+  id_ = next_univ_id;
+
+  // create a new DAGMC instance
+  dagmc_instance_ = std::make_shared<moab::DagMC>();
+
+  // --- Materials ---
+
+  // create uwuw instance
+  UWUW uwuw(filename.c_str());
+
+  // check for uwuw material definitions
+  bool using_uwuw = !uwuw.material_library.empty();
+
+  // notify user if UWUW materials are going to be used
+  if (using_uwuw) {
+    write_message("Found UWUW Materials in the DAGMC geometry file.", 6);
+  }
+
+  // load the DAGMC geometry
+  moab::ErrorCode rval = dagmc_instance_->load_file(filename.c_str());
+  MB_CHK_ERR_CONT(rval);
+
+  // initialize acceleration data structures
+  rval = dagmc_instance_->init_OBBTree();
+  MB_CHK_ERR_CONT(rval);
+
+  // parse model metadata
+  dagmcMetaData DMD(dagmc_instance_.get(), false, false);
+  DMD.load_property_data();
+
+  std::vector<std::string> keywords {"temp"};
+  std::map<std::string, std::string> dum;
+  std::string delimiters = ":/";
+  rval = dagmc_instance_->parse_properties(keywords, dum, delimiters.c_str());
+  MB_CHK_ERR_CONT(rval);
+
+  // --- Cells (Volumes) ---
+
+  // initialize cell objects
+  int n_cells = dagmc_instance_->num_entities(3);
+  moab::EntityHandle graveyard = 0;
+  for (int i = 0; i < n_cells; i++) {
+    moab::EntityHandle vol_handle = dagmc_instance_->entity_by_index(3, i + 1);
+
+    // set cell ids using global IDs
+    DAGCell* c = new DAGCell();
+    c->dag_index_ = i + 1;
+    std::cout << "Cell id: " << next_cell_id << std::endl;
+    c->id_ = next_cell_id++;
+    // c->id_ = model::DAG->id_by_index(3, c->dag_index_);
+    c->dagmc_ptr_ = dagmc_instance_;
+    c->universe_ = id_; // set to zero for now
+    c->fill_ = C_NONE; // no fill, single universe
+
+    model::cells.emplace_back(c);
+    model::cell_map[c->id_] = model::cells.size() - 1;
+    cells_.push_back(model::cell_map[c->id_]);
+
+    // Populate the Universe vector and dict
+    // auto it = model::universe_map.find(dagmc_univ_id);
+    // if (it == model::universe_map.end()) {
+    //   model::universes.push_back(std::make_unique<Universe>());
+    //   model::universes.back()->id_ = dagmc_univ_id;
+    //   model::universes.back()->cells_.push_back(i);
+    //   model::universe_map[id_] = model::universes.size() - 1;
+    // } else {
+    //   model::universes[it->second]->cells_.push_back(i);
+    // }
+
+    // MATERIALS
+
+    // determine volume material assignment
+    std::string mat_str = DMD.get_volume_property("material", vol_handle);
+
+    if (mat_str.empty()) {
+      fatal_error(fmt::format("Volume {} has no material assignment.", c->id_));
+    }
+
+    to_lower(mat_str);
+
+    if (mat_str == "graveyard") {
+      graveyard = vol_handle;
+    }
+
+    // material void checks
+    if (mat_str == "void" || mat_str == "vacuum" || mat_str == "graveyard") {
+      c->material_.push_back(MATERIAL_VOID);
+    } else {
+      if (using_uwuw) {
+        // lookup material in uwuw if present
+        std::string uwuw_mat = DMD.volume_material_property_data_eh[vol_handle];
+        if (uwuw.material_library.count(uwuw_mat) != 0) {
+          // Note: material numbers are set by UWUW
+          int mat_number = uwuw.material_library[uwuw_mat].metadata["mat_number"].asInt();
+          c->material_.push_back(mat_number);
+        } else {
+          fatal_error(fmt::format("Material with value {} not found in the "
+            "UWUW material library", mat_str));
+        }
+      } else {
+        legacy_assign_material(mat_str, c);
+      }
+    }
+
+    // check for temperature assignment
+    std::string temp_value;
+
+    // no temperature if void
+    if (c->material_[0] == MATERIAL_VOID) continue;
+
+    // assign cell temperature
+    const auto& mat = model::materials[model::material_map.at(c->material_[0])];
+    if (dagmc_instance_->has_prop(vol_handle, "temp")) {
+      rval = dagmc_instance_->prop_value(vol_handle, "temp", temp_value);
+      MB_CHK_ERR_CONT(rval);
+      double temp = std::stod(temp_value);
+      c->sqrtkT_.push_back(std::sqrt(K_BOLTZMANN * temp));
+    } else if (mat->temperature_ > 0.0) {
+      c->sqrtkT_.push_back(std::sqrt(K_BOLTZMANN * mat->temperature_));
+    } else {
+      c->sqrtkT_.push_back(std::sqrt(K_BOLTZMANN * settings::temperature_default));
+    }
+
+  }
+
+  // allocate the cell overlap count if necessary
+  if (settings::check_overlaps) {
+    model::overlap_check_count.resize(model::cells.size(), 0);
+  }
+
+  if (!graveyard) {
+    warning("No graveyard volume found in the DagMC model."
+            "This may result in lost particles and rapid simulation failure.");
+  }
+
+  // --- Surfaces ---
+
+  // initialize surface objects
+  int n_surfaces = dagmc_instance_->num_entities(2);
+
+  for (int i = 0; i < n_surfaces; i++) {
+    moab::EntityHandle surf_handle = dagmc_instance_->entity_by_index(2, i+1);
+
+    // set cell ids using global IDs
+    DAGSurface* s = new DAGSurface();
+    s->dag_index_ = i+1;
+    s->id_ = next_surf_id++;
+    s->dagmc_ptr_ = dagmc_instance_;
+
+    // set BCs
+    std::string bc_value = DMD.get_surface_property("boundary", surf_handle);
+    to_lower(bc_value);
+    if (bc_value.empty() || bc_value == "transmit" || bc_value == "transmission") {
+      // set to transmission by default
+      s->bc_ = Surface::BoundaryType::TRANSMIT;
+    } else if (bc_value == "vacuum") {
+      s->bc_ = Surface::BoundaryType::VACUUM;
+    } else if (bc_value == "reflective" || bc_value == "reflect" || bc_value == "reflecting") {
+      s->bc_ = Surface::BoundaryType::REFLECT;
+    } else if (bc_value == "periodic") {
+      fatal_error("Periodic boundary condition not supported in DAGMC.");
+    } else {
+      fatal_error(fmt::format("Unknown boundary condition \"{}\" specified "
+        "on surface {}", bc_value, s->id_));
+    }
+
+    // graveyard check
+    moab::Range parent_vols;
+    rval = dagmc_instance_->moab_instance()->get_parent_meshsets(surf_handle, parent_vols);
+    MB_CHK_ERR_CONT(rval);
+
+    // if this surface belongs to the graveyard
+    if (graveyard && parent_vols.find(graveyard) != parent_vols.end()) {
+      // set graveyard surface BC's to vacuum
+      s->bc_ = Surface::BoundaryType::VACUUM;
+    }
+
+    // add to global array and map
+    model::surfaces.emplace_back(s);
+    model::surface_map[s->id_] = model::surfaces.size() - 1;
+  }
 }
 
 int32_t create_dagmc_universe(const std::string& filename) {
@@ -204,6 +412,7 @@ int32_t create_dagmc_universe(const std::string& filename) {
     DAGCell* c = new DAGCell();
     c->dag_index_ = i+1;
     c->id_ = model::DAG->id_by_index(3, c->dag_index_);
+    std::cout << "Cell id: " << c->id_ << std::endl;
     c->dagmc_ptr_ = model::DAG;
     c->universe_ = dagmc_univ_id; // set to zero for now
     c->fill_ = C_NONE; // no fill, single universe
@@ -346,7 +555,7 @@ void read_geometry_dagmc()
   write_message("Reading DAGMC geometry...", 5);
   load_dagmc_geometry();
 
-  model::root_universe = find_root_universe();
+  model::root_universe = 0;
 }
 
 }

tests/regression_tests/dagmc/refl/inputs_true.dat

@@ -14,7 +14,7 @@
 <?xml version='1.0' encoding='utf-8'?>
 <tallies>
   <filter id="1" type="cell">
-    <bins>1</bins>
+    <bins>2</bins>
   </filter>
   <tally id="1">
     <filters>1</filters>

tests/regression_tests/dagmc/refl/test.py

@@ -30,7 +30,7 @@ class UWUWTest(PyAPITestHarness):
         # tally
         tally = openmc.Tally()
         tally.scores = ['total']
-        tally.filters = [openmc.CellFilter(1)]
+        tally.filters = [openmc.CellFilter(2)]
         model.tallies = [tally]
 
         model.tallies.export_to_xml()