From 4e3350ebcf685547795f5fa2ce5a7275d3d766f5 Mon Sep 17 00:00:00 2001 From: Patrick Shriwise Date: Mon, 5 Jul 2021 12:35:14 -0500 Subject: [PATCH 01/18] Adding a function to determine the instance of a cell at any level of Particle::coords. --- include/openmc/geometry.h | 11 +++++++++++ src/geometry.cpp | 39 +++++++++++++++++++++++++-------------- 2 files changed, 36 insertions(+), 14 deletions(-) diff --git a/include/openmc/geometry.h b/include/openmc/geometry.h index 0ba29ef48..8590e33e1 100644 --- a/include/openmc/geometry.h +++ b/include/openmc/geometry.h @@ -40,6 +40,17 @@ inline bool coincident(double d1, double d2) { bool check_cell_overlap(Particle& p, bool error=true); +//============================================================================== +//! Get the cell instance for a particle at the specified universe level +//! +//! \param p A particle for which to compute the instance using +//! it's vector of LocalCoord. +//! \param level The level (zero indexed) of the geometry where the instance should be computed. +//! \return The instance of the cell at the specified level. +//============================================================================== + +int cell_instance_at_level(const Particle& p, int level); + //============================================================================== //! Locate a particle in the geometry tree and set its geometry data fields. //! diff --git a/src/geometry.cpp b/src/geometry.cpp index c93e1940b..1cbb74935 100644 --- a/src/geometry.cpp +++ b/src/geometry.cpp @@ -65,6 +65,30 @@ bool check_cell_overlap(Particle& p, bool error) //============================================================================== +int cell_instance_at_level(const Particle& p, int level) { + if (level > p.n_coord()) { + fatal_error(fmt::format("Cell instance at level {} requested, but only {} levels exist in the model.", level, p.n_coord())); + } + // determine the cell instance + Cell& c {*model::cells[p.coord(level).cell]}; + int instance = 0; + for (int i = 0; i < level; i++) { + const auto& c_i {*model::cells[p.coord(i).cell]}; + if (c_i.type_ == Fill::UNIVERSE) { + instance += c_i.offset_[c.distribcell_index_]; + } else if (c_i.type_ == Fill::LATTICE) { + auto& lat {*model::lattices[p.coord(i + 1).lattice]}; + const auto& i_xyz {p.coord(i + 1).lattice_i}; + if (lat.are_valid_indices(i_xyz)) { + instance += lat.offset(c.distribcell_index_, i_xyz); + } + } + } + return instance; +} + +//============================================================================== + bool find_cell_inner(Particle& p, const NeighborList* neighbor_list) { @@ -130,22 +154,9 @@ find_cell_inner(Particle& p, const NeighborList* neighbor_list) // Found a material cell which means this is the lowest coord level. // Find the distribcell instance number. - int offset = 0; if (c.distribcell_index_ >= 0) { - for (int i = 0; i < p.n_coord(); i++) { - const auto& c_i {*model::cells[p.coord(i).cell]}; - if (c_i.type_ == Fill::UNIVERSE) { - offset += c_i.offset_[c.distribcell_index_]; - } else if (c_i.type_ == Fill::LATTICE) { - auto& lat {*model::lattices[p.coord(i + 1).lattice]}; - const auto& i_xyz {p.coord(i + 1).lattice_i}; - if (lat.are_valid_indices(i_xyz)) { - offset += lat.offset(c.distribcell_index_, i_xyz); - } - } - } + p.cell_instance() = cell_instance_at_level(p, p.n_coord() - 1); } - p.cell_instance() = offset; // Set the material and temperature. p.material_last() = p.material(); From 39c52dcf94639a3bb9c8ca02a3e47aac032d08d8 Mon Sep 17 00:00:00 2001 From: Patrick Shriwise Date: Wed, 28 Jul 2021 17:49:33 -0500 Subject: [PATCH 02/18] Incorporate cell instance levels into plotting. Return C_NONE if it is not a distributed cell. --- src/geometry.cpp | 5 +++++ src/plot.cpp | 2 +- 2 files changed, 6 insertions(+), 1 deletion(-) diff --git a/src/geometry.cpp b/src/geometry.cpp index 1cbb74935..563e5467b 100644 --- a/src/geometry.cpp +++ b/src/geometry.cpp @@ -66,11 +66,16 @@ bool check_cell_overlap(Particle& p, bool error) //============================================================================== int cell_instance_at_level(const Particle& p, int level) { + // throw error if the requested level is too deep for the geometry if (level > p.n_coord()) { fatal_error(fmt::format("Cell instance at level {} requested, but only {} levels exist in the model.", level, p.n_coord())); } + // determine the cell instance Cell& c {*model::cells[p.coord(level).cell]}; + + if (c.distribcell_index_ == C_NONE) return C_NONE; + int instance = 0; for (int i = 0; i < level; i++) { const auto& c_i {*model::cells[p.coord(i).cell]}; diff --git a/src/plot.cpp b/src/plot.cpp index 999cf3d29..c3ef36467 100644 --- a/src/plot.cpp +++ b/src/plot.cpp @@ -47,7 +47,7 @@ IdData::set_value(size_t y, size_t x, const Particle& p, int level) { data_(y, x, 1) = NOT_FOUND; } else { data_(y, x, 0) = model::cells.at(p.coord(level).cell)->id_; - data_(y, x, 1) = p.cell_instance(); + data_(y, x, 1) = level == p.n_coord() - 1 ? p.cell_instance() : cell_instance_at_level(p, level); } // set material data From 058e5a9094a8f7c35d40263f198bfffda9d9efc6 Mon Sep 17 00:00:00 2001 From: Patrick Shriwise Date: Sat, 17 Jul 2021 15:46:52 -0500 Subject: [PATCH 03/18] Allowing cell instance filters for non material-filled cells. --- include/openmc/tallies/filter_cell_instance.h | 5 +++++ src/tallies/filter_cell_instance.cpp | 21 ++++++++++++++++++- 2 files changed, 25 insertions(+), 1 deletion(-) diff --git a/include/openmc/tallies/filter_cell_instance.h b/include/openmc/tallies/filter_cell_instance.h index dee9b8bea..a2a204db6 100644 --- a/include/openmc/tallies/filter_cell_instance.h +++ b/include/openmc/tallies/filter_cell_instance.h @@ -46,6 +46,8 @@ public: void set_cell_instances(gsl::span instances); + void set_geom_level(int32_t level); + private: //---------------------------------------------------------------------------- // Data members @@ -55,6 +57,9 @@ private: //! A map from cell/instance indices to filter bin indices. std::unordered_map map_; + + //! Level in the geometry to check for the cell instance + int32_t geom_level_ {-1}; }; } // namespace openmc diff --git a/src/tallies/filter_cell_instance.cpp b/src/tallies/filter_cell_instance.cpp index c41077631..c5a88dc7a 100644 --- a/src/tallies/filter_cell_instance.cpp +++ b/src/tallies/filter_cell_instance.cpp @@ -7,6 +7,7 @@ #include "openmc/capi.h" #include "openmc/cell.h" #include "openmc/error.h" +#include "openmc/geometry.h" #include "openmc/xml_interface.h" namespace openmc { @@ -53,7 +54,7 @@ CellInstanceFilter::set_cell_instances(gsl::span instances) Expects(x.index_cell >= 0); Expects(x.index_cell < model::cells.size()); const auto& c {model::cells[x.index_cell]}; - if (c->type_ != Fill::MATERIAL) { + if (c->type_ != Fill::MATERIAL && geom_level_ < 0) { throw std::invalid_argument{fmt::format( "Cell {} is not filled with a material. Only material cells can be " "used in a cell instance filter.", c->id_)}; @@ -65,12 +66,30 @@ CellInstanceFilter::set_cell_instances(gsl::span instances) n_bins_ = cell_instances_.size(); } +void +CellInstanceFilter::set_geom_level(int32_t level) { + Expects(level >= 0); + geom_level_ = level; +} + void CellInstanceFilter::get_all_bins(const Particle& p, TallyEstimator estimator, FilterMatch& match) const { gsl::index index_cell = p.coord(p.n_coord() - 1).cell; gsl::index instance = p.cell_instance(); + + if (geom_level_ >= 0) { + // if the particle has fewer levels than the cell we're looking for, + // return no bins + if (p.coord(geom_level_).cell == C_NONE) return; + + // otherwise use the cell at the requested level + // and compute the cell instance for this particle's position + index_cell = p.coord(geom_level_).cell; + instance = cell_instance_at_level(p, geom_level_); + } + auto search = map_.find({index_cell, instance}); if (search != map_.end()) { int index_bin = search->second; From ca5dc796b9bf5abd298e4f29a6de863590a4fefb Mon Sep 17 00:00:00 2001 From: Patrick Shriwise Date: Thu, 29 Jul 2021 10:53:14 -0500 Subject: [PATCH 04/18] Updating error messages. Adding cell instance cells to our distribcell set. --- include/openmc/tallies/filter_cell_instance.h | 2 +- src/geometry.cpp | 4 ++-- src/geometry_aux.cpp | 4 ++++ src/tallies/filter_cell_instance.cpp | 6 +++--- 4 files changed, 10 insertions(+), 6 deletions(-) diff --git a/include/openmc/tallies/filter_cell_instance.h b/include/openmc/tallies/filter_cell_instance.h index a2a204db6..a3fbad181 100644 --- a/include/openmc/tallies/filter_cell_instance.h +++ b/include/openmc/tallies/filter_cell_instance.h @@ -59,7 +59,7 @@ private: std::unordered_map map_; //! Level in the geometry to check for the cell instance - int32_t geom_level_ {-1}; + int32_t geom_level_ {C_NONE}; }; } // namespace openmc diff --git a/src/geometry.cpp b/src/geometry.cpp index 563e5467b..7eb731660 100644 --- a/src/geometry.cpp +++ b/src/geometry.cpp @@ -67,8 +67,8 @@ bool check_cell_overlap(Particle& p, bool error) int cell_instance_at_level(const Particle& p, int level) { // throw error if the requested level is too deep for the geometry - if (level > p.n_coord()) { - fatal_error(fmt::format("Cell instance at level {} requested, but only {} levels exist in the model.", level, p.n_coord())); + if (level > model::n_coord_levels) { + fatal_error(fmt::format("Cell instance at level {} requested, but only {} levels exist in the geometry.", level, p.n_coord())); } // determine the cell instance diff --git a/src/geometry_aux.cpp b/src/geometry_aux.cpp index af6ebac9f..9d3062901 100644 --- a/src/geometry_aux.cpp +++ b/src/geometry_aux.cpp @@ -319,9 +319,13 @@ prepare_distribcell() std::unordered_set distribcells; for (auto& filt : model::tally_filters) { auto* distrib_filt = dynamic_cast(filt.get()); + auto* cell_inst_filt = dynamic_cast(filt.get()); if (distrib_filt) { distribcells.insert(distrib_filt->cell()); } + if (cell_inst_filt) { + for (const auto& c_inst : cell_inst_filt->cell_instances()) distribcells.insert(c_inst.index_cell); + } } // By default, add material cells to the list of distributed cells diff --git a/src/tallies/filter_cell_instance.cpp b/src/tallies/filter_cell_instance.cpp index c5a88dc7a..ba70b06bb 100644 --- a/src/tallies/filter_cell_instance.cpp +++ b/src/tallies/filter_cell_instance.cpp @@ -56,8 +56,8 @@ CellInstanceFilter::set_cell_instances(gsl::span instances) const auto& c {model::cells[x.index_cell]}; if (c->type_ != Fill::MATERIAL && geom_level_ < 0) { throw std::invalid_argument{fmt::format( - "Cell {} is not filled with a material. Only material cells can be " - "used in a cell instance filter.", c->id_)}; + "Cell {} is not filled with a material. A geometry level must be specified to" + "use cells filled with a universe or lattice.", c->id_)}; } cell_instances_.push_back(x); map_[x] = cell_instances_.size() - 1; @@ -82,7 +82,7 @@ CellInstanceFilter::get_all_bins(const Particle& p, TallyEstimator estimator, if (geom_level_ >= 0) { // if the particle has fewer levels than the cell we're looking for, // return no bins - if (p.coord(geom_level_).cell == C_NONE) return; + if (p.n_coord() - 1 < geom_level_ || p.coord(geom_level_).cell == C_NONE) return; // otherwise use the cell at the requested level // and compute the cell instance for this particle's position From 885ed95facc879334b170d92bf6191ec3b202e60 Mon Sep 17 00:00:00 2001 From: Patrick Shriwise Date: Thu, 29 Jul 2021 11:28:39 -0500 Subject: [PATCH 05/18] Adding geometry_level to the cell instance filter XML output. --- openmc/filter.py | 15 ++++++++++++++- 1 file changed, 14 insertions(+), 1 deletion(-) diff --git a/openmc/filter.py b/openmc/filter.py index e46effa7b..2edca37ba 100644 --- a/openmc/filter.py +++ b/openmc/filter.py @@ -550,6 +550,8 @@ class CellInstanceFilter(Filter): Unique identifier for the filter num_bins : Integral The number of filter bins + geometry_level : Integral + Geometry level of the specified cell (used to compute instances in the simulation) See Also -------- @@ -559,6 +561,7 @@ class CellInstanceFilter(Filter): def __init__(self, bins, filter_id=None): self.bins = bins self.id = filter_id + self._geom_level = None @Filter.bins.setter def bins(self, bins): @@ -570,6 +573,15 @@ class CellInstanceFilter(Filter): pairs[i, 1] = instance self._bins = pairs + @property + def geometry_level(self): + return self._geometry_level + + @geometry_level.setter + def geometry_level(self, val): + cv.check_type('geometry level', val, Integral) + self._geometry_level = val + def get_pandas_dataframe(self, data_size, stride, **kwargs): """Builds a Pandas DataFrame for the Filter's bins. @@ -618,7 +630,8 @@ class CellInstanceFilter(Filter): element = ET.Element('filter') element.set('id', str(self.id)) element.set('type', self.short_name.lower()) - + if self.geometry_level: + element.set('geometry_level', self.geometry_level) subelement = ET.SubElement(element, 'bins') subelement.text = ' '.join(str(i) for i in self.bins.ravel()) return element From 991c3d1f76d63e73b022261cf19d48e536f2c461 Mon Sep 17 00:00:00 2001 From: Patrick Shriwise Date: Thu, 29 Jul 2021 11:29:00 -0500 Subject: [PATCH 06/18] Adding read of 'geometry_level' for cell instance filters. --- src/tallies/filter_cell_instance.cpp | 4 ++++ 1 file changed, 4 insertions(+) diff --git a/src/tallies/filter_cell_instance.cpp b/src/tallies/filter_cell_instance.cpp index ba70b06bb..568a989f0 100644 --- a/src/tallies/filter_cell_instance.cpp +++ b/src/tallies/filter_cell_instance.cpp @@ -38,6 +38,10 @@ CellInstanceFilter::from_xml(pugi::xml_node node) instances.push_back({index, instance}); } + if (check_for_node(node, "geometry_level")) { + this->set_geom_level(std::stoi(get_node_value(node, "geometry_level"))); + } + this->set_cell_instances(instances); } From 4c55737a9421607b612caceb33a5d1937809c9e2 Mon Sep 17 00:00:00 2001 From: Patrick Shriwise Date: Thu, 29 Jul 2021 14:55:23 -0500 Subject: [PATCH 07/18] Better method for supporting non-material filled cells. --- include/openmc/tallies/filter_cell_instance.h | 6 +++ src/tallies/filter_cell_instance.cpp | 50 +++++++++++-------- 2 files changed, 36 insertions(+), 20 deletions(-) diff --git a/include/openmc/tallies/filter_cell_instance.h b/include/openmc/tallies/filter_cell_instance.h index a3fbad181..89f8d18af 100644 --- a/include/openmc/tallies/filter_cell_instance.h +++ b/include/openmc/tallies/filter_cell_instance.h @@ -55,11 +55,17 @@ private: //! The indices of the cells binned by this filter. vector cell_instances_; + //! The set of cells used in this filter + std::unordered_set cells_; + //! A map from cell/instance indices to filter bin indices. std::unordered_map map_; //! Level in the geometry to check for the cell instance int32_t geom_level_ {C_NONE}; + + //! Indicates if filter uses only material-filled cells + bool material_cells_only_; }; } // namespace openmc diff --git a/src/tallies/filter_cell_instance.cpp b/src/tallies/filter_cell_instance.cpp index 568a989f0..a4476d739 100644 --- a/src/tallies/filter_cell_instance.cpp +++ b/src/tallies/filter_cell_instance.cpp @@ -57,17 +57,20 @@ CellInstanceFilter::set_cell_instances(gsl::span instances) for (auto& x : instances) { Expects(x.index_cell >= 0); Expects(x.index_cell < model::cells.size()); - const auto& c {model::cells[x.index_cell]}; - if (c->type_ != Fill::MATERIAL && geom_level_ < 0) { - throw std::invalid_argument{fmt::format( - "Cell {} is not filled with a material. A geometry level must be specified to" - "use cells filled with a universe or lattice.", c->id_)}; - } cell_instances_.push_back(x); + cells_.insert(x.index_cell); map_[x] = cell_instances_.size() - 1; } n_bins_ = cell_instances_.size(); + + material_cells_only_ = true; + for (const auto& cell_inst : cell_instances_) { + const auto& c = *model::cells[cell_inst.index_cell]; + if (c.type_ == Fill::MATERIAL) continue; + material_cells_only_ = false; + break; + } } void @@ -83,22 +86,29 @@ CellInstanceFilter::get_all_bins(const Particle& p, TallyEstimator estimator, gsl::index index_cell = p.coord(p.n_coord() - 1).cell; gsl::index instance = p.cell_instance(); - if (geom_level_ >= 0) { - // if the particle has fewer levels than the cell we're looking for, - // return no bins - if (p.n_coord() - 1 < geom_level_ || p.coord(geom_level_).cell == C_NONE) return; - - // otherwise use the cell at the requested level - // and compute the cell instance for this particle's position - index_cell = p.coord(geom_level_).cell; - instance = cell_instance_at_level(p, geom_level_); + if (cells_.count(index_cell) > 0) { + auto search = map_.find({index_cell, instance}); + if (search != map_.end()) { + int index_bin = search->second; + match.bins_.push_back(index_bin); + match.weights_.push_back(1.0); + } } - auto search = map_.find({index_cell, instance}); - if (search != map_.end()) { - int index_bin = search->second; - match.bins_.push_back(index_bin); - match.weights_.push_back(1.0); + if (!material_cells_only_) { + for (int i = 0; i < p.n_coord() - 1; i++) { + gsl::index index_cell = p.coord(i).cell; + // if this cell isn't used on the filter, move on + if (cells_.count(index_cell) == 0) continue; + + // if this cell is used in the filter, check the instance as well + gsl::index instance = cell_instance_at_level(p, i); + auto search = map_.find({index_cell, instance}); + if (search != map_.end()) { + match.bins_.push_back(search->second); + match.bins_.push_back(1.0); + } + } } } From 8d4c6f87083878dcebb507f80b05f0ff0c027b60 Mon Sep 17 00:00:00 2001 From: Patrick Shriwise Date: Thu, 29 Jul 2021 14:57:58 -0500 Subject: [PATCH 08/18] Reverting changes to the cell instance filter Python object. --- openmc/filter.py | 15 +-------------- 1 file changed, 1 insertion(+), 14 deletions(-) diff --git a/openmc/filter.py b/openmc/filter.py index 2edca37ba..e46effa7b 100644 --- a/openmc/filter.py +++ b/openmc/filter.py @@ -550,8 +550,6 @@ class CellInstanceFilter(Filter): Unique identifier for the filter num_bins : Integral The number of filter bins - geometry_level : Integral - Geometry level of the specified cell (used to compute instances in the simulation) See Also -------- @@ -561,7 +559,6 @@ class CellInstanceFilter(Filter): def __init__(self, bins, filter_id=None): self.bins = bins self.id = filter_id - self._geom_level = None @Filter.bins.setter def bins(self, bins): @@ -573,15 +570,6 @@ class CellInstanceFilter(Filter): pairs[i, 1] = instance self._bins = pairs - @property - def geometry_level(self): - return self._geometry_level - - @geometry_level.setter - def geometry_level(self, val): - cv.check_type('geometry level', val, Integral) - self._geometry_level = val - def get_pandas_dataframe(self, data_size, stride, **kwargs): """Builds a Pandas DataFrame for the Filter's bins. @@ -630,8 +618,7 @@ class CellInstanceFilter(Filter): element = ET.Element('filter') element.set('id', str(self.id)) element.set('type', self.short_name.lower()) - if self.geometry_level: - element.set('geometry_level', self.geometry_level) + subelement = ET.SubElement(element, 'bins') subelement.text = ' '.join(str(i) for i in self.bins.ravel()) return element From 930fb20010d1c1b6d2a6cfc8a3f7c9ce90405cdd Mon Sep 17 00:00:00 2001 From: Patrick Shriwise Date: Thu, 29 Jul 2021 15:03:08 -0500 Subject: [PATCH 09/18] Removing geom_level form cell instance filter. --- src/tallies/filter_cell_instance.cpp | 11 +---------- 1 file changed, 1 insertion(+), 10 deletions(-) diff --git a/src/tallies/filter_cell_instance.cpp b/src/tallies/filter_cell_instance.cpp index a4476d739..0db6c7d21 100644 --- a/src/tallies/filter_cell_instance.cpp +++ b/src/tallies/filter_cell_instance.cpp @@ -38,10 +38,6 @@ CellInstanceFilter::from_xml(pugi::xml_node node) instances.push_back({index, instance}); } - if (check_for_node(node, "geometry_level")) { - this->set_geom_level(std::stoi(get_node_value(node, "geometry_level"))); - } - this->set_cell_instances(instances); } @@ -51,6 +47,7 @@ CellInstanceFilter::set_cell_instances(gsl::span instances) // Clear existing cells cell_instances_.clear(); cell_instances_.reserve(instances.size()); + cells_.clear(); map_.clear(); // Update cells and mapping @@ -73,12 +70,6 @@ CellInstanceFilter::set_cell_instances(gsl::span instances) } } -void -CellInstanceFilter::set_geom_level(int32_t level) { - Expects(level >= 0); - geom_level_ = level; -} - void CellInstanceFilter::get_all_bins(const Particle& p, TallyEstimator estimator, FilterMatch& match) const From 3114c935629fd6812bcb04cbd0115a4b776d5f49 Mon Sep 17 00:00:00 2001 From: Patrick Shriwise Date: Thu, 29 Jul 2021 20:01:21 -0500 Subject: [PATCH 10/18] Adding an extra layer to the geometry to test the cell instance filter regression test for a universe-filled cell. --- .../filter_cellinstance/inputs_true.dat | 21 ++++++++++--------- .../filter_cellinstance/test.py | 8 +++++-- 2 files changed, 17 insertions(+), 12 deletions(-) diff --git a/tests/regression_tests/filter_cellinstance/inputs_true.dat b/tests/regression_tests/filter_cellinstance/inputs_true.dat index 5abf73433..5df992627 100644 --- a/tests/regression_tests/filter_cellinstance/inputs_true.dat +++ b/tests/regression_tests/filter_cellinstance/inputs_true.dat @@ -2,18 +2,19 @@ - - - - + + + + + 2 2 4 4 -4 -4 -1 2 2 2 -2 1 2 2 -2 2 1 2 -2 2 2 1 +2 3 3 3 +3 2 3 3 +3 3 2 3 +3 3 3 2 @@ -48,10 +49,10 @@ - 3 0 3 1 3 2 3 3 3 4 3 5 3 6 3 7 3 8 3 9 3 10 3 11 2 0 2 1 2 2 2 3 + 4 0 4 1 4 2 4 3 4 4 4 5 4 6 4 7 4 8 4 9 4 10 4 11 2 0 2 1 2 2 2 3 3 0 3 1 3 2 3 3 - 2 3 2 2 2 1 2 0 3 11 3 10 3 9 3 8 3 7 3 6 3 5 3 4 3 3 3 2 3 1 3 0 + 3 3 3 2 3 1 3 0 2 3 2 2 2 1 2 0 4 11 4 10 4 9 4 8 4 7 4 6 4 5 4 4 4 3 4 2 4 1 4 0 1 diff --git a/tests/regression_tests/filter_cellinstance/test.py b/tests/regression_tests/filter_cellinstance/test.py index dc00394fc..bd708788c 100644 --- a/tests/regression_tests/filter_cellinstance/test.py +++ b/tests/regression_tests/filter_cellinstance/test.py @@ -22,7 +22,10 @@ def model(): cyl1 = openmc.ZCylinder(r=0.7) c1 = openmc.Cell(fill=m1, region=-cyl1) c2 = openmc.Cell(fill=m2, region=+cyl1) - u1 = openmc.Universe(cells=[c1, c2]) + # intermediate universe + ui = openmc.Universe(cells=[c2]) + ci = openmc.Cell(fill=ui) + u1 = openmc.Universe(cells=[c1, ci]) cyl2 = openmc.ZCylinder(r=0.5) c3 = openmc.Cell(fill=m1, region=-cyl2) @@ -50,7 +53,8 @@ def model(): model.settings.source = openmc.Source(space=openmc.stats.Point()) instances = ([(c3, i) for i in range(c3.num_instances)] + - [(c2, i) for i in range(c2.num_instances)]) + [(c2, i) for i in range(c2.num_instances)] + + [(ci, i) for i in range(ci.num_instances)]) f1 = openmc.CellInstanceFilter(instances) f2 = openmc.CellInstanceFilter(instances[::-1]) t1 = openmc.Tally() From f5651337e1275f485d5c2c9af6411c10a0d19ff0 Mon Sep 17 00:00:00 2001 From: Patrick Shriwise Date: Thu, 29 Jul 2021 20:02:15 -0500 Subject: [PATCH 11/18] Correction to match's weights vector. --- src/tallies/filter_cell_instance.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/tallies/filter_cell_instance.cpp b/src/tallies/filter_cell_instance.cpp index 0db6c7d21..d846aa294 100644 --- a/src/tallies/filter_cell_instance.cpp +++ b/src/tallies/filter_cell_instance.cpp @@ -97,7 +97,7 @@ CellInstanceFilter::get_all_bins(const Particle& p, TallyEstimator estimator, auto search = map_.find({index_cell, instance}); if (search != map_.end()) { match.bins_.push_back(search->second); - match.bins_.push_back(1.0); + match.weights_.push_back(1.0); } } } From 82719c0fdeb53f6da73b96bf10d62376c7f917e3 Mon Sep 17 00:00:00 2001 From: Patrick Shriwise Date: Fri, 30 Jul 2021 10:33:19 -0500 Subject: [PATCH 12/18] Updating cell instance regression test to include a higher level cell tally. --- .../filter_cellinstance/inputs_true.dat | 2 +- .../filter_cellinstance/results_true.dat | 16 +++++++ .../filter_cellinstance/test.py | 45 +++++++++++++------ tests/testing_harness.py | 4 ++ 4 files changed, 53 insertions(+), 14 deletions(-) diff --git a/tests/regression_tests/filter_cellinstance/inputs_true.dat b/tests/regression_tests/filter_cellinstance/inputs_true.dat index 5df992627..f9d3ca58d 100644 --- a/tests/regression_tests/filter_cellinstance/inputs_true.dat +++ b/tests/regression_tests/filter_cellinstance/inputs_true.dat @@ -1,6 +1,6 @@ - + diff --git a/tests/regression_tests/filter_cellinstance/results_true.dat b/tests/regression_tests/filter_cellinstance/results_true.dat index 053ecf576..c6832c6e1 100644 --- a/tests/regression_tests/filter_cellinstance/results_true.dat +++ b/tests/regression_tests/filter_cellinstance/results_true.dat @@ -33,6 +33,14 @@ tally 1: 2.015188E+02 9.951936E+00 2.345815E+01 +1.061333E+01 +2.508486E+01 +2.767786E+01 +1.781264E+02 +2.898957E+01 +2.015188E+02 +9.951936E+00 +2.345815E+01 tally 2: 9.951936E+00 2.345815E+01 @@ -42,6 +50,14 @@ tally 2: 1.781264E+02 1.061333E+01 2.508486E+01 +9.951936E+00 +2.345815E+01 +2.898957E+01 +2.015188E+02 +2.767786E+01 +1.781264E+02 +1.061333E+01 +2.508486E+01 6.794395E-02 1.135006E-03 1.299309E-01 diff --git a/tests/regression_tests/filter_cellinstance/test.py b/tests/regression_tests/filter_cellinstance/test.py index bd708788c..ccd547f91 100644 --- a/tests/regression_tests/filter_cellinstance/test.py +++ b/tests/regression_tests/filter_cellinstance/test.py @@ -2,8 +2,27 @@ import openmc import openmc.model import pytest +from numpy.testing import assert_array_almost_equal + from tests.testing_harness import PyAPITestHarness +class CellInstanceFilterTest(PyAPITestHarness): + + def _compare_results(self): + with openmc.StatePoint(self.statepoint_name) as sp: + # we expect the tally results for the instances of + # cells 2 and 3 to be the same as 2 is nested + # in a universe directly under 3 + t1 = sp.tallies[1] + f1 = sp.filters[1] + c2_bins = [tuple(tuple(i) for i in f1.bins if i[0] == 2)] + c2_mean = t1.get_values(filters=[openmc.CellInstanceFilter], filter_bins=c2_bins) + c3_bins = [tuple(tuple(i) for i in f1.bins if i[0] == 3)] + c3_mean = t1.get_values(filters=[openmc.CellInstanceFilter], filter_bins=c3_bins) + assert_array_almost_equal(c2_mean, c3_mean) + + return super()._compare_results() + @pytest.fixture def model(): @@ -23,23 +42,23 @@ def model(): c1 = openmc.Cell(fill=m1, region=-cyl1) c2 = openmc.Cell(fill=m2, region=+cyl1) # intermediate universe - ui = openmc.Universe(cells=[c2]) - ci = openmc.Cell(fill=ui) - u1 = openmc.Universe(cells=[c1, ci]) + u1 = openmc.Universe(cells=[c2]) + c3 = openmc.Cell(fill=u1) + u2 = openmc.Universe(cells=[c1, c3]) cyl2 = openmc.ZCylinder(r=0.5) - c3 = openmc.Cell(fill=m1, region=-cyl2) - c4 = openmc.Cell(fill=m2, region=+cyl2) - u2 = openmc.Universe(cells=[c3, c4]) + c4 = openmc.Cell(fill=m1, region=-cyl2) + c5 = openmc.Cell(fill=m2, region=+cyl2) + u3 = openmc.Universe(cells=[c4, c5]) lat = openmc.RectLattice() lat.lower_left = (-4, -4) lat.pitch = (2, 2) lat.universes = [ - [u1, u2, u2, u2], - [u2, u1, u2, u2], - [u2, u2, u1, u2], - [u2, u2, u2, u1] + [u2, u3, u3, u3], + [u3, u2, u3, u3], + [u3, u3, u2, u3], + [u3, u3, u3, u2] ] box = openmc.model.rectangular_prism(8.0, 8.0, boundary_type='reflective') main_cell = openmc.Cell(fill=lat, region=box) @@ -52,9 +71,9 @@ def model(): model.settings.particles = 1000 model.settings.source = openmc.Source(space=openmc.stats.Point()) - instances = ([(c3, i) for i in range(c3.num_instances)] + + instances = ([(c4, i) for i in range(c4.num_instances)] + [(c2, i) for i in range(c2.num_instances)] + - [(ci, i) for i in range(ci.num_instances)]) + [(c3, i) for i in range(c3.num_instances)]) f1 = openmc.CellInstanceFilter(instances) f2 = openmc.CellInstanceFilter(instances[::-1]) t1 = openmc.Tally() @@ -69,5 +88,5 @@ def model(): def test_cell_instance(model): - harness = PyAPITestHarness('statepoint.5.h5', model) + harness = CellInstanceFilterTest('statepoint.5.h5', model) harness.main() diff --git a/tests/testing_harness.py b/tests/testing_harness.py index 48f121d45..92d34ccdc 100644 --- a/tests/testing_harness.py +++ b/tests/testing_harness.py @@ -115,6 +115,10 @@ class TestHarness: return outstr + @property + def statepoint_name(self): + return self._sp_name + def _write_results(self, results_string): """Write the results to an ASCII file.""" with open('results_test.dat', 'w') as fh: From f5f087ef07ea53016e843924436aed11039051e6 Mon Sep 17 00:00:00 2001 From: Patrick Shriwise Date: Fri, 30 Jul 2021 10:41:25 -0500 Subject: [PATCH 13/18] Moving import and adding comment to test. --- tests/regression_tests/filter_cellinstance/test.py | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/tests/regression_tests/filter_cellinstance/test.py b/tests/regression_tests/filter_cellinstance/test.py index ccd547f91..3f5fa54b2 100644 --- a/tests/regression_tests/filter_cellinstance/test.py +++ b/tests/regression_tests/filter_cellinstance/test.py @@ -1,11 +1,11 @@ +from numpy.testing import assert_array_almost_equal import openmc import openmc.model import pytest -from numpy.testing import assert_array_almost_equal - from tests.testing_harness import PyAPITestHarness + class CellInstanceFilterTest(PyAPITestHarness): def _compare_results(self): @@ -41,7 +41,7 @@ def model(): cyl1 = openmc.ZCylinder(r=0.7) c1 = openmc.Cell(fill=m1, region=-cyl1) c2 = openmc.Cell(fill=m2, region=+cyl1) - # intermediate universe + # intermediate universe containing only cell 2 u1 = openmc.Universe(cells=[c2]) c3 = openmc.Cell(fill=u1) u2 = openmc.Universe(cells=[c1, c3]) From b2efde103307cba6e7b2f797a74b742a6666a690 Mon Sep 17 00:00:00 2001 From: Patrick Shriwise Date: Fri, 30 Jul 2021 10:41:53 -0500 Subject: [PATCH 14/18] Removing unused method declaration. --- include/openmc/tallies/filter_cell_instance.h | 5 ----- 1 file changed, 5 deletions(-) diff --git a/include/openmc/tallies/filter_cell_instance.h b/include/openmc/tallies/filter_cell_instance.h index 89f8d18af..789a4151e 100644 --- a/include/openmc/tallies/filter_cell_instance.h +++ b/include/openmc/tallies/filter_cell_instance.h @@ -46,8 +46,6 @@ public: void set_cell_instances(gsl::span instances); - void set_geom_level(int32_t level); - private: //---------------------------------------------------------------------------- // Data members @@ -61,9 +59,6 @@ private: //! A map from cell/instance indices to filter bin indices. std::unordered_map map_; - //! Level in the geometry to check for the cell instance - int32_t geom_level_ {C_NONE}; - //! Indicates if filter uses only material-filled cells bool material_cells_only_; }; From b42a04cc7614d2c68b6ba707df0f588984a0d981 Mon Sep 17 00:00:00 2001 From: Patrick Shriwise Date: Fri, 30 Jul 2021 10:42:11 -0500 Subject: [PATCH 15/18] Adding comment to cell_instance_at_level. --- src/geometry.cpp | 1 + 1 file changed, 1 insertion(+) diff --git a/src/geometry.cpp b/src/geometry.cpp index 7eb731660..f38f649f3 100644 --- a/src/geometry.cpp +++ b/src/geometry.cpp @@ -74,6 +74,7 @@ int cell_instance_at_level(const Particle& p, int level) { // determine the cell instance Cell& c {*model::cells[p.coord(level).cell]}; + // quick exit if this cell doesn't have distribcell instances if (c.distribcell_index_ == C_NONE) return C_NONE; int instance = 0; From 595a00d8d5a291413a00d4bd4eda04580d68cf7f Mon Sep 17 00:00:00 2001 From: Patrick Shriwise Date: Fri, 30 Jul 2021 10:42:52 -0500 Subject: [PATCH 16/18] Lowering nesting of check for non-material-filled cell instances. --- src/tallies/filter_cell_instance.cpp | 24 ++++++++++++------------ 1 file changed, 12 insertions(+), 12 deletions(-) diff --git a/src/tallies/filter_cell_instance.cpp b/src/tallies/filter_cell_instance.cpp index d846aa294..0dcc4b094 100644 --- a/src/tallies/filter_cell_instance.cpp +++ b/src/tallies/filter_cell_instance.cpp @@ -86,19 +86,19 @@ CellInstanceFilter::get_all_bins(const Particle& p, TallyEstimator estimator, } } - if (!material_cells_only_) { - for (int i = 0; i < p.n_coord() - 1; i++) { - gsl::index index_cell = p.coord(i).cell; - // if this cell isn't used on the filter, move on - if (cells_.count(index_cell) == 0) continue; + if (material_cells_only_) return; - // if this cell is used in the filter, check the instance as well - gsl::index instance = cell_instance_at_level(p, i); - auto search = map_.find({index_cell, instance}); - if (search != map_.end()) { - match.bins_.push_back(search->second); - match.weights_.push_back(1.0); - } + for (int i = 0; i < p.n_coord() - 1; i++) { + gsl::index index_cell = p.coord(i).cell; + // if this cell isn't used on the filter, move on + if (cells_.count(index_cell) == 0) continue; + + // if this cell is used in the filter, check the instance as well + gsl::index instance = cell_instance_at_level(p, i); + auto search = map_.find({index_cell, instance}); + if (search != map_.end()) { + match.bins_.push_back(search->second); + match.weights_.push_back(1.0); } } } From 822858dfbbdfee9a290150627618f54da29b2dd0 Mon Sep 17 00:00:00 2001 From: Patrick Shriwise Date: Mon, 2 Aug 2021 16:39:10 -0500 Subject: [PATCH 17/18] Update include/openmc/geometry.h Co-authored-by: Paul Romano --- include/openmc/geometry.h | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/include/openmc/geometry.h b/include/openmc/geometry.h index 8590e33e1..617921f45 100644 --- a/include/openmc/geometry.h +++ b/include/openmc/geometry.h @@ -44,7 +44,7 @@ bool check_cell_overlap(Particle& p, bool error=true); //! Get the cell instance for a particle at the specified universe level //! //! \param p A particle for which to compute the instance using -//! it's vector of LocalCoord. +//! its coordinates //! \param level The level (zero indexed) of the geometry where the instance should be computed. //! \return The instance of the cell at the specified level. //============================================================================== From 404dc1ecc9571aab924ea404212c38e2735da5a9 Mon Sep 17 00:00:00 2001 From: Patrick Shriwise Date: Mon, 2 Aug 2021 17:38:27 -0500 Subject: [PATCH 18/18] Ensuring Particle::cell_instance_ is set to zero if not a distribcell. --- src/geometry.cpp | 2 ++ 1 file changed, 2 insertions(+) diff --git a/src/geometry.cpp b/src/geometry.cpp index f38f649f3..8f5bda38f 100644 --- a/src/geometry.cpp +++ b/src/geometry.cpp @@ -77,6 +77,7 @@ int cell_instance_at_level(const Particle& p, int level) { // quick exit if this cell doesn't have distribcell instances if (c.distribcell_index_ == C_NONE) return C_NONE; + // compute the cell's instance int instance = 0; for (int i = 0; i < level; i++) { const auto& c_i {*model::cells[p.coord(i).cell]}; @@ -159,6 +160,7 @@ find_cell_inner(Particle& p, const NeighborList* neighbor_list) if (c.type_ == Fill::MATERIAL) { // Found a material cell which means this is the lowest coord level. + p.cell_instance() = 0; // Find the distribcell instance number. if (c.distribcell_index_ >= 0) { p.cell_instance() = cell_instance_at_level(p, p.n_coord() - 1);