From 964fdfde3001b7040eedfaddf1a3667dba6ca76e Mon Sep 17 00:00:00 2001 From: Paul Romano Date: Fri, 15 Feb 2019 06:35:59 -0600 Subject: [PATCH] Remove cells array on Fortran side --- CMakeLists.txt | 1 - include/openmc/cell.h | 2 + include/openmc/dagmc.h | 1 + include/openmc/geometry_aux.h | 4 +- include/openmc/lattice.h | 6 + include/openmc/surface.h | 6 + openmc/capi/cell.py | 5 +- src/api.F90 | 6 +- src/cell.cpp | 127 ++++++++-------- src/dagmc.cpp | 19 ++- src/geometry.F90 | 1 - src/geometry_aux.cpp | 47 +++++- src/geometry_header.F90 | 279 ---------------------------------- src/initialize.cpp | 1 - src/input_xml.F90 | 226 --------------------------- src/lattice.cpp | 7 +- src/source.cpp | 2 +- src/surface.cpp | 7 +- src/tallies/tally.F90 | 1 - tests/unit_tests/test_capi.py | 2 +- 20 files changed, 162 insertions(+), 588 deletions(-) delete mode 100644 src/geometry_header.F90 diff --git a/CMakeLists.txt b/CMakeLists.txt index 43e250c34..65292367d 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -310,7 +310,6 @@ add_library(libopenmc SHARED src/dict_header.F90 src/error.F90 src/geometry.F90 - src/geometry_header.F90 src/hdf5_interface.F90 src/initialize.F90 src/input_xml.F90 diff --git a/include/openmc/cell.h b/include/openmc/cell.h index ba202c343..cb77d7e41 100644 --- a/include/openmc/cell.h +++ b/include/openmc/cell.h @@ -199,6 +199,8 @@ public: // Non-member functions //============================================================================== +void read_cells(pugi::xml_node node); + #ifdef DAGMC int32_t next_cell(DAGCell* cur_cell, DAGSurface* surf_xed); #endif diff --git a/include/openmc/dagmc.h b/include/openmc/dagmc.h index 4a129ef09..5a944a42d 100644 --- a/include/openmc/dagmc.h +++ b/include/openmc/dagmc.h @@ -26,6 +26,7 @@ extern moab::DagMC* DAG; extern "C" void load_dagmc_geometry(); extern "C" void free_memory_dagmc(); +void read_geometry_dagmc(); extern "C" pugi::xml_document* read_uwuw_materials(); bool get_uwuw_materials_xml(std::string& s); diff --git a/include/openmc/geometry_aux.h b/include/openmc/geometry_aux.h index 3e5cb93e1..a6bfe06e2 100644 --- a/include/openmc/geometry_aux.h +++ b/include/openmc/geometry_aux.h @@ -10,6 +10,8 @@ namespace openmc { +void read_geometry_xml(); + //============================================================================== //! Replace Universe, Lattice, and Material IDs with indices. //============================================================================== @@ -111,7 +113,7 @@ extern "C" int maximum_levels(int32_t univ); //! Deallocates global vectors and maps for cells, universes, and lattices. //============================================================================== -extern "C" void free_memory_geometry_c(); +extern "C" void free_memory_geometry(); } // namespace openmc #endif // OPENMC_GEOMETRY_AUX_H diff --git a/include/openmc/lattice.h b/include/openmc/lattice.h index 5a5eff0b4..68728de5f 100644 --- a/include/openmc/lattice.h +++ b/include/openmc/lattice.h @@ -273,5 +273,11 @@ private: std::array pitch_; //!< Lattice tile width and height }; +//============================================================================== +// Non-member functions +//============================================================================== + +void read_lattices(pugi::xml_node node); + } // namespace openmc #endif // OPENMC_LATTICE_H diff --git a/include/openmc/surface.h b/include/openmc/surface.h index 3df791ba3..ce9f78399 100644 --- a/include/openmc/surface.h +++ b/include/openmc/surface.h @@ -395,6 +395,12 @@ public: void to_hdf5_inner(hid_t group_id) const; }; +//============================================================================== +// Non-member functions +//============================================================================== + +void read_surfaces(pugi::xml_node node); + //============================================================================== // Fortran compatibility functions //============================================================================== diff --git a/openmc/capi/cell.py b/openmc/capi/cell.py index 1a4c8996f..894e4c79e 100644 --- a/openmc/capi/cell.py +++ b/openmc/capi/cell.py @@ -38,6 +38,7 @@ _dll.openmc_cell_set_temperature.errcheck = _error_handler _dll.openmc_get_cell_index.argtypes = [c_int32, POINTER(c_int32)] _dll.openmc_get_cell_index.restype = c_int _dll.openmc_get_cell_index.errcheck = _error_handler +_dll.cells_size.restype = c_int class Cell(_FortranObjectWithID): @@ -156,10 +157,10 @@ class _CellMapping(Mapping): def __iter__(self): for i in range(len(self)): - yield Cell(index=i + 1).id + yield Cell(index=i).id def __len__(self): - return c_int32.in_dll(_dll, 'n_cells').value + return _dll.cells_size() def __repr__(self): return repr(dict(self)) diff --git a/src/api.F90 b/src/api.F90 index 5200927cf..cb3856861 100644 --- a/src/api.F90 +++ b/src/api.F90 @@ -39,7 +39,7 @@ contains err = E_UNASSIGNED if (found) then - index = p % coord(p % n_coord) % cell + 1 + index = p % coord(p % n_coord) % cell instance = p % cell_instance err = 0 else @@ -59,7 +59,6 @@ contains subroutine free_memory() bind(C) use bank_header - use geometry_header use material_header use nuclide_header use photon_header @@ -93,6 +92,9 @@ contains subroutine free_memory_surfaces() bind(C) end subroutine + subroutine free_memory_geometry() bind(C) + end subroutine + subroutine sab_clear() bind(C) end subroutine end interface diff --git a/src/cell.cpp b/src/cell.cpp index 56635309f..2a2bd82d6 100644 --- a/src/cell.cpp +++ b/src/cell.cpp @@ -23,8 +23,6 @@ namespace openmc { namespace model { -int32_t n_cells {0}; - std::vector cells; std::unordered_map cell_map; @@ -643,18 +641,18 @@ void DAGCell::to_hdf5(hid_t group_id) const { return; } // Non-method functions //============================================================================== -extern "C" void -read_cells(pugi::xml_node* node) +void read_cells(pugi::xml_node node) { // Count the number of cells. - for (pugi::xml_node cell_node: node->children("cell")) {model::n_cells++;} - if (model::n_cells == 0) { + int n_cells = 0; + for (pugi::xml_node cell_node: node.children("cell")) {n_cells++;} + if (n_cells == 0) { fatal_error("No cells found in geometry.xml!"); } // Loop over XML cell elements and populate the array. - model::cells.reserve(model::n_cells); - for (pugi::xml_node cell_node : node->children("cell")) { + model::cells.reserve(n_cells); + for (pugi::xml_node cell_node : node.children("cell")) { model::cells.push_back(new CSGCell(cell_node)); } @@ -699,9 +697,8 @@ read_cells(pugi::xml_node* node) extern "C" int openmc_cell_get_fill(int32_t index, int* type, int32_t** indices, int32_t* n) { - if (index >= 1 && index <= model::cells.size()) { - //TODO: off-by-one - Cell& c {*model::cells[index - 1]}; + if (index >= 0 && index < model::cells.size()) { + Cell& c {*model::cells[index]}; *type = c.type_; if (c.type_ == FILL_MATERIAL) { *indices = c.material_.data(); @@ -721,9 +718,8 @@ extern "C" int openmc_cell_set_fill(int32_t index, int type, int32_t n, const int32_t* indices) { - if (index >= 1 && index <= model::cells.size()) { - //TODO: off-by-one - Cell& c {*model::cells[index - 1]}; + if (index >= 0 && index < model::cells.size()) { + Cell& c {*model::cells[index]}; if (type == FILL_MATERIAL) { c.type_ = FILL_MATERIAL; c.material_.clear(); @@ -756,9 +752,9 @@ openmc_cell_set_fill(int32_t index, int type, int32_t n, extern "C" int openmc_cell_set_temperature(int32_t index, double T, const int32_t* instance) { - if (index >= 1 && index <= model::cells.size()) { + if (index >= 0 && index < model::cells.size()) { //TODO: off-by-one - Cell& c {*model::cells[index - 1]}; + Cell& c {*model::cells[index]}; if (instance) { if (*instance >= 0 && *instance < c.sqrtkT_.size()) { @@ -781,28 +777,66 @@ openmc_cell_set_temperature(int32_t index, double T, const int32_t* instance) return 0; } +//! Return the index in the cells array of a cell with a given ID +extern "C" int +openmc_get_cell_index(int32_t id, int32_t* index) +{ + auto it = model::cell_map.find(id); + if (it != model::cell_map.end()) { + *index = it->second; + return 0; + } else { + set_errmsg("No cell exists with ID=" + std::to_string(id) + "."); + return OPENMC_E_INVALID_ID; + } +} + +//! Return the ID of a cell +extern "C" int +openmc_cell_get_id(int32_t index, int32_t* id) +{ + if (index >= 0 && index < model::cells.size()) { + *id = model::cells[index]->id_; + return 0; + } else { + set_errmsg("Index in cells array is out of bounds."); + return OPENMC_E_OUT_OF_BOUNDS; + } +} + +//! Set the ID of a cell +extern "C" int +openmc_cell_set_id(int32_t index, int32_t id) +{ + if (index >= 0 && index < model::cells.size()) { + model::cells[index]->id_ = id; + model::cell_map[id] = index; + return 0; + } else { + set_errmsg("Index in cells array is out of bounds."); + return OPENMC_E_OUT_OF_BOUNDS; + } +} + +//! Extend the cells array by n elements +extern "C" int +openmc_extend_cells(int32_t n, int32_t* index_start, int32_t* index_end) +{ + if (index_start) *index_start = model::cells.size(); + if (index_end) *index_end = model::cells.size() + n - 1; + for (int32_t i = 0; i < n; i++) { + model::cells.push_back(new CSGCell()); + } + return 0; +} + + //============================================================================== // Fortran compatibility functions //============================================================================== extern "C" { - Cell* cell_pointer(int32_t cell_ind) {return model::cells[cell_ind];} - - int32_t cell_id(Cell* c) {return c->id_;} - - void - cell_set_id(Cell* c, int32_t id) - { - c->id_ = id; - - // Find the index of this cell and update the cell map. - for (int i = 0; i < model::cells.size(); i++) { - if (model::cells[i] == c) { - model::cell_map[id] = i; - break; - } - } - } + int cells_size() { return model::cells.size(); } #ifdef DAGMC int32_t next_cell(DAGCell* cur_cell, DAGSurface* surf_xed) @@ -818,33 +852,6 @@ extern "C" { return cur_cell->dagmc_ptr_->index_by_handle(new_vol); } #endif - - int32_t cell_universe(Cell* c) {return c->universe_;} - - int32_t cell_fill(Cell* c) {return c->fill_;} - - int cell_material_size(Cell* c) {return c->material_.size();} - - //TODO: off-by-one - int32_t cell_material(Cell* c, int i) - { - int32_t mat = c->material_[i-1]; - if (mat == MATERIAL_VOID) return MATERIAL_VOID; - return mat + 1; - } - - int cell_sqrtkT_size(Cell* c) {return c->sqrtkT_.size();} - - double cell_sqrtkT(Cell* c, int i) {return c->sqrtkT_[i];} - - void extend_cells_c(int32_t n) - { - model::cells.reserve(model::cells.size() + n); - for (int32_t i = 0; i < n; i++) { - model::cells.push_back(new CSGCell()); - } - model::n_cells = model::cells.size(); - } } diff --git a/src/dagmc.cpp b/src/dagmc.cpp index 22f9d392a..3af23bc46 100644 --- a/src/dagmc.cpp +++ b/src/dagmc.cpp @@ -4,6 +4,7 @@ #include "openmc/constants.h" #include "openmc/error.h" #include "openmc/file_utils.h" +#include "openmc/geometry_aux.h" #include "openmc/string_utils.h" #include "openmc/settings.h" #include "openmc/geometry.h" @@ -133,9 +134,9 @@ void load_dagmc_geometry() /// Cells (Volumes) \\\ // initialize cell objects - model::n_cells = model::DAG->num_entities(3); + int n_cells = model::DAG->num_entities(3); moab::EntityHandle graveyard = 0; - for (int i = 0; i < model::n_cells; i++) { + for (int i = 0; i < n_cells; i++) { moab::EntityHandle vol_handle = model::DAG->entity_by_index(3, i+1); // set cell ids using global IDs @@ -312,6 +313,20 @@ void load_dagmc_geometry() return; } +void read_geometry_dagmc() +{ + // Check if dagmc.h5m exists + std::string filename = settings::path_input + "dagmc.h5m" + if (!file_exists(filename)) { + fatal_error("Geometry DAGMC file '" + filename + "' does not exist!"); + } + + write_message("Reading DAGMC geometry...", 5); + load_dagmc_geometry(); + + model::root_universe = find_root_universe() +} + void free_memory_dagmc() { delete model::DAG; diff --git a/src/geometry.F90 b/src/geometry.F90 index 0abe96ee4..1b7ead90a 100644 --- a/src/geometry.F90 +++ b/src/geometry.F90 @@ -1,7 +1,6 @@ module geometry use constants - use geometry_header use particle_header use simulation_header use settings diff --git a/src/geometry_aux.cpp b/src/geometry_aux.cpp index 9073c511d..f69a63d9b 100644 --- a/src/geometry_aux.cpp +++ b/src/geometry_aux.cpp @@ -4,10 +4,14 @@ #include #include +#include "pugixml.hpp" + #include "openmc/cell.h" #include "openmc/constants.h" #include "openmc/container_util.h" +#include "openmc/dagmc.h" #include "openmc/error.h" +#include "openmc/file_utils.h" #include "openmc/geometry.h" #include "openmc/lattice.h" #include "openmc/material.h" @@ -19,6 +23,46 @@ namespace openmc { +void read_geometry_xml() +{ +#ifdef DAGMC + if (settings::dagmc) { + read_geometry_dagmc(); + return; + } +#endif + + // Display output message + write_message("Reading geometry XML file...", 5); + + // Check if geometry.xml exists + std::string filename = settings::path_input + "geometry.xml"; + if (!file_exists(filename)) { + fatal_error("Geometry XML file '" + filename + "' does not exist!"); + } + + // Parse settings.xml file + pugi::xml_document doc; + auto result = doc.load_file(filename.c_str()); + if (!result) { + fatal_error("Error processing geometry.xml file."); + } + + // Get root element + pugi::xml_node root = doc.document_element(); + + // Read surfaces, cells, lattice + read_surfaces(root); + read_cells(root); + read_lattices(root); + + // ========================================================================== + // SETUP UNIVERSES + + // Allocate universes, universe cell arrays, and assign base universe + model::root_universe = find_root_universe(); +} + //============================================================================== void @@ -473,12 +517,11 @@ maximum_levels(int32_t univ) //============================================================================== void -free_memory_geometry_c() +free_memory_geometry() { for (Cell* c : model::cells) {delete c;} model::cells.clear(); model::cell_map.clear(); - model::n_cells = 0; for (Universe* u : model::universes) {delete u;} model::universes.clear(); diff --git a/src/geometry_header.F90 b/src/geometry_header.F90 deleted file mode 100644 index afb9e2c0d..000000000 --- a/src/geometry_header.F90 +++ /dev/null @@ -1,279 +0,0 @@ -module geometry_header - - use, intrinsic :: ISO_C_BINDING - - use dict_header, only: DictIntInt - use error - use string, only: to_str - - implicit none - - interface - function cell_pointer(cell_ind) bind(C) result(ptr) - import C_PTR, C_INT32_T - integer(C_INT32_T), intent(in), value :: cell_ind - type(C_PTR) :: ptr - end function cell_pointer - - function cell_id_c(cell_ptr) bind(C, name='cell_id') result(id) - import C_PTR, C_INT32_T - type(C_PTR), intent(in), value :: cell_ptr - integer(C_INT32_T) :: id - end function cell_id_c - - subroutine cell_set_id_c(cell_ptr, id) bind(C, name='cell_set_id') - import C_PTR, C_INT32_T - type(C_PTR), intent(in), value :: cell_ptr - integer(C_INT32_T), intent(in), value :: id - end subroutine cell_set_id_c - - function cell_universe_c(cell_ptr) bind(C, name='cell_universe') & - result(universe) - import C_PTR, C_INT32_T - type(C_PTR), intent(in), value :: cell_ptr - integer(C_INT32_T) :: universe - end function cell_universe_c - - function cell_fill_c(cell_ptr) bind(C, name="cell_fill") result(fill) - import C_PTR, C_INT32_T - type(C_PTR), intent(in), value :: cell_ptr - integer(C_INT32_T) :: fill - end function cell_fill_c - - function cell_material_size_c(cell_ptr) bind(C, name='cell_material_size') & - result(n) - import C_PTR, C_INT - type(C_PTR), intent(in), value :: cell_ptr - integer(C_INT) :: n - end function cell_material_size_c - - function cell_material_c(cell_ptr, i) bind(C, name='cell_material') & - result(mat) - import C_PTR, C_INT, C_INT32_T - type(C_PTR), intent(in), value :: cell_ptr - integer(C_INT), intent(in), value :: i - integer(C_INT32_T) :: mat - end function cell_material_c - - function cell_sqrtkT_size_c(cell_ptr) bind(C, name='cell_sqrtkT_size') & - result(n) - import C_PTR, C_INT - type(C_PTR), intent(in), value :: cell_ptr - integer(C_INT) :: n - end function cell_sqrtkT_size_c - - function cell_sqrtkT_c(cell_ptr, i) bind(C, name='cell_sqrtkT') & - result(sqrtkT) - import C_PTR, C_INT, C_DOUBLE - type(C_PTR), intent(in), value :: cell_ptr - integer(C_INT), intent(in), value :: i - real(C_DOUBLE) :: sqrtkT - end function cell_sqrtkT_c - - subroutine extend_cells_c(n) bind(C) - import C_INT32_t - integer(C_INT32_T), intent(in), value :: n - end subroutine extend_cells_c - end interface - -!=============================================================================== -! CELL defines a closed volume by its bounding surfaces -!=============================================================================== - - type Cell - type(C_PTR) :: ptr - contains - - procedure :: id => cell_id - procedure :: set_id => cell_set_id - procedure :: universe => cell_universe - procedure :: fill => cell_fill - procedure :: material_size => cell_material_size - procedure :: material => cell_material - procedure :: sqrtkT_size => cell_sqrtkT_size - procedure :: sqrtkT => cell_sqrtkT - - end type Cell - - ! array index of the root universe - integer(C_INT), bind(C) :: root_universe - - integer(C_INT32_T), bind(C) :: n_cells ! # of cells - integer(C_INT32_T), bind(C) :: n_universes ! # of universes - - type(Cell), allocatable, target :: cells(:) - - ! Dictionaries which map user IDs to indices in the global arrays - type(DictIntInt) :: cell_dict - -contains - -!=============================================================================== - - function cell_id(this) result(id) - class(Cell), intent(in) :: this - integer(C_INT32_T) :: id - id = cell_id_c(this % ptr) - end function cell_id - - subroutine cell_set_id(this, id) - class(Cell), intent(in) :: this - integer(C_INT32_T), intent(in) :: id - call cell_set_id_c(this % ptr, id) - end subroutine cell_set_id - - function cell_universe(this) result(universe) - class(Cell), intent(in) :: this - integer(C_INT32_T) :: universe - universe = cell_universe_c(this % ptr) - end function cell_universe - - function cell_fill(this) result(fill) - class(Cell), intent(in) :: this - integer(C_INT32_T) :: fill - fill = cell_fill_c(this % ptr) - end function cell_fill - - function cell_material_size(this) result(n) - class(Cell), intent(in) :: this - integer(C_INT) :: n - n = cell_material_size_c(this % ptr) - end function cell_material_size - - function cell_material(this, i) result(mat) - class(Cell), intent(in) :: this - integer, intent(in) :: i - integer(C_INT32_T) :: mat - mat = cell_material_c(this % ptr, i) - end function cell_material - - function cell_sqrtkT_size(this) result(n) - class(Cell), intent(in) :: this - integer :: n - n = cell_sqrtkT_size_c(this % ptr) - end function cell_sqrtkT_size - - function cell_sqrtkT(this, i) result(sqrtkT) - class(Cell), intent(in) :: this - integer, intent(in) :: i - real(C_DOUBLE) :: sqrtkT - sqrtkT = cell_sqrtkT_c(this % ptr, i) - end function cell_sqrtkT - -!=============================================================================== -! FREE_MEMORY_GEOMETRY deallocates global arrays defined in this module -!=============================================================================== - - subroutine free_memory_geometry() - interface - subroutine free_memory_geometry_c() bind(C) - end subroutine free_memory_geometry_c - end interface - - call free_memory_geometry_c() - - n_cells = 0 - n_universes = 0 - - if (allocated(cells)) deallocate(cells) - - call cell_dict % clear() - - end subroutine free_memory_geometry - -!=============================================================================== -! C API FUNCTIONS -!=============================================================================== - - function openmc_extend_cells(n, index_start, index_end) result(err) bind(C) - ! Extend the cells array by n elements - integer(C_INT32_T), value, intent(in) :: n - integer(C_INT32_T), optional, intent(out) :: index_start - integer(C_INT32_T), optional, intent(out) :: index_end - integer(C_INT32_T) :: i - integer(C_INT) :: err - - type(Cell), allocatable :: temp(:) ! temporary cells array - - if (n_cells == 0) then - ! Allocate cells array - allocate(cells(n)) - else - ! Allocate cells array with increased size - allocate(temp(n_cells + n)) - - ! Copy original cells to temporary array - temp(1:n_cells) = cells - - ! Move allocation from temporary array - call move_alloc(FROM=temp, TO=cells) - end if - - ! Return indices in cells array - if (present(index_start)) index_start = n_cells + 1 - if (present(index_end)) index_end = n_cells + n - - ! Extend the C++ cells array and get pointers to the C++ objects - call extend_cells_c(n) - do i = n_cells - n, n_cells - cells(i) % ptr = cell_pointer(i - 1) - end do - - err = 0 - end function openmc_extend_cells - - - function openmc_get_cell_index(id, index) result(err) bind(C) - ! Return the index in the cells array of a cell with a given ID - integer(C_INT32_T), value :: id - integer(C_INT32_T), intent(out) :: index - integer(C_INT) :: err - - if (allocated(cells)) then - if (cell_dict % has(id)) then - index = cell_dict % get(id) - err = 0 - else - err = E_INVALID_ID - call set_errmsg("No cell exists with ID=" // trim(to_str(id)) // ".") - end if - else - err = E_ALLOCATE - call set_errmsg("Memory has not been allocated for cells.") - end if - end function openmc_get_cell_index - - - function openmc_cell_get_id(index, id) result(err) bind(C) - ! Return the ID of a cell - integer(C_INT32_T), value :: index - integer(C_INT32_T), intent(out) :: id - integer(C_INT) :: err - - if (index >= 1 .and. index <= size(cells)) then - id = cells(index) % id() - err = 0 - else - err = E_OUT_OF_BOUNDS - call set_errmsg("Index in cells array is out of bounds.") - end if - end function openmc_cell_get_id - - - function openmc_cell_set_id(index, id) result(err) bind(C) - ! Set the ID of a cell - integer(C_INT32_T), value, intent(in) :: index - integer(C_INT32_T), value, intent(in) :: id - integer(C_INT) :: err - - if (index >= 1 .and. index <= n_cells) then - call cells(index) % set_id(id) - call cell_dict % set(id, index) - err = 0 - else - err = E_OUT_OF_BOUNDS - call set_errmsg("Index in cells array is out of bounds.") - end if - end function openmc_cell_set_id - -end module geometry_header diff --git a/src/initialize.cpp b/src/initialize.cpp index e454917d0..59662d577 100644 --- a/src/initialize.cpp +++ b/src/initialize.cpp @@ -31,7 +31,6 @@ // data/functions from Fortran side extern "C" void read_command_line(); -extern "C" void read_geometry_xml(); extern "C" void read_materials_xml(); extern "C" void read_plots_xml(); extern "C" void read_tallies_xml(); diff --git a/src/input_xml.F90 b/src/input_xml.F90 index 60c6485ac..8da108e35 100644 --- a/src/input_xml.F90 +++ b/src/input_xml.F90 @@ -3,9 +3,7 @@ module input_xml use, intrinsic :: ISO_C_BINDING use constants - use dict_header, only: DictIntInt, DictCharInt, DictEntryCI use error, only: fatal_error, warning, write_message, openmc_err_msg - use geometry_header #ifdef DAGMC use dagmc_header #endif @@ -30,42 +28,11 @@ module input_xml save interface - subroutine count_cell_instances(univ_indx) bind(C) - import C_INT32_T - integer(C_INT32_T), intent(in), value :: univ_indx - end subroutine count_cell_instances - - subroutine read_surfaces(node_ptr) bind(C) - import C_PTR - type(C_PTR) :: node_ptr - end subroutine read_surfaces - - subroutine read_cells(node_ptr) bind(C) - import C_PTR - type(C_PTR) :: node_ptr - end subroutine read_cells - - subroutine read_lattices(node_ptr) bind(C) - import C_PTR - type(C_PTR) :: node_ptr - end subroutine read_lattices - subroutine read_materials(node_ptr) bind(C) import C_PTR type(C_PTR) :: node_ptr end subroutine read_materials - function find_root_universe() bind(C) result(root) - import C_INT32_T - integer(C_INT32_T) :: root - end function find_root_universe - - function maximum_levels(univ) bind(C) result(n) - import C_INT32_T, C_INT - integer(C_INT32_T), intent(in), value :: univ - integer(C_INT) :: n - end function maximum_levels - subroutine read_plots(node_ptr) bind(C) import C_PTR type(C_PTR) :: node_ptr @@ -81,199 +48,6 @@ module input_xml contains -#ifdef DAGMC - -!=============================================================================== -! READ_GEOMETRY_DAGMC reads data from a DAGMC .h5m file, checking -! for material properties and surface boundary conditions -! some universe information is spoofed for now -!=============================================================================== - - subroutine read_geometry_dagmc() - - integer :: i - integer :: univ_id - integer :: n_cells_in_univ - logical :: file_exists - character(MAX_LINE_LEN) :: filename - type(Cell), pointer :: c - type(VectorInt) :: univ_ids ! List of all universe IDs - type(DictIntInt) :: cells_in_univ_dict ! Used to count how many cells each - ! universe contains - - ! Check if dagmc.h5m exists - filename = trim(path_input) // "dagmc.h5m" - inquire(FILE=filename, EXIST=file_exists) - if (.not. file_exists) then - call fatal_error("Geometry DAGMC file '" // trim(filename) // "' does not & - &exist!") - end if - - call write_message("Reading DAGMC geometry...", 5) - call load_dagmc_geometry() - call allocate_cells() - - ! setup universe data structs - do i = 1, n_cells - c => cells(i) - ! additional metadata spoofing - univ_id = c % universe() - - if (.not. cells_in_univ_dict % has(univ_id)) then - n_universes = n_universes + 1 - n_cells_in_univ = 1 - call univ_ids % push_back(univ_id) - else - n_cells_in_univ = 1 + cells_in_univ_dict % get(univ_id) - end if - call cells_in_univ_dict % set(univ_id, n_cells_in_univ) - end do - - root_universe = find_root_universe() - - end subroutine read_geometry_dagmc - -#endif - -!=============================================================================== -! READ_GEOMETRY_XML reads data from a geometry.xml file and parses it, checking -! for errors and placing properly-formatted data in the right data structures -!=============================================================================== - - subroutine read_geometry_xml() bind(C) - - integer :: i, n - integer :: univ_id - integer :: n_cells_in_univ - real(8) :: phi, theta, psi - logical :: file_exists - character(MAX_LINE_LEN) :: filename - type(Cell), pointer :: c - type(XMLDocument) :: doc - type(XMLNode) :: root - type(XMLNode) :: node_cell - type(XMLNode), allocatable :: node_cell_list(:) - type(VectorInt) :: univ_ids ! List of all universe IDs - type(DictIntInt) :: cells_in_univ_dict ! Used to count how many cells each - ! universe contains -#ifdef DAGMC - if (dagmc) then - call read_geometry_dagmc() - return - end if -#endif - - ! Display output message - call write_message("Reading geometry XML file...", 5) - - ! Check if geometry.xml exists - filename = trim(path_input) // "geometry.xml" - inquire(FILE=filename, EXIST=file_exists) - if (.not. file_exists) then - call fatal_error("Geometry XML file '" // trim(filename) // "' does not & - &exist!") - end if - - ! Parse geometry.xml file - call doc % load_file(filename) - root = doc % document_element() - - ! ========================================================================== - ! READ SURFACES FROM GEOMETRY.XML - - call read_surfaces(root % ptr) - - ! ========================================================================== - ! READ CELLS FROM GEOMETRY.XML - - call read_cells(root % ptr) - - ! Get pointer to list of XML - call get_node_list(root, "cell", node_cell_list) - - ! Get number of tags - n_cells = size(node_cell_list) - - ! Check for no cells - if (n_cells == 0) then - call fatal_error("No cells found in geometry.xml!") - end if - - ! Allocate cells array - allocate(cells(n_cells)) - - n_universes = 0 - do i = 1, n_cells - c => cells(i) - - c % ptr = cell_pointer(i - 1) - - ! Get pointer to i-th cell node - node_cell = node_cell_list(i) - - ! Check to make sure 'id' hasn't been used - if (cell_dict % has(c % id())) then - call fatal_error("Two or more cells use the same unique ID: " & - // to_str(c % id())) - end if - - ! Add cell to dictionary - call cell_dict % set(c % id(), i) - - ! For cells, we also need to check if there's a new universe -- - ! also for every cell add 1 to the count of cells for the - ! specified universe - univ_id = c % universe() - if (.not. cells_in_univ_dict % has(univ_id)) then - n_universes = n_universes + 1 - n_cells_in_univ = 1 - call univ_ids % push_back(univ_id) - else - n_cells_in_univ = 1 + cells_in_univ_dict % get(univ_id) - end if - call cells_in_univ_dict % set(univ_id, n_cells_in_univ) - - end do - - ! ========================================================================== - ! READ LATTICES FROM GEOMETRY.XML - - call read_lattices(root % ptr) - - ! ========================================================================== - ! SETUP UNIVERSES - - ! Allocate universes, universe cell arrays, and assign base universe - root_universe = find_root_universe() - - ! Clear dictionary - call cells_in_univ_dict%clear() - - ! Close geometry XML file - call doc % clear() - - end subroutine read_geometry_xml - - subroutine allocate_cells() - integer :: i - type(Cell), pointer :: c - - ! Allocate cells array - allocate(cells(n_cells)) - - do i = 1, n_cells - c => cells(i) - c % ptr = cell_pointer(i - 1) - ! Check to make sure 'id' hasn't been used - if (cell_dict % has(c % id())) then - call fatal_error("Two or more cells use the same unique ID: " & - // to_str(c % id())) - end if - ! Add cell to dictionary - call cell_dict % set(c % id(), i) - end do - end subroutine allocate_cells - subroutine read_materials_xml() bind(C) logical :: file_exists ! does materials.xml exist? character(MAX_LINE_LEN) :: filename ! absolute path to materials.xml diff --git a/src/lattice.cpp b/src/lattice.cpp index d42440439..d8185de62 100644 --- a/src/lattice.cpp +++ b/src/lattice.cpp @@ -864,13 +864,12 @@ HexLattice::to_hdf5_inner(hid_t lat_group) const // Non-method functions //============================================================================== -extern "C" void -read_lattices(pugi::xml_node *node) +void read_lattices(pugi::xml_node node) { - for (pugi::xml_node lat_node : node->children("lattice")) { + for (pugi::xml_node lat_node : node.children("lattice")) { model::lattices.push_back(new RectLattice(lat_node)); } - for (pugi::xml_node lat_node : node->children("hex_lattice")) { + for (pugi::xml_node lat_node : node.children("hex_lattice")) { model::lattices.push_back(new HexLattice(lat_node)); } diff --git a/src/source.cpp b/src/source.cpp index 0e23bb3dd..231a95e21 100644 --- a/src/source.cpp +++ b/src/source.cpp @@ -172,7 +172,7 @@ Bank SourceDistribution::sample() const if (space_box) { if (space_box->only_fissionable()) { // Determine material - auto c = model::cells[cell_index - 1]; + auto c = model::cells[cell_index]; auto mat_index = c->material_.size() == 1 ? c->material_[0] : c->material_[instance]; diff --git a/src/surface.cpp b/src/surface.cpp index 0cceee9f4..cd1f4fc2b 100644 --- a/src/surface.cpp +++ b/src/surface.cpp @@ -1061,12 +1061,11 @@ void SurfaceQuadric::to_hdf5_inner(hid_t group_id) const //============================================================================== -extern "C" void -read_surfaces(pugi::xml_node* node) +void read_surfaces(pugi::xml_node node) { // Count the number of surfaces. int n_surfaces = 0; - for (pugi::xml_node surf_node : node->children("surface")) {n_surfaces++;} + for (pugi::xml_node surf_node : node.children("surface")) {n_surfaces++;} if (n_surfaces == 0) { fatal_error("No surfaces found in geometry.xml!"); } @@ -1076,7 +1075,7 @@ read_surfaces(pugi::xml_node* node) { pugi::xml_node surf_node; int i_surf; - for (surf_node = node->child("surface"), i_surf = 0; surf_node; + for (surf_node = node.child("surface"), i_surf = 0; surf_node; surf_node = surf_node.next_sibling("surface"), i_surf++) { std::string surf_type = get_node_value(surf_node, "type", true, true); diff --git a/src/tallies/tally.F90 b/src/tallies/tally.F90 index e7cdcbc1d..823ef1acc 100644 --- a/src/tallies/tally.F90 +++ b/src/tallies/tally.F90 @@ -4,7 +4,6 @@ module tally use bank_header use constants - use geometry_header use material_header use message_passing use mgxs_interface diff --git a/tests/unit_tests/test_capi.py b/tests/unit_tests/test_capi.py index 89cf5d708..d5943b7be 100644 --- a/tests/unit_tests/test_capi.py +++ b/tests/unit_tests/test_capi.py @@ -72,7 +72,7 @@ def test_cell(capi_init): cell = openmc.capi.cells[1] assert isinstance(cell.fill, openmc.capi.Material) cell.fill = openmc.capi.materials[1] - assert str(cell) == 'Cell[1]' + assert str(cell) == 'Cell[0]' def test_new_cell(capi_init):