diff --git a/docs/source/pythonapi/base.rst b/docs/source/pythonapi/base.rst index d0e7d2439b8..ff2c7e4ba0d 100644 --- a/docs/source/pythonapi/base.rst +++ b/docs/source/pythonapi/base.rst @@ -133,6 +133,7 @@ Constructing Tallies openmc.MeshSurfaceFilter openmc.EnergyFilter openmc.EnergyoutFilter + openmc.SecondaryEnergyFilter openmc.MuFilter openmc.MuSurfaceFilter openmc.PolarFilter diff --git a/include/openmc/constants.h b/include/openmc/constants.h index 252194528c7..c326786c09e 100644 --- a/include/openmc/constants.h +++ b/include/openmc/constants.h @@ -290,7 +290,7 @@ enum class TallyEstimator { ANALOG, TRACKLENGTH, COLLISION }; enum class TallyEvent { SURFACE, LATTICE, KILL, SCATTER, ABSORB }; // Tally score type -- if you change these, make sure you also update the -// _SCORES dictionary in openmc/capi/tally.py +// _SCORES dictionary in openmc/lib/tally.py // // These are kept as a normal enum and made negative, since variables which // store one of these enum values usually also may be responsible for storing diff --git a/include/openmc/particle_data.h b/include/openmc/particle_data.h index a07a99ffbed..9ce5e47b91a 100644 --- a/include/openmc/particle_data.h +++ b/include/openmc/particle_data.h @@ -466,6 +466,11 @@ class ParticleData : public GeometryState { vector secondary_bank_; + // Keeps track of how many secondary particles were created + // in the collision the particle most recently experienced. + // This is used to tally secondary photon energies. + int secondaries_this_collision_ {0}; + int64_t current_work_; vector flux_derivs_; @@ -615,7 +620,17 @@ class ParticleData : public GeometryState { // secondary particle bank SourceSite& secondary_bank(int i) { return secondary_bank_[i]; } + const SourceSite& secondary_bank(int i) const { return secondary_bank_[i]; } decltype(secondary_bank_)& secondary_bank() { return secondary_bank_; } + decltype(secondary_bank_) const& secondary_bank() const + { + return secondary_bank_; + } + const int& secondaries_this_collision() const + { + return secondaries_this_collision_; + } + int& secondaries_this_collision() { return secondaries_this_collision_; } // Current simulation work index int64_t& current_work() { return current_work_; } diff --git a/include/openmc/tallies/filter.h b/include/openmc/tallies/filter.h index 65098597a5d..6ec733516d9 100644 --- a/include/openmc/tallies/filter.h +++ b/include/openmc/tallies/filter.h @@ -28,6 +28,7 @@ enum class FilterType { ENERGY_FUNCTION, ENERGY, ENERGY_OUT, + SECONDARY_ENERGY, LEGENDRE, MATERIAL, MATERIALFROM, diff --git a/include/openmc/tallies/filter_energy.h b/include/openmc/tallies/filter_energy.h index cf8a8aa0f58..1e40d7e262c 100644 --- a/include/openmc/tallies/filter_energy.h +++ b/include/openmc/tallies/filter_energy.h @@ -1,6 +1,7 @@ #ifndef OPENMC_TALLIES_FILTER_ENERGY_H #define OPENMC_TALLIES_FILTER_ENERGY_H +#include "openmc/particle.h" #include "openmc/span.h" #include "openmc/tallies/filter.h" #include "openmc/vector.h" @@ -72,5 +73,38 @@ class EnergyoutFilter : public EnergyFilter { std::string text_label(int bin) const override; }; +//============================================================================== +//! Bins the outgoing energy of secondary particles +//! +//! This is used to get the photon production matrix for multigroup photon +//! transport. It could also be used to find the energy distribution of +//! neutron secondaries or others, for example. +//! +//! Using anything other than analog estimators here would be complicated +//============================================================================== + +class SecondaryEnergyFilter : public EnergyFilter { +public: + //---------------------------------------------------------------------------- + // Methods + + std::string type_str() const override { return "secondaryenergy"; } + FilterType type() const override { return FilterType::SECONDARY_ENERGY; } + + void get_all_bins(const Particle& p, TallyEstimator estimator, + FilterMatch& match) const override; + + std::string text_label(int bin) const override; + + void from_xml(pugi::xml_node node) override; + +protected: + // This filter could simultaneously use different particle types, for + // example creating one set of energy bins for photons and another for + // electrons. However, for typical use only photons are of interest. + // Covering multiple particles can be done by defining separate tallies. + ParticleType secondary_type_ {ParticleType::photon}; +}; + } // namespace openmc #endif // OPENMC_TALLIES_FILTER_ENERGY_H diff --git a/openmc/filter.py b/openmc/filter.py index f29c540857c..5bd6a0612bd 100644 --- a/openmc/filter.py +++ b/openmc/filter.py @@ -127,9 +127,9 @@ def __eq__(self, other): def __gt__(self, other): if type(self) is not type(other): if self.short_name in _FILTER_TYPES and \ - other.short_name in _FILTER_TYPES: + other.short_name in _FILTER_TYPES: delta = _FILTER_TYPES.index(self.short_name) - \ - _FILTER_TYPES.index(other.short_name) + _FILTER_TYPES.index(other.short_name) return delta > 0 else: return False @@ -278,7 +278,6 @@ def from_xml_element(cls, elem, **kwargs): if filter_type == subclass.short_name.lower(): return subclass.from_xml_element(elem, **kwargs) - def can_merge(self, other): """Determine if filter can be merged with another. @@ -429,6 +428,7 @@ def get_pandas_dataframe(self, data_size, stride, **kwargs): class WithIDFilter(Filter): """Abstract parent for filters of types with IDs (Cell, Material, etc.).""" + def __init__(self, bins, filter_id=None): bins = np.atleast_1d(bins) @@ -632,6 +632,7 @@ class CellInstanceFilter(Filter): DistribcellFilter """ + def __init__(self, bins, filter_id=None): self.bins = bins self.id = filter_id @@ -753,6 +754,7 @@ class ParticleFilter(Filter): The number of filter bins """ + def __eq__(self, other): if type(self) is not type(other): return False @@ -1059,6 +1061,7 @@ class MeshMaterialFilter(MeshFilter): Unique identifier for the filter """ + def __init__(self, mesh: openmc.MeshBase, bins, filter_id=None): self.mesh = mesh self.bins = bins @@ -1388,6 +1391,7 @@ class RealFilter(Filter): The number of filter bins """ + def __init__(self, values, filter_id=None): self.values = np.asarray(values) self.bins = np.vstack((self.values[:-1], self.values[1:])).T @@ -1663,7 +1667,8 @@ def from_group_structure(cls, group_structure): """ - cv.check_value('group_structure', group_structure, openmc.mgxs.GROUP_STRUCTURES.keys()) + cv.check_value('group_structure', group_structure, + openmc.mgxs.GROUP_STRUCTURES.keys()) return cls(openmc.mgxs.GROUP_STRUCTURES[group_structure.upper()]) @@ -1693,6 +1698,66 @@ class EnergyoutFilter(EnergyFilter): """ + +class SecondaryEnergyFilter(EnergyFilter): + """Bins tally events based on energy of secondary particles. + + .. versionadded:: 0.15.3 + + This filter collects the energies of secondary particles (e.g., photons + or electrons) produced in a reaction. This is useful for constructing + production matrices or analyzing secondary particle spectra. + + The primary particle type should be filtered by the usual ParticleFilter. + + Parameters + ---------- + values : Iterable of Real + A list of energy boundaries in [eV]; each successive pair defines a bin. + filter_id : int, optional + Unique identifier for the filter + particle : str + Type of secondary particle to tally ('photon', 'neutron', etc.) + + Attributes + ---------- + values : numpy.ndarray + Energy boundaries in [eV] + bins : numpy.ndarray + Array of (low, high) energy bin pairs + num_bins : int + Number of filter bins + particle : str + The secondary particle type this filter applies to + """ + + def __init__(self, values, filter_id=None, particle='photon'): + self.particle = particle + super().__init__(values, filter_id) + + @property + def particle(self): + return self._particle + + @particle.setter + def particle(self, particle): + cv.check_value('particle', particle, _PARTICLES) + self._particle = particle + + def to_xml_element(self): + element = super().to_xml_element() + subelement = ET.SubElement(element, 'particle') + subelement.text = self.particle + return element + + @classmethod + def from_xml_element(cls, elem, **kwargs): + filter_id = int(elem.get('id')) + values = [float(x) for x in get_text(elem, 'bins').split()] + particle = get_text(elem, 'particle') + return cls(values, filter_id=filter_id, particle=particle) + + class TimeFilter(RealFilter): """Bins tally events based on the particle's time. @@ -1851,7 +1916,7 @@ def bins(self, bins): # Make sure there is only 1 bin. if not len(bins) == 1: msg = (f'Unable to add bins "{bins}" to a DistribcellFilter since ' - 'only a single distribcell can be used per tally') + 'only a single distribcell can be used per tally') raise ValueError(msg) # Check the type and extract the id, if necessary. @@ -2004,7 +2069,7 @@ def get_pandas_dataframe(self, data_size, stride, **kwargs): # requests Summary geometric information filter_bins = _repeat_and_tile( np.arange(self.num_bins), stride, data_size) - df = pd.DataFrame({self.short_name.lower() : filter_bins}) + df = pd.DataFrame({self.short_name.lower(): filter_bins}) # Concatenate with DataFrame of distribcell instance IDs if level_df is not None: @@ -2043,6 +2108,7 @@ class MuFilter(RealFilter): The number of filter bins """ + def __init__(self, values, filter_id=None): if isinstance(values, Integral): values = np.linspace(-1., 1., values + 1) @@ -2208,6 +2274,7 @@ class DelayedGroupFilter(Filter): The number of filter bins """ + def check_bins(self, bins): # Check the bin values. for g in bins: @@ -2276,9 +2343,9 @@ def __eq__(self, other): def __gt__(self, other): if type(self) is not type(other): if self.short_name in _FILTER_TYPES and \ - other.short_name in _FILTER_TYPES: + other.short_name in _FILTER_TYPES: delta = _FILTER_TYPES.index(self.short_name) - \ - _FILTER_TYPES.index(other.short_name) + _FILTER_TYPES.index(other.short_name) return delta > 0 else: return False @@ -2288,9 +2355,9 @@ def __gt__(self, other): def __lt__(self, other): if type(self) is not type(other): if self.short_name in _FILTER_TYPES and \ - other.short_name in _FILTER_TYPES: + other.short_name in _FILTER_TYPES: delta = _FILTER_TYPES.index(self.short_name) - \ - _FILTER_TYPES.index(other.short_name) + _FILTER_TYPES.index(other.short_name) return delta < 0 else: return False @@ -2301,14 +2368,16 @@ def __hash__(self): string = type(self).__name__ + '\n' string += '{: <16}=\t{}\n'.format('\tEnergy', self.energy) string += '{: <16}=\t{}\n'.format('\tInterpolant', self.y) - string += '{: <16}=\t{}\n'.format('\tInterpolation', self.interpolation) + string += '{: <16}=\t{}\n'.format('\tInterpolation', + self.interpolation) return hash(string) def __repr__(self): string = type(self).__name__ + '\n' string += '{: <16}=\t{}\n'.format('\tEnergy', self.energy) string += '{: <16}=\t{}\n'.format('\tInterpolant', self.y) - string += '{: <16}=\t{}\n'.format('\tInterpolation', self.interpolation) + string += '{: <16}=\t{}\n'.format('\tInterpolation', + self.interpolation) string += '{: <16}=\t{}\n'.format('\tID', self.id) return string @@ -2394,10 +2463,12 @@ def interpolation(self): @interpolation.setter def interpolation(self, val): cv.check_type('interpolation', val, str) - cv.check_value('interpolation', val, self.INTERPOLATION_SCHEMES.values()) + cv.check_value('interpolation', val, + self.INTERPOLATION_SCHEMES.values()) if val == 'quadratic' and len(self.energy) < 3: - raise ValueError('Quadratic interpolation requires 3 or more values.') + raise ValueError( + 'Quadratic interpolation requires 3 or more values.') if val == 'cubic' and len(self.energy) < 4: raise ValueError('Cubic interpolation requires 3 or more values.') diff --git a/src/particle.cpp b/src/particle.cpp index 324276d15f3..000fa13441f 100644 --- a/src/particle.cpp +++ b/src/particle.cpp @@ -84,6 +84,11 @@ bool Particle::create_secondary( return false; } + // This is used to count backward in the secondary source bank for tallying + // outgoing photon energies from a neutron in the SecondaryPhotonEnergy + // filter. + secondaries_this_collision()++; + auto& bank = secondary_bank().emplace_back(); bank.particle = type; bank.wgt = wgt; @@ -320,6 +325,7 @@ void Particle::event_cross_surface() void Particle::event_collide() { + // Score collision estimate of keff if (settings::run_mode == RunMode::EIGENVALUE && type() == ParticleType::neutron) { @@ -364,6 +370,11 @@ void Particle::event_collide() n_bank() = 0; bank_second_E() = 0.0; wgt_bank() = 0.0; + + // Clear number of secondaries in this collision. This is + // distinct from the number of created neutrons n_bank() above! + secondaries_this_collision() = 0; + zero_delayed_bank(); // Reset fission logical diff --git a/src/tallies/filter.cpp b/src/tallies/filter.cpp index 79817981db0..62559b9144e 100644 --- a/src/tallies/filter.cpp +++ b/src/tallies/filter.cpp @@ -124,6 +124,8 @@ Filter* Filter::create(const std::string& type, int32_t id) return Filter::create(id); } else if (type == "energyout") { return Filter::create(id); + } else if (type == "secondaryenergy") { + return Filter::create(id); } else if (type == "legendre") { return Filter::create(id); } else if (type == "material") { diff --git a/src/tallies/filter_delayedgroup.cpp b/src/tallies/filter_delayedgroup.cpp index 01e39e554a9..46adf529b95 100644 --- a/src/tallies/filter_delayedgroup.cpp +++ b/src/tallies/filter_delayedgroup.cpp @@ -27,7 +27,7 @@ void DelayedGroupFilter::set_groups(span groups) } else if (group > MAX_DELAYED_GROUPS) { throw std::invalid_argument { "Encountered delayedgroup bin with index " + std::to_string(group) + - " which is greater than MAX_DELATED_GROUPS (" + + " which is greater than MAX_DELAYED_GROUPS (" + std::to_string(MAX_DELAYED_GROUPS) + ")"}; } groups_.push_back(group); @@ -39,6 +39,10 @@ void DelayedGroupFilter::set_groups(span groups) void DelayedGroupFilter::get_all_bins( const Particle& p, TallyEstimator estimator, FilterMatch& match) const { + // Note that the bin is set to zero here, but bins outside zero are + // tallied to regardless. This is because that logic has to be handled + // in the scoring code instead where looping over the delayed + // group takes place (tally_scoring.cpp). match.bins_.push_back(0); match.weights_.push_back(1.0); } diff --git a/src/tallies/filter_energy.cpp b/src/tallies/filter_energy.cpp index 0b954cce3a3..71a21aced32 100644 --- a/src/tallies/filter_energy.cpp +++ b/src/tallies/filter_energy.cpp @@ -117,6 +117,47 @@ std::string EnergyoutFilter::text_label(int bin) const "Outgoing Energy [{}, {})", bins_.at(bin), bins_.at(bin + 1)); } +//============================================================================== +// SecondaryEnergyoutFilter implementation +//============================================================================== + +void SecondaryEnergyFilter::get_all_bins( + const Particle& p, TallyEstimator estimator, FilterMatch& match) const +{ + assert(p.secondary_bank().size() >= p.secondaries_this_collision()); + + // Loop over secondary bank entries from latest to earliest + for (int secondary_idx = 0; secondary_idx < p.secondaries_this_collision(); + secondary_idx++) { + + int bank_idx = p.secondary_bank().size() - 1 - secondary_idx; + + // Check if this is the correct type of secondary, then + // match its energy if it's the right type + if (p.secondary_bank(bank_idx).particle == secondary_type_) { + const double E = p.secondary_bank(bank_idx).E; + if (E >= bins_.front() && E <= bins_.back()) { + auto bin = lower_bound_index(bins_.begin(), bins_.end(), E); + match.bins_.push_back(bin); + match.weights_.push_back(1.0); + } + } + } +} + +std::string SecondaryEnergyFilter::text_label(int bin) const +{ + return fmt::format( + "Secondary outgoing Energy [{}, {})", bins_.at(bin), bins_.at(bin + 1)); +} + +void SecondaryEnergyFilter::from_xml(pugi::xml_node node) +{ + EnergyFilter::from_xml(node); + std::string p = get_node_value(node, "particle"); + secondary_type_ = str_to_particle_type(p); +} + //============================================================================== // C-API functions //============================================================================== diff --git a/src/tallies/tally.cpp b/src/tallies/tally.cpp index 54840b5188e..7bbe58628e9 100644 --- a/src/tallies/tally.cpp +++ b/src/tallies/tally.cpp @@ -163,6 +163,8 @@ Tally::Tally(pugi::xml_node node) filt_type == FilterType::ZERNIKE || filt_type == FilterType::ZERNIKE_RADIAL) { estimator_ = TallyEstimator::COLLISION; + } else if (filt_type == FilterType::SECONDARY_ENERGY) { + estimator_ = TallyEstimator::ANALOG; } } diff --git a/tests/regression_tests/model_xml/photon_production_inputs_true.dat b/tests/regression_tests/model_xml/photon_production_inputs_true.dat index 10f0bad9841..13db17263eb 100644 --- a/tests/regression_tests/model_xml/photon_production_inputs_true.dat +++ b/tests/regression_tests/model_xml/photon_production_inputs_true.dat @@ -41,6 +41,16 @@ neutron photon electron positron + + neutron + + + 0.0 20000000.0 + + + 0.0 100000.0 300000.0 500000.0 2000000.0 20000000.0 + photon + 1 2 current @@ -63,5 +73,10 @@ total heating (n,gamma) analog + + 5 3 4 + total + analog + diff --git a/tests/regression_tests/photon_production/inputs_true.dat b/tests/regression_tests/photon_production/inputs_true.dat index 10f0bad9841..13db17263eb 100644 --- a/tests/regression_tests/photon_production/inputs_true.dat +++ b/tests/regression_tests/photon_production/inputs_true.dat @@ -41,6 +41,16 @@ neutron photon electron positron + + neutron + + + 0.0 20000000.0 + + + 0.0 100000.0 300000.0 500000.0 2000000.0 20000000.0 + photon + 1 2 current @@ -63,5 +73,10 @@ total heating (n,gamma) analog + + 5 3 4 + total + analog + diff --git a/tests/regression_tests/photon_production/results_true.dat b/tests/regression_tests/photon_production/results_true.dat index 413f6f0ca45..ec319aa0eda 100644 --- a/tests/regression_tests/photon_production/results_true.dat +++ b/tests/regression_tests/photon_production/results_true.dat @@ -138,3 +138,14 @@ tally 4: 2.173936E+08 0.000000E+00 0.000000E+00 +tally 5: +2.000000E-02 +4.000000E-04 +8.200000E-03 +6.724000E-05 +5.280000E-02 +2.787840E-03 +3.546000E-01 +1.257412E-01 +5.063000E-01 +2.563397E-01 diff --git a/tests/regression_tests/photon_production/test.py b/tests/regression_tests/photon_production/test.py index 150448a12f4..579319e9a31 100644 --- a/tests/regression_tests/photon_production/test.py +++ b/tests/regression_tests/photon_production/test.py @@ -25,7 +25,8 @@ def model(): inner_cyl_right.region = -cyl & +x_plane_center & -x_plane_right outer_cyl.region = ~(-cyl & +x_plane_left & -x_plane_right) inner_cyl_right.fill = mat - model.geometry = openmc.Geometry([inner_cyl_left, inner_cyl_right, outer_cyl]) + model.geometry = openmc.Geometry( + [inner_cyl_left, inner_cyl_right, outer_cyl]) source = openmc.IndependentSource() source.space = openmc.stats.Point((0, 0, 0)) @@ -38,17 +39,19 @@ def model(): model.settings.batches = 1 model.settings.photon_transport = True model.settings.electron_treatment = 'ttb' - model.settings.cutoff = {'energy_photon' : 1000.0} + model.settings.cutoff = {'energy_photon': 1000.0} model.settings.source = source surface_filter = openmc.SurfaceFilter(cyl) - particle_filter = openmc.ParticleFilter(['neutron', 'photon', 'electron', 'positron']) + particle_filter = openmc.ParticleFilter( + ['neutron', 'photon', 'electron', 'positron']) current_tally = openmc.Tally() current_tally.filters = [surface_filter, particle_filter] current_tally.scores = ['current'] tally_tracklength = openmc.Tally() tally_tracklength.filters = [particle_filter] - tally_tracklength.scores = ['total', '(n,gamma)'] # heating doesn't work with tracklength + # heating doesn't work with tracklength + tally_tracklength.scores = ['total', '(n,gamma)'] tally_tracklength.nuclides = ['Al27', 'total'] tally_tracklength.estimator = 'tracklength' tally_collision = openmc.Tally() @@ -61,8 +64,26 @@ def model(): tally_analog.scores = ['total', 'heating', '(n,gamma)'] tally_analog.nuclides = ['Al27', 'total'] tally_analog.estimator = 'analog' + + # This is an analog tally tracking the energy distribution of photons + # generated by neutrons. The sum of the tally should give the total + # number of photons generated per source neutron. + ene_filter = openmc.EnergyFilter([0.0, 20e6]) # incident neutron energy + + # Track source energies of secondary gammas + ene2_filter = openmc.SecondaryEnergyFilter( + [0.0, 100e3, 300e3, 500e3, 2e6, 20e6]) + ene2_filter.particle = 'photon' + + neutron_only = openmc.ParticleFilter(['neutron']) + tally_gam_ene = openmc.Tally() + tally_gam_ene.filters = [neutron_only, ene_filter, ene2_filter] + tally_gam_ene.scores = ['total'] + tally_gam_ene.estimator = 'analog' + model.tallies.extend([current_tally, tally_tracklength, - tally_collision, tally_analog]) + tally_collision, tally_analog, + tally_gam_ene]) return model diff --git a/tests/unit_tests/test_filters.py b/tests/unit_tests/test_filters.py index 8c56a310e17..9af22122ae5 100644 --- a/tests/unit_tests/test_filters.py +++ b/tests/unit_tests/test_filters.py @@ -17,14 +17,16 @@ def box_model(): model.settings.particles = 100 model.settings.batches = 10 model.settings.inactive = 0 - model.settings.source = openmc.IndependentSource(space=openmc.stats.Point()) + model.settings.source = openmc.IndependentSource( + space=openmc.stats.Point()) return model def test_cell_instance(): c1 = openmc.Cell() c2 = openmc.Cell() - f = openmc.CellInstanceFilter([(c1, 0), (c1, 1), (c1, 2), (c2, 0), (c2, 1)]) + f = openmc.CellInstanceFilter( + [(c1, 0), (c1, 1), (c1, 2), (c2, 0), (c2, 1)]) # Make sure __repr__ works repr(f) @@ -234,7 +236,7 @@ def test_first_moment(run_in_tmpdir, box_model): flux, scatter = sp.tallies[plain_tally.id].mean.ravel() # Check that first moment matches - first_score = lambda t: sp.tallies[t.id].mean.ravel()[0] + def first_score(t): return sp.tallies[t.id].mean.ravel()[0] assert first_score(leg_tally) == scatter assert first_score(leg_sptl_tally) == scatter assert first_score(sph_scat_tally) == scatter @@ -258,7 +260,8 @@ def test_lethargy_bin_width(): assert len(f.lethargy_bin_width) == 175 energy_bins = openmc.mgxs.GROUP_STRUCTURES['VITAMIN-J-175'] assert f.lethargy_bin_width[0] == np.log10(energy_bins[1]/energy_bins[0]) - assert f.lethargy_bin_width[-1] == np.log10(energy_bins[-1]/energy_bins[-2]) + assert f.lethargy_bin_width[-1] == np.log10( + energy_bins[-1]/energy_bins[-2]) def test_energyfunc(): @@ -292,7 +295,8 @@ def test_tabular_from_energyfilter(): # 'histogram' is the default assert tab.interpolation == 'histogram' - tab = efilter.get_tabular(values=np.array([10, 10, 5]), interpolation='linear-linear') + tab = efilter.get_tabular(values=np.array( + [10, 10, 5]), interpolation='linear-linear') assert tab.interpolation == 'linear-linear' @@ -314,6 +318,38 @@ def test_energy_filter(): openmc.EnergyFilter([-1.2, 0.25, 0.5]) +def test_secondary_energy_filter(): + energy_bins = [1e3, 1e4, 1e5, 1e6] + f = openmc.SecondaryEnergyFilter(energy_bins, particle='photon') + + assert f.particle == 'photon' + assert f.num_bins == 3 + assert f.bins.shape == (3, 2) + assert np.allclose(f.values, energy_bins) + + # __repr__ check + repr(f) + + # to_xml_element() + elem = f.to_xml_element() + assert elem.tag == 'filter' + assert elem.attrib['type'] == 'secondaryenergy' + assert elem.find('particle').text == 'photon' + assert elem.find('bins').text.split()[0] == str(energy_bins[0]) + + # from_xml_element() + new_f = openmc.Filter.from_xml_element(elem) + assert new_f.id == f.id + assert new_f.particle == f.particle + assert np.allclose(new_f.bins, f.bins) + + # pandas output + df = f.get_pandas_dataframe(data_size=3, stride=1) + assert df.shape[0] == 3 + assert "secondaryenergy low [eV]" in df.columns + assert "secondaryenergy high [eV]" in df.columns + + def test_weight(): f = openmc.WeightFilter([0.01, 0.1, 1.0, 10.0]) expected_bins = [[0.01, 0.1], [0.1, 1.0], [1.0, 10.0]]