From 13c2c837b1795e1c392671cb6a4620450e2aee90 Mon Sep 17 00:00:00 2001 From: Paul Romano Date: Fri, 28 Jun 2019 14:43:17 -0500 Subject: [PATCH] Major cleanup of tally filter interfaces --- include/openmc/capi.h | 8 +-- include/openmc/tallies/filter.h | 34 +++++++++-- include/openmc/tallies/filter_azimuthal.h | 4 ++ include/openmc/tallies/filter_cell.h | 6 +- include/openmc/tallies/filter_delayedgroup.h | 4 ++ include/openmc/tallies/filter_distribcell.h | 4 +- include/openmc/tallies/filter_energy.h | 4 ++ include/openmc/tallies/filter_legendre.h | 2 + include/openmc/tallies/filter_material.h | 6 +- include/openmc/tallies/filter_mu.h | 4 ++ include/openmc/tallies/filter_particle.h | 2 + include/openmc/tallies/filter_polar.h | 4 ++ include/openmc/tallies/filter_sph_harm.h | 6 ++ include/openmc/tallies/filter_sptl_legendre.h | 6 ++ include/openmc/tallies/filter_surface.h | 6 +- include/openmc/tallies/filter_universe.h | 6 +- openmc/capi/filter.py | 12 ++-- src/tallies/filter.cpp | 58 +++++++++++++++++-- src/tallies/filter_azimuthal.cpp | 33 +++++++---- src/tallies/filter_cell.cpp | 44 ++++++++------ src/tallies/filter_delayedgroup.cpp | 25 +++++--- src/tallies/filter_distribcell.cpp | 27 +++++---- src/tallies/filter_energy.cpp | 31 +++++++--- src/tallies/filter_legendre.cpp | 11 +++- src/tallies/filter_material.cpp | 53 +++++++++-------- src/tallies/filter_mu.cpp | 34 +++++++---- src/tallies/filter_particle.cpp | 19 +++++- src/tallies/filter_polar.cpp | 34 +++++++---- src/tallies/filter_sph_harm.cpp | 54 ++++++++++------- src/tallies/filter_sptl_legendre.cpp | 53 ++++++++++++----- src/tallies/filter_surface.cpp | 41 ++++++++----- src/tallies/filter_universe.cpp | 41 ++++++++----- src/tallies/tally.cpp | 31 +--------- 33 files changed, 477 insertions(+), 230 deletions(-) diff --git a/include/openmc/capi.h b/include/openmc/capi.h index 9f05a56dcc..afca455ddf 100644 --- a/include/openmc/capi.h +++ b/include/openmc/capi.h @@ -17,8 +17,8 @@ extern "C" { int openmc_cell_set_fill(int32_t index, int type, int32_t n, const int32_t* indices); int openmc_cell_set_id(int32_t index, int32_t id); int openmc_cell_set_temperature(int32_t index, double T, const int32_t* instance); - int openmc_energy_filter_get_bins(int32_t index, double** energies, int32_t* n); - int openmc_energy_filter_set_bins(int32_t index, int32_t n, const double* energies); + int openmc_energy_filter_get_bins(int32_t index, double** energies, size_t* n); + int openmc_energy_filter_set_bins(int32_t index, size_t n, const double* energies); int openmc_extend_cells(int32_t n, int32_t* index_start, int32_t* index_end); int openmc_extend_filters(int32_t n, int32_t* index_start, int32_t* index_end); int openmc_extend_materials(int32_t n, int32_t* index_start, int32_t* index_end); @@ -57,8 +57,8 @@ extern "C" { int openmc_material_set_densities(int32_t index, int n, const char** name, const double* density); int openmc_material_set_id(int32_t index, int32_t id); int openmc_material_set_volume(int32_t index, double volume); - int openmc_material_filter_get_bins(int32_t index, int32_t** bins, int32_t* n); - int openmc_material_filter_set_bins(int32_t index, int32_t n, const int32_t* bins); + int openmc_material_filter_get_bins(int32_t index, int32_t** bins, size_t* n); + int openmc_material_filter_set_bins(int32_t index, size_t n, const int32_t* bins); int openmc_mesh_filter_get_mesh(int32_t index, int32_t* index_mesh); int openmc_mesh_filter_set_mesh(int32_t index, int32_t index_mesh); int openmc_mesh_get_id(int32_t index, int32_t* id); diff --git a/include/openmc/tallies/filter.h b/include/openmc/tallies/filter.h index c2fd4dfb37..2c0dc2aa2e 100644 --- a/include/openmc/tallies/filter.h +++ b/include/openmc/tallies/filter.h @@ -7,6 +7,8 @@ #include #include +#include + #include "openmc/hdf5_interface.h" #include "openmc/particle.h" #include "pugixml.hpp" @@ -43,6 +45,23 @@ namespace openmc { class Filter { public: + // Default constructor + Filter(); + + //! Create a new tally filter + // + //! \param[in] type Type of the filter + //! \param[in] id Unique ID for the filter. If none is passed, an ID is + //! automatically assigned + //! \return Pointer to the new filter object + static Filter* create(const std::string& type, int32_t id = -1); + + //! Create a new tally filter from an XML node + // + //! \param[in] node XML node + //! \return Pointer to the new filter object + static Filter* create(pugi::xml_node node); + virtual ~Filter() = default; virtual std::string type() const = 0; @@ -57,6 +76,13 @@ public: virtual void get_all_bins(const Particle* p, int estimator, FilterMatch& match) const = 0; + //! Assign a unique ID to the filter + // + //! \param[in] Unique ID to assign + void set_id(int32_t id); + + gsl::index index() const { return index_; } + //! Writes data describing this filter to an HDF5 statepoint group. virtual void to_statepoint(hid_t filter_group) const @@ -71,11 +97,11 @@ public: //! "Incoming Energy [0.625E-6, 20.0)". virtual std::string text_label(int bin) const = 0; - virtual void initialize() {} - - int32_t id_; + int32_t id_ {-1}; int n_bins_; +private: + gsl::index index_; }; //============================================================================== @@ -99,8 +125,6 @@ namespace model { // Non-member functions //============================================================================== -Filter* allocate_filter(const std::string& type); - //! Make sure index corresponds to a valid filter int verify_filter(int32_t index); diff --git a/include/openmc/tallies/filter_azimuthal.h b/include/openmc/tallies/filter_azimuthal.h index b896400315..f753bb2806 100644 --- a/include/openmc/tallies/filter_azimuthal.h +++ b/include/openmc/tallies/filter_azimuthal.h @@ -4,6 +4,8 @@ #include #include +#include + #include "openmc/tallies/filter.h" namespace openmc { @@ -24,6 +26,8 @@ public: void get_all_bins(const Particle* p, int estimator, FilterMatch& match) const override; + void set_bins(gsl::span bins); + void to_statepoint(hid_t filter_group) const override; std::string text_label(int bin) const override; diff --git a/include/openmc/tallies/filter_cell.h b/include/openmc/tallies/filter_cell.h index 9c654e82ed..eb3ff52367 100644 --- a/include/openmc/tallies/filter_cell.h +++ b/include/openmc/tallies/filter_cell.h @@ -5,6 +5,8 @@ #include #include +#include + #include "openmc/tallies/filter.h" namespace openmc { @@ -22,11 +24,11 @@ public: void from_xml(pugi::xml_node node) override; - void initialize() override; - void get_all_bins(const Particle* p, int estimator, FilterMatch& match) const override; + void set_cells(gsl::span cells); + void to_statepoint(hid_t filter_group) const override; std::string text_label(int bin) const override; diff --git a/include/openmc/tallies/filter_delayedgroup.h b/include/openmc/tallies/filter_delayedgroup.h index 2f65a50593..ae2ad0cb44 100644 --- a/include/openmc/tallies/filter_delayedgroup.h +++ b/include/openmc/tallies/filter_delayedgroup.h @@ -3,6 +3,8 @@ #include +#include + #include "openmc/tallies/filter.h" namespace openmc { @@ -26,6 +28,8 @@ public: void get_all_bins(const Particle* p, int estimator, FilterMatch& match) const override; + void set_groups(gsl::span groups); + void to_statepoint(hid_t filter_group) const override; std::string text_label(int bin) const override; diff --git a/include/openmc/tallies/filter_distribcell.h b/include/openmc/tallies/filter_distribcell.h index 73857e7cb0..db4cd4b1a3 100644 --- a/include/openmc/tallies/filter_distribcell.h +++ b/include/openmc/tallies/filter_distribcell.h @@ -20,11 +20,11 @@ public: void from_xml(pugi::xml_node node) override; - void initialize() override; - void get_all_bins(const Particle* p, int estimator, FilterMatch& match) const override; + void set_cell(int32_t cell); + void to_statepoint(hid_t filter_group) const override; std::string text_label(int bin) const override; diff --git a/include/openmc/tallies/filter_energy.h b/include/openmc/tallies/filter_energy.h index 025a77c624..b2123d7e7a 100644 --- a/include/openmc/tallies/filter_energy.h +++ b/include/openmc/tallies/filter_energy.h @@ -3,6 +3,8 @@ #include +#include + #include "openmc/tallies/filter.h" namespace openmc { @@ -23,6 +25,8 @@ public: void get_all_bins(const Particle* p, int estimator, FilterMatch& match) const override; + void set_bins(gsl::span bins); + void to_statepoint(hid_t filter_group) const override; std::string text_label(int bin) const override; diff --git a/include/openmc/tallies/filter_legendre.h b/include/openmc/tallies/filter_legendre.h index 054ba14e79..25cf99e156 100644 --- a/include/openmc/tallies/filter_legendre.h +++ b/include/openmc/tallies/filter_legendre.h @@ -23,6 +23,8 @@ public: void get_all_bins(const Particle* p, int estimator, FilterMatch& match) const override; + void set_order(int order); + void to_statepoint(hid_t filter_group) const override; std::string text_label(int bin) const override; diff --git a/include/openmc/tallies/filter_material.h b/include/openmc/tallies/filter_material.h index d272451db5..34dffa2b83 100644 --- a/include/openmc/tallies/filter_material.h +++ b/include/openmc/tallies/filter_material.h @@ -5,6 +5,8 @@ #include #include +#include + #include "openmc/tallies/filter.h" namespace openmc { @@ -22,11 +24,11 @@ public: void from_xml(pugi::xml_node node) override; - void initialize() override; - void get_all_bins(const Particle* p, int estimator, FilterMatch& match) const override; + void set_materials(gsl::span materials); + void to_statepoint(hid_t filter_group) const override; std::string text_label(int bin) const override; diff --git a/include/openmc/tallies/filter_mu.h b/include/openmc/tallies/filter_mu.h index 68b2f7b023..06044e760b 100644 --- a/include/openmc/tallies/filter_mu.h +++ b/include/openmc/tallies/filter_mu.h @@ -3,6 +3,8 @@ #include +#include + #include "openmc/tallies/filter.h" namespace openmc { @@ -24,6 +26,8 @@ public: void get_all_bins(const Particle* p, int estimator, FilterMatch& match) const override; + void set_bins(gsl::span bins); + void to_statepoint(hid_t filter_group) const override; std::string text_label(int bin) const override; diff --git a/include/openmc/tallies/filter_particle.h b/include/openmc/tallies/filter_particle.h index 268c19c6de..05e3d0505e 100644 --- a/include/openmc/tallies/filter_particle.h +++ b/include/openmc/tallies/filter_particle.h @@ -24,6 +24,8 @@ public: void get_all_bins(const Particle* p, int estimator, FilterMatch& match) const override; + void set_particles(gsl::span particles); + void to_statepoint(hid_t filter_group) const override; std::string text_label(int bin) const override; diff --git a/include/openmc/tallies/filter_polar.h b/include/openmc/tallies/filter_polar.h index 2bcdcd07f4..93488e5b09 100644 --- a/include/openmc/tallies/filter_polar.h +++ b/include/openmc/tallies/filter_polar.h @@ -4,6 +4,8 @@ #include #include +#include + #include "openmc/tallies/filter.h" namespace openmc { @@ -24,6 +26,8 @@ public: void get_all_bins(const Particle* p, int estimator, FilterMatch& match) const override; + void set_bins(gsl::span bins); + void to_statepoint(hid_t filter_group) const override; std::string text_label(int bin) const override; diff --git a/include/openmc/tallies/filter_sph_harm.h b/include/openmc/tallies/filter_sph_harm.h index b53b624417..522e89e786 100644 --- a/include/openmc/tallies/filter_sph_harm.h +++ b/include/openmc/tallies/filter_sph_harm.h @@ -3,6 +3,8 @@ #include +#include + #include "openmc/tallies/filter.h" namespace openmc { @@ -27,6 +29,10 @@ public: void get_all_bins(const Particle* p, int estimator, FilterMatch& match) const override; + void set_order(int order); + + void set_cosine(gsl::cstring_span cosine); + void to_statepoint(hid_t filter_group) const override; std::string text_label(int bin) const override; diff --git a/include/openmc/tallies/filter_sptl_legendre.h b/include/openmc/tallies/filter_sptl_legendre.h index 995bc73608..c65bb0cfc4 100644 --- a/include/openmc/tallies/filter_sptl_legendre.h +++ b/include/openmc/tallies/filter_sptl_legendre.h @@ -27,6 +27,12 @@ public: void get_all_bins(const Particle* p, int estimator, FilterMatch& match) const override; + void set_order(int order); + + void set_axis(LegendreAxis axis); + + void set_minmax(double min, double max); + void to_statepoint(hid_t filter_group) const override; std::string text_label(int bin) const override; diff --git a/include/openmc/tallies/filter_surface.h b/include/openmc/tallies/filter_surface.h index ea4dc5b990..61b320ad63 100644 --- a/include/openmc/tallies/filter_surface.h +++ b/include/openmc/tallies/filter_surface.h @@ -5,6 +5,8 @@ #include #include +#include + #include "openmc/tallies/filter.h" namespace openmc { @@ -22,11 +24,11 @@ public: void from_xml(pugi::xml_node node) override; - void initialize() override; - void get_all_bins(const Particle* p, int estimator, FilterMatch& match) const override; + void set_surfaces(gsl::span surfaces); + void to_statepoint(hid_t filter_group) const override; std::string text_label(int bin) const override; diff --git a/include/openmc/tallies/filter_universe.h b/include/openmc/tallies/filter_universe.h index 3ae0093cbd..0b43bb85cc 100644 --- a/include/openmc/tallies/filter_universe.h +++ b/include/openmc/tallies/filter_universe.h @@ -5,6 +5,8 @@ #include #include +#include + #include "openmc/tallies/filter.h" namespace openmc { @@ -22,11 +24,11 @@ public: void from_xml(pugi::xml_node node) override; - void initialize() override; - void get_all_bins(const Particle* p, int estimator, FilterMatch& match) const override; + void set_universes(gsl::span universes); + void to_statepoint(hid_t filter_group) const override; std::string text_label(int bin) const override; diff --git a/openmc/capi/filter.py b/openmc/capi/filter.py index c0716bcef6..404cd3898e 100644 --- a/openmc/capi/filter.py +++ b/openmc/capi/filter.py @@ -28,10 +28,10 @@ _dll.openmc_cell_filter_get_bins.argtypes = [ _dll.openmc_cell_filter_get_bins.restype = c_int _dll.openmc_cell_filter_get_bins.errcheck = _error_handler _dll.openmc_energy_filter_get_bins.argtypes = [ - c_int32, POINTER(POINTER(c_double)), POINTER(c_int32)] + c_int32, POINTER(POINTER(c_double)), POINTER(c_size_t)] _dll.openmc_energy_filter_get_bins.restype = c_int _dll.openmc_energy_filter_get_bins.errcheck = _error_handler -_dll.openmc_energy_filter_set_bins.argtypes = [c_int32, c_int32, POINTER(c_double)] +_dll.openmc_energy_filter_set_bins.argtypes = [c_int32, c_size_t, POINTER(c_double)] _dll.openmc_energy_filter_set_bins.restype = c_int _dll.openmc_energy_filter_set_bins.errcheck = _error_handler _dll.openmc_filter_get_id.argtypes = [c_int32, POINTER(c_int32)] @@ -53,10 +53,10 @@ _dll.openmc_legendre_filter_set_order.argtypes = [c_int32, c_int] _dll.openmc_legendre_filter_set_order.restype = c_int _dll.openmc_legendre_filter_set_order.errcheck = _error_handler _dll.openmc_material_filter_get_bins.argtypes = [ - c_int32, POINTER(POINTER(c_int32)), POINTER(c_int32)] + c_int32, POINTER(POINTER(c_int32)), POINTER(c_size_t)] _dll.openmc_material_filter_get_bins.restype = c_int _dll.openmc_material_filter_get_bins.errcheck = _error_handler -_dll.openmc_material_filter_set_bins.argtypes = [c_int32, c_int32, POINTER(c_int32)] +_dll.openmc_material_filter_set_bins.argtypes = [c_int32, c_size_t, POINTER(c_int32)] _dll.openmc_material_filter_set_bins.restype = c_int _dll.openmc_material_filter_set_bins.errcheck = _error_handler _dll.openmc_mesh_filter_get_mesh.argtypes = [c_int32, POINTER(c_int32)] @@ -148,7 +148,7 @@ class EnergyFilter(Filter): @property def bins(self): energies = POINTER(c_double)() - n = c_int32() + n = c_size_t() _dll.openmc_energy_filter_get_bins(self._index, energies, n) return as_array(energies, (n.value,)) @@ -231,7 +231,7 @@ class MaterialFilter(Filter): @property def bins(self): materials = POINTER(c_int32)() - n = c_int32() + n = c_size_t() _dll.openmc_material_filter_get_bins(self._index, materials, n) return [Material(index=materials[i]) for i in range(n.value)] diff --git a/src/tallies/filter.cpp b/src/tallies/filter.cpp index e9aae53a7a..36a36044ab 100644 --- a/src/tallies/filter.cpp +++ b/src/tallies/filter.cpp @@ -7,6 +7,7 @@ #include "openmc/capi.h" #include "openmc/constants.h" // for MAX_LINE_LEN; #include "openmc/error.h" +#include "openmc/xml_interface.h" #include "openmc/tallies/filter_azimuthal.h" #include "openmc/tallies/filter_cell.h" #include "openmc/tallies/filter_cellborn.h" @@ -50,8 +51,41 @@ namespace model { // Non-member functions //============================================================================== -Filter* -allocate_filter(const std::string& type) +extern "C" size_t tally_filters_size() +{ + return model::tally_filters.size(); +} + +//============================================================================== +// Filter implementation +//============================================================================== + +Filter::Filter() : index_{model::tally_filters.size()} +{ } + +Filter* Filter::create(pugi::xml_node node) +{ + // Copy filter id + if (!check_for_node(node, "id")) { + fatal_error("Must specify id for filter in tally XML file."); + } + int filter_id = std::stoi(get_node_value(node, "id")); + + // Convert filter type to lower case + std::string s; + if (check_for_node(node, "type")) { + s = get_node_value(node, "type", true); + } + + // Allocate according to the filter type + auto f = Filter::create(s, filter_id); + + // Read filter data from XML + f->from_xml(node); + return f; +} + +Filter* Filter::create(const std::string& type, int32_t id) { if (type == "azimuthal") { model::tally_filters.push_back(std::make_unique()); @@ -100,12 +134,26 @@ allocate_filter(const std::string& type) } else { throw std::runtime_error{"Unknown filter type: " + type}; } + + // Assign ID + model::tally_filters.back()->set_id(id); + return model::tally_filters.back().get(); } -extern "C" size_t tally_filters_size() +void Filter::set_id(int32_t id) { - return model::tally_filters.size(); + Expects(id >= 0); + if (model::filter_map.find(id) != model::filter_map.end()) { + throw std::runtime_error{"Two filters have the same ID: " + std::to_string(id)}; + } + + // Clear entry in filter map if an ID was already assigned before + if (id_ != -1) model::filter_map.erase(id_); + + // Update ID and entry in filter map + id_ = id; + model::filter_map[id] = index_; } //============================================================================== @@ -181,7 +229,7 @@ extern "C" int openmc_new_filter(const char* type, int32_t* index) { *index = model::tally_filters.size(); - allocate_filter(type); + Filter::create(type); return 0; } diff --git a/src/tallies/filter_azimuthal.cpp b/src/tallies/filter_azimuthal.cpp index 1a33a6bfe3..eb02319c20 100644 --- a/src/tallies/filter_azimuthal.cpp +++ b/src/tallies/filter_azimuthal.cpp @@ -15,22 +15,35 @@ AzimuthalFilter::from_xml(pugi::xml_node node) { auto bins = get_node_array(node, "bins"); - if (bins.size() > 1) { - bins_ = bins; - - } else { + if (bins.size() == 1) { // Allow a user to input a lone number which will mean that you subdivide // [-pi,pi) evenly with the input being the number of bins int n_angle = bins[0]; - - if (n_angle <= 1) fatal_error("Number of bins for azimuthal filter must " - "be greater than 1."); + if (n_angle <= 1) throw std::runtime_error{ + "Number of bins for azimuthal filter must be greater than 1."}; double d_angle = 2.0 * PI / n_angle; - bins_.resize(n_angle + 1); - for (int i = 0; i < n_angle; i++) bins_[i] = -PI + i * d_angle; - bins_[n_angle] = PI; + bins.resize(n_angle + 1); + for (int i = 0; i < n_angle; i++) bins[i] = -PI + i * d_angle; + bins[n_angle] = PI; + } + + this->set_bins(bins); +} + +void AzimuthalFilter::set_bins(gsl::span bins) +{ + // Clear existing bins + bins_.clear(); + bins_.reserve(bins.size()); + + // Copy bins, ensuring they are valid + for (gsl::index i = 0; i < bins.size(); ++i) { + if (i > 0 && bins[i] <= bins[i-1]) { + throw std::runtime_error{"Azimuthal bins must be monotonically increasing."}; + } + bins_.push_back(bins[i]); } n_bins_ = bins_.size() - 1; diff --git a/src/tallies/filter_cell.cpp b/src/tallies/filter_cell.cpp index 8286a69760..685d8ff745 100644 --- a/src/tallies/filter_cell.cpp +++ b/src/tallies/filter_cell.cpp @@ -12,29 +12,39 @@ namespace openmc { void CellFilter::from_xml(pugi::xml_node node) { - cells_ = get_node_array(node, "bins"); - n_bins_ = cells_.size(); + // Get cell IDs and convert into indices into the global cells vector + auto cells = get_node_array(node, "bins"); + for (auto& c : cells) { + auto search = model::cell_map.find(c); + if (search == model::cell_map.end()) { + std::stringstream err_msg; + err_msg << "Could not find cell " << c + << " specified on tally filter."; + throw std::runtime_error{err_msg.str()}; + } + c = search->second; + } + + this->set_cells(cells); } void -CellFilter::initialize() +CellFilter::set_cells(gsl::span cells) { - // Convert cell IDs to indices of the global array. - for (auto& c : cells_) { - auto search = model::cell_map.find(c); - if (search != model::cell_map.end()) { - c = search->second; - } else { - std::stringstream err_msg; - err_msg << "Could not find cell " << c << " specified on tally filter."; - fatal_error(err_msg); - } + // Clear existing cells + cells_.clear(); + cells_.reserve(cells.size()); + map_.clear(); + + // Update cells and mapping + for (auto& index : cells) { + Expects(index >= 0); + Expects(index < model::cells.size()); + cells_.push_back(index); + map_[index] = cells_.size() - 1; } - // Populate the index->bin map. - for (int i = 0; i < cells_.size(); i++) { - map_[cells_[i]] = i; - } + n_bins_ = cells_.size(); } void diff --git a/src/tallies/filter_delayedgroup.cpp b/src/tallies/filter_delayedgroup.cpp index 4f448358e3..2b69783081 100644 --- a/src/tallies/filter_delayedgroup.cpp +++ b/src/tallies/filter_delayedgroup.cpp @@ -8,21 +8,32 @@ namespace openmc { void DelayedGroupFilter::from_xml(pugi::xml_node node) { - groups_ = get_node_array(node, "bins"); - n_bins_ = groups_.size(); + auto groups = get_node_array(node, "bins"); + this->set_groups(groups); +} + +void +DelayedGroupFilter::set_groups(gsl::span groups) +{ + // Clear existing groups + groups_.clear(); + groups_.reserve(groups.size()); // Make sure all the group index values are valid. // TODO: do these need to be decremented for zero-based indexing? - for (auto group : groups_) { + for (auto group : groups) { if (group < 1) { - fatal_error("Encountered delayedgroup bin with index " - + std::to_string(group) + " which is less than 1"); + throw std::invalid_argument{"Encountered delayedgroup bin with index " + + std::to_string(group) + " which is less than 1"}; } else if (group > MAX_DELAYED_GROUPS) { - fatal_error("Encountered delayedgroup bin with index " + throw std::invalid_argument{"Encountered delayedgroup bin with index " + std::to_string(group) + " which is greater than MAX_DELATED_GROUPS (" - + std::to_string(MAX_DELAYED_GROUPS) + ")"); + + std::to_string(MAX_DELAYED_GROUPS) + ")"}; } + groups_.push_back(group); } + + n_bins_ = groups_.size(); } void diff --git a/src/tallies/filter_distribcell.cpp b/src/tallies/filter_distribcell.cpp index 94f23621b1..f9c7ab4d66 100644 --- a/src/tallies/filter_distribcell.cpp +++ b/src/tallies/filter_distribcell.cpp @@ -15,23 +15,26 @@ DistribcellFilter::from_xml(pugi::xml_node node) if (cells.size() != 1) { fatal_error("Only one cell can be specified per distribcell filter."); } - cell_ = cells[0]; -} -void -DistribcellFilter::initialize() -{ - // Convert the cell ID to an index of the global array. - auto search = model::cell_map.find(cell_); - if (search != model::cell_map.end()) { - cell_ = search->second; - n_bins_ = model::cells[cell_]->n_instances_; - } else { + // Find index in global cells vector corresponding to cell ID + auto search = model::cell_map.find(cells[0]); + if (search == model::cell_map.end()) { std::stringstream err_msg; err_msg << "Could not find cell " << cell_ << " specified on tally filter."; - fatal_error(err_msg); + throw std::runtime_error{err_msg.str()}; } + + this->set_cell(search->second); +} + +void +DistribcellFilter::set_cell(int32_t cell) +{ + Expects(cell >= 0); + Expects(cell < model::cells.size()); + cell_ = cell; + n_bins_ = model::cells[cell]->n_instances_; } void diff --git a/src/tallies/filter_energy.cpp b/src/tallies/filter_energy.cpp index e11ab30803..d4e26a1f0d 100644 --- a/src/tallies/filter_energy.cpp +++ b/src/tallies/filter_energy.cpp @@ -16,7 +16,25 @@ namespace openmc { void EnergyFilter::from_xml(pugi::xml_node node) { - bins_ = get_node_array(node, "bins"); + auto bins = get_node_array(node, "bins"); + this->set_bins(bins); +} + +void +EnergyFilter::set_bins(gsl::span bins) +{ + // Clear existing bins + bins_.clear(); + bins_.reserve(bins.size()); + + // Copy bins, ensuring they are valid + for (gsl::index i = 0; i < bins.size(); ++i) { + if (i > 0 && bins[i] <= bins[i-1]) { + throw std::runtime_error{"Energy bins must be monotonically increasing."}; + } + bins_.push_back(bins[i]); + } + n_bins_ = bins_.size() - 1; // In MG mode, check if the filter bins match the transport bins. @@ -27,7 +45,7 @@ EnergyFilter::from_xml(pugi::xml_node node) if (!settings::run_CE) { if (n_bins_ == data::num_energy_groups) { matches_transport_groups_ = true; - for (auto i = 0; i < n_bins_ + 1; i++) { + for (gsl::index i = 0; i < n_bins_ + 1; ++i) { if (data::rev_energy_bins[i] != bins_[i]) { matches_transport_groups_ = false; break; @@ -111,7 +129,7 @@ EnergyoutFilter::text_label(int bin) const //============================================================================== extern"C" int -openmc_energy_filter_get_bins(int32_t index, double** energies, int32_t* n) +openmc_energy_filter_get_bins(int32_t index, double** energies, size_t* n) { // Make sure this is a valid index to an allocated filter. if (int err = verify_filter(index)) return err; @@ -133,7 +151,7 @@ openmc_energy_filter_get_bins(int32_t index, double** energies, int32_t* n) } extern "C" int -openmc_energy_filter_set_bins(int32_t index, int32_t n, const double* energies) +openmc_energy_filter_set_bins(int32_t index, size_t n, const double* energies) { // Make sure this is a valid index to an allocated filter. if (int err = verify_filter(index)) return err; @@ -149,10 +167,7 @@ openmc_energy_filter_set_bins(int32_t index, int32_t n, const double* energies) } // Update the filter. - filt->bins_.clear(); - filt->bins_.resize(n); - for (int i = 0; i < n; i++) filt->bins_[i] = energies[i]; - filt->n_bins_ = n - 1; + filt->set_bins({energies, n}); return 0; } diff --git a/src/tallies/filter_legendre.cpp b/src/tallies/filter_legendre.cpp index 6fa041c76b..21939b7d45 100644 --- a/src/tallies/filter_legendre.cpp +++ b/src/tallies/filter_legendre.cpp @@ -10,7 +10,13 @@ namespace openmc { void LegendreFilter::from_xml(pugi::xml_node node) { - order_ = std::stoi(get_node_value(node, "order")); + this->set_order(std::stoi(get_node_value(node, "order"))); +} + +void +LegendreFilter::set_order(int order) +{ + order_ = order; n_bins_ = order_ + 1; } @@ -81,8 +87,7 @@ openmc_legendre_filter_set_order(int32_t index, int order) } // Update the filter. - filt->order_ = order; - filt->n_bins_ = order + 1; + filt->set_order(order); return 0; } diff --git a/src/tallies/filter_material.cpp b/src/tallies/filter_material.cpp index fd9e637bde..c9fd644adc 100644 --- a/src/tallies/filter_material.cpp +++ b/src/tallies/filter_material.cpp @@ -3,7 +3,6 @@ #include #include "openmc/capi.h" -#include "openmc/error.h" #include "openmc/material.h" #include "openmc/xml_interface.h" @@ -12,30 +11,39 @@ namespace openmc { void MaterialFilter::from_xml(pugi::xml_node node) { - materials_ = get_node_array(node, "bins"); - n_bins_ = materials_.size(); -} - -void -MaterialFilter::initialize() -{ - // Convert material IDs to indices of the global array. - for (auto& m : materials_) { + // Get material IDs and convert to indices in the global materials vector + auto mats = get_node_array(node, "bins"); + for (auto& m : mats) { auto search = model::material_map.find(m); - if (search != model::material_map.end()) { - m = search->second; - } else { + if (search == model::material_map.end()) { std::stringstream err_msg; err_msg << "Could not find material " << m << " specified on tally filter."; - fatal_error(err_msg); + throw std::runtime_error{err_msg.str()}; } + m = search->second; } - // Populate the index->bin map. - for (int i = 0; i < materials_.size(); i++) { - map_[materials_[i]] = i; + this->set_materials(mats); +} + +void +MaterialFilter::set_materials(gsl::span materials) +{ + // Clear existing materials + materials_.clear(); + materials_.reserve(materials.size()); + map_.clear(); + + // Update materials and mapping + for (auto& index : materials) { + Expects(index >= 0); + Expects(index < model::materials.size()); + materials_.push_back(index); + map_[index] = materials_.size() - 1; } + + n_bins_ = materials_.size(); } void @@ -69,7 +77,7 @@ MaterialFilter::text_label(int bin) const //============================================================================== extern "C" int -openmc_material_filter_get_bins(int32_t index, int32_t** bins, int32_t* n) +openmc_material_filter_get_bins(int32_t index, int32_t** bins, size_t* n) { // Make sure this is a valid index to an allocated filter. if (int err = verify_filter(index)) return err; @@ -91,7 +99,7 @@ openmc_material_filter_get_bins(int32_t index, int32_t** bins, int32_t* n) } extern "C" int -openmc_material_filter_set_bins(int32_t index, int32_t n, const int32_t* bins) +openmc_material_filter_set_bins(int32_t index, size_t n, const int32_t* bins) { // Make sure this is a valid index to an allocated filter. if (int err = verify_filter(index)) return err; @@ -107,12 +115,7 @@ openmc_material_filter_set_bins(int32_t index, int32_t n, const int32_t* bins) } // Update the filter. - filt->materials_.clear(); - filt->materials_.resize(n); - for (int i = 0; i < n; i++) filt->materials_[i] = bins[i]; - filt->n_bins_ = filt->materials_.size(); - filt->map_.clear(); - for (int i = 0; i < n; i++) filt->map_[filt->materials_[i]] = i; + filt->set_materials({bins, n}); return 0; } diff --git a/src/tallies/filter_mu.cpp b/src/tallies/filter_mu.cpp index aed8371b80..9aaebea8cf 100644 --- a/src/tallies/filter_mu.cpp +++ b/src/tallies/filter_mu.cpp @@ -13,22 +13,36 @@ MuFilter::from_xml(pugi::xml_node node) { auto bins = get_node_array(node, "bins"); - if (bins.size() > 1) { - bins_ = bins; - - } else { + if (bins.size() == 1) { // Allow a user to input a lone number which will mean that you subdivide // [-1,1) evenly with the input being the number of bins int n_angle = bins[0]; - - if (n_angle <= 1) fatal_error("Number of bins for mu filter must " - "be greater than 1."); + if (n_angle <= 1) throw std::runtime_error{ + "Number of bins for mu filter must be greater than 1."}; double d_angle = 2.0 / n_angle; - bins_.resize(n_angle + 1); - for (int i = 0; i < n_angle; i++) bins_[i] = -1 + i * d_angle; - bins_[n_angle] = 1; + bins.resize(n_angle + 1); + for (int i = 0; i < n_angle; i++) bins[i] = -1 + i * d_angle; + bins[n_angle] = 1; + } + + this->set_bins(bins); +} + +void +MuFilter::set_bins(gsl::span bins) +{ + // Clear existing bins + bins_.clear(); + bins_.reserve(bins.size()); + + // Copy bins, ensuring they are valid + for (gsl::index i = 0; i < bins.size(); ++i) { + if (i > 0 && bins[i] <= bins[i-1]) { + throw std::runtime_error{"Mu bins must be monotonically increasing."}; + } + bins_.push_back(bins[i]); } n_bins_ = bins_.size() - 1; diff --git a/src/tallies/filter_particle.cpp b/src/tallies/filter_particle.cpp index 23a305b8dd..de142489ff 100644 --- a/src/tallies/filter_particle.cpp +++ b/src/tallies/filter_particle.cpp @@ -8,8 +8,25 @@ void ParticleFilter::from_xml(pugi::xml_node node) { auto particles = get_node_array(node, "bins"); + + // Convert to vector of Particle::Type + std::vector types; for (auto& p : particles) { - particles_.push_back(static_cast(p - 1)); + types.push_back(static_cast(p - 1)); + } + this->set_particles(types); +} + +void +ParticleFilter::set_particles(gsl::span particles) +{ + // Clear existing particles + particles_.clear(); + particles_.reserve(particles.size()); + + // Set particles and number of bins + for (auto p : particles) { + particles_.push_back(p); } n_bins_ = particles_.size(); } diff --git a/src/tallies/filter_polar.cpp b/src/tallies/filter_polar.cpp index 9730bbb584..edbbc89229 100644 --- a/src/tallies/filter_polar.cpp +++ b/src/tallies/filter_polar.cpp @@ -14,22 +14,36 @@ PolarFilter::from_xml(pugi::xml_node node) { auto bins = get_node_array(node, "bins"); - if (bins.size() > 1) { - bins_ = bins; - - } else { + if (bins.size() == 1) { // Allow a user to input a lone number which will mean that you subdivide // [0,pi] evenly with the input being the number of bins int n_angle = bins[0]; - - if (n_angle <= 1) fatal_error("Number of bins for polar filter must " - "be greater than 1."); + if (n_angle <= 1) throw std::runtime_error{ + "Number of bins for polar filter must be greater than 1."}; double d_angle = PI / n_angle; - bins_.resize(n_angle + 1); - for (int i = 0; i < n_angle; i++) bins_[i] = i * d_angle; - bins_[n_angle] = PI; + bins.resize(n_angle + 1); + for (int i = 0; i < n_angle; i++) bins[i] = i * d_angle; + bins[n_angle] = PI; + } + + this->set_bins(bins); +} + +void +PolarFilter::set_bins(gsl::span bins) +{ + // Clear existing bins + bins_.clear(); + bins_.reserve(bins.size()); + + // Copy bins, ensuring they are valid + for (gsl::index i = 0; i < bins.size(); ++i) { + if (i > 0 && bins[i] <= bins[i-1]) { + throw std::runtime_error{"Polar bins must be monotonically increasing."}; + } + bins_.push_back(bins[i]); } n_bins_ = bins_.size() - 1; diff --git a/src/tallies/filter_sph_harm.cpp b/src/tallies/filter_sph_harm.cpp index 0b3a371187..a01c2e87f2 100644 --- a/src/tallies/filter_sph_harm.cpp +++ b/src/tallies/filter_sph_harm.cpp @@ -12,21 +12,34 @@ namespace openmc { void SphericalHarmonicsFilter::from_xml(pugi::xml_node node) { - order_ = std::stoi(get_node_value(node, "order")); - n_bins_ = (order_ + 1) * (order_ + 1); - + this->set_order(std::stoi(get_node_value(node, "order"))); if (check_for_node(node, "cosine")) { - auto cos = get_node_value(node, "cosine", true); - if (cos == "scatter") { - cosine_ = SphericalHarmonicsCosine::scatter; - } else if (cos == "particle") { - cosine_ = SphericalHarmonicsCosine::particle; - } else { - std::stringstream err_msg; - err_msg << "Unrecognized cosine type, \"" << cos - << "\" in spherical harmonics filter"; - fatal_error(err_msg); - } + this->set_cosine(get_node_value(node, "cosine", true)); + } +} + +void +SphericalHarmonicsFilter::set_order(int order) +{ + if (order < 0) { + throw std::invalid_argument{"Spherical harmonics order must be non-negative."}; + } + order_ = order; + n_bins_ = (order_ + 1) * (order_ + 1); +} + +void +SphericalHarmonicsFilter::set_cosine(gsl::cstring_span cosine) +{ + if (cosine == "scatter") { + cosine_ = SphericalHarmonicsCosine::scatter; + } else if (cosine == "particle") { + cosine_ = SphericalHarmonicsCosine::particle; + } else { + std::stringstream err_msg; + err_msg << "Unrecognized cosine type, \"" << cos + << "\" in spherical harmonics filter"; + throw std::invalid_argument{err_msg.str()}; } } @@ -153,8 +166,7 @@ openmc_sphharm_filter_set_order(int32_t index, int order) if (err) return err; // Update the filter. - filt->order_ = order; - filt->n_bins_ = (order + 1) * (order + 1); + filt->set_order(order); return 0; } @@ -168,12 +180,10 @@ openmc_sphharm_filter_set_cosine(int32_t index, const char cosine[]) if (err) return err; // Update the filter. - if (strcmp(cosine, "scatter") == 0) { - filt->cosine_ = SphericalHarmonicsCosine::scatter; - } else if (strcmp(cosine, "particle") == 0) { - filt->cosine_ = SphericalHarmonicsCosine::particle; - } else { - set_errmsg("Invalid spherical harmonics cosine."); + try { + filt->set_cosine(cosine); + } catch (const std::invalid_argument& e) { + set_errmsg(e.what()); return OPENMC_E_INVALID_ARGUMENT; } return 0; diff --git a/src/tallies/filter_sptl_legendre.cpp b/src/tallies/filter_sptl_legendre.cpp index 6562ff01c7..e0aeff9c17 100644 --- a/src/tallies/filter_sptl_legendre.cpp +++ b/src/tallies/filter_sptl_legendre.cpp @@ -12,25 +12,51 @@ namespace openmc { void SpatialLegendreFilter::from_xml(pugi::xml_node node) { - order_ = std::stoi(get_node_value(node, "order")); + this->set_order(std::stoi(get_node_value(node, "order"))); auto axis = get_node_value(node, "axis"); - if (axis == "x") { - axis_ = LegendreAxis::x; - } else if (axis == "y") { - axis_ = LegendreAxis::y; - } else if (axis == "z") { - axis_ = LegendreAxis::z; - } else { - fatal_error("Unrecognized axis on SpatialLegendreFilter"); + switch (axis[0]) { + case 'x': + this->set_axis(LegendreAxis::x); + break; + case 'y': + this->set_axis(LegendreAxis::y); + break; + case 'z': + this->set_axis(LegendreAxis::z); + break; + default: + throw std::runtime_error{"Axis for SpatialLegendreFilter must be 'x', 'y', or 'z'"}; } - min_ = std::stod(get_node_value(node, "min")); - max_ = std::stod(get_node_value(node, "max")); + double min = std::stod(get_node_value(node, "min")); + double max = std::stod(get_node_value(node, "max")); + this->set_minmax(min, max); +} +void +SpatialLegendreFilter::set_order(int order) +{ + order_ = order; n_bins_ = order_ + 1; } +void +SpatialLegendreFilter::set_axis(LegendreAxis axis) +{ + axis_ = axis; +} + +void +SpatialLegendreFilter::set_minmax(double min, double max) +{ + if (max < min) { + throw std::invalid_argument{"Maximum value must be greater than minimum value"}; + } + min_ = min; + max_ = max; +} + void SpatialLegendreFilter::get_all_bins(const Particle* p, int estimator, FilterMatch& match) const @@ -157,8 +183,7 @@ openmc_spatial_legendre_filter_set_order(int32_t index, int order) if (err) return err; // Update the filter. - filt->order_ = order; - filt->n_bins_ = order + 1; + filt->set_order(order); return 0; } @@ -173,7 +198,7 @@ openmc_spatial_legendre_filter_set_params(int32_t index, const int* axis, if (err) return err; // Update the filter. - if (axis) filt->axis_ = static_cast(*axis); + if (axis) filt->set_axis(static_cast(*axis)); if (min) filt->min_ = *min; if (max) filt->max_ = *max; return 0; diff --git a/src/tallies/filter_surface.cpp b/src/tallies/filter_surface.cpp index faa8aafbf7..72b356f216 100644 --- a/src/tallies/filter_surface.cpp +++ b/src/tallies/filter_surface.cpp @@ -11,30 +11,41 @@ namespace openmc { void SurfaceFilter::from_xml(pugi::xml_node node) { - surfaces_ = get_node_array(node, "bins"); - n_bins_ = surfaces_.size(); -} + auto surfaces = get_node_array(node, "bins"); -void -SurfaceFilter::initialize() -{ - // Convert surface IDs to indices of the global array. - for (auto& s : surfaces_) { + // Convert surface IDs to indices of the global surfaces vector. + for (auto& s : surfaces) { auto search = model::surface_map.find(s); - if (search != model::surface_map.end()) { - s = search->second; - } else { + if (search == model::surface_map.end()) { std::stringstream err_msg; err_msg << "Could not find surface " << s << " specified on tally filter."; - fatal_error(err_msg); + throw std::runtime_error{err_msg.str()}; } + + s = search->second; } - // Populate the index->bin map. - for (int i = 0; i < surfaces_.size(); i++) { - map_[surfaces_[i]] = i; + this->set_surfaces(surfaces); +} + +void +SurfaceFilter::set_surfaces(gsl::span surfaces) +{ + // Clear existing surfaces + surfaces_.clear(); + surfaces_.reserve(surfaces.size()); + map_.clear(); + + // Update surfaces and mapping + for (auto& index : surfaces) { + Expects(index >= 0); + Expects(index < model::surfaces.size()); + surfaces_.push_back(index); + map_[index] = surfaces_.size() - 1; } + + n_bins_ = surfaces_.size(); } void diff --git a/src/tallies/filter_universe.cpp b/src/tallies/filter_universe.cpp index 50c058c05b..dffdee621f 100644 --- a/src/tallies/filter_universe.cpp +++ b/src/tallies/filter_universe.cpp @@ -11,30 +11,39 @@ namespace openmc { void UniverseFilter::from_xml(pugi::xml_node node) { - universes_ = get_node_array(node, "bins"); - n_bins_ = universes_.size(); -} - -void -UniverseFilter::initialize() -{ - // Convert universe IDs to indices of the global array. - for (auto& u : universes_) { + // Get material IDs and convert to indices in the global materials vector + auto universes = get_node_array(node, "bins"); + for (auto& u : universes) { auto search = model::universe_map.find(u); - if (search != model::universe_map.end()) { - u = search->second; - } else { + if (search == model::universe_map.end()) { std::stringstream err_msg; err_msg << "Could not find universe " << u << " specified on tally filter."; - fatal_error(err_msg); + throw std::runtime_error{err_msg.str()}; } + u = search->second; } - // Populate the index->bin map. - for (int i = 0; i < universes_.size(); i++) { - map_[universes_[i]] = i; + this->set_universes(universes); +} + +void +UniverseFilter::set_universes(gsl::span universes) +{ + // Clear existing universes + universes_.clear(); + universes_.reserve(universes.size()); + map_.clear(); + + // Update universes and mapping + for (auto& index : universes) { + Expects(index >= 0); + Expects(index < model::universes.size()); + universes_.push_back(index); + map_[index] = universes_.size() - 1; } + + n_bins_ = universes_.size(); } void diff --git a/src/tallies/tally.cpp b/src/tallies/tally.cpp index 2cb0945ab3..a4cb041694 100644 --- a/src/tallies/tally.cpp +++ b/src/tallies/tally.cpp @@ -614,36 +614,7 @@ void read_tallies_xml() // Check for user filters and allocate for (auto node_filt : root.children("filter")) { - // Copy filter id - if (!check_for_node(node_filt, "id")) { - fatal_error("Must specify id for filter in tally XML file."); - } - int filter_id = std::stoi(get_node_value(node_filt, "id")); - - // Check to make sure 'id' hasn't been used - if (model::filter_map.find(filter_id) != model::filter_map.end()) { - fatal_error("Two or more filters use the same unique ID: " - + std::to_string(filter_id)); - } - - // Convert filter type to lower case - std::string s; - if (check_for_node(node_filt, "type")) { - s = get_node_value(node_filt, "type", true); - } - - // Allocate according to the filter type - Filter* f = allocate_filter(s); - - // Read filter data from XML - f->from_xml(node_filt); - - // Set filter id - f->id_ = filter_id; - model::filter_map[filter_id] = model::tally_filters.size() - 1; - - // Initialize filter - f->initialize(); + auto f = Filter::create(node_filt); } // ==========================================================================