diff --git a/include/openmc/cell.h b/include/openmc/cell.h index d975750a16..01f0375865 100644 --- a/include/openmc/cell.h +++ b/include/openmc/cell.h @@ -300,5 +300,8 @@ struct CellInstanceHash { void read_cells(pugi::xml_node node); +//! Add cells to universes +void populate_universes(); + } // namespace openmc #endif // OPENMC_CELL_H diff --git a/include/openmc/dagmc.h b/include/openmc/dagmc.h index 64bf4d2091..9930763474 100644 --- a/include/openmc/dagmc.h +++ b/include/openmc/dagmc.h @@ -22,6 +22,7 @@ void check_dagmc_root_univ(); #ifdef DAGMC #include "DagMC.hpp" +#include "dagmcmetadata.hpp" #include "openmc/cell.h" #include "openmc/particle.h" @@ -87,6 +88,11 @@ public: explicit DAGUniverse(const std::string& filename, bool auto_geom_ids = false, bool auto_mat_ids = false); + //! Alternative DAGMC universe constructor for external DAGMC instance + explicit DAGUniverse(std::shared_ptr external_dagmc_ptr, + const std::string& filename = "", bool auto_geom_ids = false, + bool auto_mat_ids = false); + //! Initialize the DAGMC accel. data structures, indices, material //! assignments, etc. void initialize(); @@ -139,10 +145,16 @@ public: bool has_graveyard() const { return has_graveyard_; } private: + void set_id(); //!< Deduce the universe id from model::universes + void init_dagmc(); //!< Create and initialise DAGMC pointer + void init_metadata(); //!< Create and initialise dagmcMetaData pointer + void init_geometry(); //!< Create cells and surfaces from DAGMC entities + std::string filename_; //!< Name of the DAGMC file used to create this universe std::shared_ptr - uwuw_; //!< Pointer to the UWUW instance for this universe + uwuw_; //!< Pointer to the UWUW instance for this universe + std::unique_ptr dmd_ptr; //! Pointer to DAGMC metadata object bool adjust_geometry_ids_; //!< Indicates whether or not to automatically //!< generate new cell and surface IDs for the //!< universe diff --git a/include/openmc/material.h b/include/openmc/material.h index 709d205738..b251a3ca85 100644 --- a/include/openmc/material.h +++ b/include/openmc/material.h @@ -113,6 +113,9 @@ public: //! \param[in] units Units of density void set_density(double density, gsl::cstring_span units); + //! Set temperature of the material + void set_temperature(double temperature) { temperature_ = temperature; }; + //! Get nuclides in material //! \return Indices into the global nuclides vector gsl::span nuclides() const diff --git a/src/cell.cpp b/src/cell.cpp index 28be2c2cff..f9786ce3eb 100644 --- a/src/cell.cpp +++ b/src/cell.cpp @@ -846,6 +846,20 @@ void read_cells(pugi::xml_node node) read_dagmc_universes(node); + populate_universes(); + + // Allocate the cell overlap count if necessary. + if (settings::check_overlaps) { + model::overlap_check_count.resize(model::cells.size(), 0); + } + + if (model::cells.size() == 0) { + fatal_error("No cells were found in the geometry.xml file"); + } +} + +void populate_universes() +{ // Populate the Universe vector and map. for (int i = 0; i < model::cells.size(); i++) { int32_t uid = model::cells[i]->universe_; @@ -860,15 +874,6 @@ void read_cells(pugi::xml_node node) } } model::universes.shrink_to_fit(); - - // Allocate the cell overlap count if necessary. - if (settings::check_overlaps) { - model::overlap_check_count.resize(model::cells.size(), 0); - } - - if (model::cells.size() == 0) { - fatal_error("No cells were found in the geometry.xml file"); - } } //============================================================================== diff --git a/src/dagmc.cpp b/src/dagmc.cpp index 8131c99297..e780d2e77f 100644 --- a/src/dagmc.cpp +++ b/src/dagmc.cpp @@ -12,7 +12,6 @@ #include "openmc/string_utils.h" #ifdef DAGMC -#include "dagmcmetadata.hpp" #include "uwuw.hpp" #endif #include @@ -74,6 +73,23 @@ DAGUniverse::DAGUniverse( const std::string& filename, bool auto_geom_ids, bool auto_mat_ids) : filename_(filename), adjust_geometry_ids_(auto_geom_ids), adjust_material_ids_(auto_mat_ids) +{ + set_id(); + initialize(); +} + +DAGUniverse::DAGUniverse(std::shared_ptr dagmc_ptr, + const std::string& filename, bool auto_geom_ids, bool auto_mat_ids) + : dagmc_instance_(dagmc_ptr), filename_(filename), + adjust_geometry_ids_(auto_geom_ids), adjust_material_ids_(auto_mat_ids) +{ + set_id(); + init_metadata(); + read_uwuw_materials(); + init_geometry(); +} + +void DAGUniverse::set_id() { // determine the next universe id int32_t next_univ_id = 0; @@ -85,48 +101,27 @@ DAGUniverse::DAGUniverse( // set the universe id id_ = next_univ_id; - - initialize(); } void DAGUniverse::initialize() { geom_type() = GeometryType::DAG; - // 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_; - } - cell_idx_offset_ = model::cells.size(); - next_cell_id++; + init_dagmc(); - // 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_; - } - surf_idx_offset_ = model::surfaces.size(); - next_surf_id++; + init_metadata(); + + read_uwuw_materials(); + + init_geometry(); +} + +void DAGUniverse::init_dagmc() +{ // create a new DAGMC instance dagmc_instance_ = std::make_shared(); - // --- Materials --- - - // read any UWUW materials from the file - read_uwuw_materials(); - - // check for uwuw material definitions - bool using_uwuw = uses_uwuw(); - - // 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 filename_ = settings::path_input + filename_; if (!file_exists(filename_)) { @@ -138,18 +133,35 @@ void DAGUniverse::initialize() // initialize acceleration data structures rval = dagmc_instance_->init_OBBTree(); MB_CHK_ERR_CONT(rval); +} +void DAGUniverse::init_metadata() +{ // parse model metadata - dagmcMetaData DMD(dagmc_instance_.get(), false, false); - DMD.load_property_data(); + dmd_ptr = + std::make_unique(dagmc_instance_.get(), false, false); + dmd_ptr->load_property_data(); std::vector keywords {"temp"}; std::map dum; std::string delimiters = ":/"; + moab::ErrorCode rval; rval = dagmc_instance_->parse_properties(keywords, dum, delimiters.c_str()); MB_CHK_ERR_CONT(rval); +} - // --- Cells (Volumes) --- +void DAGUniverse::init_geometry() +{ + moab::ErrorCode rval; + + // 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_; + } + cell_idx_offset_ = model::cells.size(); + next_cell_id++; // initialize cell objects int n_cells = dagmc_instance_->num_entities(3); @@ -178,7 +190,7 @@ void DAGUniverse::initialize() // --- Materials --- // determine volume material assignment - std::string mat_str = DMD.get_volume_property("material", vol_handle); + std::string mat_str = dmd_ptr->get_volume_property("material", vol_handle); if (mat_str.empty()) { fatal_error(fmt::format("Volume {} has no material assignment.", c->id_)); @@ -194,9 +206,10 @@ void DAGUniverse::initialize() if (mat_str == "void" || mat_str == "vacuum" || mat_str == "graveyard") { c->material_.push_back(MATERIAL_VOID); } else { - if (using_uwuw) { + if (uses_uwuw()) { // lookup material in uwuw if present - std::string uwuw_mat = DMD.volume_material_property_data_eh[vol_handle]; + std::string uwuw_mat = + dmd_ptr->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.get_material(uwuw_mat) @@ -246,7 +259,14 @@ void DAGUniverse::initialize() has_graveyard_ = graveyard; - // --- Surfaces --- + // 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_; + } + surf_idx_offset_ = model::surfaces.size(); + next_surf_id++; // initialize surface objects int n_surfaces = dagmc_instance_->num_entities(2); @@ -259,7 +279,8 @@ void DAGUniverse::initialize() : dagmc_instance_->id_by_index(2, i + 1); // set BCs - std::string bc_value = DMD.get_surface_property("boundary", surf_handle); + std::string bc_value = + dmd_ptr->get_surface_property("boundary", surf_handle); to_lower(bc_value); if (bc_value.empty() || bc_value == "transmit" || bc_value == "transmission") { @@ -398,7 +419,7 @@ void DAGUniverse::to_hdf5(hid_t universes_group) const bool DAGUniverse::uses_uwuw() const { - return !uwuw_->material_library.empty(); + return uwuw_ && !uwuw_->material_library.empty(); } std::string DAGUniverse::get_uwuw_materials_xml() const @@ -487,33 +508,32 @@ void DAGUniverse::legacy_assign_material( void DAGUniverse::read_uwuw_materials() { - - int32_t next_material_id = 0; - for (const auto& m : model::materials) { - next_material_id = std::max(m->id_, next_material_id); - } - next_material_id++; + // If no filename was provided, don't read UWUW materials + if (filename_ == "") + return; uwuw_ = std::make_shared(filename_.c_str()); - const auto& mat_lib = uwuw_->material_library; - if (mat_lib.size() == 0) + + if (!uses_uwuw()) return; + // Notify user if UWUW materials are going to be used + write_message("Found UWUW Materials in the DAGMC geometry file.", 6); + // if we're using automatic IDs, update the UWUW material metadata if (adjust_material_ids_) { + int32_t next_material_id = 0; + for (const auto& m : model::materials) { + next_material_id = std::max(m->id_, next_material_id); + } + next_material_id++; + for (auto& mat : uwuw_->material_library) { mat.second->metadata["mat_number"] = next_material_id++; } } - std::stringstream ss; - ss << "\n"; - ss << "\n"; - for (auto mat : mat_lib) { - ss << mat.second->openmc("atom"); - } - ss << ""; - std::string mat_xml_string = ss.str(); + std::string mat_xml_string = get_uwuw_materials_xml(); // create a pugi XML document from this string pugi::xml_document doc; diff --git a/tests/regression_tests/dagmc/external/__init__.py b/tests/regression_tests/dagmc/external/__init__.py new file mode 100644 index 0000000000..e69de29bb2 diff --git a/tests/regression_tests/dagmc/external/dagmc.h5m b/tests/regression_tests/dagmc/external/dagmc.h5m new file mode 120000 index 0000000000..92c41719c5 --- /dev/null +++ b/tests/regression_tests/dagmc/external/dagmc.h5m @@ -0,0 +1 @@ +../legacy/dagmc.h5m \ No newline at end of file diff --git a/tests/regression_tests/dagmc/external/inputs_true.dat b/tests/regression_tests/dagmc/external/inputs_true.dat new file mode 100644 index 0000000000..dd75bbdbd9 --- /dev/null +++ b/tests/regression_tests/dagmc/external/inputs_true.dat @@ -0,0 +1,40 @@ + + + + + + + + + + + + + + + + + + + + eigenvalue + 100 + 5 + 0 + + + -4 -4 -4 4 4 4 + + + 293 + + + + + 1 + + + 1 + total + + diff --git a/tests/regression_tests/dagmc/external/main.cpp b/tests/regression_tests/dagmc/external/main.cpp new file mode 100644 index 0000000000..42f7d97c90 --- /dev/null +++ b/tests/regression_tests/dagmc/external/main.cpp @@ -0,0 +1,94 @@ +#include "openmc/capi.h" +#include "openmc/cross_sections.h" +#include "openmc/dagmc.h" +#include "openmc/error.h" +#include "openmc/geometry.h" +#include "openmc/geometry_aux.h" +#include "openmc/material.h" +#include "openmc/nuclide.h" +#include + +int main(int argc, char* argv[]) +{ + using namespace openmc; + int openmc_err; + + // Initialise OpenMC + openmc_err = openmc_init(argc, argv, nullptr); + if (openmc_err == -1) { + // This happens for the -h and -v flags + return EXIT_SUCCESS; + } else if (openmc_err) { + fatal_error(openmc_err_msg); + } + + // Create DAGMC ptr + std::string filename = "dagmc.h5m"; + std::shared_ptr dag_ptr = std::make_shared(); + moab::ErrorCode rval = dag_ptr->load_file(filename.c_str()); + if (rval != moab::MB_SUCCESS) { + fatal_error("Failed to load file"); + } + + // Initialize acceleration data structures + rval = dag_ptr->init_OBBTree(); + if (rval != moab::MB_SUCCESS) { + fatal_error("Failed to initialise OBB tree"); + } + + // Get rid of existing geometry + std::unordered_map nuclide_map_copy = + openmc::data::nuclide_map; + openmc::data::nuclides.clear(); + openmc::data::nuclide_map = nuclide_map_copy; + openmc::model::surfaces.clear(); + openmc::model::surface_map.clear(); + openmc::model::cells.clear(); + openmc::model::cell_map.clear(); + openmc::model::universes.clear(); + openmc::model::universe_map.clear(); + + // Update materials (emulate an external program) + for (auto& mat_ptr : openmc::model::materials) { + mat_ptr->set_temperature(300); + } + + // Create new DAGMC universe + openmc::model::universes.push_back( + std::make_unique(dag_ptr)); + model::universe_map[model::universes.back()->id_] = + model::universes.size() - 1; + + // Add cells to universes + openmc::populate_universes(); + + // Set root universe + openmc::model::root_universe = openmc::find_root_universe(); + openmc::check_dagmc_root_univ(); + + // Final geometry setup and assign temperatures + openmc::finalize_geometry(); + + // Finalize cross sections having assigned temperatures + openmc::finalize_cross_sections(); + + // Check that we correctly assigned cell temperatures with non-void fill + for (auto& cell_ptr : openmc::model::cells) { + if (cell_ptr->material_.front() != openmc::C_NONE && + cell_ptr->temperature() != 300) { + fatal_error("Failed to set cell temperature"); + } + } + + // Run OpenMC + openmc_err = openmc_run(); + if (openmc_err) + fatal_error(openmc_err_msg); + + // Deallocate memory + openmc_err = openmc_finalize(); + if (openmc_err) + fatal_error(openmc_err_msg); + + return EXIT_SUCCESS; +} diff --git a/tests/regression_tests/dagmc/external/results_true.dat b/tests/regression_tests/dagmc/external/results_true.dat new file mode 100644 index 0000000000..e5127c8272 --- /dev/null +++ b/tests/regression_tests/dagmc/external/results_true.dat @@ -0,0 +1,5 @@ +k-combined: +8.426936E-01 5.715847E-02 +tally 1: +8.093843E+00 +1.328829E+01 diff --git a/tests/regression_tests/dagmc/external/test.py b/tests/regression_tests/dagmc/external/test.py new file mode 100644 index 0000000000..93123ae1ed --- /dev/null +++ b/tests/regression_tests/dagmc/external/test.py @@ -0,0 +1,128 @@ +from pathlib import Path +import os +import shutil +import subprocess +import textwrap + +import openmc +import openmc.lib +import pytest + +from tests.regression_tests import config +from tests.testing_harness import PyAPITestHarness + +pytestmark = pytest.mark.skipif( + not openmc.lib._dagmc_enabled(), + reason="DAGMC is not enabled.") + +# Test that an external DAGMC instance can be passed in through the C API + +@pytest.fixture +def cpp_driver(request): + """Compile the external source""" + + # Get build directory and write CMakeLists.txt file + openmc_dir = Path(str(request.config.rootdir)) / 'build' + with open('CMakeLists.txt', 'w') as f: + f.write(textwrap.dedent(""" + cmake_minimum_required(VERSION 3.3 FATAL_ERROR) + project(openmc_cpp_driver CXX) + add_executable(main main.cpp) + find_package(OpenMC REQUIRED HINTS {}) + target_link_libraries(main OpenMC::libopenmc) + target_compile_features(main PUBLIC cxx_std_14) + set(CMAKE_CXX_FLAGS "-pedantic-errors") + add_compile_definitions(DAGMC=1) + """.format(openmc_dir))) + + # Create temporary build directory and change to there + local_builddir = Path('build') + local_builddir.mkdir(exist_ok=True) + os.chdir(local_builddir) + + mpi_arg = "On" if config['mpi'] else "Off" + + try: + # Run cmake/make to build the shared libary + subprocess.run(['cmake', os.path.pardir, f'-DOPENMC_USE_MPI={mpi_arg}'], check=True) + subprocess.run(['make'], check=True) + os.chdir(os.path.pardir) + + yield "./build/main" + + finally: + # Remove local build directory when test is complete + shutil.rmtree('build') + os.remove('CMakeLists.txt') + +@pytest.fixture +def model(): + model = openmc.model.Model() + + # Settings + model.settings.batches = 5 + model.settings.inactive = 0 + model.settings.particles = 100 + source_box = openmc.stats.Box([-4, -4, -4], + [ 4, 4, 4]) + source = openmc.Source(space=source_box) + model.settings.source = source + model.settings.temperature['default'] = 293 + + # Geometry + dag_univ = openmc.DAGMCUniverse("dagmc.h5m") + model.geometry = openmc.Geometry(dag_univ) + + # Tallies + tally = openmc.Tally() + tally.scores = ['total'] + tally.filters = [openmc.CellFilter(1)] + model.tallies = [tally] + + # Materials + u235 = openmc.Material(name="no-void fuel") + u235.add_nuclide('U235', 1.0, 'ao') + u235.set_density('g/cc', 11) + u235.id = 40 + water = openmc.Material(name="water") + water.add_nuclide('H1', 2.0, 'ao') + water.add_nuclide('O16', 1.0, 'ao') + water.set_density('g/cc', 1.0) + water.add_s_alpha_beta('c_H_in_H2O') + water.id = 41 + mats = openmc.Materials([u235, water]) + model.materials = mats + + return model + +class ExternalDAGMCTest(PyAPITestHarness): + def __init__(self, executable, statepoint_name, model): + super().__init__(statepoint_name, model) + self.executable = executable + + def _run_openmc(self): + """ + Just test if results generated with the external C++ API are + self-consistent with the internal python API. + We generate the "truth" results with the python API but + the main test produces the results file by running the + executable compiled from main.cpp. This future-proofs + the test - we only care that the two routes are equivalent. + """ + if config['update']: + # Generate the results file with internal python API + openmc.run(openmc_exec=config['exe'], event_based=config['event']) + elif config['mpi']: + mpi_args = [config['mpiexec'], '-n', config['mpi_np']] + # Run main cpp executable with MPI + openmc.run(openmc_exec=self.executable, + mpi_args=mpi_args, + event_based=config['event']) + else: + # Run main cpp executable + openmc.run(openmc_exec=self.executable, + event_based=config['event']) + +def test_external_dagmc(cpp_driver, model): + harness = ExternalDAGMCTest(cpp_driver, 'statepoint.5.h5', model) + harness.main()