From 146d7002224cfdf676168ca55f1caee89299190a Mon Sep 17 00:00:00 2001 From: shriwise Date: Tue, 4 Sep 2018 21:49:58 -0500 Subject: [PATCH] Some updates to get the CAD build working again. --- include/openmc/settings.h | 2 +- include/openmc/surface.h | 2 ++ src/cad.cpp | 42 ++++++++++++++-------------- src/cad.h | 4 +-- src/cell.cpp | 8 +++--- src/input_xml.F90 | 59 ++++++++++++++++----------------------- src/settings.F90 | 2 +- src/settings.cpp | 14 +++++++++- src/summary.cpp | 4 +-- src/surface.cpp | 8 ++++-- 10 files changed, 75 insertions(+), 70 deletions(-) diff --git a/include/openmc/settings.h b/include/openmc/settings.h index e5684c7ae..6565d3e46 100644 --- a/include/openmc/settings.h +++ b/include/openmc/settings.h @@ -45,6 +45,7 @@ extern "C" bool ufs_on; //!< uniform fission site method on? extern "C" bool urr_ptables_on; //!< use unresolved resonance prob. tables? extern "C" bool write_all_tracks; //!< write track files for every particle? extern "C" bool write_initial_source; //!< write out initial source file? +extern "C" bool dagmc; //!< indicator of DAGMC geometry // Paths to various files extern std::string path_cross_sections; //!< path to cross_sections.xml @@ -86,7 +87,6 @@ extern "C" int trigger_batch_interval; //!< Batch interval for triggers extern "C" int verbosity; //!< How verbose to make output extern "C" double weight_cutoff; //!< Weight cutoff for Russian roulette extern "C" double weight_survive; //!< Survival weight after Russian roulette - } // namespace settings //! Read settings from XML file diff --git a/include/openmc/surface.h b/include/openmc/surface.h index d1e3f2b49..f5b30e413 100644 --- a/include/openmc/surface.h +++ b/include/openmc/surface.h @@ -147,6 +147,8 @@ class CADSurface : public Surface Direction normal(Position p) const; //! Get the bounding box of this surface. BoundingBox bounding_box() const; + + void to_hdf5(hid_t group_id) const; }; #endif //============================================================================== diff --git a/src/cad.cpp b/src/cad.cpp index a3d1c7648..ef937a346 100644 --- a/src/cad.cpp +++ b/src/cad.cpp @@ -2,7 +2,7 @@ #ifdef CAD #include "cad.h" -#include "error.h" +#include "openmc/error.h" #include #include @@ -59,25 +59,25 @@ void load_cad_geometry_c() c->fill_ = openmc::C_NONE; - openmc::global_cells.push_back(c); + openmc::cells.push_back(c); openmc::cell_map[c->id_] = c->id_; // Populate the Universe vector and dict auto it = openmc::universe_map.find(cad_univ_id); if (it == openmc::universe_map.end()) { - openmc::global_universes.push_back(new openmc::Universe()); - openmc::global_universes.back()-> id = cad_univ_id; - openmc::global_universes.back()->cells.push_back(i); - openmc::universe_map[cad_univ_id] = openmc::global_universes.size() - 1; + openmc::universes.push_back(new openmc::Universe()); + openmc::universes.back()-> id_ = cad_univ_id; + openmc::universes.back()->cells_.push_back(i); + openmc::universe_map[cad_univ_id] = openmc::universes.size() - 1; } else { - openmc::global_universes[it->second]->cells.push_back(i); + openmc::universes[it->second]->cells_.push_back(i); } if(DAGMC->is_implicit_complement(vol_handle)) { // assuming implicit complement is void for now - c->material.push_back(openmc::MATERIAL_VOID); + c->material_.push_back(openmc::MATERIAL_VOID); continue; } @@ -88,21 +88,21 @@ void load_cad_geometry_c() TOLOWER(mat_value); if(mat_value == "void") { - c->material.push_back(openmc::MATERIAL_VOID); + c->material_.push_back(openmc::MATERIAL_VOID); } else { - c->material.push_back(std::stoi(mat_value)); + c->material_.push_back(std::stoi(mat_value)); } } else { std::cout << "Warning: volume without material found!" << std::endl; - c->material.push_back(openmc::MATERIAL_VOID); + c->material_.push_back(openmc::MATERIAL_VOID); } } // initialize surface objects openmc::n_surfaces = DAGMC->num_entities(2); - openmc::global_surfaces.resize(openmc::n_surfaces); + openmc::surfaces.resize(openmc::n_surfaces); for(int i = 0; i < openmc::n_surfaces; i++) { @@ -110,7 +110,7 @@ void load_cad_geometry_c() // set cell ids using global IDs openmc::CADSurface* s = new openmc::CADSurface(); - s->id = DAGMC->id_by_index(2, i+1); + s->id_ = DAGMC->id_by_index(2, i+1); s->dagmc_ptr = DAGMC; if(DAGMC->has_prop(surf_handle, "boundary")) { @@ -120,32 +120,32 @@ void load_cad_geometry_c() TOLOWER(bc_value); if(bc_value == "transmit" || bc_value == "transmission") { - s->bc = openmc::BC_TRANSMIT; + s->bc_ = openmc::BC_TRANSMIT; } else if(bc_value == "vacuum") { - s->bc = openmc::BC_VACUUM; + s->bc_ = openmc::BC_VACUUM; } else if(bc_value == "reflective" || bc_value == "reflect" || bc_value == "reflecting") { - s->bc = openmc::BC_REFLECT; + s->bc_ = openmc::BC_REFLECT; } else if(bc_value == "periodic") { openmc::fatal_error("Periodic boundary condition not supported in CAD."); } else { std::stringstream err_msg; - err_msg << "Unknown boundary condition \"" << s->bc - << "\" specified on surface " << s->id; + err_msg << "Unknown boundary condition \"" << s->bc_ + << "\" specified on surface " << s->id_; openmc::fatal_error(err_msg); } } // if no BC property is found, set to transmit else { - s->bc = openmc::BC_TRANSMIT; + s->bc_ = openmc::BC_TRANSMIT; } // add to global array and map - openmc::global_surfaces[i] = s; - openmc::surface_map[s->id] = s->id; + openmc::surfaces[i] = s; + openmc::surface_map[s->id_] = s->id_; } return; diff --git a/src/cad.h b/src/cad.h index 8f14d2530..1a1663cdb 100644 --- a/src/cad.h +++ b/src/cad.h @@ -5,8 +5,8 @@ #define CAD_H #include "DagMC.hpp" -#include "cell.h" -#include "surface.h" +#include "openmc/cell.h" +#include "openmc/surface.h" extern moab::DagMC* DAGMC; diff --git a/src/cell.cpp b/src/cell.cpp index 01f4eb6d4..047798a24 100644 --- a/src/cell.cpp +++ b/src/cell.cpp @@ -599,7 +599,7 @@ CADCell::CADCell() : Cell{} {}; std::pair CADCell::distance(Position p, Direction u, int32_t on_surface) const { - moab::EntityHandle vol = dagmc_ptr->entity_by_id(3, id); + moab::EntityHandle vol = dagmc_ptr->entity_by_id(3, id_); moab::EntityHandle hit_surf; double dist; dagmc_ptr->ray_fire(vol, p.xyz, u.xyz, hit_surf, dist); @@ -620,7 +620,7 @@ std::pair CADCell::distance(Position p, Direction u, int32_t on bool CADCell::contains(Position p, Direction u, int32_t on_surface) const { - moab::EntityHandle vol = dagmc_ptr->entity_by_id(3, id); + moab::EntityHandle vol = dagmc_ptr->entity_by_id(3, id_); int result = 0; dagmc_ptr->point_in_volume(vol, p.xyz, result, u.xyz); @@ -772,8 +772,8 @@ extern "C" { #ifdef CAD int32_t next_cell(CADCell* cur_cell, CADSurface *surf_xed ) { - moab::EntityHandle surf = surf_xed->dagmc_ptr->entity_by_id(2,surf_xed->id); - moab::EntityHandle vol = cur_cell->dagmc_ptr->entity_by_id(3,cur_cell->id); + moab::EntityHandle surf = surf_xed->dagmc_ptr->entity_by_id(2,surf_xed->id_); + moab::EntityHandle vol = cur_cell->dagmc_ptr->entity_by_id(3,cur_cell->id_); moab::EntityHandle new_vol; cur_cell->dagmc_ptr->next_vol(surf, vol, new_vol); diff --git a/src/input_xml.F90 b/src/input_xml.F90 index 3875bc80a..d3531cbe5 100644 --- a/src/input_xml.F90 +++ b/src/input_xml.F90 @@ -219,17 +219,6 @@ contains ! Get proper XMLNode type given pointer root % ptr = root_ptr - ! Check for use of CAD geometry - if (check_for_node(root, "dagmc")) then -#ifdef CAD - call get_node_value_bool(root, "dagmc", dagmc) -#else - if (dagmc) then - call fatal_error("CAD mode unsupported for this build of OpenMC") - end if -#endif - end if - if (run_mode == MODE_EIGENVALUE) then ! Preallocate space for keff and entropy by generation call k_generation % reserve(n_max_batches*gen_per_batch) @@ -597,10 +586,10 @@ contains do i = 1, n_cells c => cells(i) ! additional metadata spoofing - allocate(c % sqrtKT(1)) - c % sqrtkT(1) = -1.0 ! rely on temps from elsewhere +! allocate(c % sqrtKT(1)) +! c % sqrtkT(1) = -1.0 ! rely on temps from elsewhere univ_id = c % universe() - if (.not. allocated(c%region)) allocate(c%region(0)) +! if (.not. allocated(c%region)) allocate(c%region(0)) if (.not. cells_in_univ_dict % has(univ_id)) then n_universes = n_universes + 1 @@ -615,30 +604,30 @@ contains ! ========================================================================== ! SETUP UNIVERSES - + ! Allocate universes, universe cell arrays, and assign base universe - allocate(universes(n_universes)) - do i = 1, n_universes - associate (u => universes(i)) - u % id = univ_ids % data(i) - ! Allocate cell list - n_cells_in_univ = cells_in_univ_dict % get(u % id) - allocate(u % cells(n_cells_in_univ)) - u % cells(:) = 0 - end associate - end do + ! allocate(universes(n_universes)) + ! do i = 1, n_universes + ! associate (u => universes(i)) + ! u % id = univ_ids % data(i) + ! ! Allocate cell list + ! n_cells_in_univ = cells_in_univ_dict % get(u % id) + ! allocate(u % cells(n_cells_in_univ)) + ! u % cells(:) = 0 + ! end associate + ! end do - root_universe = 0 + 1 + root_universe = find_root_universe() - do i = 1, n_cells - ! Get index in universes array - j = universe_dict % get(cells(i) % universe()) - ! Set the first zero entry in the universe cells array to the index in the - ! global cells array - associate (u => universes(j)) - u % cells(find(u % cells, 0)) = i - end associate - end do + ! do i = 1, n_cells + ! ! Get index in universes array + ! j = universe_dict % get(cells(i) % universe()) + ! ! Set the first zero entry in the universe cells array to the index in the + ! ! global cells array + ! associate (u => universes(j)) + ! u % cells(find(u % cells, 0)) = i + ! end associate + ! end do ! Clear dictionary call cells_in_univ_dict%clear() diff --git a/src/settings.F90 b/src/settings.F90 index 66af618f7..6625ccef0 100644 --- a/src/settings.F90 +++ b/src/settings.F90 @@ -81,7 +81,7 @@ module settings integer(C_INT), bind(C) :: run_mode ! flag for use of CAD geometry - logical, bind(C) :: dagmc = .false. + logical(C_BOOL), bind(C) :: dagmc ! Restart run logical(C_BOOL), bind(C) :: restart_run diff --git a/src/settings.cpp b/src/settings.cpp index 007c515b0..60a805a01 100644 --- a/src/settings.cpp +++ b/src/settings.cpp @@ -56,7 +56,8 @@ bool ufs_on {false}; bool urr_ptables_on {true}; bool write_all_tracks {false}; bool write_initial_source {false}; - +bool dagmc {false}; + std::string path_cross_sections; std::string path_input; std::string path_multipole; @@ -207,6 +208,17 @@ void read_settings_xml() verbosity = std::stoi(get_node_value(root, "verbosity")); } + // DAGMC geometry check + if (check_for_node(root, "dagmc")) { + dagmc = get_node_value_bool(root, "dagmc"); + } + +#ifndef CAD + if (dagmc) { + fatal_error("CAD mode unsupported for this build of OpenMC"); + } +#endif + // To this point, we haven't displayed any output since we didn't know what // the verbosity is. Now that we checked for it, show the title if necessary if (openmc_master) { diff --git a/src/summary.cpp b/src/summary.cpp index c371eb56c..67a9a7f07 100644 --- a/src/summary.cpp +++ b/src/summary.cpp @@ -2,7 +2,7 @@ #include "openmc/hdf5_interface.h" #include "openmc/lattice.h" #include "openmc/surface.h" - +#include "openmc/settings.h" namespace openmc { @@ -10,7 +10,7 @@ extern "C" void write_geometry(hid_t file_id) { #ifdef CAD - if (dagmc) { + if (settings::dagmc) { return; } #endif diff --git a/src/surface.cpp b/src/surface.cpp index aee21287e..b95339386 100644 --- a/src/surface.cpp +++ b/src/surface.cpp @@ -257,7 +257,7 @@ double CADSurface::evaluate(Position r) const } double CADSurface::distance(Position p, Direction u, bool coincident) const { - moab::EntityHandle surf = dagmc_ptr->entity_by_id(2, id); + moab::EntityHandle surf = dagmc_ptr->entity_by_id(2, id_); moab::EntityHandle hit_surf; double dist; dagmc_ptr->ray_fire(surf, p.xyz, u.xyz, hit_surf, dist, NULL, 0, 0); @@ -267,17 +267,19 @@ double CADSurface::distance(Position p, Direction u, bool coincident) const { Direction CADSurface::normal(Position p) const { Direction u; - moab::EntityHandle surf = dagmc_ptr->entity_by_id(2, id); + moab::EntityHandle surf = dagmc_ptr->entity_by_id(2, id_); dagmc_ptr->get_angle(surf, p.xyz, u.xyz); return u; } BoundingBox CADSurface::bounding_box() const { - moab::EntityHandle surf = dagmc_ptr->entity_by_id(2, id); + moab::EntityHandle surf = dagmc_ptr->entity_by_id(2, id_); double min[3], max[3]; dagmc_ptr->getobb(surf, min, max); return BoundingBox(min[0], max[0], min[1], max[1], min[2], max[2]); } + +void CADSurface::to_hdf5(hid_t group_id) const {} #endif //============================================================================== // PeriodicSurface implementation