From ef8bf2b8204df6f2cf76861ce9d924f64346b2fa Mon Sep 17 00:00:00 2001 From: Sterling Harper Date: Sun, 19 Aug 2018 18:39:28 -0400 Subject: [PATCH] Move Cell.distribcell_index to C++ --- src/cell.cpp | 3 +- src/cell.h | 4 + src/geometry.F90 | 2 +- src/geometry_aux.cpp | 100 +++++++++++++++++------ src/geometry_aux.h | 13 +-- src/geometry_header.F90 | 16 +++- src/input_xml.F90 | 81 ++++-------------- src/tallies/tally_filter_distribcell.F90 | 4 +- 8 files changed, 114 insertions(+), 109 deletions(-) diff --git a/src/cell.cpp b/src/cell.cpp index da5998b90..a68ab9d0b 100644 --- a/src/cell.cpp +++ b/src/cell.cpp @@ -4,7 +4,6 @@ #include #include -#include "constants.h" #include "error.h" #include "geometry.h" #include "hdf5_interface.h" @@ -617,6 +616,8 @@ extern "C" { int32_t cell_n_instances(Cell* c) {return c->n_instances;} + int cell_distribcell_index(Cell* c) {return c->distribcell_index;} + int cell_material_size(Cell* c) {return c->material.size();} //TODO: off-by-one diff --git a/src/cell.h b/src/cell.h index 42b34ae40..8b7c3a78b 100644 --- a/src/cell.h +++ b/src/cell.h @@ -7,6 +7,7 @@ #include #include +#include "constants.h" #include "hdf5.h" #include "pugixml.hpp" @@ -64,6 +65,9 @@ public: int32_t fill; //!< Universe # filling this cell int32_t n_instances{0}; //!< Number of instances of this cell + //! \brief Index corresponding to this cell in distribcell arrays + int distribcell_index{C_NONE}; + //! \brief Material(s) within this cell. //! //! May be multiple materials for distribcell. diff --git a/src/geometry.F90 b/src/geometry.F90 index 60c4948bb..f717810bb 100644 --- a/src/geometry.F90 +++ b/src/geometry.F90 @@ -98,7 +98,7 @@ contains ! Distributed instances of this cell have different ! materials/temperatures. Determine which instance this is for ! assigning the matching material/temperature. - distribcell_index = c % distribcell_index + distribcell_index = c % distribcell_index() offset = 0 do k = 1, p % n_coord if (cells(p % coord(k) % cell) % type() == FILL_UNIVERSE) then diff --git a/src/geometry_aux.cpp b/src/geometry_aux.cpp index 93c08fbe1..6d7d5badb 100644 --- a/src/geometry_aux.cpp +++ b/src/geometry_aux.cpp @@ -153,17 +153,90 @@ find_root_universe() //============================================================================== void -allocate_offset_tables(int n_maps) +prepare_distribcell(int32_t* filter_cell_list, int n) { + // Read the list of cells contained in distribcell filters from Fortran. + std::unordered_set distribcells; + for (int i = 0; i < n; i++) { + distribcells.insert(filter_cell_list[i]); + } + + // Find all cells with distributed materials or temperatures. Make sure that + // the number of materials/temperatures matches the number of cell instances. + for (int i = 0; i < global_cells.size(); i++) { + Cell& c {*global_cells[i]}; + + if (c.material.size() > 1) { + if (c.material.size() != c.n_instances) { + std::stringstream err_msg; + err_msg << "Cell " << c.id << " was specified with " + << c.material.size() << " materials but has " << c.n_instances + << " distributed instances. The number of materials must equal " + "one or the number of instances."; + fatal_error(err_msg); + } + distribcells.insert(i); + } + + if (c.sqrtkT.size() > 1) { + if (c.sqrtkT.size() != c.n_instances) { + std::stringstream err_msg; + err_msg << "Cell " << c.id << " was specified with " + << c.sqrtkT.size() << " temperatures but has " << c.n_instances + << " distributed instances. The number of temperatures must equal " + "one or the number of instances."; + fatal_error(err_msg); + } + distribcells.insert(i); + } + } + + // Search through universes for distributed cells and assign each one a + // unique distribcell array index. + //TODO: off-by-one + int distribcell_index = 1; + std::vector target_univ_ids; + for (Universe* u : global_universes) { + for (auto cell_indx : u->cells) { + if (distribcells.find(cell_indx) != distribcells.end()) { + global_cells[cell_indx]->distribcell_index = distribcell_index; + target_univ_ids.push_back(u->id); + ++distribcell_index; + } + } + } + + // Allocate the cell and lattice offset tables. + int n_maps = target_univ_ids.size(); for (Cell* c : global_cells) { if (c->type != FILL_MATERIAL) { c->offset.resize(n_maps, C_NONE); } } - for (Lattice* lat : lattices_c) { lat->allocate_offset_table(n_maps); } + + // Fill the cell and lattice offset tables. + for (int map = 0; map < target_univ_ids.size(); map++) { + auto target_univ_id = target_univ_ids[map]; + for (Universe* univ : global_universes) { + int32_t offset {0}; // TODO: is this a bug? It matches F90 implementation. + for (int32_t cell_indx : univ->cells) { + Cell& c = *global_cells[cell_indx]; + + if (c.type == FILL_UNIVERSE) { + c.offset[map] = offset; + int32_t search_univ = c.fill; + offset += count_universe_instances(search_univ, target_univ_id); + + } else if (c.type == FILL_LATTICE) { + Lattice& lat = *lattices_c[c.fill]; + offset = lat.fill_offset_table(offset, target_univ_id, map); + } + } + } + } } //============================================================================== @@ -221,29 +294,6 @@ count_universe_instances(int32_t search_univ, int32_t target_univ_id) //============================================================================== -void -fill_offset_tables(int32_t target_univ_id, int map) -{ - for (Universe* univ : global_universes) { - int32_t offset {0}; // TODO: is this a bug? It matches F90 implementation. - for (int32_t cell_indx : univ->cells) { - Cell& c = *global_cells[cell_indx]; - - if (c.type == FILL_UNIVERSE) { - c.offset[map] = offset; - int32_t search_univ = c.fill; - offset += count_universe_instances(search_univ, target_univ_id); - - } else if (c.type == FILL_LATTICE) { - Lattice& lat = *lattices_c[c.fill]; - offset = lat.fill_offset_table(offset, target_univ_id, map); - } - } - } -} - -//============================================================================== - std::string distribcell_path_inner(int32_t target_cell, int32_t map, int32_t target_offset, const Universe& search_univ, int32_t offset) diff --git a/src/geometry_aux.h b/src/geometry_aux.h index fdc6a16b7..e3fd034c6 100644 --- a/src/geometry_aux.h +++ b/src/geometry_aux.h @@ -32,10 +32,10 @@ extern "C" void assign_temperatures(); extern "C" int32_t find_root_universe(); //============================================================================== -//! Allocate storage in Lattice and Cell objects for distribcell offset tables. +//! Populate all data structures needed for distribcells. //============================================================================== -extern "C" void allocate_offset_tables(int n_maps); +extern "C" void prepare_distribcell(int32_t* filter_cell_list, int n); //============================================================================== //! Recursively search through the geometry and count cell instances. @@ -59,15 +59,6 @@ extern "C" void count_cell_instances(int32_t univ_indx); extern "C" int count_universe_instances(int32_t search_univ, int32_t target_univ_id); -//============================================================================== -//! Populate Cell and Lattice distribcell offset tables. -//! @param target_univ_id The ID of the universe to be counted. -//! @param map The index of the distribcell map that defines the offsets for the -//! target universe. -//============================================================================== - -extern "C" void fill_offset_tables(int32_t target_univ_id, int map); - //============================================================================== //! Find the length necessary for a string to contain a distribcell path. //! @param target_cell The index of the Cell in the global Cell array. diff --git a/src/geometry_header.F90 b/src/geometry_header.F90 index f7f533bd7..272966b33 100644 --- a/src/geometry_header.F90 +++ b/src/geometry_header.F90 @@ -60,6 +60,13 @@ module geometry_header integer(C_INT32_T) :: n_instances end function cell_n_instances_c + function cell_distribcell_index_c(cell_ptr) & + bind(C, name='cell_distribcell_index') result(distribcell_index) + import C_PTR, C_INT + type(C_PTR), intent(in), value :: cell_ptr + integer(C_INT) :: distribcell_index + end function cell_distribcell_index_c + function cell_material_size_c(cell_ptr) bind(C, name='cell_material_size') & result(n) import C_PTR, C_INT @@ -266,8 +273,6 @@ module geometry_header integer, allocatable :: region(:) ! Definition of spatial region as ! Boolean expression of half-spaces - integer :: distribcell_index ! Index corresponding to this cell in - ! distribcell arrays ! Rotation matrix and translation vector real(8), allocatable :: translation(:) @@ -282,6 +287,7 @@ module geometry_header procedure :: universe => cell_universe procedure :: fill => cell_fill procedure :: n_instances => cell_n_instances + procedure :: distribcell_index => cell_distribcell_index procedure :: material_size => cell_material_size procedure :: material => cell_material procedure :: sqrtkT_size => cell_sqrtkT_size @@ -414,6 +420,12 @@ contains n_instances = cell_n_instances_c(this % ptr) end function cell_n_instances + function cell_distribcell_index(this) result(distribcell_index) + class(Cell), intent(in) :: this + integer(C_INT) :: distribcell_index + distribcell_index = cell_distribcell_index_c(this % ptr) + end function cell_distribcell_index + function cell_material_size(this) result(n) class(Cell), intent(in) :: this integer(C_INT) :: n diff --git a/src/input_xml.F90 b/src/input_xml.F90 index 6c26da0b0..7d8a8393d 100644 --- a/src/input_xml.F90 +++ b/src/input_xml.F90 @@ -51,25 +51,21 @@ module input_xml subroutine adjust_indices() bind(C) end subroutine adjust_indices - subroutine allocate_offset_tables(n_maps) bind(C) - import C_INT - integer(C_INT), intent(in), value :: n_maps - end subroutine allocate_offset_tables - subroutine assign_temperatures() bind(C) end subroutine assign_temperatures - subroutine fill_offset_tables(target_univ_id, map) bind(C) - import C_INT32_T, C_INT - integer(C_INT32_T), intent(in), value :: target_univ_id - integer(C_INT), intent(in), value :: map - end subroutine fill_offset_tables - 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 prepare_distribcell_c(cell_list, n) & + bind(C, name="prepare_distribcell") + import C_INT32_T, C_INT + integer(C_INT), intent(in), value :: n + integer(C_INT32_T), intent(in) :: cell_list(n) + end subroutine prepare_distribcell_c + subroutine read_surfaces(node_ptr) bind(C) import C_PTR type(C_PTR) :: node_ptr @@ -1007,7 +1003,7 @@ contains subroutine read_geometry_xml() integer :: i, j, k - integer :: n, n_mats, n_rlats, n_hlats + integer :: n, n_rlats, n_hlats integer :: id integer :: univ_id integer :: n_cells_in_univ @@ -1015,7 +1011,6 @@ contains logical :: file_exists logical :: boundary_exists character(MAX_LINE_LEN) :: filename - character(MAX_WORD_LEN), allocatable :: sarray(:) character(:), allocatable :: region_spec type(Cell), pointer :: c class(Lattice), pointer :: lat @@ -1100,9 +1095,6 @@ contains c % ptr = cell_pointer(i - 1) - ! Initialize distribcell instances and distribcell index - c % distribcell_index = NONE - ! Get pointer to i-th cell node node_cell = node_cell_list(i) @@ -3821,9 +3813,9 @@ contains subroutine prepare_distribcell() - integer :: i, j, k + integer :: i, j type(SetInt) :: cell_list ! distribcells to track - type(ListInt) :: univ_list ! universes containing distribcells + integer(C_INT32_T), allocatable :: cell_list_c(:) ! Find all cells listed in a distribcell filter. do i = 1, n_tallies @@ -3835,56 +3827,11 @@ contains end do end do - ! Find all cells with multiple (distributed) materials or temperatures. - do i = 1, n_cells - if (cells(i) % material_size() > 1 .or. cells(i) % sqrtkT_size() > 1) then - call cell_list % add(i) - end if - end do - - ! Make sure the number of distributed materials and temperatures matches the - ! number of respective cell instances. - do i = 1, n_cells - associate (c => cells(i)) - if (c % material_size() > 1) then - if (c % material_size() /= c % n_instances()) then - call fatal_error("Cell " // trim(to_str(c % id())) // " was & - &specified with " // trim(to_str(c % material_size())) & - // " materials but has " // trim(to_str(c % n_instances())) & - // " distributed instances. The number of materials must & - &equal one or the number of instances.") - end if - end if - if (c % sqrtkT_size() > 1) then - if (c % sqrtkT_size() /= c % n_instances()) then - call fatal_error("Cell " // trim(to_str(c % id())) // " was & - &specified with " // trim(to_str(c % sqrtkT_size())) & - // " temperatures but has " // trim(to_str(c % n_instances()))& - // " distributed instances. The number of temperatures must & - &equal one or the number of instances.") - end if - end if - end associate - end do - - ! Search through universes for distributed cells and assign each one a - ! unique distribcell array index. - k = 1 - do i = 1, n_universes - do j = 1, size(universes(i) % cells) - if (cell_list % contains(universes(i) % cells(j))) then - cells(universes(i) % cells(j)) % distribcell_index = k - call univ_list % append(universes(i) % id) - k = k + 1 - end if - end do - end do - - ! Allocate and fill cell and lattice offset tables. - call allocate_offset_tables(univ_list % size()) - do i = 1, univ_list % size() - call fill_offset_tables(univ_list % get_item(i), i-1) + allocate(cell_list_c(cell_list % size())) + do i = 1, cell_list % size() + cell_list_c(i) = cell_list % get_item(i) - 1 end do + call prepare_distribcell_c(cell_list_c, cell_list % size()) end subroutine prepare_distribcell diff --git a/src/tallies/tally_filter_distribcell.F90 b/src/tallies/tally_filter_distribcell.F90 index c9564083d..9ef282126 100644 --- a/src/tallies/tally_filter_distribcell.F90 +++ b/src/tallies/tally_filter_distribcell.F90 @@ -54,7 +54,7 @@ contains integer :: distribcell_index, offset, i - distribcell_index = cells(this % cell) % distribcell_index + distribcell_index = cells(this % cell) % distribcell_index() offset = 0 do i = 1, p % n_coord if (cells(p % coord(i) % cell) % type() == FILL_UNIVERSE) then @@ -154,7 +154,7 @@ contains end interface ! Get the distribcell index for this cell - map = cells(i_cell) % distribcell_index + map = cells(i_cell) % distribcell_index() path_len = distribcell_path_len(i_cell-1, map-1, target_offset, & root_universe-1)