From eed5dae4f32e325ba22bedf6069a6d4b23299a83 Mon Sep 17 00:00:00 2001 From: Paul Romano Date: Mon, 3 Aug 2026 08:55:14 -0500 Subject: [PATCH 1/6] Initial support for photon MGXS --- docs/source/pythonapi/mgxs.rst | 1 + openmc/mgxs/library.py | 104 ++++++++++++- openmc/mgxs/mgxs.py | 137 +++++++++++++++++- openmc/mgxs_library.py | 73 +++++++++- .../mgxs_library_condense/test.py | 5 +- .../mgxs_library_distribcell/test.py | 5 +- .../mgxs_library_hdf5/test.py | 5 +- .../mgxs_library_mesh/test.py | 5 +- .../mgxs_library_no_nuclides/test.py | 5 +- .../mgxs_library_nuclides/test.py | 5 +- .../mgxs_library_specific_nuclides/test.py | 5 +- tests/unit_tests/test_photon_mgxs.py | 115 +++++++++++++++ 12 files changed, 443 insertions(+), 22 deletions(-) create mode 100644 tests/unit_tests/test_photon_mgxs.py diff --git a/docs/source/pythonapi/mgxs.rst b/docs/source/pythonapi/mgxs.rst index 4141aa0a09a..bda07adf35d 100644 --- a/docs/source/pythonapi/mgxs.rst +++ b/docs/source/pythonapi/mgxs.rst @@ -51,6 +51,7 @@ Multi-group Cross Sections openmc.mgxs.KappaFissionXS openmc.mgxs.MultiplicityMatrixXS openmc.mgxs.NuFissionMatrixXS + openmc.mgxs.PhotonProductionMatrixXS openmc.mgxs.ReducedAbsorptionXS openmc.mgxs.ScatterXS openmc.mgxs.ScatterMatrixXS diff --git a/openmc/mgxs/library.py b/openmc/mgxs/library.py index faa83c0481f..19da8d7cad5 100644 --- a/openmc/mgxs/library.py +++ b/openmc/mgxs/library.py @@ -35,6 +35,9 @@ class Library: A geometry which has been initialized with a root universe by_nuclide : bool If true, computes cross sections for each nuclide in each domain + particle_type : {'neutron', 'photon'}, optional + Particle type for which cross sections are computed. If not specified, + tallies are not filtered by particle type. mgxs_types : Iterable of str The types of cross sections in the library (e.g., ['total', 'scatter']) name : str, optional @@ -47,6 +50,8 @@ class Library: An geometry which has been initialized with a root universe by_nuclide : bool If true, computes cross sections for each nuclide in each domain + particle_type : openmc.ParticleType or None + Particle type for which cross sections are computed mgxs_types : Iterable of str The types of cross sections in the library (e.g., ['total', 'scatter']) domain_type : {'material', 'cell', 'distribcell', 'universe', 'mesh'} @@ -102,11 +107,12 @@ class Library: """ def __init__(self, geometry, by_nuclide=False, - mgxs_types=None, name=''): + mgxs_types=None, name='', particle_type=None): self._name = '' self._geometry = None self._by_nuclide = None + self._particle_type = None self._mgxs_types = [] self._domain_type = None self._domains = 'all' @@ -129,6 +135,8 @@ def __init__(self, geometry, by_nuclide=False, self.name = name self.geometry = geometry self.by_nuclide = by_nuclide + if particle_type is not None: + self.particle_type = particle_type if mgxs_types is not None: self.mgxs_types = mgxs_types @@ -142,6 +150,7 @@ def __deepcopy__(self, memo): clone._name = self.name clone._geometry = self.geometry clone._by_nuclide = self.by_nuclide + clone._particle_type = self.particle_type clone._mgxs_types = self.mgxs_types clone._domain_type = self.domain_type clone._domains = copy.deepcopy(self.domains) @@ -203,7 +212,13 @@ def mgxs_types(self, mgxs_types): openmc.mgxs.ARBITRARY_VECTOR_TYPES + \ openmc.mgxs.ARBITRARY_MATRIX_TYPES if mgxs_types == 'all': - self._mgxs_types = all_mgxs_types + if self.particle_type == openmc.ParticleType.PHOTON: + self._mgxs_types = ( + 'total', 'absorption', 'photon-production matrix') + else: + self._mgxs_types = tuple( + xs for xs in all_mgxs_types + if xs != 'photon-production matrix') else: cv.check_iterable_type('mgxs_types', mgxs_types, str) for mgxs_type in mgxs_types: @@ -224,6 +239,18 @@ def by_nuclide(self, by_nuclide): self._by_nuclide = by_nuclide + @property + def particle_type(self): + return self._particle_type + + @particle_type.setter + def particle_type(self, particle_type): + particle_type = openmc.ParticleType(particle_type) + cv.check_value('particle type', particle_type, + (openmc.ParticleType.NEUTRON, + openmc.ParticleType.PHOTON)) + self._particle_type = particle_type + @property def domain_type(self): return self._domain_type @@ -513,13 +540,15 @@ def build_library(self): else: mgxs = openmc.mgxs.MGXS.get_mgxs( mgxs_type, name=self.name, num_polar=self.num_polar, - num_azimuthal=self.num_azimuthal) + num_azimuthal=self.num_azimuthal, + particle_type=self.particle_type) mgxs.domain = domain mgxs.domain_type = self.domain_type mgxs.energy_groups = self.energy_groups mgxs.by_nuclide = self.by_nuclide - if self.estimator is not None: + if self.estimator is not None and not isinstance( + mgxs, openmc.mgxs.PhotonProductionMatrixXS): mgxs.estimator = self.estimator if mgxs_type in openmc.mgxs.MDGXS_TYPES: @@ -1194,14 +1223,27 @@ def get_xsdata(self, domain, xsdata_name, nuclide='total', xs_type='macro', xsdata.set_decay_rate_mgxs(mymgxs, temperature=temperature, xs_type=xs_type, nuclide=[nuclide], subdomain=subdomain) + # Photon production includes both the surviving primary photon and all + # banked secondary photons. Its multiplicity is therefore already + # folded into the matrix values. + if 'photon-production matrix' in self.mgxs_types: + production = self.get_mgxs( + domain, 'photon-production matrix') + xsdata.set_photon_production_mgxs( + production, temperature=temperature, xs_type=xs_type, + subdomain=subdomain) + using_photon_production = True + using_multiplicity = False + # If multiplicity matrix is available, prefer that - if 'multiplicity matrix' in self.mgxs_types: + elif 'multiplicity matrix' in self.mgxs_types: mymgxs = self.get_mgxs(domain, 'multiplicity matrix') xsdata.set_multiplicity_matrix_mgxs(mymgxs, temperature=temperature, xs_type=xs_type, nuclide=[nuclide], subdomain=subdomain) using_multiplicity = True + using_photon_production = False # multiplicity will fall back to using scatter and nu-scatter elif 'scatter matrix' in self.mgxs_types and \ @@ -1214,6 +1256,7 @@ def get_xsdata(self, domain, xsdata_name, nuclide='total', xs_type='macro', nuclide=[nuclide], subdomain=subdomain) using_multiplicity = True + using_photon_production = False # multiplicity will fall back to using scatter and nu-scatter elif 'consistent scatter matrix' in self.mgxs_types and \ @@ -1227,11 +1270,15 @@ def get_xsdata(self, domain, xsdata_name, nuclide='total', xs_type='macro', nuclide=[nuclide], subdomain=subdomain) using_multiplicity = True + using_photon_production = False else: using_multiplicity = False + using_photon_production = False - if using_multiplicity: + if using_photon_production: + pass + elif using_multiplicity: if 'nu-scatter matrix' in self.mgxs_types: nuscatt_mgxs = self.get_mgxs(domain, 'nu-scatter matrix') else: @@ -1369,7 +1416,8 @@ def create_mg_library(self, xs_type='macro', xsdata_names=None, # Initialize file mgxs_file = openmc.MGXSLibrary( - self.energy_groups, num_delayed_groups=self.num_delayed_groups) + self.energy_groups, num_delayed_groups=self.num_delayed_groups, + particle_type=self.particle_type) if self.domain_type == 'mesh': # Create the xsdata objects and add to the mgxs_file @@ -1585,6 +1633,48 @@ def check_library_for_openmc_mgxs(self): error_flag = False + if 'photon-production matrix' in self.mgxs_types: + photon_mgxs_types = { + 'total', 'absorption', 'photon-production matrix'} + unsupported = set(self.mgxs_types) - photon_mgxs_types + if unsupported: + warn('Photon MGXS libraries do not support the following ' + f'MGXS types: {sorted(unsupported)}.') + error_flag = True + if self.particle_type != openmc.ParticleType.PHOTON: + warn('A photon production matrix requires particle_type to ' + 'be "photon".') + error_flag = True + if self.by_nuclide: + warn('Photon production cannot be tallied by nuclide.') + error_flag = True + if self.correction is not None: + warn('Photon MGXS libraries do not support a transport ' + 'correction; correction must be None.') + error_flag = True + if self.num_polar != 1 or self.num_azimuthal != 1: + warn('Photon production only supports an isotropic ' + 'representation.') + error_flag = True + if self.estimator not in (None, 'tracklength'): + warn('Photon MGXS libraries require the default tracklength ' + 'estimator for flux-weighted cross sections.') + error_flag = True + for mgxs_type in ('total', 'absorption'): + if mgxs_type not in self.mgxs_types: + warn(f'A "{mgxs_type}" MGXS type is required for a ' + 'photon MGXS library.') + error_flag = True + if error_flag: + raise ValueError('Invalid photon MGXS configuration ' + 'encountered.') + return + + if self.particle_type == openmc.ParticleType.PHOTON: + warn('A "photon-production matrix" MGXS type is required for a ' + 'photon MGXS library.') + raise ValueError('Invalid photon MGXS configuration encountered.') + # if correction is 'P0', then transport must be provided # otherwise total must be provided if self.correction == 'P0': diff --git a/openmc/mgxs/mgxs.py b/openmc/mgxs/mgxs.py index 7f0200c0976..d415e2939a0 100644 --- a/openmc/mgxs/mgxs.py +++ b/openmc/mgxs/mgxs.py @@ -30,6 +30,7 @@ 'nu-scatter', 'scatter matrix', 'nu-scatter matrix', + 'photon-production matrix', 'multiplicity matrix', 'nu-fission matrix', 'scatter probability matrix', @@ -177,6 +178,9 @@ class MGXS: The energy group structure for energy condensation by_nuclide : bool If true, computes cross sections for each nuclide in domain + particle_type : {'neutron', 'photon'}, optional + Particle type for which cross sections are computed. If not specified, + tallies are not filtered by particle type. name : str, optional Name of the multi-group cross section. Used as a label to identify tallies in OpenMC 'tallies.xml' file. @@ -195,6 +199,8 @@ class MGXS: Reaction type (e.g., 'total', 'nu-fission', etc.) by_nuclide : bool If true, computes cross sections for each nuclide in domain + particle_type : openmc.ParticleType or None + Particle type for which cross sections are computed domain : openmc.Material or openmc.Cell or openmc.Universe or openmc.RegularMesh Domain for spatial homogenization domain_type : {'material', 'cell', 'distribcell', 'universe', 'mesh'} @@ -261,10 +267,11 @@ class MGXS: def __init__(self, domain=None, domain_type=None, energy_groups=None, by_nuclide=False, name='', num_polar=1, - num_azimuthal=1): + num_azimuthal=1, particle_type=None): self._name = '' self._rxn_type = None self._by_nuclide = None + self._particle_type = None self._nuclides = None self._estimator = 'tracklength' self._domain = None @@ -284,6 +291,8 @@ def __init__(self, domain=None, domain_type=None, self.name = name self.by_nuclide = by_nuclide + if particle_type is not None: + self.particle_type = particle_type if domain_type is not None: self.domain_type = domain_type @@ -306,6 +315,7 @@ def __deepcopy__(self, memo): clone._name = self.name clone._rxn_type = self.rxn_type clone._by_nuclide = self.by_nuclide + clone._particle_type = self.particle_type clone._nuclides = copy.deepcopy(self._nuclides, memo) clone._domain = self.domain clone._domain_type = self.domain_type @@ -472,6 +482,18 @@ def by_nuclide(self, by_nuclide): cv.check_type('by_nuclide', by_nuclide, bool) self._by_nuclide = by_nuclide + @property + def particle_type(self): + return self._particle_type + + @particle_type.setter + def particle_type(self, particle_type): + particle_type = openmc.ParticleType(particle_type) + cv.check_value('particle type', particle_type, + (openmc.ParticleType.NEUTRON, + openmc.ParticleType.PHOTON)) + self._particle_type = particle_type + @property def domain(self): return self._domain @@ -611,6 +633,10 @@ def tallies(self): for add_filter in filters: self._tallies[key].filters.append(add_filter) + if self.particle_type is not None: + self._tallies[key].filters.append( + openmc.ParticleFilter(self.particle_type)) + # If this is a by-nuclide cross-section, add nuclides to Tally if self.by_nuclide and score != 'flux': self._tallies[key].nuclides += self.get_nuclides() @@ -719,7 +745,7 @@ def mgxs_type(self): @staticmethod def get_mgxs(mgxs_type, domain=None, domain_type=None, energy_groups=None, by_nuclide=False, name='', num_polar=1, - num_azimuthal=1): + num_azimuthal=1, particle_type=None): """Return a MGXS subclass object for some energy group structure within some spatial domain for some reaction type. @@ -753,6 +779,9 @@ def get_mgxs(mgxs_type, domain=None, domain_type=None, num_azimuthal : Integral, optional Number of equi-width azimuthal angles for angle discretization; defaults to no discretization + particle_type : {'neutron', 'photon'}, optional + Particle type for which cross sections are computed. If not + specified, tallies are not filtered by particle type. Returns ------- @@ -792,6 +821,10 @@ def get_mgxs(mgxs_type, domain=None, domain_type=None, mgxs = ScatterMatrixXS(domain, domain_type, energy_groups) elif mgxs_type == 'nu-scatter matrix': mgxs = ScatterMatrixXS(domain, domain_type, energy_groups, nu=True) + elif mgxs_type == 'photon-production matrix': + mgxs = PhotonProductionMatrixXS( + domain, domain_type, energy_groups, by_nuclide, name, + num_polar, num_azimuthal) elif mgxs_type == 'multiplicity matrix': mgxs = MultiplicityMatrixXS(domain, domain_type, energy_groups) elif mgxs_type == 'scatter probability matrix': @@ -835,6 +868,8 @@ def get_mgxs(mgxs_type, domain=None, domain_type=None, mgxs.name = name mgxs.num_polar = num_polar mgxs.num_azimuthal = num_azimuthal + if particle_type is not None: + mgxs.particle_type = particle_type return mgxs def get_nuclides(self): @@ -4923,6 +4958,104 @@ def print_groups_and_histogram(avg_xs, err_xs, num_groups, print(string) +@add_params +class PhotonProductionMatrixXS(MatrixMGXS): + r"""A photon production matrix multigroup cross section. + + This matrix includes both the photon that survives a coherent or incoherent + scattering event and secondary photons banked during photon interactions. + The latter include photons from atomic relaxation, thick-target + bremsstrahlung, and positron annihilation. Since each banked photon is + scored with its statistical weight, photon multiplicity is included directly + in the production matrix. + + Photon production is only available as a macroscopic, isotropic cross + section. Per-nuclide production cannot be determined because secondary + photons are tallied from the collision bank rather than a reaction score. + """ + + def __init__(self, domain=None, domain_type=None, energy_groups=None, + by_nuclide=False, name='', num_polar=1, num_azimuthal=1): + if by_nuclide: + raise ValueError('Photon production cannot be tallied by nuclide') + if num_polar != 1 or num_azimuthal != 1: + raise ValueError('Photon production only supports an isotropic ' + 'representation') + super().__init__(domain, domain_type, energy_groups, False, name) + self._rxn_type = 'photon-production' + self._mgxs_type = 'photon-production matrix' + self._particle_type = openmc.ParticleType.PHOTON + self._valid_estimators = ['analog'] + + @property + def by_nuclide(self): + return self._by_nuclide + + @by_nuclide.setter + def by_nuclide(self, by_nuclide): + cv.check_type('by_nuclide', by_nuclide, bool) + if by_nuclide: + raise ValueError('Photon production cannot be tallied by nuclide') + self._by_nuclide = False + + @property + def particle_type(self): + return self._particle_type + + @particle_type.setter + def particle_type(self, particle_type): + particle_type = openmc.ParticleType(particle_type) + if particle_type != openmc.ParticleType.PHOTON: + raise ValueError('Photon production requires photon tallies') + self._particle_type = particle_type + + @property + def scores(self): + return ['flux', 'scatter', 'events'] + + @property + def tally_keys(self): + return ['flux', 'primary photon production', 'secondary photon production'] + + @property + def filters(self): + group_edges = self.energy_groups.group_edges + energy = openmc.EnergyFilter(group_edges) + energyout = openmc.EnergyoutFilter(group_edges) + production = openmc.ParticleProductionFilter( + 'photon', group_edges) + return [[energy], [energy, energyout], [energy, production]] + + @property + def estimator(self): + return ['tracklength', 'analog', 'analog'] + + @property + def rxn_rate_tally(self): + if self._rxn_rate_tally is None: + primary = copy.deepcopy( + self.tallies['primary photon production']) + secondary = copy.deepcopy( + self.tallies['secondary photon production']) + + # Give the secondary tally the same outgoing-energy filter as the + # primary tally so that the two production contributions can be + # combined with standard tally arithmetic. + energyout = copy.deepcopy( + primary.find_filter(openmc.EnergyoutFilter)) + secondary.filters = [ + energyout + if isinstance(f, openmc.ParticleProductionFilter) else f + for f in secondary.filters + ] + primary._scores = ['photon-production'] + secondary._scores = ['photon-production'] + self._rxn_rate_tally = primary + secondary + self._rxn_rate_tally.sparse = self.sparse + + return self._rxn_rate_tally + + @add_params class MultiplicityMatrixXS(MatrixMGXS): r"""The scattering multiplicity matrix. diff --git a/openmc/mgxs_library.py b/openmc/mgxs_library.py index ab9b58b7a3a..776aa8b5bec 100644 --- a/openmc/mgxs_library.py +++ b/openmc/mgxs_library.py @@ -1634,6 +1634,49 @@ def set_scatter_matrix_mgxs(self, scatter, temperature=ROOM_TEMPERATURE_KELVIN, scatter.get_xs(nuclides=nuclide, xs_type=xs_type, subdomains=subdomain) + def set_photon_production_mgxs( + self, production, temperature=ROOM_TEMPERATURE_KELVIN, + xs_type='macro', subdomain=None): + """Set the scattering data from a photon production matrix MGXS. + + The production matrix is written as the zeroth-order scattering + matrix. Its multiplicity is already included in the matrix values, so + no separate multiplicity matrix is written. + + Parameters + ---------- + production : openmc.mgxs.PhotonProductionMatrixXS + Photon production matrix to use + temperature : float, optional + Temperature of the cross section data in kelvin + xs_type : {'macro'}, optional + Cross section type. Photon production is only available as a + macroscopic cross section. + subdomain : iterable of int, optional + Mesh subdomain for which data are requested + """ + + check_type('production', production, + openmc.mgxs.PhotonProductionMatrixXS) + check_value('energy_groups', production.energy_groups, + [self.energy_groups]) + check_value('domain_type', production.domain_type, + openmc.mgxs.DOMAIN_TYPES) + check_value('xs_type', xs_type, ['macro']) + self._check_temperature(temperature) + + if self.representation != REPRESENTATION_ISOTROPIC: + raise ValueError('Photon production only supports an isotropic ' + 'representation') + + self.scatter_format = SCATTER_LEGENDRE + self.order = 0 + i = self._temperature_index(temperature) + self._scatter_matrix[i] = np.zeros( + self.xs_shapes["[G][G'][Order]"]) + self._scatter_matrix[i][:, :, 0] = production.get_xs( + xs_type=xs_type, subdomains=subdomain) + def set_multiplicity_matrix_mgxs(self, nuscatter, scatter=None, temperature=ROOM_TEMPERATURE_KELVIN, nuclide='total', xs_type='macro', subdomain=None): @@ -2364,6 +2407,8 @@ class MGXSLibrary: Energy group structure num_delayed_groups : int Num delayed groups + particle_type : {'neutron', 'photon'}, optional + Particle type represented by the library Attributes ---------- @@ -2371,13 +2416,19 @@ class MGXSLibrary: Energy group structure. num_delayed_groups : int Num delayed groups + particle_type : openmc.ParticleType or None + Particle type represented by the library xsdatas : Iterable of openmc.XSdata Iterable of multi-Group cross section data objects """ - def __init__(self, energy_groups, num_delayed_groups=0): + def __init__(self, energy_groups, num_delayed_groups=0, + particle_type=None): self.energy_groups = energy_groups self.num_delayed_groups = num_delayed_groups + self._particle_type = None + if particle_type is not None: + self.particle_type = particle_type self._xsdatas = [] def __deepcopy__(self, memo): @@ -2388,6 +2439,7 @@ def __deepcopy__(self, memo): clone = type(self).__new__(type(self)) clone._energy_groups = copy.deepcopy(self.energy_groups, memo) clone._num_delayed_groups = self.num_delayed_groups + clone._particle_type = self.particle_type clone._xsdatas = copy.deepcopy(self.xsdatas, memo) memo[id(self)] = clone @@ -2420,6 +2472,18 @@ def num_delayed_groups(self, num_delayed_groups): openmc.mgxs.MAX_DELAYED_GROUPS, equality=True) self._num_delayed_groups = num_delayed_groups + @property + def particle_type(self): + return self._particle_type + + @particle_type.setter + def particle_type(self, particle_type): + particle_type = openmc.ParticleType(particle_type) + check_value('particle type', particle_type, + (openmc.ParticleType.NEUTRON, + openmc.ParticleType.PHOTON)) + self._particle_type = particle_type + @property def xsdatas(self): return self._xsdatas @@ -2592,6 +2656,8 @@ def export_to_hdf5(self, filename='mgxs.h5', libver='earliest'): file.attrs['energy_groups'] = self.energy_groups.num_groups file.attrs['delayed_groups'] = self.num_delayed_groups file.attrs['group structure'] = self.energy_groups.group_edges + if self.particle_type is not None: + file.attrs['particle_type'] = np.bytes_(str(self.particle_type)) for xsdata in self._xsdatas: xsdata.to_hdf5(file) @@ -2633,7 +2699,10 @@ def from_hdf5(cls, filename=None): group_structure = file.attrs['group structure'] num_delayed_groups = file.attrs['delayed_groups'] energy_groups = openmc.mgxs.EnergyGroups(group_structure) - data = cls(energy_groups, num_delayed_groups) + particle_type = file.attrs.get('particle_type') + if isinstance(particle_type, bytes): + particle_type = particle_type.decode() + data = cls(energy_groups, num_delayed_groups, particle_type) for group_name, group in file.items(): data.add_xsdata(openmc.XSdata.from_hdf5(group, group_name, diff --git a/tests/regression_tests/mgxs_library_condense/test.py b/tests/regression_tests/mgxs_library_condense/test.py index bbc4c11bfa9..c9eae94e6c5 100644 --- a/tests/regression_tests/mgxs_library_condense/test.py +++ b/tests/regression_tests/mgxs_library_condense/test.py @@ -19,7 +19,10 @@ def __init__(self, *args, **kwargs): self.mgxs_lib.by_nuclide = False # Test all MGXS types - self.mgxs_lib.mgxs_types = openmc.mgxs.MGXS_TYPES + \ + neutron_mgxs_types = tuple( + xs for xs in openmc.mgxs.MGXS_TYPES + if xs != 'photon-production matrix') + self.mgxs_lib.mgxs_types = neutron_mgxs_types + \ openmc.mgxs.MDGXS_TYPES self.mgxs_lib.energy_groups = energy_groups self.mgxs_lib.num_delayed_groups = 6 diff --git a/tests/regression_tests/mgxs_library_distribcell/test.py b/tests/regression_tests/mgxs_library_distribcell/test.py index 464b309c002..c0df9fce8d9 100644 --- a/tests/regression_tests/mgxs_library_distribcell/test.py +++ b/tests/regression_tests/mgxs_library_distribcell/test.py @@ -23,8 +23,9 @@ def __init__(self, *args, **kwargs): self.mgxs_lib.by_nuclide = False # Test all relevant MGXS types - relevant_MGXS_TYPES = [item for item in openmc.mgxs.MGXS_TYPES - if item != 'current'] + relevant_MGXS_TYPES = [ + item for item in openmc.mgxs.MGXS_TYPES + if item not in ('current', 'photon-production matrix')] self.mgxs_lib.mgxs_types = tuple(relevant_MGXS_TYPES) + \ openmc.mgxs.MDGXS_TYPES self.mgxs_lib.energy_groups = energy_groups diff --git a/tests/regression_tests/mgxs_library_hdf5/test.py b/tests/regression_tests/mgxs_library_hdf5/test.py index 4fb4bf09369..709512cce2c 100644 --- a/tests/regression_tests/mgxs_library_hdf5/test.py +++ b/tests/regression_tests/mgxs_library_hdf5/test.py @@ -23,7 +23,10 @@ def __init__(self, *args, **kwargs): self.mgxs_lib.by_nuclide = False # Test all MGXS types - self.mgxs_lib.mgxs_types = openmc.mgxs.MGXS_TYPES + \ + neutron_mgxs_types = tuple( + xs for xs in openmc.mgxs.MGXS_TYPES + if xs != 'photon-production matrix') + self.mgxs_lib.mgxs_types = neutron_mgxs_types + \ openmc.mgxs.MDGXS_TYPES self.mgxs_lib.energy_groups = energy_groups self.mgxs_lib.num_delayed_groups = 6 diff --git a/tests/regression_tests/mgxs_library_mesh/test.py b/tests/regression_tests/mgxs_library_mesh/test.py index c1a5980b5d1..3d220ee8c54 100644 --- a/tests/regression_tests/mgxs_library_mesh/test.py +++ b/tests/regression_tests/mgxs_library_mesh/test.py @@ -40,7 +40,10 @@ def model(): model.mgxs_lib.by_nuclide = False # Test all MGXS types - model.mgxs_lib.mgxs_types = openmc.mgxs.MGXS_TYPES + openmc.mgxs.MDGXS_TYPES + neutron_mgxs_types = tuple( + xs for xs in openmc.mgxs.MGXS_TYPES + if xs != 'photon-production matrix') + model.mgxs_lib.mgxs_types = neutron_mgxs_types + openmc.mgxs.MDGXS_TYPES model.mgxs_lib.energy_groups = energy_groups model.mgxs_lib.num_delayed_groups = 6 model.mgxs_lib.correction = None # Avoid warning about P0 correction diff --git a/tests/regression_tests/mgxs_library_no_nuclides/test.py b/tests/regression_tests/mgxs_library_no_nuclides/test.py index a02086af3ec..f7588bafc78 100644 --- a/tests/regression_tests/mgxs_library_no_nuclides/test.py +++ b/tests/regression_tests/mgxs_library_no_nuclides/test.py @@ -20,8 +20,9 @@ def __init__(self, *args, **kwargs): self.mgxs_lib.by_nuclide = False # Test relevant MGXS types - relevant_MGXS_TYPES = [item for item in openmc.mgxs.MGXS_TYPES - if item != 'current'] + relevant_MGXS_TYPES = [ + item for item in openmc.mgxs.MGXS_TYPES + if item not in ('current', 'photon-production matrix')] # Add in a subset of openmc.mgxs.ARBITRARY_VECTOR_TYPES and # openmc.mgxs.ARBITRARY_MATRIX_TYPES so we can see the code works, # but not use too much resources diff --git a/tests/regression_tests/mgxs_library_nuclides/test.py b/tests/regression_tests/mgxs_library_nuclides/test.py index a10070358ad..e11bb995e27 100644 --- a/tests/regression_tests/mgxs_library_nuclides/test.py +++ b/tests/regression_tests/mgxs_library_nuclides/test.py @@ -19,8 +19,9 @@ def __init__(self, *args, **kwargs): self.mgxs_lib.by_nuclide = True # Test relevant MGXS types - relevant_MGXS_TYPES = [item for item in openmc.mgxs.MGXS_TYPES - if item != 'current'] + relevant_MGXS_TYPES = [ + item for item in openmc.mgxs.MGXS_TYPES + if item not in ('current', 'photon-production matrix')] # Add in a subset of openmc.mgxs.ARBITRARY_VECTOR_TYPES and # openmc.mgxs.ARBITRARY_MATRIX_TYPES so we can see the code works, # but not use too much resources diff --git a/tests/regression_tests/mgxs_library_specific_nuclides/test.py b/tests/regression_tests/mgxs_library_specific_nuclides/test.py index 0ccbb83bdbe..688c0660c84 100644 --- a/tests/regression_tests/mgxs_library_specific_nuclides/test.py +++ b/tests/regression_tests/mgxs_library_specific_nuclides/test.py @@ -19,8 +19,9 @@ def __init__(self, *args, **kwargs): self.mgxs_lib.by_nuclide = True # Test relevant MGXS types - relevant_MGXS_TYPES = [item for item in openmc.mgxs.MGXS_TYPES - if item != 'current'] + relevant_MGXS_TYPES = [ + item for item in openmc.mgxs.MGXS_TYPES + if item not in ('current', 'photon-production matrix')] # Add in a subset of openmc.mgxs.ARBITRARY_VECTOR_TYPES and # openmc.mgxs.ARBITRARY_MATRIX_TYPES so we can see the code works, # but not use too much resources diff --git a/tests/unit_tests/test_photon_mgxs.py b/tests/unit_tests/test_photon_mgxs.py new file mode 100644 index 00000000000..08a8eccc07d --- /dev/null +++ b/tests/unit_tests/test_photon_mgxs.py @@ -0,0 +1,115 @@ +import h5py +import numpy as np +import pytest + +import openmc + + +def test_photon_production_matrix_tallies(): + material = openmc.Material() + groups = openmc.mgxs.EnergyGroups([1.0, 10.0, 100.0]) + production = openmc.mgxs.PhotonProductionMatrixXS( + material, 'material', groups) + + assert production.particle_type == openmc.ParticleType.PHOTON + assert production.estimator == ['tracklength', 'analog', 'analog'] + assert production.tally_keys == [ + 'flux', 'primary photon production', 'secondary photon production'] + + tallies = production.tallies + assert tallies['primary photon production'].scores == ['scatter'] + assert tallies['secondary photon production'].scores == ['events'] + assert tallies['secondary photon production'].contains_filter( + openmc.ParticleProductionFilter) + for tally in tallies.values(): + particle_filter = tally.find_filter(openmc.ParticleFilter) + assert particle_filter.bins == ['photon'] + + +def test_photon_production_matrix_constraints(): + with pytest.raises(ValueError, match='by nuclide'): + openmc.mgxs.MGXS.get_mgxs( + 'photon-production matrix', by_nuclide=True) + with pytest.raises(ValueError, match='isotropic'): + openmc.mgxs.MGXS.get_mgxs( + 'photon-production matrix', num_polar=2) + with pytest.raises(ValueError, match='photon tallies'): + openmc.mgxs.MGXS.get_mgxs( + 'photon-production matrix', particle_type='neutron') + + +def test_photon_production_matrix_combines_primary_and_secondary(): + material = openmc.Material() + groups = openmc.mgxs.EnergyGroups([1.0, 10.0, 100.0]) + production = openmc.mgxs.PhotonProductionMatrixXS( + material, 'material', groups) + + for tally in production.tallies.values(): + values = np.arange(1, tally.num_filter_bins + 1, dtype=float) + values.shape = (tally.num_filter_bins, 1, 1) + tally._sum = 2.0 * values + tally._sum_sq = 2.0 * values**2 + tally._num_realizations = 2 + tally._sp_filename = 'statepoint.h5' + tally._results_read = True + + np.testing.assert_allclose( + production.rxn_rate_tally.mean.ravel(), [2.0, 4.0, 6.0, 8.0]) + + +def test_set_photon_production_mgxs(monkeypatch): + material = openmc.Material() + groups = openmc.mgxs.EnergyGroups([1.0, 10.0, 100.0]) + production = openmc.mgxs.PhotonProductionMatrixXS( + material, 'material', groups) + values = np.array([[0.1, 0.2], [0.3, 0.4]]) + monkeypatch.setattr(production, 'get_xs', lambda **kwargs: values) + + xsdata = openmc.XSdata('photon', groups) + xsdata.order = 0 + absorption = np.array([0.5, 0.6]) + xsdata.set_absorption(absorption) + xsdata.set_photon_production_mgxs(production) + + assert xsdata.scatter_format == 'legendre' + assert xsdata.order == 0 + np.testing.assert_allclose(xsdata._scatter_matrix[0][:, :, 0], values) + np.testing.assert_allclose(xsdata._absorption[0], absorption) + assert xsdata._multiplicity_matrix[0] is None + + +def test_photon_mgxs_library(tmp_path): + material = openmc.Material() + geometry = openmc.Geometry([openmc.Cell(fill=material)]) + groups = openmc.mgxs.EnergyGroups([1.0, 10.0, 100.0]) + library = openmc.mgxs.Library( + geometry, mgxs_types=[ + 'total', 'absorption', 'photon-production matrix'], + particle_type='photon') + library.domain_type = 'material' + library.energy_groups = groups + library.correction = None + library.build_library() + library.check_library_for_openmc_mgxs() + + for mgxs in library.all_mgxs[material.id].values(): + assert mgxs.particle_type == openmc.ParticleType.PHOTON + + mg_library = openmc.MGXSLibrary(groups, particle_type='photon') + path = tmp_path / 'mgxs.h5' + mg_library.export_to_hdf5(path) + with h5py.File(path) as h5file: + assert h5file.attrs['particle_type'] == b'photon' + assert openmc.MGXSLibrary.from_hdf5(path).particle_type == \ + openmc.ParticleType.PHOTON + + +def test_all_mgxs_types_respects_particle_type(): + geometry = openmc.Geometry([openmc.Cell()]) + neutron_library = openmc.mgxs.Library(geometry, mgxs_types='all') + photon_library = openmc.mgxs.Library( + geometry, mgxs_types='all', particle_type='photon') + + assert 'photon-production matrix' not in neutron_library.mgxs_types + assert photon_library.mgxs_types == ( + 'total', 'absorption', 'photon-production matrix') From 39e4fb26061f592e1830a7a93f3b5286e15d1f82 Mon Sep 17 00:00:00 2001 From: Paul Romano Date: Sat, 8 Aug 2026 21:35:18 -0500 Subject: [PATCH 2/6] Add photon MGXS regression test --- .../regression_tests/mgxs_photon/__init__.py | 0 .../mgxs_photon/inputs_true.dat | 95 ++++++++++++++++ .../mgxs_photon/results_true.dat | 8 ++ tests/regression_tests/mgxs_photon/test.py | 102 ++++++++++++++++++ 4 files changed, 205 insertions(+) create mode 100644 tests/regression_tests/mgxs_photon/__init__.py create mode 100644 tests/regression_tests/mgxs_photon/inputs_true.dat create mode 100644 tests/regression_tests/mgxs_photon/results_true.dat create mode 100644 tests/regression_tests/mgxs_photon/test.py diff --git a/tests/regression_tests/mgxs_photon/__init__.py b/tests/regression_tests/mgxs_photon/__init__.py new file mode 100644 index 00000000000..e69de29bb2d diff --git a/tests/regression_tests/mgxs_photon/inputs_true.dat b/tests/regression_tests/mgxs_photon/inputs_true.dat new file mode 100644 index 00000000000..aac04060864 --- /dev/null +++ b/tests/regression_tests/mgxs_photon/inputs_true.dat @@ -0,0 +1,95 @@ + + + + + + + + + + + + + + + + + fixed source + 2000 + 2 + + + 0.0 0.0 0.0 + + + 1000000.0 1.0 + + + ttb + true + true + + 1000.0 + + + + + 1 + + + 1000.0 100000.0 500000.0 1100000.0 + + + photon + + + 1000.0 100000.0 500000.0 1100000.0 + + + photon + 1000.0 100000.0 500000.0 1100000.0 + + + 1 2 3 + total + flux + tracklength + + + 1 2 3 + total + total + tracklength + + + 1 2 3 + total + flux + tracklength + + + 1 2 3 + total + absorption + tracklength + + + 1 2 3 + total + flux + tracklength + + + 1 2 11 3 + total + scatter + analog + + + 1 2 12 3 + total + events + analog + + + diff --git a/tests/regression_tests/mgxs_photon/results_true.dat b/tests/regression_tests/mgxs_photon/results_true.dat new file mode 100644 index 00000000000..0a84849600c --- /dev/null +++ b/tests/regression_tests/mgxs_photon/results_true.dat @@ -0,0 +1,8 @@ +primary photons: 6.55500000e-01 +secondary photons: 1.25150000e+00 +production matrix: +[[ 0.33790683 0.34765112 0.79557413] + [ 0. 1.17578035 5.49497349] + [ 0. 0. 19.5842782 ]] +absorption: +[ 0.23475329 3.62317561 76.43950604] diff --git a/tests/regression_tests/mgxs_photon/test.py b/tests/regression_tests/mgxs_photon/test.py new file mode 100644 index 00000000000..b7ea928c7c7 --- /dev/null +++ b/tests/regression_tests/mgxs_photon/test.py @@ -0,0 +1,102 @@ +import hashlib + +import numpy as np +import openmc + +from tests.testing_harness import PyAPITestHarness + + +class PhotonMGXSTestHarness(PyAPITestHarness): + """Run photon transport and verify photon MGXS post-processing.""" + + def __init__(self, *args, **kwargs): + super().__init__(*args, **kwargs) + + self.material = openmc.Material() + self.material.set_density('g/cm3', 11.35) + self.material.add_element('Pb', 1.0) + self._model.materials = openmc.Materials([self.material]) + + sphere = openmc.Sphere(r=1.0, boundary_type='vacuum') + cell = openmc.Cell(fill=self.material, region=-sphere) + self._model.geometry = openmc.Geometry([cell]) + + self._model.settings.run_mode = 'fixed source' + self._model.settings.particles = 2000 + self._model.settings.batches = 2 + self._model.settings.photon_transport = True + self._model.settings.atomic_relaxation = True + self._model.settings.electron_treatment = 'ttb' + self._model.settings.cutoff = {'energy_photon': 1000.0} + self._model.settings.source = openmc.IndependentSource( + particle='photon', + space=openmc.stats.Point((0.0, 0.0, 0.0)), + energy=openmc.stats.Discrete([1.0e6], [1.0])) + + groups = openmc.mgxs.EnergyGroups( + group_edges=[1.0e3, 1.0e5, 5.0e5, 1.1e6]) + self.mgxs_lib = openmc.mgxs.Library( + self._model.geometry, + mgxs_types=['total', 'absorption', 'photon-production matrix'], + particle_type='photon') + self.mgxs_lib.energy_groups = groups + self.mgxs_lib.correction = None + self.mgxs_lib.domain_type = 'material' + self.mgxs_lib.build_library() + self.mgxs_lib.add_to_tallies(self._model.tallies, merge=False) + + def _get_results(self, hash_output=False): + with openmc.StatePoint(self._sp_name) as statepoint: + self.mgxs_lib.load_from_statepoint(statepoint) + + production = self.mgxs_lib.get_mgxs( + self.material, 'photon-production matrix') + absorption = self.mgxs_lib.get_mgxs(self.material, 'absorption') + + primary = production.tallies[ + 'primary photon production'].mean.sum() + secondary = production.tallies[ + 'secondary photon production'].mean.sum() + if secondary <= 0.0: + raise AssertionError( + 'Photon transport did not score any secondary photons') + + production_xs = production.get_xs() + absorption_xs = absorption.get_xs() + + # Exercise the XSdata conversion used by subsequent MG/RR + # workflows. + mg_library = self.mgxs_lib.create_mg_library() + xsdata = mg_library.xsdatas[0] + if xsdata.multiplicity_matrix[0] is not None: + raise AssertionError( + 'Photon production should not create a multiplicity ' + 'matrix') + if not np.allclose( + xsdata.scatter_matrix[0][:, :, 0], production_xs): + raise AssertionError('Photon production matrix was not ' + 'exported') + if not np.allclose(xsdata.absorption[0], absorption_xs): + raise AssertionError('Photon absorption was changed during ' + 'export') + + output = [ + f'primary photons: {primary:.8e}', + f'secondary photons: {secondary:.8e}', + 'production matrix:', + np.array2string(production_xs, precision=8), + 'absorption:', + np.array2string(absorption_xs, precision=8), + ] + output = '\n'.join(output) + '\n' + + if hash_output: + digest = hashlib.sha512(output.encode('utf-8')) + output = digest.hexdigest() + return output + + +def test_photon_mgxs(): + harness = PhotonMGXSTestHarness( + 'statepoint.2.h5', model=openmc.Model()) + harness.main() From 8d695ce382d990df631112a3a1f6104beb72f8ff Mon Sep 17 00:00:00 2001 From: Paul Romano Date: Sun, 9 Aug 2026 16:47:34 -0500 Subject: [PATCH 3/6] Use nu-scatter formalism for photon production --- openmc/mgxs/library.py | 43 ++++---------- openmc/mgxs/mgxs.py | 22 +++++-- openmc/mgxs_library.py | 46 +-------------- .../mgxs_library_condense/test.py | 5 +- .../mgxs_library_distribcell/test.py | 5 +- .../mgxs_library_hdf5/test.py | 5 +- .../mgxs_library_mesh/test.py | 5 +- .../mgxs_library_no_nuclides/test.py | 5 +- .../mgxs_library_nuclides/test.py | 5 +- .../mgxs_library_specific_nuclides/test.py | 5 +- tests/regression_tests/mgxs_photon/test.py | 59 +++++++++++-------- tests/unit_tests/test_photon_mgxs.py | 20 +++---- 12 files changed, 86 insertions(+), 139 deletions(-) diff --git a/openmc/mgxs/library.py b/openmc/mgxs/library.py index 19da8d7cad5..6338c7d732d 100644 --- a/openmc/mgxs/library.py +++ b/openmc/mgxs/library.py @@ -214,11 +214,9 @@ def mgxs_types(self, mgxs_types): if mgxs_types == 'all': if self.particle_type == openmc.ParticleType.PHOTON: self._mgxs_types = ( - 'total', 'absorption', 'photon-production matrix') + 'total', 'absorption', 'nu-scatter matrix') else: - self._mgxs_types = tuple( - xs for xs in all_mgxs_types - if xs != 'photon-production matrix') + self._mgxs_types = all_mgxs_types else: cv.check_iterable_type('mgxs_types', mgxs_types, str) for mgxs_type in mgxs_types: @@ -1223,16 +1221,14 @@ def get_xsdata(self, domain, xsdata_name, nuclide='total', xs_type='macro', xsdata.set_decay_rate_mgxs(mymgxs, temperature=temperature, xs_type=xs_type, nuclide=[nuclide], subdomain=subdomain) - # Photon production includes both the surviving primary photon and all - # banked secondary photons. Its multiplicity is therefore already - # folded into the matrix values. - if 'photon-production matrix' in self.mgxs_types: - production = self.get_mgxs( - domain, 'photon-production matrix') - xsdata.set_photon_production_mgxs( + # Photon nu-scatter includes both the surviving primary photon and all + # banked secondary photons, so its multiplicity is already folded into + # the matrix values. + if self.particle_type == openmc.ParticleType.PHOTON: + production = self.get_mgxs(domain, 'nu-scatter matrix') + xsdata.set_scatter_matrix_mgxs( production, temperature=temperature, xs_type=xs_type, - subdomain=subdomain) - using_photon_production = True + nuclide=[nuclide], subdomain=subdomain) using_multiplicity = False # If multiplicity matrix is available, prefer that @@ -1243,7 +1239,6 @@ def get_xsdata(self, domain, xsdata_name, nuclide='total', xs_type='macro', nuclide=[nuclide], subdomain=subdomain) using_multiplicity = True - using_photon_production = False # multiplicity will fall back to using scatter and nu-scatter elif 'scatter matrix' in self.mgxs_types and \ @@ -1256,7 +1251,6 @@ def get_xsdata(self, domain, xsdata_name, nuclide='total', xs_type='macro', nuclide=[nuclide], subdomain=subdomain) using_multiplicity = True - using_photon_production = False # multiplicity will fall back to using scatter and nu-scatter elif 'consistent scatter matrix' in self.mgxs_types and \ @@ -1270,13 +1264,11 @@ def get_xsdata(self, domain, xsdata_name, nuclide='total', xs_type='macro', nuclide=[nuclide], subdomain=subdomain) using_multiplicity = True - using_photon_production = False else: using_multiplicity = False - using_photon_production = False - if using_photon_production: + if self.particle_type == openmc.ParticleType.PHOTON: pass elif using_multiplicity: if 'nu-scatter matrix' in self.mgxs_types: @@ -1633,18 +1625,14 @@ def check_library_for_openmc_mgxs(self): error_flag = False - if 'photon-production matrix' in self.mgxs_types: + if self.particle_type == openmc.ParticleType.PHOTON: photon_mgxs_types = { - 'total', 'absorption', 'photon-production matrix'} + 'total', 'absorption', 'nu-scatter matrix'} unsupported = set(self.mgxs_types) - photon_mgxs_types if unsupported: warn('Photon MGXS libraries do not support the following ' f'MGXS types: {sorted(unsupported)}.') error_flag = True - if self.particle_type != openmc.ParticleType.PHOTON: - warn('A photon production matrix requires particle_type to ' - 'be "photon".') - error_flag = True if self.by_nuclide: warn('Photon production cannot be tallied by nuclide.') error_flag = True @@ -1660,7 +1648,7 @@ def check_library_for_openmc_mgxs(self): warn('Photon MGXS libraries require the default tracklength ' 'estimator for flux-weighted cross sections.') error_flag = True - for mgxs_type in ('total', 'absorption'): + for mgxs_type in ('total', 'absorption', 'nu-scatter matrix'): if mgxs_type not in self.mgxs_types: warn(f'A "{mgxs_type}" MGXS type is required for a ' 'photon MGXS library.') @@ -1670,11 +1658,6 @@ def check_library_for_openmc_mgxs(self): 'encountered.') return - if self.particle_type == openmc.ParticleType.PHOTON: - warn('A "photon-production matrix" MGXS type is required for a ' - 'photon MGXS library.') - raise ValueError('Invalid photon MGXS configuration encountered.') - # if correction is 'P0', then transport must be provided # otherwise total must be provided if self.correction == 'P0': diff --git a/openmc/mgxs/mgxs.py b/openmc/mgxs/mgxs.py index d415e2939a0..b4eeae63045 100644 --- a/openmc/mgxs/mgxs.py +++ b/openmc/mgxs/mgxs.py @@ -30,7 +30,6 @@ 'nu-scatter', 'scatter matrix', 'nu-scatter matrix', - 'photon-production matrix', 'multiplicity matrix', 'nu-fission matrix', 'scatter probability matrix', @@ -819,12 +818,15 @@ def get_mgxs(mgxs_type, domain=None, domain_type=None, mgxs = ScatterXS(domain, domain_type, energy_groups, nu=True) elif mgxs_type == 'scatter matrix': mgxs = ScatterMatrixXS(domain, domain_type, energy_groups) - elif mgxs_type == 'nu-scatter matrix': - mgxs = ScatterMatrixXS(domain, domain_type, energy_groups, nu=True) - elif mgxs_type == 'photon-production matrix': + elif mgxs_type == 'nu-scatter matrix' and \ + particle_type is not None and \ + openmc.ParticleType(particle_type) == \ + openmc.ParticleType.PHOTON: mgxs = PhotonProductionMatrixXS( domain, domain_type, energy_groups, by_nuclide, name, num_polar, num_azimuthal) + elif mgxs_type == 'nu-scatter matrix': + mgxs = ScatterMatrixXS(domain, domain_type, energy_groups, nu=True) elif mgxs_type == 'multiplicity matrix': mgxs = MultiplicityMatrixXS(domain, domain_type, energy_groups) elif mgxs_type == 'scatter probability matrix': @@ -4982,8 +4984,8 @@ def __init__(self, domain=None, domain_type=None, energy_groups=None, raise ValueError('Photon production only supports an isotropic ' 'representation') super().__init__(domain, domain_type, energy_groups, False, name) - self._rxn_type = 'photon-production' - self._mgxs_type = 'photon-production matrix' + self._rxn_type = 'nu-scatter' + self._mgxs_type = 'nu-scatter matrix' self._particle_type = openmc.ParticleType.PHOTON self._valid_estimators = ['analog'] @@ -5030,6 +5032,14 @@ def filters(self): def estimator(self): return ['tracklength', 'analog', 'analog'] + @property + def scatter_format(self): + return SCATTER_LEGENDRE + + @property + def legendre_order(self): + return 0 + @property def rxn_rate_tally(self): if self._rxn_rate_tally is None: diff --git a/openmc/mgxs_library.py b/openmc/mgxs_library.py index 776aa8b5bec..7d3c3d4151b 100644 --- a/openmc/mgxs_library.py +++ b/openmc/mgxs_library.py @@ -1586,7 +1586,8 @@ def set_scatter_matrix_mgxs(self, scatter, temperature=ROOM_TEMPERATURE_KELVIN, """ - check_type('scatter', scatter, openmc.mgxs.ScatterMatrixXS) + check_type('scatter', scatter, (openmc.mgxs.ScatterMatrixXS, + openmc.mgxs.PhotonProductionMatrixXS)) check_value('energy_groups', scatter.energy_groups, [self.energy_groups]) check_value('domain_type', scatter.domain_type, @@ -1634,49 +1635,6 @@ def set_scatter_matrix_mgxs(self, scatter, temperature=ROOM_TEMPERATURE_KELVIN, scatter.get_xs(nuclides=nuclide, xs_type=xs_type, subdomains=subdomain) - def set_photon_production_mgxs( - self, production, temperature=ROOM_TEMPERATURE_KELVIN, - xs_type='macro', subdomain=None): - """Set the scattering data from a photon production matrix MGXS. - - The production matrix is written as the zeroth-order scattering - matrix. Its multiplicity is already included in the matrix values, so - no separate multiplicity matrix is written. - - Parameters - ---------- - production : openmc.mgxs.PhotonProductionMatrixXS - Photon production matrix to use - temperature : float, optional - Temperature of the cross section data in kelvin - xs_type : {'macro'}, optional - Cross section type. Photon production is only available as a - macroscopic cross section. - subdomain : iterable of int, optional - Mesh subdomain for which data are requested - """ - - check_type('production', production, - openmc.mgxs.PhotonProductionMatrixXS) - check_value('energy_groups', production.energy_groups, - [self.energy_groups]) - check_value('domain_type', production.domain_type, - openmc.mgxs.DOMAIN_TYPES) - check_value('xs_type', xs_type, ['macro']) - self._check_temperature(temperature) - - if self.representation != REPRESENTATION_ISOTROPIC: - raise ValueError('Photon production only supports an isotropic ' - 'representation') - - self.scatter_format = SCATTER_LEGENDRE - self.order = 0 - i = self._temperature_index(temperature) - self._scatter_matrix[i] = np.zeros( - self.xs_shapes["[G][G'][Order]"]) - self._scatter_matrix[i][:, :, 0] = production.get_xs( - xs_type=xs_type, subdomains=subdomain) - def set_multiplicity_matrix_mgxs(self, nuscatter, scatter=None, temperature=ROOM_TEMPERATURE_KELVIN, nuclide='total', xs_type='macro', subdomain=None): diff --git a/tests/regression_tests/mgxs_library_condense/test.py b/tests/regression_tests/mgxs_library_condense/test.py index c9eae94e6c5..bbc4c11bfa9 100644 --- a/tests/regression_tests/mgxs_library_condense/test.py +++ b/tests/regression_tests/mgxs_library_condense/test.py @@ -19,10 +19,7 @@ def __init__(self, *args, **kwargs): self.mgxs_lib.by_nuclide = False # Test all MGXS types - neutron_mgxs_types = tuple( - xs for xs in openmc.mgxs.MGXS_TYPES - if xs != 'photon-production matrix') - self.mgxs_lib.mgxs_types = neutron_mgxs_types + \ + self.mgxs_lib.mgxs_types = openmc.mgxs.MGXS_TYPES + \ openmc.mgxs.MDGXS_TYPES self.mgxs_lib.energy_groups = energy_groups self.mgxs_lib.num_delayed_groups = 6 diff --git a/tests/regression_tests/mgxs_library_distribcell/test.py b/tests/regression_tests/mgxs_library_distribcell/test.py index c0df9fce8d9..464b309c002 100644 --- a/tests/regression_tests/mgxs_library_distribcell/test.py +++ b/tests/regression_tests/mgxs_library_distribcell/test.py @@ -23,9 +23,8 @@ def __init__(self, *args, **kwargs): self.mgxs_lib.by_nuclide = False # Test all relevant MGXS types - relevant_MGXS_TYPES = [ - item for item in openmc.mgxs.MGXS_TYPES - if item not in ('current', 'photon-production matrix')] + relevant_MGXS_TYPES = [item for item in openmc.mgxs.MGXS_TYPES + if item != 'current'] self.mgxs_lib.mgxs_types = tuple(relevant_MGXS_TYPES) + \ openmc.mgxs.MDGXS_TYPES self.mgxs_lib.energy_groups = energy_groups diff --git a/tests/regression_tests/mgxs_library_hdf5/test.py b/tests/regression_tests/mgxs_library_hdf5/test.py index 709512cce2c..4fb4bf09369 100644 --- a/tests/regression_tests/mgxs_library_hdf5/test.py +++ b/tests/regression_tests/mgxs_library_hdf5/test.py @@ -23,10 +23,7 @@ def __init__(self, *args, **kwargs): self.mgxs_lib.by_nuclide = False # Test all MGXS types - neutron_mgxs_types = tuple( - xs for xs in openmc.mgxs.MGXS_TYPES - if xs != 'photon-production matrix') - self.mgxs_lib.mgxs_types = neutron_mgxs_types + \ + self.mgxs_lib.mgxs_types = openmc.mgxs.MGXS_TYPES + \ openmc.mgxs.MDGXS_TYPES self.mgxs_lib.energy_groups = energy_groups self.mgxs_lib.num_delayed_groups = 6 diff --git a/tests/regression_tests/mgxs_library_mesh/test.py b/tests/regression_tests/mgxs_library_mesh/test.py index 3d220ee8c54..c1a5980b5d1 100644 --- a/tests/regression_tests/mgxs_library_mesh/test.py +++ b/tests/regression_tests/mgxs_library_mesh/test.py @@ -40,10 +40,7 @@ def model(): model.mgxs_lib.by_nuclide = False # Test all MGXS types - neutron_mgxs_types = tuple( - xs for xs in openmc.mgxs.MGXS_TYPES - if xs != 'photon-production matrix') - model.mgxs_lib.mgxs_types = neutron_mgxs_types + openmc.mgxs.MDGXS_TYPES + model.mgxs_lib.mgxs_types = openmc.mgxs.MGXS_TYPES + openmc.mgxs.MDGXS_TYPES model.mgxs_lib.energy_groups = energy_groups model.mgxs_lib.num_delayed_groups = 6 model.mgxs_lib.correction = None # Avoid warning about P0 correction diff --git a/tests/regression_tests/mgxs_library_no_nuclides/test.py b/tests/regression_tests/mgxs_library_no_nuclides/test.py index f7588bafc78..a02086af3ec 100644 --- a/tests/regression_tests/mgxs_library_no_nuclides/test.py +++ b/tests/regression_tests/mgxs_library_no_nuclides/test.py @@ -20,9 +20,8 @@ def __init__(self, *args, **kwargs): self.mgxs_lib.by_nuclide = False # Test relevant MGXS types - relevant_MGXS_TYPES = [ - item for item in openmc.mgxs.MGXS_TYPES - if item not in ('current', 'photon-production matrix')] + relevant_MGXS_TYPES = [item for item in openmc.mgxs.MGXS_TYPES + if item != 'current'] # Add in a subset of openmc.mgxs.ARBITRARY_VECTOR_TYPES and # openmc.mgxs.ARBITRARY_MATRIX_TYPES so we can see the code works, # but not use too much resources diff --git a/tests/regression_tests/mgxs_library_nuclides/test.py b/tests/regression_tests/mgxs_library_nuclides/test.py index e11bb995e27..a10070358ad 100644 --- a/tests/regression_tests/mgxs_library_nuclides/test.py +++ b/tests/regression_tests/mgxs_library_nuclides/test.py @@ -19,9 +19,8 @@ def __init__(self, *args, **kwargs): self.mgxs_lib.by_nuclide = True # Test relevant MGXS types - relevant_MGXS_TYPES = [ - item for item in openmc.mgxs.MGXS_TYPES - if item not in ('current', 'photon-production matrix')] + relevant_MGXS_TYPES = [item for item in openmc.mgxs.MGXS_TYPES + if item != 'current'] # Add in a subset of openmc.mgxs.ARBITRARY_VECTOR_TYPES and # openmc.mgxs.ARBITRARY_MATRIX_TYPES so we can see the code works, # but not use too much resources diff --git a/tests/regression_tests/mgxs_library_specific_nuclides/test.py b/tests/regression_tests/mgxs_library_specific_nuclides/test.py index 688c0660c84..0ccbb83bdbe 100644 --- a/tests/regression_tests/mgxs_library_specific_nuclides/test.py +++ b/tests/regression_tests/mgxs_library_specific_nuclides/test.py @@ -19,9 +19,8 @@ def __init__(self, *args, **kwargs): self.mgxs_lib.by_nuclide = True # Test relevant MGXS types - relevant_MGXS_TYPES = [ - item for item in openmc.mgxs.MGXS_TYPES - if item not in ('current', 'photon-production matrix')] + relevant_MGXS_TYPES = [item for item in openmc.mgxs.MGXS_TYPES + if item != 'current'] # Add in a subset of openmc.mgxs.ARBITRARY_VECTOR_TYPES and # openmc.mgxs.ARBITRARY_MATRIX_TYPES so we can see the code works, # but not use too much resources diff --git a/tests/regression_tests/mgxs_photon/test.py b/tests/regression_tests/mgxs_photon/test.py index b7ea928c7c7..d52dc1e9723 100644 --- a/tests/regression_tests/mgxs_photon/test.py +++ b/tests/regression_tests/mgxs_photon/test.py @@ -2,42 +2,52 @@ import numpy as np import openmc +import pytest from tests.testing_harness import PyAPITestHarness +@pytest.fixture +def model(): + model = openmc.Model() + + material = openmc.Material() + material.set_density('g/cm3', 11.35) + material.add_element('Pb', 1.0) + model.materials.append(material) + + sphere = openmc.Sphere(r=1.0, boundary_type='vacuum') + cell = openmc.Cell(fill=material, region=-sphere) + model.geometry = openmc.Geometry([cell]) + + model.settings.run_mode = 'fixed source' + model.settings.particles = 2000 + model.settings.batches = 2 + model.settings.photon_transport = True + model.settings.atomic_relaxation = True + model.settings.electron_treatment = 'ttb' + model.settings.cutoff = {'energy_photon': 1000.0} + model.settings.source = openmc.IndependentSource( + particle='photon', + space=openmc.stats.Point((0.0, 0.0, 0.0)), + energy=openmc.stats.Discrete([1.0e6], [1.0])) + + return model + + class PhotonMGXSTestHarness(PyAPITestHarness): """Run photon transport and verify photon MGXS post-processing.""" def __init__(self, *args, **kwargs): super().__init__(*args, **kwargs) - self.material = openmc.Material() - self.material.set_density('g/cm3', 11.35) - self.material.add_element('Pb', 1.0) - self._model.materials = openmc.Materials([self.material]) - - sphere = openmc.Sphere(r=1.0, boundary_type='vacuum') - cell = openmc.Cell(fill=self.material, region=-sphere) - self._model.geometry = openmc.Geometry([cell]) - - self._model.settings.run_mode = 'fixed source' - self._model.settings.particles = 2000 - self._model.settings.batches = 2 - self._model.settings.photon_transport = True - self._model.settings.atomic_relaxation = True - self._model.settings.electron_treatment = 'ttb' - self._model.settings.cutoff = {'energy_photon': 1000.0} - self._model.settings.source = openmc.IndependentSource( - particle='photon', - space=openmc.stats.Point((0.0, 0.0, 0.0)), - energy=openmc.stats.Discrete([1.0e6], [1.0])) + self.material = self._model.materials[0] groups = openmc.mgxs.EnergyGroups( group_edges=[1.0e3, 1.0e5, 5.0e5, 1.1e6]) self.mgxs_lib = openmc.mgxs.Library( self._model.geometry, - mgxs_types=['total', 'absorption', 'photon-production matrix'], + mgxs_types=['total', 'absorption', 'nu-scatter matrix'], particle_type='photon') self.mgxs_lib.energy_groups = groups self.mgxs_lib.correction = None @@ -50,7 +60,7 @@ def _get_results(self, hash_output=False): self.mgxs_lib.load_from_statepoint(statepoint) production = self.mgxs_lib.get_mgxs( - self.material, 'photon-production matrix') + self.material, 'nu-scatter matrix') absorption = self.mgxs_lib.get_mgxs(self.material, 'absorption') primary = production.tallies[ @@ -96,7 +106,6 @@ def _get_results(self, hash_output=False): return output -def test_photon_mgxs(): - harness = PhotonMGXSTestHarness( - 'statepoint.2.h5', model=openmc.Model()) +def test_photon_mgxs(model): + harness = PhotonMGXSTestHarness('statepoint.2.h5', model) harness.main() diff --git a/tests/unit_tests/test_photon_mgxs.py b/tests/unit_tests/test_photon_mgxs.py index 08a8eccc07d..175ac245422 100644 --- a/tests/unit_tests/test_photon_mgxs.py +++ b/tests/unit_tests/test_photon_mgxs.py @@ -29,13 +29,14 @@ def test_photon_production_matrix_tallies(): def test_photon_production_matrix_constraints(): with pytest.raises(ValueError, match='by nuclide'): openmc.mgxs.MGXS.get_mgxs( - 'photon-production matrix', by_nuclide=True) + 'nu-scatter matrix', by_nuclide=True, particle_type='photon') with pytest.raises(ValueError, match='isotropic'): openmc.mgxs.MGXS.get_mgxs( - 'photon-production matrix', num_polar=2) - with pytest.raises(ValueError, match='photon tallies'): - openmc.mgxs.MGXS.get_mgxs( - 'photon-production matrix', particle_type='neutron') + 'nu-scatter matrix', num_polar=2, particle_type='photon') + + neutron_production = openmc.mgxs.MGXS.get_mgxs( + 'nu-scatter matrix', particle_type='neutron') + assert isinstance(neutron_production, openmc.mgxs.ScatterMatrixXS) def test_photon_production_matrix_combines_primary_and_secondary(): @@ -57,7 +58,7 @@ def test_photon_production_matrix_combines_primary_and_secondary(): production.rxn_rate_tally.mean.ravel(), [2.0, 4.0, 6.0, 8.0]) -def test_set_photon_production_mgxs(monkeypatch): +def test_set_photon_nu_scatter_mgxs(monkeypatch): material = openmc.Material() groups = openmc.mgxs.EnergyGroups([1.0, 10.0, 100.0]) production = openmc.mgxs.PhotonProductionMatrixXS( @@ -69,7 +70,7 @@ def test_set_photon_production_mgxs(monkeypatch): xsdata.order = 0 absorption = np.array([0.5, 0.6]) xsdata.set_absorption(absorption) - xsdata.set_photon_production_mgxs(production) + xsdata.set_scatter_matrix_mgxs(production) assert xsdata.scatter_format == 'legendre' assert xsdata.order == 0 @@ -84,7 +85,7 @@ def test_photon_mgxs_library(tmp_path): groups = openmc.mgxs.EnergyGroups([1.0, 10.0, 100.0]) library = openmc.mgxs.Library( geometry, mgxs_types=[ - 'total', 'absorption', 'photon-production matrix'], + 'total', 'absorption', 'nu-scatter matrix'], particle_type='photon') library.domain_type = 'material' library.energy_groups = groups @@ -110,6 +111,5 @@ def test_all_mgxs_types_respects_particle_type(): photon_library = openmc.mgxs.Library( geometry, mgxs_types='all', particle_type='photon') - assert 'photon-production matrix' not in neutron_library.mgxs_types assert photon_library.mgxs_types == ( - 'total', 'absorption', 'photon-production matrix') + 'total', 'absorption', 'nu-scatter matrix') From e9dd574b3ce70efc625bc9ec944840480e418e55 Mon Sep 17 00:00:00 2001 From: Paul Romano Date: Sun, 16 Aug 2026 16:12:12 -0500 Subject: [PATCH 4/6] Simplify mgxs_photon test --- tests/regression_tests/mgxs_photon/inputs_true.dat | 8 -------- tests/regression_tests/mgxs_photon/test.py | 6 +----- 2 files changed, 1 insertion(+), 13 deletions(-) diff --git a/tests/regression_tests/mgxs_photon/inputs_true.dat b/tests/regression_tests/mgxs_photon/inputs_true.dat index aac04060864..130deaa0698 100644 --- a/tests/regression_tests/mgxs_photon/inputs_true.dat +++ b/tests/regression_tests/mgxs_photon/inputs_true.dat @@ -18,19 +18,11 @@ 2000 2 - - 0.0 0.0 0.0 - 1000000.0 1.0 - ttb - true true - - 1000.0 - diff --git a/tests/regression_tests/mgxs_photon/test.py b/tests/regression_tests/mgxs_photon/test.py index d52dc1e9723..e2d93a7e86c 100644 --- a/tests/regression_tests/mgxs_photon/test.py +++ b/tests/regression_tests/mgxs_photon/test.py @@ -24,13 +24,9 @@ def model(): model.settings.particles = 2000 model.settings.batches = 2 model.settings.photon_transport = True - model.settings.atomic_relaxation = True - model.settings.electron_treatment = 'ttb' - model.settings.cutoff = {'energy_photon': 1000.0} model.settings.source = openmc.IndependentSource( particle='photon', - space=openmc.stats.Point((0.0, 0.0, 0.0)), - energy=openmc.stats.Discrete([1.0e6], [1.0])) + energy=openmc.stats.delta_function(1.0e6)) return model From 2f954c1201c0782ce56ee94c432a09e41426e978 Mon Sep 17 00:00:00 2001 From: Paul Romano Date: Sun, 16 Aug 2026 23:49:46 -0500 Subject: [PATCH 5/6] Simplify unit tests --- tests/unit_tests/test_photon_mgxs.py | 60 ---------------------------- 1 file changed, 60 deletions(-) diff --git a/tests/unit_tests/test_photon_mgxs.py b/tests/unit_tests/test_photon_mgxs.py index 175ac245422..3bdd79c764f 100644 --- a/tests/unit_tests/test_photon_mgxs.py +++ b/tests/unit_tests/test_photon_mgxs.py @@ -1,7 +1,3 @@ -import h5py -import numpy as np -import pytest - import openmc @@ -26,59 +22,6 @@ def test_photon_production_matrix_tallies(): assert particle_filter.bins == ['photon'] -def test_photon_production_matrix_constraints(): - with pytest.raises(ValueError, match='by nuclide'): - openmc.mgxs.MGXS.get_mgxs( - 'nu-scatter matrix', by_nuclide=True, particle_type='photon') - with pytest.raises(ValueError, match='isotropic'): - openmc.mgxs.MGXS.get_mgxs( - 'nu-scatter matrix', num_polar=2, particle_type='photon') - - neutron_production = openmc.mgxs.MGXS.get_mgxs( - 'nu-scatter matrix', particle_type='neutron') - assert isinstance(neutron_production, openmc.mgxs.ScatterMatrixXS) - - -def test_photon_production_matrix_combines_primary_and_secondary(): - material = openmc.Material() - groups = openmc.mgxs.EnergyGroups([1.0, 10.0, 100.0]) - production = openmc.mgxs.PhotonProductionMatrixXS( - material, 'material', groups) - - for tally in production.tallies.values(): - values = np.arange(1, tally.num_filter_bins + 1, dtype=float) - values.shape = (tally.num_filter_bins, 1, 1) - tally._sum = 2.0 * values - tally._sum_sq = 2.0 * values**2 - tally._num_realizations = 2 - tally._sp_filename = 'statepoint.h5' - tally._results_read = True - - np.testing.assert_allclose( - production.rxn_rate_tally.mean.ravel(), [2.0, 4.0, 6.0, 8.0]) - - -def test_set_photon_nu_scatter_mgxs(monkeypatch): - material = openmc.Material() - groups = openmc.mgxs.EnergyGroups([1.0, 10.0, 100.0]) - production = openmc.mgxs.PhotonProductionMatrixXS( - material, 'material', groups) - values = np.array([[0.1, 0.2], [0.3, 0.4]]) - monkeypatch.setattr(production, 'get_xs', lambda **kwargs: values) - - xsdata = openmc.XSdata('photon', groups) - xsdata.order = 0 - absorption = np.array([0.5, 0.6]) - xsdata.set_absorption(absorption) - xsdata.set_scatter_matrix_mgxs(production) - - assert xsdata.scatter_format == 'legendre' - assert xsdata.order == 0 - np.testing.assert_allclose(xsdata._scatter_matrix[0][:, :, 0], values) - np.testing.assert_allclose(xsdata._absorption[0], absorption) - assert xsdata._multiplicity_matrix[0] is None - - def test_photon_mgxs_library(tmp_path): material = openmc.Material() geometry = openmc.Geometry([openmc.Cell(fill=material)]) @@ -99,15 +42,12 @@ def test_photon_mgxs_library(tmp_path): mg_library = openmc.MGXSLibrary(groups, particle_type='photon') path = tmp_path / 'mgxs.h5' mg_library.export_to_hdf5(path) - with h5py.File(path) as h5file: - assert h5file.attrs['particle_type'] == b'photon' assert openmc.MGXSLibrary.from_hdf5(path).particle_type == \ openmc.ParticleType.PHOTON def test_all_mgxs_types_respects_particle_type(): geometry = openmc.Geometry([openmc.Cell()]) - neutron_library = openmc.mgxs.Library(geometry, mgxs_types='all') photon_library = openmc.mgxs.Library( geometry, mgxs_types='all', particle_type='photon') From 33c379303e1a0cb2b2e91e42b0608afa1846e153 Mon Sep 17 00:00:00 2001 From: Paul Romano Date: Sun, 23 Aug 2026 09:24:59 -0500 Subject: [PATCH 6/6] Fix slicing and condensation on photon production MGXS --- openmc/mgxs/library.py | 3 ++ openmc/mgxs/mgxs.py | 56 ++++++++++++++++++++++++++++ openmc/tallies.py | 2 +- tests/unit_tests/test_photon_mgxs.py | 49 +++++++++++++++++++++++- 4 files changed, 108 insertions(+), 2 deletions(-) diff --git a/openmc/mgxs/library.py b/openmc/mgxs/library.py index 6338c7d732d..08cf5350d90 100644 --- a/openmc/mgxs/library.py +++ b/openmc/mgxs/library.py @@ -168,6 +168,7 @@ def __deepcopy__(self, memo): clone._sp_filename = self._sp_filename clone._keff = self._keff clone._sparse = self.sparse + clone._estimator = self.estimator clone._all_mgxs = {} for domain in self.domains: @@ -248,6 +249,8 @@ def particle_type(self, particle_type): (openmc.ParticleType.NEUTRON, openmc.ParticleType.PHOTON)) self._particle_type = particle_type + if particle_type == openmc.ParticleType.PHOTON: + self._correction = None @property def domain_type(self): diff --git a/openmc/mgxs/mgxs.py b/openmc/mgxs/mgxs.py index b4eeae63045..e7ba529ecf7 100644 --- a/openmc/mgxs/mgxs.py +++ b/openmc/mgxs/mgxs.py @@ -5040,6 +5040,62 @@ def scatter_format(self): def legendre_order(self): return 0 + def _with_energyout_filter(self): + """Return a copy using an outgoing-energy production filter.""" + production = copy.deepcopy(self) + secondary = production.tallies['secondary photon production'] + if not secondary.contains_filter(openmc.ParticleProductionFilter): + return production + particle_production = secondary.find_filter( + openmc.ParticleProductionFilter) + energyout = openmc.EnergyoutFilter(particle_production.energies) + secondary.filters = [ + energyout + if isinstance(f, openmc.ParticleProductionFilter) else f + for f in secondary.filters + ] + return production + + def get_condensed_xs(self, coarse_groups): + """Construct an energy-condensed photon production matrix. + + Parameters + ---------- + coarse_groups : openmc.mgxs.EnergyGroups + Coarse energy group structure + + Returns + ------- + openmc.mgxs.PhotonProductionMatrixXS + Photon production matrix condensed to the coarse group structure + + """ + production = self._with_energyout_filter() + return MGXS.get_condensed_xs(production, coarse_groups) + + def get_slice(self, nuclides=[], in_groups=[], out_groups=[]): + """Build a sliced photon production matrix. + + Parameters + ---------- + nuclides : list of str + Nuclides to include; photon production only supports macroscopic + data + in_groups : list of int + Incoming energy groups to include + out_groups : list of int + Outgoing energy groups to include + + Returns + ------- + openmc.mgxs.PhotonProductionMatrixXS + Sliced photon production matrix + + """ + production = self._with_energyout_filter() + return MatrixMGXS.get_slice( + production, nuclides, in_groups, out_groups) + @property def rxn_rate_tally(self): if self._rxn_rate_tally is None: diff --git a/openmc/tallies.py b/openmc/tallies.py index a503d04bc91..4393dc90b71 100644 --- a/openmc/tallies.py +++ b/openmc/tallies.py @@ -3383,7 +3383,7 @@ def get_slice(self, scores=[], filters=[], filter_bins=[], nuclides=[], # Replace existing filter with new one for j, test_filter in enumerate(new_tally.filters): - if isinstance(test_filter, filter_type): + if type(test_filter) is filter_type: new_tally.filters[j] = new_filter # If original tally was sparse, sparsify the sliced tally diff --git a/tests/unit_tests/test_photon_mgxs.py b/tests/unit_tests/test_photon_mgxs.py index 3bdd79c764f..a1f1e3eb0cc 100644 --- a/tests/unit_tests/test_photon_mgxs.py +++ b/tests/unit_tests/test_photon_mgxs.py @@ -1,6 +1,16 @@ +import copy + +import numpy as np import openmc +def _add_tally_data(mgxs): + for tally in mgxs.tallies.values(): + tally._derived = True + tally._mean = np.ones(tally.shape) + tally._std_dev = np.zeros(tally.shape) + + def test_photon_production_matrix_tallies(): material = openmc.Material() groups = openmc.mgxs.EnergyGroups([1.0, 10.0, 100.0]) @@ -32,12 +42,23 @@ def test_photon_mgxs_library(tmp_path): particle_type='photon') library.domain_type = 'material' library.energy_groups = groups - library.correction = None library.build_library() library.check_library_for_openmc_mgxs() + assert library.correction is None + + clone = copy.deepcopy(library) + clone.check_library_for_openmc_mgxs() + assert clone.estimator is None for mgxs in library.all_mgxs[material.id].values(): assert mgxs.particle_type == openmc.ParticleType.PHOTON + _add_tally_data(mgxs) + + library._sp_filename = 'statepoint.h5' + coarse_groups = openmc.mgxs.EnergyGroups([1.0, 100.0]) + condensed = library.get_condensed_library(coarse_groups) + xsdata = condensed.create_mg_library().xsdatas[0] + assert xsdata.scatter_matrix[0].shape == (1, 1, 1) mg_library = openmc.MGXSLibrary(groups, particle_type='photon') path = tmp_path / 'mgxs.h5' @@ -53,3 +74,29 @@ def test_all_mgxs_types_respects_particle_type(): assert photon_library.mgxs_types == ( 'total', 'absorption', 'nu-scatter matrix') + assert photon_library.correction is None + + +def test_photon_production_matrix_group_transformations(): + material = openmc.Material() + groups = openmc.mgxs.EnergyGroups([1.0, 10.0, 100.0, 1000.0]) + production = openmc.mgxs.PhotonProductionMatrixXS( + material, 'material', groups) + _add_tally_data(production) + + coarse_groups = openmc.mgxs.EnergyGroups([1.0, 100.0, 1000.0]) + condensed = production.get_condensed_xs(coarse_groups) + assert condensed.get_xs().shape == (2, 2) + assert condensed.tallies['secondary photon production'].contains_filter( + openmc.EnergyoutFilter) + + sliced = production.get_slice(in_groups=[1, 2], out_groups=[1, 2]) + assert sliced.get_xs().shape == (2, 2) + assert sliced.tallies['secondary photon production'].contains_filter( + openmc.EnergyoutFilter) + assert production.tallies[ + 'secondary photon production'].contains_filter( + openmc.ParticleProductionFilter) + + resliced = sliced.get_slice(in_groups=[1], out_groups=[1]) + assert resliced.get_xs().shape == (1, 1)