Generalize get decay particle energy

docs/source/io_formats/depletion_chain.rst

@@ -60,16 +60,17 @@ attributes:
 ``<source>`` Element
 --------------------
 
-The ``<source>`` element represents photon and electron sources associated with
-the decay of a nuclide and contains information to construct an
+The ``<source>`` element represents decay particle sources associated with the
+decay of a nuclide and contains information to construct an
 :class:`openmc.stats.Univariate` object that represents this emission as an
 energy distribution. This element has the following attributes:
 
   :type:
     The type of :class:`openmc.stats.Univariate` source term.
 
   :particle:
-    The type of particle emitted, e.g., 'photon' or 'electron'
+    The type of particle emitted, e.g., ``photon``, ``electron``,
+    ``positron``, ``alpha``, ``neutron``, ``proton``, or ``fragment``
 
   :parameters:
     The parameters of the source term, e.g., for a

docs/source/pythonapi/data.rst

@@ -64,6 +64,7 @@ Core Functions
     combine_distributions
     decay_constant
     decay_energy
+    decay_particle_energy
     decay_photon_energy
     dose_coefficients
     gnds_name

docs/source/usersguide/decay_sources.rst

@@ -39,9 +39,13 @@ Generally, it involves the following steps:
    with the timestep index returns a :class:`~openmc.deplete.StepResult` object,
    which itself has a :meth:`~openmc.deplete.StepResult.get_material` method
    returning an activated material.
-4. The :meth:`openmc.Material.get_decay_photon_energy` method is used to obtain
-   the energy spectrum of the decay photon source. The integral of the spectrum
-   also indicates the intensity of the source in units of [Bq].
+4. The :meth:`openmc.Material.get_decay_particle_energy` method can be used to
+   obtain the energy spectrum of decay particles from the activated material.
+   For R2S, the photon spectrum is needed and can be obtained either with
+   :meth:`openmc.Material.get_decay_particle_energy` using
+   ``particle='photon'`` or with the wrapper
+   :meth:`openmc.Material.get_decay_photon_energy`. The integral of the
+   spectrum indicates the intensity of the source in units of [Bq].
 5. A new photon source is defined using one of OpenMC's source classes with the
    energy distribution set equal to the object returned by the
    :meth:`openmc.Material.get_decay_photon_energy` method. The source is then
@@ -76,9 +80,10 @@ Altogether, the workflow looks as follows::
     model.settings.source = photon_source
     model.run()
 
-Note that by default, the :meth:`~openmc.Material.get_decay_photon_energy`
-method will eliminate spectral lines with very low intensity, but this behavior
-can be configured with the ``clip_tolerance`` argument.
+Note that :meth:`~openmc.Material.get_decay_particle_energy` and
+:meth:`~openmc.Material.get_decay_photon_energy` will, by default, eliminate
+spectral lines with very low intensity, but this behavior can be configured
+with the ``clip_tolerance`` argument.
 
 Cell-based R2S
 --------------
@@ -236,4 +241,3 @@ relevant tallies. This can be done with the aid of the
 
     # Apply time correction factors
     tally = d1s.apply_time_correction(dose_tally, factors, time_index)
-

openmc/config.py

@@ -20,7 +20,7 @@
 from typing import Any, Dict, Iterator
 
 from openmc.data import DataLibrary
-from openmc.data.decay import _DECAY_ENERGY, _DECAY_PHOTON_ENERGY
+from openmc.data.decay import _DECAY_ENERGY, _DECAY_PARTICLE_ENERGY
 
 __all__ = ["config"]
 
@@ -77,7 +77,7 @@ def __delitem__(self, key: str):
             if env_var in os.environ:
                 del os.environ[env_var]
         if key == 'chain_file':
-            _DECAY_PHOTON_ENERGY.clear()
+            _DECAY_PARTICLE_ENERGY.clear()
             _DECAY_ENERGY.clear()
 
     def __setitem__(self, key: str, value: Any):
@@ -103,7 +103,7 @@ def __setitem__(self, key: str, value: Any):
             os.environ[self._PATH_KEYS[key]] = str(stored_path)
 
             if key == 'chain_file':
-                _DECAY_PHOTON_ENERGY.clear()
+                _DECAY_PARTICLE_ENERGY.clear()
                 _DECAY_ENERGY.clear()
 
             if not stored_path.exists():

openmc/data/decay.py

@@ -561,7 +561,7 @@ def sources(self):
         return merged_sources
 
 
-_DECAY_PHOTON_ENERGY = {}
+_DECAY_PARTICLE_ENERGY = {}
 
 
 def decay_photon_energy(nuclide: str) -> Univariate | None:
@@ -571,7 +571,8 @@ def decay_photon_energy(nuclide: str) -> Univariate | None:
     for the first time, you need to ensure that a depletion chain has been
     specified in openmc.config['chain_file'].
 
-    .. versionadded:: 0.13.2
+    This is a backward-compatible convenience wrapper around
+    :func:`decay_particle_energy` with ``particle='photon'``.
 
     Parameters
     ----------
@@ -585,7 +586,34 @@ def decay_photon_energy(nuclide: str) -> Univariate | None:
         if no photon source exists. Note that the probabilities represent
         intensities, given as [Bq/atom] (in other words, decay constants).
     """
-    if not _DECAY_PHOTON_ENERGY:
+    return decay_particle_energy(nuclide, 'photon')
+
+
+def decay_particle_energy(nuclide: str, particle: str) -> Univariate | None:
+    """Get a decay particle energy distribution resulting from the decay of a nuclide
+
+    This function relies on data stored in a depletion chain. Before calling it
+    for the first time, you need to ensure that a depletion chain has been
+    specified in openmc.config['chain_file'].
+
+    Parameters
+    ----------
+    nuclide : str
+        Name of nuclide, e.g., 'Co58'
+    particle : str
+        Name of the decay particle, e.g., ``photon``, ``electron``,
+        ``positron``, ``alpha``, ``neutron``, ``proton``, or ``fragment``.
+
+    Returns
+    -------
+    openmc.stats.Univariate or None
+        Distribution of energies in [eV] of decay particles emitted from decay,
+        or None if no decay source of that particle exists. Note that the
+        probabilities represent intensities, given as [Bq/atom] (in other
+        words, decay constants).
+    """
+
+    if not _DECAY_PARTICLE_ENERGY:
         chain_file = openmc.config.get('chain_file')
         if chain_file is None:
             raise DataError(
@@ -596,15 +624,15 @@ def decay_photon_energy(nuclide: str) -> Univariate | None:
         from openmc.deplete import Chain
         chain = Chain.from_xml(chain_file)
         for nuc in chain.nuclides:
-            if 'photon' in nuc.sources:
-                _DECAY_PHOTON_ENERGY[nuc.name] = nuc.sources['photon']
+            for part, dist in nuc.sources.items():
+                _DECAY_PARTICLE_ENERGY[(nuc.name, part)] = dist
 
         # If the chain file contained no sources at all, warn the user
-        if not _DECAY_PHOTON_ENERGY:
-            warn(f"Chain file '{chain_file}' does not have any decay photon "
-                 "sources listed.")
+        if not _DECAY_PARTICLE_ENERGY:
+            warn(f"Chain file '{chain_file}' does not have any decay particle "
+                 "source listed.")
 
-    return _DECAY_PHOTON_ENERGY.get(nuclide)
+    return _DECAY_PARTICLE_ENERGY.get((nuclide, particle))
 
 
 _DECAY_ENERGY = {}
@@ -649,5 +677,3 @@ def decay_energy(nuclide: str):
             warn(f"Chain file '{chain_file}' does not have any decay energy.")
 
     return _DECAY_ENERGY.get(nuclide, 0.0)
-
-

openmc/material.py

@@ -332,20 +332,23 @@ def decay_photon_energy(self) -> Univariate | None:
             "version.", FutureWarning)
         return self.get_decay_photon_energy(0.0)
 
-    def get_decay_photon_energy(
+    def get_decay_particle_energy(
         self,
+        particle: str = 'photon',
         clip_tolerance: float = 1e-6,
         units: str = 'Bq',
         volume: float | None = None,
         exclude_nuclides: list[str] | None = None,
         include_nuclides: list[str] | None = None
     ) -> Univariate | None:
-        r"""Return energy distribution of decay photons from unstable nuclides.
-
-        .. versionadded:: 0.14.0
+        r"""Return energy distribution of decay particles from unstable nuclides.
 
         Parameters
         ----------
+        particle : string
+            Particle type to return the energy distribution for. Supported
+            values are ``'photon'``, ``'electron'``, ``'positron'``,
+            ``'alpha'``, ``'neutron'``, ``'proton'``, and ``'fragment'``.
         clip_tolerance : float
             Maximum fraction of :math:`\sum_i x_i p_i` for discrete distributions
             that will be discarded.
@@ -355,19 +358,21 @@ def get_decay_photon_energy(
             Volume of the material. If not passed, defaults to using the
             :attr:`Material.volume` attribute.
         exclude_nuclides : list of str, optional
-            Nuclides to exclude from the photon source calculation.
+            Nuclides to exclude from the particle source calculation.
         include_nuclides : list of str, optional
-            Nuclides to include in the photon source calculation. If specified,
+            Nuclides to include in the particle source calculation. If specified,
             only these nuclides are used.
 
         Returns
         -------
         Univariate or None
-            Decay photon energy distribution. The integral of this distribution is
-            the total intensity of the photon source in the requested units.
+            Decay particle energy distribution. The integral of this distribution is
+            the total intensity of the particle source in the requested units.
 
         """
         cv.check_value('units', units, {'Bq', 'Bq/g', 'Bq/kg', 'Bq/cm3', 'Bq/m3'})
+        cv.check_type('particle', particle, str)
+        cv.check_value('particle', particle, {'alpha', 'fragment', 'neutron','photon', 'electron', 'positron','proton'})
 
         if exclude_nuclides is not None and include_nuclides is not None:
             raise ValueError("Cannot specify both exclude_nuclides and include_nuclides")
@@ -393,12 +398,12 @@ def get_decay_photon_energy(
             if include_nuclides is not None and nuc not in include_nuclides:
                 continue
 
-            source_per_atom = openmc.data.decay_photon_energy(nuc)
+            source_per_atom = openmc.data.decay_particle_energy(nuc,particle)
             if source_per_atom is not None and atoms_per_bcm > 0.0:
                 dists.append(source_per_atom)
                 probs.append(1e24 * atoms_per_bcm * multiplier)
 
-        # If no photon sources, exit early
+        # If no particle sources, exit early
         if not dists:
             return None
 
@@ -414,6 +419,53 @@ def get_decay_photon_energy(
 
         return combined
 
+    def get_decay_photon_energy(
+        self,
+        clip_tolerance: float = 1e-6,
+        units: str = 'Bq',
+        volume: float | None = None,
+        exclude_nuclides: list[str] | None = None,
+        include_nuclides: list[str] | None = None
+    ) -> Univariate | None:
+        r"""Return energy distribution of decay photons from unstable nuclides.
+
+        This is a backward-compatible convenience wrapper around
+        :meth:`Material.get_decay_particle_energy` with
+        ``particle='photon'``.
+
+        Parameters
+        ----------
+        clip_tolerance : float
+            Maximum fraction of :math:`\sum_i x_i p_i` for discrete distributions
+            that will be discarded.
+        units : {'Bq', 'Bq/g', 'Bq/kg', 'Bq/cm3', 'Bq/m3'}
+            Specifies the units on the integral of the distribution.
+        volume : float, optional
+            Volume of the material. If not passed, defaults to using the
+            :attr:`Material.volume` attribute.
+        exclude_nuclides : list of str, optional
+            Nuclides to exclude from the photon source calculation.
+        include_nuclides : list of str, optional
+            Nuclides to include in the photon source calculation. If specified,
+            only these nuclides are used.
+
+        Returns
+        -------
+        Univariate or None
+            Decay photon energy distribution. The integral of this distribution is
+            the total intensity of the photon source in the requested units.
+
+        """
+
+        return self.get_decay_particle_energy(
+            particle='photon',
+            clip_tolerance=clip_tolerance,
+            units=units,
+            volume=volume,
+            exclude_nuclides=exclude_nuclides,
+            include_nuclides=include_nuclides,
+        )
+
     def get_photon_contact_dose_rate(
         self,
         dose_quantity: str = "absorbed-air",
@@ -540,7 +592,7 @@ def get_photon_contact_dose_rate(
                           * sv_per_psv * 1e24)
 
         for nuc, nuc_atoms_per_bcm in self.get_nuclide_atom_densities().items():
-            photon_source_per_atom = openmc.data.decay_photon_energy(nuc)
+            photon_source_per_atom = openmc.data.decay_particle_energy(nuc,'photon')
 
             # nuclides with no contribution
             if photon_source_per_atom is None or nuc_atoms_per_bcm <= 0.0:

tests/unit_tests/test_config.py

@@ -99,6 +99,6 @@ def test_config_chain_side_effect(tmp_path):
     """Test that modifying chain_file clears decay data caches."""
     chain_file = tmp_path / "chain.xml"; chain_file.touch()
     decay._DECAY_ENERGY['U235'] = (1.0, 2.0)
-    decay._DECAY_PHOTON_ENERGY['PU239'] = {}
+    decay._DECAY_PARTICLE_ENERGY['PU239','photon'] = {}
     openmc.config['chain_file'] = chain_file
-    assert not decay._DECAY_ENERGY and not decay._DECAY_PHOTON_ENERGY
+    assert not decay._DECAY_ENERGY and not decay._DECAY_PARTICLE_ENERGY

tests/unit_tests/test_data_decay.py

@@ -146,3 +146,23 @@ def test_decay_photon_energy():
     src = openmc.data.decay_photon_energy('Xe135')
     assert isinstance(src, openmc.stats.Tabular)
     assert src.integral() == pytest.approx(2.076506258964966e-05)
+
+
+def test_decay_particle_energy():
+    # If chain file is not set, we should get a data error
+    if 'chain_file' in openmc.config:
+        del openmc.config['chain_file']
+    with pytest.raises(DataError):
+        openmc.data.decay_particle_energy('Fe55', 'electron')
+
+    # Set chain file to Ni chain and check electron/positron sources
+    with openmc.config.patch('chain_file', Path(__file__).parents[1] / "chain_ni.xml"):
+        src = openmc.data.decay_particle_energy('Fe55', 'electron')
+        assert isinstance(src, openmc.stats.Discrete)
+        assert src.integral() == pytest.approx(4.225437849490689e-08)
+        assert 6377.571 in src.x
+
+        src = openmc.data.decay_particle_energy('Fe55', 'positron')
+        assert isinstance(src, openmc.stats.Discrete)
+        assert src.integral() == pytest.approx(8.004558990612365e-09)
+        assert 231210.0 in src.x