From a7594ec16911e6256ad18874d4a8f25ca4a963bf Mon Sep 17 00:00:00 2001 From: Guillaume Giudicelli Date: Fri, 28 Aug 2026 16:52:21 -0600 Subject: [PATCH 1/3] Add per-domain transport correction ratios to mgxs.Library Add a transport_correction_ratios attribute to openmc.mgxs.Library that lets users supply a per-domain, per-group transport correction to apply when the library is written to an openmc.MGXSLibrary. The ratios are a nested dict keyed by domain type and then domain ID, with one ratio r_g = sigma_tr,g / sigma_t,g per energy group. For each listed domain the total is replaced by the transport-corrected value r_g * sigma_t,g and the same correction (1 - r_g) * sigma_t,g is subtracted from the in-group P0 element of the scattering matrix, which preserves the absorption balance. This is an alternative to the tally-based 'P0' correction, so it requires correction=None, a 'legendre' scatter_format, and a 'total' MGXS type; domains without an entry are left uncorrected. Both isotropic and 'angle' representations are supported. Add unit tests covering the setter validation/normalization and the application logic (total, scatter diagonal, absorption balance, angle representation, no-op without an entry, and mutual exclusivity with the 'P0' correction). Co-Authored-By: Claude --- openmc/mgxs/library.py | 210 +++++++++++++++++- .../test_mgxs_transport_correction.py | 185 +++++++++++++++ 2 files changed, 394 insertions(+), 1 deletion(-) create mode 100644 tests/unit_tests/test_mgxs_transport_correction.py diff --git a/openmc/mgxs/library.py b/openmc/mgxs/library.py index faa83c0481f..497b880e2dd 100644 --- a/openmc/mgxs/library.py +++ b/openmc/mgxs/library.py @@ -1,4 +1,4 @@ -from collections.abc import Iterable +from collections.abc import Iterable, Mapping import copy from numbers import Integral import os @@ -55,6 +55,25 @@ class Library: The spatial domain(s) for which MGXS in the Library are computed correction : {'P0', None} Apply the P0 correction to scattering matrices if set to 'P0' + transport_correction_ratios : dict or None + An optional, user-supplied transport correction to apply when the + library is written to an :class:`openmc.MGXSLibrary`. This is a nested + dictionary keyed first by domain type (e.g., ``'material'``) and then + by domain ID, whose values are iterables of per-group transport + correction ratios :math:`r_g = \\sigma_{tr,g} / \\sigma_{t,g}` (one + ratio per energy group, ordered from fast to thermal to match the + :class:`openmc.XSdata` group indexing). For each listed domain the + transport-corrected total cross section is computed as + :math:`\\sigma_{tr,g} = r_g \\sigma_{t,g}` and the same correction + :math:`\\Delta_g = (1 - r_g)\\sigma_{t,g}` is subtracted from the + in-group (diagonal) :math:`P_0` element of the scattering matrix, + preserving the absorption balance. This is an alternative to the + tally-based ``'P0'`` correction and therefore requires + :attr:`correction` to be ``None``, :attr:`scatter_format` to be + ``'legendre'``, and a ``'total'`` MGXS type to be present. Domains + without an entry are left uncorrected. Defaults to ``None``. + + .. versionadded:: 0.16.1 scatter_format : {'legendre', 'histogram'} Representation of the angular scattering distribution (default is 'legendre') @@ -116,6 +135,7 @@ def __init__(self, geometry, by_nuclide=False, self._nuclides = None self._num_delayed_groups = 0 self._correction = 'P0' + self._transport_correction_ratios = None self._scatter_format = 'legendre' self._legendre_order = 0 self._histogram_bins = 16 @@ -146,6 +166,8 @@ def __deepcopy__(self, memo): clone._domain_type = self.domain_type clone._domains = copy.deepcopy(self.domains) clone._correction = self.correction + clone._transport_correction_ratios = \ + copy.deepcopy(self._transport_correction_ratios) clone._scatter_format = self.scatter_format clone._legendre_order = self.legendre_order clone._histogram_bins = self.histogram_bins @@ -367,6 +389,58 @@ def correction(self, correction): self._correction = correction + @property + def transport_correction_ratios(self): + return self._transport_correction_ratios + + @transport_correction_ratios.setter + def transport_correction_ratios(self, ratios): + if ratios is None: + self._transport_correction_ratios = None + return + + cv.check_type('transport_correction_ratios', ratios, Mapping) + + normalized = {} + for domain_type, domain_ratios in ratios.items(): + cv.check_value('transport_correction_ratios domain type', + domain_type, openmc.mgxs.DOMAIN_TYPES) + cv.check_type(f'transport_correction_ratios["{domain_type}"]', + domain_ratios, Mapping) + + normalized[domain_type] = {} + for domain_id, group_ratios in domain_ratios.items(): + cv.check_type('transport correction ratio domain ID', + domain_id, Integral) + + try: + arr = np.asarray(group_ratios, dtype=float) + except (ValueError, TypeError): + raise ValueError( + 'Transport correction ratios for domain ' + f'{domain_id} must be a 1-D iterable of real numbers.') + + if arr.ndim != 1: + raise ValueError( + 'Transport correction ratios for domain ' + f'{domain_id} must be a 1-D iterable (one ratio per ' + 'energy group).') + + if np.any(arr <= 0.0): + raise ValueError( + 'Transport correction ratios for domain ' + f'{domain_id} must be positive.') + + # Validate the number of groups when it is already known + if self._energy_groups is not None: + cv.check_length( + 'transport correction ratios for domain ' + f'{domain_id}', arr, self.num_groups, self.num_groups) + + normalized[domain_type][int(domain_id)] = arr + + self._transport_correction_ratios = normalized + @property def scatter_format(self): return self._scatter_format @@ -954,6 +1028,102 @@ def load_from_file(filename='mgxs', directory='mgxs'): with open(full_filename, 'rb') as f: return pickle.load(f) + def _get_transport_correction_ratios(self, domain): + """Return the per-group transport correction ratios for a domain. + + Parameters + ---------- + domain : openmc.Material or openmc.Cell or openmc.Universe or openmc.RegularMesh + The domain of interest + + Returns + ------- + numpy.ndarray or None + The transport correction ratios (one per energy group) for the + domain, or None if none were provided. + + """ + + if not self._transport_correction_ratios: + return None + + domain_ratios = self._transport_correction_ratios.get(self.domain_type) + if not domain_ratios: + return None + + return domain_ratios.get(domain.id) + + def _apply_transport_correction_ratios(self, xsdata, domain, temperature): + """Apply user-supplied transport correction ratios to an XSdata object. + + For each energy group ``g``, the transport correction + :math:`\\Delta_g = (1 - r_g)\\sigma_{t,g}` is computed from the plain + total cross section :math:`\\sigma_{t,g}` and the user-supplied ratio + :math:`r_g`. The total cross section is replaced by the + transport-corrected value :math:`r_g \\sigma_{t,g}`, and the same + :math:`\\Delta_g` is subtracted from the in-group (diagonal) + :math:`P_0` element of the scattering matrix. Subtracting the same + correction from both quantities leaves the absorption balance + unchanged. + + Parameters + ---------- + xsdata : openmc.XSdata + The dataset to correct in place. Its total cross section is + expected to be the plain (un-corrected) total. + domain : openmc.Material or openmc.Cell or openmc.Universe or openmc.RegularMesh + The domain the dataset describes + temperature : float + Temperature (in Kelvin) of the data to correct + + """ + + ratios = self._get_transport_correction_ratios(domain) + if ratios is None: + return + + # The ratio-based correction takes the place of the tally-based 'P0' + # correction, so the total and scattering matrix must be un-corrected. + if self.correction is not None: + raise ValueError( + 'The "correction" parameter must be None when ' + 'transport_correction_ratios are provided, otherwise the ' + 'transport correction would be applied twice.') + + if self.scatter_format != 'legendre': + raise ValueError( + 'transport_correction_ratios require a "legendre" ' + 'scatter_format.') + + if len(ratios) != self.num_groups: + raise ValueError( + f'Expected {self.num_groups} transport correction ratios for ' + f'domain {domain.id} but got {len(ratios)}.') + + i = xsdata._temperature_index(temperature) + + total = xsdata._total[i] + if total is None: + raise ValueError( + 'A "total" MGXS type is required to apply ' + 'transport_correction_ratios.') + + # Compute the correction from the plain total cross section. ratios has + # shape (G,) and broadcasts against the trailing group axis of the + # total cross section for both isotropic and angle representations. + delta = (1.0 - ratios) * total + xsdata._total[i] = total - delta + + # Subtract the same correction from the in-group (diagonal) P0 element + # of the scattering matrix to preserve the absorption balance. + scatter_matrix = xsdata._scatter_matrix[i] + if scatter_matrix is not None: + groups = np.arange(self.num_groups) + if xsdata.representation == 'angle': + scatter_matrix[:, :, groups, groups, 0] -= delta + else: + scatter_matrix[groups, groups, 0] -= delta + def get_xsdata(self, domain, xsdata_name, nuclide='total', xs_type='macro', subdomain=None, apply_domain_chi=False, temperature=ROOM_TEMPERATURE_KELVIN): """Generates an openmc.XSdata object describing a multi-group cross section @@ -1303,6 +1473,10 @@ def get_xsdata(self, domain, xsdata_name, nuclide='total', xs_type='macro', nuclide=[nuclide], subdomain=subdomain) + # Apply any user-supplied transport correction ratios to the total + # cross section and the scattering matrix diagonal + self._apply_transport_correction_ratios(xsdata, domain, temperature) + return xsdata def create_mg_library(self, xs_type='macro', xsdata_names=None, @@ -1639,5 +1813,39 @@ def check_library_for_openmc_mgxs(self): error_flag = True warn('An "absorption" MGXS type is required but not provided.') + # Validate user-supplied transport correction ratios + if self._transport_correction_ratios: + if self.correction is not None: + error_flag = True + warn('The "correction" parameter must be None when ' + 'transport_correction_ratios are provided.') + if self.scatter_format != 'legendre': + error_flag = True + warn('transport_correction_ratios require a "legendre" ' + 'scatter_format.') + if 'total' not in self.mgxs_types: + error_flag = True + warn('A "total" MGXS type is required when ' + 'transport_correction_ratios are provided.') + + domain_ratios = \ + self._transport_correction_ratios.get(self.domain_type) + if not domain_ratios: + warn('The transport_correction_ratios do not contain any ' + f'entries for the "{self.domain_type}" domain type, so ' + 'no transport correction will be applied.') + else: + domain_ids = [domain.id for domain in self.domains] + for domain_id, ratios in domain_ratios.items(): + if domain_id not in domain_ids: + warn(f'Domain {domain_id} in ' + 'transport_correction_ratios is not in the ' + 'Library and will be ignored.') + elif len(ratios) != self.num_groups: + error_flag = True + warn(f'Expected {self.num_groups} transport ' + f'correction ratios for domain {domain_id} but ' + f'got {len(ratios)}.') + if error_flag: raise ValueError('Invalid MGXS configuration encountered.') diff --git a/tests/unit_tests/test_mgxs_transport_correction.py b/tests/unit_tests/test_mgxs_transport_correction.py new file mode 100644 index 00000000000..9704ec01872 --- /dev/null +++ b/tests/unit_tests/test_mgxs_transport_correction.py @@ -0,0 +1,185 @@ +"""Tests for user-supplied transport correction ratios in openmc.mgxs.Library.""" + +import numpy as np +import pytest + +import openmc +import openmc.mgxs + + +@pytest.fixture +def simple_geometry(): + openmc.reset_auto_ids() + mat = openmc.Material(material_id=1) + mat.add_nuclide('U235', 1.0) + mat.set_density('g/cm3', 10.0) + sph = openmc.Sphere(r=1.0, boundary_type='vacuum') + cell = openmc.Cell(fill=mat, region=-sph) + return openmc.Geometry([cell]), mat + + +@pytest.fixture +def library(simple_geometry): + geometry, _ = simple_geometry + groups = openmc.mgxs.EnergyGroups(group_edges=[0.0, 0.625, 2.0e7]) + lib = openmc.mgxs.Library(geometry) + lib.energy_groups = groups + lib.domain_type = 'material' + lib.correction = None + lib.scatter_format = 'legendre' + return lib + + +def test_setter_accepts_and_normalizes(library): + library.transport_correction_ratios = {'material': {1: [0.9, 0.8]}} + stored = library.transport_correction_ratios + assert set(stored) == {'material'} + assert list(stored['material']) == [1] + assert isinstance(stored['material'][1], np.ndarray) + np.testing.assert_allclose(stored['material'][1], [0.9, 0.8]) + + +def test_setter_none_clears(library): + library.transport_correction_ratios = {'material': {1: [0.9, 0.8]}} + library.transport_correction_ratios = None + assert library.transport_correction_ratios is None + + +def test_setter_validation(library): + # Not a mapping + with pytest.raises(TypeError): + library.transport_correction_ratios = [0.9, 0.8] + + # Invalid domain type + with pytest.raises(ValueError): + library.transport_correction_ratios = {'banana': {1: [0.9, 0.8]}} + + # Non-integer domain ID + with pytest.raises(TypeError): + library.transport_correction_ratios = {'material': {'1': [0.9, 0.8]}} + + # Wrong number of groups + with pytest.raises(ValueError): + library.transport_correction_ratios = {'material': {1: [0.9, 0.8, 0.7]}} + + # Non-positive ratio + with pytest.raises(ValueError): + library.transport_correction_ratios = {'material': {1: [0.9, -0.1]}} + + +def _make_xsdata(groups, sigma_t, scatter, representation='isotropic', + num_polar=1, num_azimuthal=1): + xsdata = openmc.XSdata('set1', groups, representation=representation) + xsdata.order = 0 + if representation == 'angle': + xsdata.num_polar = num_polar + xsdata.num_azimuthal = num_azimuthal + xsdata.set_total(sigma_t) + xsdata.set_scatter_matrix(scatter) + return xsdata + + +def test_apply_isotropic(library, simple_geometry): + _, mat = simple_geometry + library.transport_correction_ratios = {'material': {1: [0.9, 0.8]}} + + sigma_t = np.array([2.0, 3.0]) + scatter = np.array([[[0.5], [0.3]], [[0.1], [1.2]]]) + scatter_orig = scatter.copy() + absorption = sigma_t - scatter_orig[:, :, 0].sum(axis=1) + + xsdata = _make_xsdata(library.energy_groups, sigma_t, scatter) + library._apply_transport_correction_ratios(xsdata, mat, 294.0) + + ratios = np.array([0.9, 0.8]) + delta = (1.0 - ratios) * sigma_t + + # Total is transport-corrected + np.testing.assert_allclose(xsdata._total[0], ratios * sigma_t) + + # In-group P0 diagonal reduced by delta, off-diagonal unchanged + sm = xsdata._scatter_matrix[0] + np.testing.assert_allclose(sm[0, 0, 0], scatter_orig[0, 0, 0] - delta[0]) + np.testing.assert_allclose(sm[1, 1, 0], scatter_orig[1, 1, 0] - delta[1]) + np.testing.assert_allclose(sm[0, 1, 0], scatter_orig[0, 1, 0]) + np.testing.assert_allclose(sm[1, 0, 0], scatter_orig[1, 0, 0]) + + # Absorption balance (total - out-scatter) is unchanged + new_absorption = xsdata._total[0] - sm[:, :, 0].sum(axis=1) + np.testing.assert_allclose(new_absorption, absorption) + + +def test_apply_angle(simple_geometry): + geometry, mat = simple_geometry + groups = openmc.mgxs.EnergyGroups(group_edges=[0.0, 0.625, 2.0e7]) + lib = openmc.mgxs.Library(geometry) + lib.energy_groups = groups + lib.domain_type = 'material' + lib.correction = None + lib.num_polar = 2 + lib.num_azimuthal = 2 + lib.transport_correction_ratios = {'material': {1: [0.9, 0.8]}} + + sigma_t = np.array([2.0, 3.0]) + scatter2d = np.array([[[0.5], [0.3]], [[0.1], [1.2]]]) + total = np.empty((2, 2, 2)) + total[...] = sigma_t + scatter = np.zeros((2, 2, 2, 2, 1)) + scatter[...] = scatter2d + scatter_orig = scatter2d.copy() + + xsdata = _make_xsdata(groups, total, scatter, representation='angle', + num_polar=2, num_azimuthal=2) + lib._apply_transport_correction_ratios(xsdata, mat, 294.0) + + ratios = np.array([0.9, 0.8]) + delta = (1.0 - ratios) * sigma_t + sm = xsdata._scatter_matrix[0] + np.testing.assert_allclose(xsdata._total[0], np.broadcast_to( + ratios * sigma_t, (2, 2, 2))) + np.testing.assert_allclose(sm[0, 0, 0, 0, 0], scatter_orig[0, 0, 0] - delta[0]) + np.testing.assert_allclose(sm[1, 1, 1, 1, 0], scatter_orig[1, 1, 0] - delta[1]) + np.testing.assert_allclose(sm[:, :, 0, 1, 0], scatter_orig[0, 1, 0]) + + +def test_apply_noop_without_entry(library, simple_geometry): + _, mat = simple_geometry + library.transport_correction_ratios = {'material': {999: [0.9, 0.8]}} + + sigma_t = np.array([2.0, 3.0]) + scatter = np.array([[[0.5], [0.3]], [[0.1], [1.2]]]) + xsdata = _make_xsdata(library.energy_groups, sigma_t, scatter) + library._apply_transport_correction_ratios(xsdata, mat, 294.0) + + # Domain 1 has no entry, so nothing is changed + np.testing.assert_allclose(xsdata._total[0], sigma_t) + + +def test_apply_requires_correction_none(library, simple_geometry): + _, mat = simple_geometry + library.correction = 'P0' + library.transport_correction_ratios = {'material': {1: [0.9, 0.8]}} + + sigma_t = np.array([2.0, 3.0]) + scatter = np.array([[[0.5], [0.3]], [[0.1], [1.2]]]) + xsdata = _make_xsdata(library.energy_groups, sigma_t, scatter) + with pytest.raises(ValueError, match='correction'): + library._apply_transport_correction_ratios(xsdata, mat, 294.0) + + +def test_check_library_rejects_p0(library): + library.mgxs_types = ['total', 'absorption', 'nu-scatter matrix', + 'scatter matrix'] + library.correction = 'P0' + library.transport_correction_ratios = {'material': {1: [0.9, 0.8]}} + with pytest.raises(ValueError, match='Invalid MGXS configuration'): + library.check_library_for_openmc_mgxs() + + +def test_check_library_warns_missing_domain_type(library): + library.mgxs_types = ['total', 'absorption', 'nu-scatter matrix', + 'scatter matrix'] + # Ratios provided for a domain type the library does not use + library.transport_correction_ratios = {'cell': {1: [0.9, 0.8]}} + with pytest.warns(UserWarning, match='do not contain any entries'): + library.check_library_for_openmc_mgxs() From fb0733740c0cbe81f69a84e983bd978091782132 Mon Sep 17 00:00:00 2001 From: Guillaume Giudicelli Date: Sun, 30 Aug 2026 20:29:41 -0600 Subject: [PATCH 2/3] Populate transport correction ratios during tally processing When a 'transport' or 'nu-transport' MGXS type has been tallied, Library.load_from_statepoint now computes the per-group transport correction ratios r_g = sigma_tr,g / sigma_t,g for each domain and stores them in the transport_correction_ratios attribute. The transport-corrected total is the transport MGXS itself; the plain total is recovered from the same MGXS's total and flux tallies. User-supplied entries are preserved and never overwritten, and only isotropic, single-subdomain data is handled. Because the ratios are now populated automatically alongside the standard 'P0' workflow, applying them is gated on correction being None: when correction is 'P0' the tally-based correction is already reflected in the total and scattering matrix, so the stored ratios are recorded but not applied (rather than raising). check_library_for_openmc_mgxs validates the ratios only in the case where they are applied. Co-Authored-By: Claude --- openmc/mgxs/library.py | 172 +++++++++++++++--- .../test_mgxs_transport_correction.py | 59 +++++- 2 files changed, 196 insertions(+), 35 deletions(-) diff --git a/openmc/mgxs/library.py b/openmc/mgxs/library.py index 497b880e2dd..e69b8f67d3b 100644 --- a/openmc/mgxs/library.py +++ b/openmc/mgxs/library.py @@ -56,22 +56,28 @@ class Library: correction : {'P0', None} Apply the P0 correction to scattering matrices if set to 'P0' transport_correction_ratios : dict or None - An optional, user-supplied transport correction to apply when the - library is written to an :class:`openmc.MGXSLibrary`. This is a nested - dictionary keyed first by domain type (e.g., ``'material'``) and then - by domain ID, whose values are iterables of per-group transport - correction ratios :math:`r_g = \\sigma_{tr,g} / \\sigma_{t,g}` (one - ratio per energy group, ordered from fast to thermal to match the - :class:`openmc.XSdata` group indexing). For each listed domain the - transport-corrected total cross section is computed as + The per-group transport correction ratios + :math:`r_g = \\sigma_{tr,g} / \\sigma_{t,g}` used when the library is + written to an :class:`openmc.MGXSLibrary`. This is a nested dictionary + keyed first by domain type (e.g., ``'material'``) and then by domain + ID, whose values are iterables of per-group ratios (one ratio per + energy group, ordered from fast to thermal to match the + :class:`openmc.XSdata` group indexing). When a ``'transport'`` or + ``'nu-transport'`` MGXS type has been tallied, + :meth:`load_from_statepoint` fills this in automatically from the + tallied data; entries supplied by the user beforehand are preserved and + never overwritten. Explicit assignment may also be used to provide or + edit the ratios directly. For each listed domain the transport- + corrected total cross section is computed as :math:`\\sigma_{tr,g} = r_g \\sigma_{t,g}` and the same correction :math:`\\Delta_g = (1 - r_g)\\sigma_{t,g}` is subtracted from the in-group (diagonal) :math:`P_0` element of the scattering matrix, - preserving the absorption balance. This is an alternative to the - tally-based ``'P0'`` correction and therefore requires - :attr:`correction` to be ``None``, :attr:`scatter_format` to be - ``'legendre'``, and a ``'total'`` MGXS type to be present. Domains - without an entry are left uncorrected. Defaults to ``None``. + preserving the absorption balance. The ratios are only applied as an + alternative to the tally-based ``'P0'`` correction, i.e. when + :attr:`correction` is ``None``, :attr:`scatter_format` is + ``'legendre'``, and a ``'total'`` MGXS type is present; otherwise they + are recorded but not applied. Domains without an entry are left + uncorrected. Defaults to ``None``. .. versionadded:: 0.16.1 scatter_format : {'legendre', 'histogram'} @@ -715,6 +721,11 @@ def load_from_statepoint(self, statepoint): mgxs.load_from_statepoint(statepoint) mgxs.sparse = self.sparse + # Record the transport correction ratios implied by any tallied + # transport-corrected total cross section (see the + # transport_correction_ratios attribute) + self._store_computed_transport_correction_ratios() + def get_mgxs(self, domain, mgxs_type): """Return the MGXS object for some domain and reaction rate type. @@ -1028,6 +1039,109 @@ def load_from_file(filename='mgxs', directory='mgxs'): with open(full_filename, 'rb') as f: return pickle.load(f) + def _uncorrected_total_xs(self, transport_mgxs, nuclides): + """Recover the plain total cross section from a transport MGXS. + + A :class:`openmc.mgxs.TransportXS` tallies both the total reaction rate + and the flux, so the un-corrected total cross section can be recovered + as ``tallies['total'] / tallies['flux (tracklength)']``. The transport + MGXS's own :meth:`~openmc.mgxs.MGXS.get_xs` machinery is reused (by + temporarily overriding its cached cross-section tally) so that the + energy-group ordering and nuclide handling match the transport- + corrected cross section exactly. + + Parameters + ---------- + transport_mgxs : openmc.mgxs.TransportXS + The transport MGXS whose plain total cross section is wanted + nuclides : str + The ``nuclides`` argument to forward to + :meth:`~openmc.mgxs.MGXS.get_xs` (e.g., ``'total'`` or ``'sum'``) + + Returns + ------- + numpy.ndarray + The plain total cross section (one value per energy group) + + """ + + saved_xs_tally = transport_mgxs._xs_tally + saved_rxn_rate_tally = transport_mgxs._rxn_rate_tally + try: + transport_mgxs._xs_tally = ( + transport_mgxs.tallies['total'] / + transport_mgxs.tallies['flux (tracklength)']) + transport_mgxs._compute_xs() + return transport_mgxs.get_xs(nuclides=nuclides, xs_type='macro') + finally: + transport_mgxs._xs_tally = saved_xs_tally + transport_mgxs._rxn_rate_tally = saved_rxn_rate_tally + + def _store_computed_transport_correction_ratios(self): + """Populate transport_correction_ratios from tallied transport data. + + When a transport-type MGXS (``'transport'`` or ``'nu-transport'``) has + been tallied, the per-group transport correction ratio + :math:`r_g = \\sigma_{tr,g} / \\sigma_{t,g}` is computed for each domain + and stored in :attr:`transport_correction_ratios`. The transport- + corrected total :math:`\\sigma_{tr}` is the transport MGXS itself, while + the plain total :math:`\\sigma_t` is recovered from the same MGXS's flux + and total tallies (see :meth:`_uncorrected_total_xs`). + + Ratios that the user already supplied for a domain are never + overwritten. Only isotropic, single-subdomain data (one ratio per + energy group) is handled; angle-dependent or multi-subdomain domains + are skipped. + + """ + + # One ratio per group is only well defined for isotropic data + if self.num_polar > 1 or self.num_azimuthal > 1: + return + + # Use the transport-corrected total that matches the scattering + # multiplicity treatment when both flavors are available + if 'nu-transport' in self.mgxs_types: + transport_type = 'nu-transport' + elif 'transport' in self.mgxs_types: + transport_type = 'transport' + else: + return + + nuclides = 'sum' if self.by_nuclide else 'total' + + # Start from any existing (e.g., user-provided) ratios so they are + # preserved, and never overwrite an entry the user already set + ratios = copy.deepcopy(self._transport_correction_ratios) or {} + domain_ratios = ratios.get(self.domain_type, {}) + + for domain in self.domains: + if domain.id in domain_ratios: + continue + + transport_mgxs = self.get_mgxs(domain, transport_type) + sigma_tr = np.asarray( + transport_mgxs.get_xs(nuclides=nuclides, xs_type='macro'), + dtype=float) + sigma_t = np.asarray( + self._uncorrected_total_xs(transport_mgxs, nuclides), + dtype=float) + + # Only isotropic, single-subdomain data yields one ratio per group + if sigma_tr.shape != (self.num_groups,) or \ + sigma_t.shape != (self.num_groups,): + continue + + # Default to unity (no correction) where the total vanishes + ratio = np.ones(self.num_groups) + nonzero = sigma_t > 0.0 + ratio[nonzero] = sigma_tr[nonzero] / sigma_t[nonzero] + domain_ratios[domain.id] = ratio + + if domain_ratios: + ratios[self.domain_type] = domain_ratios + self._transport_correction_ratios = ratios + def _get_transport_correction_ratios(self, domain): """Return the per-group transport correction ratios for a domain. @@ -1054,11 +1168,11 @@ def _get_transport_correction_ratios(self, domain): return domain_ratios.get(domain.id) def _apply_transport_correction_ratios(self, xsdata, domain, temperature): - """Apply user-supplied transport correction ratios to an XSdata object. + """Apply stored transport correction ratios to an XSdata object. For each energy group ``g``, the transport correction :math:`\\Delta_g = (1 - r_g)\\sigma_{t,g}` is computed from the plain - total cross section :math:`\\sigma_{t,g}` and the user-supplied ratio + total cross section :math:`\\sigma_{t,g}` and the stored ratio :math:`r_g`. The total cross section is replaced by the transport-corrected value :math:`r_g \\sigma_{t,g}`, and the same :math:`\\Delta_g` is subtracted from the in-group (diagonal) @@ -1066,6 +1180,12 @@ def _apply_transport_correction_ratios(self, xsdata, domain, temperature): correction from both quantities leaves the absorption balance unchanged. + The ratios are only applied when :attr:`correction` is ``None``, so + that the total and scattering matrix in ``xsdata`` are the plain, + un-corrected quantities. When :attr:`correction` is ``'P0'`` the + tally-based correction is already reflected in those quantities, so the + stored ratios are left as a record and this method makes no change. + Parameters ---------- xsdata : openmc.XSdata @@ -1083,12 +1203,11 @@ def _apply_transport_correction_ratios(self, xsdata, domain, temperature): return # The ratio-based correction takes the place of the tally-based 'P0' - # correction, so the total and scattering matrix must be un-corrected. + # correction, so it is only applied when the total and scattering + # matrix are the plain (un-corrected) quantities. When a 'P0' + # correction was tallied, the stored ratios merely record it. if self.correction is not None: - raise ValueError( - 'The "correction" parameter must be None when ' - 'transport_correction_ratios are provided, otherwise the ' - 'transport correction would be applied twice.') + return if self.scatter_format != 'legendre': raise ValueError( @@ -1813,12 +1932,11 @@ def check_library_for_openmc_mgxs(self): error_flag = True warn('An "absorption" MGXS type is required but not provided.') - # Validate user-supplied transport correction ratios - if self._transport_correction_ratios: - if self.correction is not None: - error_flag = True - warn('The "correction" parameter must be None when ' - 'transport_correction_ratios are provided.') + # Validate the transport correction ratios that will actually be + # applied. They are only applied when correction is None (otherwise the + # tally-based 'P0' correction is used and the stored ratios simply + # record it), so only validate them in that case. + if self._transport_correction_ratios and self.correction is None: if self.scatter_format != 'legendre': error_flag = True warn('transport_correction_ratios require a "legendre" ' @@ -1826,7 +1944,7 @@ def check_library_for_openmc_mgxs(self): if 'total' not in self.mgxs_types: error_flag = True warn('A "total" MGXS type is required when ' - 'transport_correction_ratios are provided.') + 'transport_correction_ratios are applied.') domain_ratios = \ self._transport_correction_ratios.get(self.domain_type) diff --git a/tests/unit_tests/test_mgxs_transport_correction.py b/tests/unit_tests/test_mgxs_transport_correction.py index 9704ec01872..d4859dd7a00 100644 --- a/tests/unit_tests/test_mgxs_transport_correction.py +++ b/tests/unit_tests/test_mgxs_transport_correction.py @@ -1,4 +1,9 @@ -"""Tests for user-supplied transport correction ratios in openmc.mgxs.Library.""" +"""Tests for transport correction ratios in openmc.mgxs.Library. + +These cover both user-supplied ratios and the ratios that +:meth:`openmc.mgxs.Library.load_from_statepoint` computes automatically from +tallied transport data. +""" import numpy as np import pytest @@ -155,25 +160,32 @@ def test_apply_noop_without_entry(library, simple_geometry): np.testing.assert_allclose(xsdata._total[0], sigma_t) -def test_apply_requires_correction_none(library, simple_geometry): +def test_apply_skips_when_p0_correction(library, simple_geometry): _, mat = simple_geometry library.correction = 'P0' library.transport_correction_ratios = {'material': {1: [0.9, 0.8]}} sigma_t = np.array([2.0, 3.0]) scatter = np.array([[[0.5], [0.3]], [[0.1], [1.2]]]) + scatter_orig = scatter.copy() xsdata = _make_xsdata(library.energy_groups, sigma_t, scatter) - with pytest.raises(ValueError, match='correction'): - library._apply_transport_correction_ratios(xsdata, mat, 294.0) + library._apply_transport_correction_ratios(xsdata, mat, 294.0) + + # With correction='P0' the tally-based correction is already applied, so + # the stored ratios are recorded but leave the dataset unchanged. + np.testing.assert_allclose(xsdata._total[0], sigma_t) + np.testing.assert_allclose(xsdata._scatter_matrix[0][:, :, 0], + scatter_orig[:, :, 0]) -def test_check_library_rejects_p0(library): - library.mgxs_types = ['total', 'absorption', 'nu-scatter matrix', +def test_check_library_allows_ratios_with_p0(library): + # With correction='P0', stored ratios are recorded but not applied, so a + # valid P0 configuration must still pass validation without error. + library.mgxs_types = ['transport', 'absorption', 'nu-scatter matrix', 'scatter matrix'] library.correction = 'P0' library.transport_correction_ratios = {'material': {1: [0.9, 0.8]}} - with pytest.raises(ValueError, match='Invalid MGXS configuration'): - library.check_library_for_openmc_mgxs() + library.check_library_for_openmc_mgxs() def test_check_library_warns_missing_domain_type(library): @@ -183,3 +195,34 @@ def test_check_library_warns_missing_domain_type(library): library.transport_correction_ratios = {'cell': {1: [0.9, 0.8]}} with pytest.warns(UserWarning, match='do not contain any entries'): library.check_library_for_openmc_mgxs() + + +def test_store_ratios_noop_without_transport(library): + # Without a 'transport' or 'nu-transport' MGXS type there is nothing to + # compute, so no ratios are stored (keeping the standard total-based + # workflow unaffected). + library.mgxs_types = ['total', 'absorption', 'scatter matrix'] + library._store_computed_transport_correction_ratios() + assert library.transport_correction_ratios is None + + +def test_store_ratios_skips_angle(library): + # Angle-dependent data does not yield a single ratio per group, so the + # helper returns before touching any tallies. + library.mgxs_types = ['transport', 'absorption', 'scatter matrix'] + library.num_polar = 2 + library.num_azimuthal = 2 + library._store_computed_transport_correction_ratios() + assert library.transport_correction_ratios is None + + +def test_store_ratios_preserves_user_entries(library): + # A user-provided entry must never be overwritten by the automatic + # computation, even when a transport MGXS type is present. + library.mgxs_types = ['transport', 'absorption', 'scatter matrix'] + library.num_polar = 2 # skip the tally-based computation for this test + library.num_azimuthal = 2 + library.transport_correction_ratios = {'material': {1: [0.9, 0.8]}} + library._store_computed_transport_correction_ratios() + np.testing.assert_allclose( + library.transport_correction_ratios['material'][1], [0.9, 0.8]) From 59a85fbe214c5d6ac362d690d749cde20e9baa58 Mon Sep 17 00:00:00 2001 From: Guillaume Giudicelli Date: Sun, 30 Aug 2026 20:46:11 -0600 Subject: [PATCH 3/3] Apply stored transport correction ratios during data extraction get_xsdata now derives the transport correction from the ratios stored in transport_correction_ratios rather than relying on the tally-based total and P0 scattering-matrix correction directly. For each domain the exported total is set to r_g * sigma_t,g and the in-group P0 scattering diagonal is shifted by the same amount, so the absorption balance is preserved and the stored ratios are the single source of truth. Editing an entry after load_from_statepoint therefore changes the correction that is applied. The plain sigma_t needed for the correction is taken from the 'total' MGXS when correction is None, or recovered from the 'transport'/'nu-transport' MGXS (via _uncorrected_total_xs) when correction is 'P0'; the latter recovery is only performed for isotropic data. When the stored ratios equal the tally-derived sigma_tr / sigma_t (as they are after automatic population), re-applying them reproduces the tally-based data unchanged. Validation in check_library_for_openmc_mgxs now covers both correction modes. Co-Authored-By: Claude --- openmc/mgxs/library.py | 189 +++++++++++------- .../test_mgxs_transport_correction.py | 75 ++++++- 2 files changed, 181 insertions(+), 83 deletions(-) diff --git a/openmc/mgxs/library.py b/openmc/mgxs/library.py index e69b8f67d3b..7883132149f 100644 --- a/openmc/mgxs/library.py +++ b/openmc/mgxs/library.py @@ -57,7 +57,7 @@ class Library: Apply the P0 correction to scattering matrices if set to 'P0' transport_correction_ratios : dict or None The per-group transport correction ratios - :math:`r_g = \\sigma_{tr,g} / \\sigma_{t,g}` used when the library is + :math:`r_g = \\sigma_{tr,g} / \\sigma_{t,g}` applied when the library is written to an :class:`openmc.MGXSLibrary`. This is a nested dictionary keyed first by domain type (e.g., ``'material'``) and then by domain ID, whose values are iterables of per-group ratios (one ratio per @@ -67,16 +67,19 @@ class Library: :meth:`load_from_statepoint` fills this in automatically from the tallied data; entries supplied by the user beforehand are preserved and never overwritten. Explicit assignment may also be used to provide or - edit the ratios directly. For each listed domain the transport- - corrected total cross section is computed as - :math:`\\sigma_{tr,g} = r_g \\sigma_{t,g}` and the same correction - :math:`\\Delta_g = (1 - r_g)\\sigma_{t,g}` is subtracted from the - in-group (diagonal) :math:`P_0` element of the scattering matrix, - preserving the absorption balance. The ratios are only applied as an - alternative to the tally-based ``'P0'`` correction, i.e. when - :attr:`correction` is ``None``, :attr:`scatter_format` is - ``'legendre'``, and a ``'total'`` MGXS type is present; otherwise they - are recorded but not applied. Domains without an entry are left + edit the ratios directly. For each listed domain the transport-corrected + total cross section is set to :math:`\\sigma_{tr,g} = r_g \\sigma_{t,g}` + and the same change is applied to the in-group (diagonal) :math:`P_0` + element of the scattering matrix, preserving the absorption balance. + These stored ratios, rather than a correction derived directly from the + tallies during data extraction, are what determine the transport + correction, so editing an entry after :meth:`load_from_statepoint` + changes the correction that is applied. When :attr:`correction` is + ``None`` the plain :math:`\\sigma_{t,g}` is taken from the ``'total'`` + MGXS already written to the dataset; when :attr:`correction` is + ``'P0'`` it is recovered from the transport MGXS. A ``'legendre'`` + :attr:`scatter_format` is required, and the ``'P0'`` recovery is only + performed for isotropic data. Domains without an entry are left uncorrected. Defaults to ``None``. .. versionadded:: 0.16.1 @@ -1039,7 +1042,8 @@ def load_from_file(filename='mgxs', directory='mgxs'): with open(full_filename, 'rb') as f: return pickle.load(f) - def _uncorrected_total_xs(self, transport_mgxs, nuclides): + def _uncorrected_total_xs(self, transport_mgxs, nuclides, xs_type='macro', + subdomains='all'): """Recover the plain total cross section from a transport MGXS. A :class:`openmc.mgxs.TransportXS` tallies both the total reaction rate @@ -1054,9 +1058,16 @@ def _uncorrected_total_xs(self, transport_mgxs, nuclides): ---------- transport_mgxs : openmc.mgxs.TransportXS The transport MGXS whose plain total cross section is wanted - nuclides : str + nuclides : str or Iterable of str The ``nuclides`` argument to forward to - :meth:`~openmc.mgxs.MGXS.get_xs` (e.g., ``'total'`` or ``'sum'``) + :meth:`~openmc.mgxs.MGXS.get_xs` (e.g., ``'total'``, ``'sum'`` or a + list with a single nuclide name) + xs_type : {'macro', 'micro'} + The ``xs_type`` argument to forward to + :meth:`~openmc.mgxs.MGXS.get_xs` + subdomains : Iterable of int or 'all' + The ``subdomains`` argument to forward to + :meth:`~openmc.mgxs.MGXS.get_xs` Returns ------- @@ -1072,7 +1083,8 @@ def _uncorrected_total_xs(self, transport_mgxs, nuclides): transport_mgxs.tallies['total'] / transport_mgxs.tallies['flux (tracklength)']) transport_mgxs._compute_xs() - return transport_mgxs.get_xs(nuclides=nuclides, xs_type='macro') + return transport_mgxs.get_xs(nuclides=nuclides, xs_type=xs_type, + subdomains=subdomains) finally: transport_mgxs._xs_tally = saved_xs_tally transport_mgxs._rxn_rate_tally = saved_rxn_rate_tally @@ -1167,34 +1179,49 @@ def _get_transport_correction_ratios(self, domain): return domain_ratios.get(domain.id) - def _apply_transport_correction_ratios(self, xsdata, domain, temperature): - """Apply stored transport correction ratios to an XSdata object. - - For each energy group ``g``, the transport correction - :math:`\\Delta_g = (1 - r_g)\\sigma_{t,g}` is computed from the plain - total cross section :math:`\\sigma_{t,g}` and the stored ratio - :math:`r_g`. The total cross section is replaced by the - transport-corrected value :math:`r_g \\sigma_{t,g}`, and the same - :math:`\\Delta_g` is subtracted from the in-group (diagonal) - :math:`P_0` element of the scattering matrix. Subtracting the same - correction from both quantities leaves the absorption balance - unchanged. - - The ratios are only applied when :attr:`correction` is ``None``, so - that the total and scattering matrix in ``xsdata`` are the plain, - un-corrected quantities. When :attr:`correction` is ``'P0'`` the - tally-based correction is already reflected in those quantities, so the - stored ratios are left as a record and this method makes no change. + def _apply_transport_correction_ratios(self, xsdata, domain, temperature, + nuclide='total', xs_type='macro', + subdomain='all'): + """Apply the stored transport correction ratios to an XSdata object. + + For each energy group ``g`` the transport-corrected total cross section + is set to :math:`r_g \\sigma_{t,g}`, where :math:`\\sigma_{t,g}` is the + plain total cross section and :math:`r_g` is the stored ratio. The same + change is applied to the in-group (diagonal) :math:`P_0` element of the + scattering matrix, so the absorption balance (total minus out-scatter) + is preserved. + + Because the total and the scattering diagonal are shifted by the same + amount, the stored ratios rather than a correction derived directly + from the tallies during data extraction are what determine the + transport correction. Editing an entry in + :attr:`transport_correction_ratios` after + :meth:`load_from_statepoint` therefore changes the correction applied + here. + + When :attr:`correction` is ``None`` the total already in ``xsdata`` is + the plain :math:`\\sigma_t`. When :attr:`correction` is ``'P0'`` the + total is the tally-based transport-corrected :math:`\\sigma_{tr}`, so + the plain :math:`\\sigma_t` is recovered from the transport MGXS (see + :meth:`_uncorrected_total_xs`); this recovery is only performed for + isotropic data, matching the ratios computed by + :meth:`_store_computed_transport_correction_ratios`. Parameters ---------- xsdata : openmc.XSdata - The dataset to correct in place. Its total cross section is - expected to be the plain (un-corrected) total. + The dataset to correct in place domain : openmc.Material or openmc.Cell or openmc.Universe or openmc.RegularMesh The domain the dataset describes temperature : float Temperature (in Kelvin) of the data to correct + nuclide : str + The nuclide the dataset describes (or ``'total'`` for material-wise + data); used to recover the plain total cross section + xs_type : {'macro', 'micro'} + Whether the dataset holds macroscopic or microscopic cross sections + subdomain : Iterable of int or 'all' + The subdomain the dataset describes (for mesh domains) """ @@ -1202,46 +1229,57 @@ def _apply_transport_correction_ratios(self, xsdata, domain, temperature): if ratios is None: return - # The ratio-based correction takes the place of the tally-based 'P0' - # correction, so it is only applied when the total and scattering - # matrix are the plain (un-corrected) quantities. When a 'P0' - # correction was tallied, the stored ratios merely record it. - if self.correction is not None: - return - + # A legendre scattering matrix is needed to correct the P0 diagonal, + # and exactly one ratio per energy group is required. if self.scatter_format != 'legendre': - raise ValueError( - 'transport_correction_ratios require a "legendre" ' - 'scatter_format.') - + return if len(ratios) != self.num_groups: - raise ValueError( - f'Expected {self.num_groups} transport correction ratios for ' - f'domain {domain.id} but got {len(ratios)}.') + return i = xsdata._temperature_index(temperature) + current_total = xsdata._total[i] + if current_total is None: + return - total = xsdata._total[i] - if total is None: - raise ValueError( - 'A "total" MGXS type is required to apply ' - 'transport_correction_ratios.') - - # Compute the correction from the plain total cross section. ratios has - # shape (G,) and broadcasts against the trailing group axis of the - # total cross section for both isotropic and angle representations. - delta = (1.0 - ratios) * total - xsdata._total[i] = total - delta + # Determine the plain total cross section sigma_t. + if self.correction is None: + # The total already in xsdata is the plain total. + sigma_t = current_total + else: + # The total in xsdata is the tally-corrected transport total, so + # recover the plain total from the transport MGXS. The stored + # ratios are only computed for isotropic data. + if xsdata.representation == 'angle': + return + if 'nu-transport' in self.mgxs_types: + transport_type = 'nu-transport' + elif 'transport' in self.mgxs_types: + transport_type = 'transport' + else: + return + transport_mgxs = self.get_mgxs(domain, transport_type) + sigma_t = np.asarray( + self._uncorrected_total_xs(transport_mgxs, [nuclide], xs_type, + subdomain), dtype=float) + if sigma_t.size != np.size(current_total): + return + sigma_t = sigma_t.reshape(np.shape(current_total)) + + # Set the transport-corrected total to r * sigma_t and shift the P0 + # scattering diagonal by the same amount. ratios has shape (G,) and + # broadcasts against the trailing group axis for both the isotropic and + # angle representations. + target_total = ratios * sigma_t + delta = target_total - current_total + xsdata._total[i] = target_total - # Subtract the same correction from the in-group (diagonal) P0 element - # of the scattering matrix to preserve the absorption balance. scatter_matrix = xsdata._scatter_matrix[i] if scatter_matrix is not None: groups = np.arange(self.num_groups) if xsdata.representation == 'angle': - scatter_matrix[:, :, groups, groups, 0] -= delta + scatter_matrix[:, :, groups, groups, 0] += delta else: - scatter_matrix[groups, groups, 0] -= delta + scatter_matrix[groups, groups, 0] += delta def get_xsdata(self, domain, xsdata_name, nuclide='total', xs_type='macro', subdomain=None, apply_domain_chi=False, temperature=ROOM_TEMPERATURE_KELVIN): @@ -1592,9 +1630,11 @@ def get_xsdata(self, domain, xsdata_name, nuclide='total', xs_type='macro', nuclide=[nuclide], subdomain=subdomain) - # Apply any user-supplied transport correction ratios to the total - # cross section and the scattering matrix diagonal - self._apply_transport_correction_ratios(xsdata, domain, temperature) + # Apply the stored transport correction ratios to the total cross + # section and the scattering matrix diagonal + self._apply_transport_correction_ratios( + xsdata, domain, temperature, nuclide=nuclide, xs_type=xs_type, + subdomain=subdomain) return xsdata @@ -1932,19 +1972,16 @@ def check_library_for_openmc_mgxs(self): error_flag = True warn('An "absorption" MGXS type is required but not provided.') - # Validate the transport correction ratios that will actually be - # applied. They are only applied when correction is None (otherwise the - # tally-based 'P0' correction is used and the stored ratios simply - # record it), so only validate them in that case. - if self._transport_correction_ratios and self.correction is None: + # Validate the transport correction ratios that will be applied to the + # exported data. The plain total cross section they require is already + # guaranteed by the total/transport checks above (a 'total' MGXS when + # correction is None, a 'transport'/'nu-transport' MGXS otherwise), so + # only the ratio-specific requirements are checked here. + if self._transport_correction_ratios: if self.scatter_format != 'legendre': error_flag = True warn('transport_correction_ratios require a "legendre" ' 'scatter_format.') - if 'total' not in self.mgxs_types: - error_flag = True - warn('A "total" MGXS type is required when ' - 'transport_correction_ratios are applied.') domain_ratios = \ self._transport_correction_ratios.get(self.domain_type) diff --git a/tests/unit_tests/test_mgxs_transport_correction.py b/tests/unit_tests/test_mgxs_transport_correction.py index d4859dd7a00..025ce5548cc 100644 --- a/tests/unit_tests/test_mgxs_transport_correction.py +++ b/tests/unit_tests/test_mgxs_transport_correction.py @@ -160,27 +160,78 @@ def test_apply_noop_without_entry(library, simple_geometry): np.testing.assert_allclose(xsdata._total[0], sigma_t) -def test_apply_skips_when_p0_correction(library, simple_geometry): +def test_apply_uses_stored_ratio_with_p0(library, simple_geometry, monkeypatch): + # With correction='P0' the dataset arrives with the tally-based transport + # total, and the plain total is recovered from the transport MGXS. The + # stored ratio (as if edited by the user) then determines the correction. _, mat = simple_geometry library.correction = 'P0' + library.mgxs_types = ['transport', 'nu-scatter matrix', 'scatter matrix'] library.transport_correction_ratios = {'material': {1: [0.9, 0.8]}} + # sigma_t recovered from the transport MGXS sigma_t = np.array([2.0, 3.0]) + monkeypatch.setattr(library, 'get_mgxs', + lambda domain, mgxs_type: object()) + monkeypatch.setattr( + library, '_uncorrected_total_xs', + lambda tm, nuclides, xs_type, subdomains: sigma_t) + + # The tally-based transport total in the dataset differs from + # ratios * sigma_t, so the stored ratio must visibly take effect. + corrected_total = np.array([1.5, 2.0]) scatter = np.array([[[0.5], [0.3]], [[0.1], [1.2]]]) scatter_orig = scatter.copy() - xsdata = _make_xsdata(library.energy_groups, sigma_t, scatter) + xsdata = _make_xsdata(library.energy_groups, corrected_total, scatter) + library._apply_transport_correction_ratios(xsdata, mat, 294.0) - # With correction='P0' the tally-based correction is already applied, so - # the stored ratios are recorded but leave the dataset unchanged. - np.testing.assert_allclose(xsdata._total[0], sigma_t) + ratios = np.array([0.9, 0.8]) + target_total = ratios * sigma_t + delta = target_total - corrected_total + + np.testing.assert_allclose(xsdata._total[0], target_total) + sm = xsdata._scatter_matrix[0] + np.testing.assert_allclose(sm[0, 0, 0], scatter_orig[0, 0, 0] + delta[0]) + np.testing.assert_allclose(sm[1, 1, 0], scatter_orig[1, 1, 0] + delta[1]) + np.testing.assert_allclose(sm[0, 1, 0], scatter_orig[0, 1, 0]) + np.testing.assert_allclose(sm[1, 0, 0], scatter_orig[1, 0, 0]) + + +def test_apply_p0_ratio_matches_tally_is_noop(library, simple_geometry, + monkeypatch): + # When the stored ratio equals the tally-derived sigma_tr / sigma_t (as it + # is after automatic population), re-applying it reproduces the tally-based + # total and leaves the scattering matrix unchanged (no double correction). + _, mat = simple_geometry + library.correction = 'P0' + library.mgxs_types = ['transport', 'nu-scatter matrix', 'scatter matrix'] + + sigma_t = np.array([2.0, 3.0]) + corrected_total = np.array([1.8, 2.4]) # sigma_tr from the tallies + ratios = corrected_total / sigma_t + library.transport_correction_ratios = {'material': {1: list(ratios)}} + + monkeypatch.setattr(library, 'get_mgxs', + lambda domain, mgxs_type: object()) + monkeypatch.setattr( + library, '_uncorrected_total_xs', + lambda tm, nuclides, xs_type, subdomains: sigma_t) + + scatter = np.array([[[0.5], [0.3]], [[0.1], [1.2]]]) + scatter_orig = scatter.copy() + xsdata = _make_xsdata(library.energy_groups, corrected_total, scatter) + library._apply_transport_correction_ratios(xsdata, mat, 294.0) + + np.testing.assert_allclose(xsdata._total[0], corrected_total) np.testing.assert_allclose(xsdata._scatter_matrix[0][:, :, 0], scatter_orig[:, :, 0]) def test_check_library_allows_ratios_with_p0(library): - # With correction='P0', stored ratios are recorded but not applied, so a - # valid P0 configuration must still pass validation without error. + # With correction='P0' the stored ratios are applied by recovering the + # plain total from the transport MGXS, so a valid P0 configuration with a + # 'transport' type must pass validation without error. library.mgxs_types = ['transport', 'absorption', 'nu-scatter matrix', 'scatter matrix'] library.correction = 'P0' @@ -188,6 +239,16 @@ def test_check_library_allows_ratios_with_p0(library): library.check_library_for_openmc_mgxs() +def test_check_library_p0_requires_transport(library): + # With correction='P0' and no 'transport'/'nu-transport' type there is no + # way to recover the plain total, so validation must fail. + library.mgxs_types = ['absorption', 'nu-scatter matrix', 'scatter matrix'] + library.correction = 'P0' + library.transport_correction_ratios = {'material': {1: [0.9, 0.8]}} + with pytest.raises(ValueError, match='Invalid MGXS configuration'): + library.check_library_for_openmc_mgxs() + + def test_check_library_warns_missing_domain_type(library): library.mgxs_types = ['total', 'absorption', 'nu-scatter matrix', 'scatter matrix']