crp/OpenMC
1
0
Fork 0
mirror of https://github.com/openmc-dev/openmc.git synced 2026-07-27 21:55:41 -04:00
Code Issues Projects Releases Packages Wiki Activity Actions 5

Add ConstantFissionYieldHelper concrete class

Browse source 
Given a requested energy, will take fission yields from that
energy on all nuclides and hold them as constant throughout
the simulation. If the requested energy is not found for a specific
nuclide, then the closest energy is used. The default is to
take the thermal 0.0253 eV yields.

This is now the helper attached to the Operator. The user can
select what energy of yields they would like to use.
This commit is contained in:
Andrew Johnson 2019-08-13 15:39:13 -05:00
parent 2cfbb21858
commit d8b9f5db54
No known key found for this signature in database
GPG key ID: 253418E91B7F6FEB
5 changed files with 99 additions and 225 deletions
Download patch file Download diff file Expand all files Collapse all files

1
docs/source/pythonapi/deplete.rst
View file

@ -78,6 +78,7 @@ data, such as number densities and reaction rates for each material.
AtomNumber
ChainFissionHelper
ConstantFissionYieldHelper
DirectReactionRateHelper
OperatorResult
ReactionRates

173
openmc/deplete/helpers.py
View file

@ -2,12 +2,19 @@
Class for normalizing fission energy deposition
"""
from itertools import product
from numbers import Real
from numpy import dot, zeros, newaxis, asarray, empty_like, where
from openmc.checkvalue import check_type, check_greater_than
from openmc.capi import (
Tally, MaterialFilter, EnergyFilter)
from .abc import ReactionRateHelper, EnergyHelper
from .abc import (
ReactionRateHelper, EnergyHelper, FissionYieldHelper)
__all__ = (
"DirectReactionRateHelper", "ChainFissionHelper",
"ConstantFissionYieldHelper")
# -------------------------------------
# Helpers for generating reaction rates
@ -150,152 +157,68 @@ class ChainFissionHelper(EnergyHelper):
# ------------------------------------
class FissionYieldHelper(object):
"""Class for using energy-dependent fission yields in depletion chain
Creates a tally across all burnable materials to score the fission
rate in nuclides with yield data. An energy filter is used to
compute this rates in a group structure corresponding to the
fission yield data. This tally data is used to compute the
relative number of fission events in each energy region,
which serve as the weights for each energy-dependent fission
yield distribution.
class ConstantFissionYieldHelper(FissionYieldHelper):
"""Class that uses a single set of fission yields on each isotope
Parameters
----------
chain_nuclides : iterable of openmc.deplete.Nuclide
Nuclides tracked in the depletion chain. Not necessary
that all have yield data.
n_bmats : int
Number of burnable materials tracked in the problem
energy : float, optional
Key in :attr:`openmc.deplete.Nuclide.yield_data` corresponding
to the desired set of fission yield data. Typically one of
``{0.0253, 500000, 14000000}`` corresponding to 0.0253 eV,
500 keV, and 14 MeV yield libraries. If the specific key is not
found, will fall back to closest energy present.
Default: 0.0253 eV for thermal yields
Attributes
----------
energy_bounds : tuple of float
Sorted energy bounds from the tally filter
results : numpy.ndarray
Array of tally results for this process with shape
``(n_local_mat, n_energy, n_nucs)``
constant_yields : dict of str to :class:`openmc.deplete.FissionYield`
Fission yields for all nuclides that only have one set of
fission yield data. Can be accessed as ``{parent: {product: yield}}``
energy : float
Energy of fission yield libraries.
"""
def __init__(self, chain_nuclides, n_bmats):
self._chain_nuclides = {}
self._chain_set = set()
self._tally_map = {}
# TODO Support user-requested minimum energy?
yield_energies = {0.0}
# Get nuclides with fission yield data, names
# and all energy points
# Names are provided from operator tally nuclides
for nuc in chain_nuclides:
if len(nuc.yield_data) == 0:
def __init__(self, chain_nuclides, energy=0.0253):
check_type("energy", energy, Real)
check_greater_than("energy", energy, 0.0, equality=True)
self._energy = energy
super().__init__(chain_nuclides)
# Iterate over all nuclides with > 1 set of yields
for name, nuc in self._chain_nuclides.items():
yield_data = nuc.yield_data.get(energy)
if yield_data is not None:
self._constant_yields[name] = yield_data
continue
self._chain_nuclides[nuc.name] = nuc
self._chain_set.add(nuc.name)
yield_energies.update(nuc.yield_energies)
# Create energy grid
self._energy_bounds = tuple(sorted(yield_energies))
self.n_bmats = n_bmats
self._reaction_tally = None
self.results = None
self.local_indexes = None
# Specific energy not found, use closest energy
distances = [abs(energy - ene) for ene in nuc.yield_energies]
min_index = min(
range(len(nuc.yield_energies)), key=distances.__getitem__)
self._constant_yields[name] = (
nuc.yield_data[nuc.yield_energies[min_index]])
@property
def energy_bounds(self):
return self._energy_bounds
def energy(self):
return self._energy
def generate_tallies(self, materials, mat_indexes):
"""Construct the fission rate tally
def weighted_yields(self, _local_mat_index=None):
"""Return fission yields for all nuclides requested
Parameters
----------
materials : iterable of C-API materials
Materials to be used in :class:`openmc.capi.MaterialFilter`
mat_indexes : iterable of int
Indexes for materials in ``materials`` tracked on this
process
"""
# Tally group-wise fission reaction rates
self._reaction_tally = Tally()
self._reaction_tally.scores = ['fission']
# Tally energy-weighted group-wise fission reaction rate
# Used to evaluated linear interpolation between fission yield points
self._weighted_reaction_tally = Tally()
self._weighted_reaction_tally.scores = ['fission']
filters = [
MaterialFilter(materials), EnergyFilter(self._energy_bounds)]
self._reaction_tally.filters = filters
self.local_indexes = asarray(mat_indexes)
def set_fissionable_nuclides(self, nuclides):
"""List of string of nuclides with data to be tallied
Parameters
----------
nuclides : iterable of str
Nuclides with non-zero densities that are candidates
for the fission tally. Not necessary that all are nuclides
with fission yields, but at least one must be
_local_mat_index : int, optional
Current material index. Not used since all yields are
constant
Returns
-------
nuclides : tuple of str
Nuclides ordered as they appear in the tally and in
the nuclide column of :attr:`results`
Raises
------
ValueError
If no nuclides in ``nuclides`` are tracked on this
object
library : dict
Dictionary of ``{parent: {product: fyield}}``
"""
# Set of all nuclides with positive density
# and fission yield data
nuc_set = self._chain_set & set(nuclides)
if len(nuc_set) == 0:
raise ValueError(
"No overlap between chain nuclides with fission yields and "
"requested tally nuclides")
nuclides = tuple(sorted(nuc_set))
self._tally_index = [self._chain_nuclides[n] for n in nuclides]
self._reaction_tally.nuclides = nuclides
return nuclides
def unpack(self):
"""Unpack fission rate tallies to produce :attr:`results`
"""
# if this process is not responsible for depleting anything
# [more processes than burnable materials]
# don't do anything
if self.local_indexes.size == 0:
return
# get view into tally results
# new shape: [material, energy, parent nuclide]
result_view = self._reaction_tally.results[..., 1].reshape(
self.n_bmats, len(self._energy_bounds) - 1,
len(self._reaction_tally.nuclides))
# Get results specific to this process
fission_rates = result_view[self.local_indexes]
self.results = empty_like(fission_rates)
# scale group fission rates proportional to total fission rate
fission_total = fission_rates.sum(axis=1)
nz_mat, nz_nuc = fission_total.nonzero()
self.results[nz_mat, :, nz_nuc] = (
fission_rates[nz_mat, :, nz_nuc]
/ fission_total[nz_mat, newaxis, nz_nuc])
# directly set values to zero where total fission rate is zero
z_mat, z_nuc = where(fission_total == 0.0)
self.results[z_mat, :, z_nuc] = 0.0
return self.constant_yields
def compute_yields(self, local_mat_index):
"""Compute single fission yields using :attr:`results`

12
openmc/deplete/operator.py
View file

@ -26,7 +26,7 @@ from .atom_number import AtomNumber
from .reaction_rates import ReactionRates
from .results_list import ResultsList
from .helpers import (
DirectReactionRateHelper, ChainFissionHelper, FissionYieldHelper)
DirectReactionRateHelper, ChainFissionHelper, ConstantFissionYieldHelper)
def _distribute(items):
@ -85,6 +85,10 @@ class Operator(TransportOperator):
in initial condition to ensure they exist in the decay chain.
Only done for nuclides with reaction rates.
Defaults to 1.0e3.
fission_yield_energy : float, optional
Energy [eV] to pull fission product yields from. Passed to the
:class:`openmc.deplete.ConstantFissionYieldHelper`. Default:
0.0253 eV
Attributes
----------
@ -123,7 +127,7 @@ class Operator(TransportOperator):
"""
def __init__(self, geometry, settings, chain_file=None, prev_results=None,
diff_burnable_mats=False, fission_q=None,
dilute_initial=1.0e3):
dilute_initial=1.0e3, fission_yield_energy=0.0253):
super().__init__(chain_file, fission_q, dilute_initial, prev_results)
self.round_number = False
self.prev_res = None
@ -178,8 +182,8 @@ class Operator(TransportOperator):
self._rate_helper = DirectReactionRateHelper(
self.reaction_rates.n_nuc, self.reaction_rates.n_react)
self._energy_helper = ChainFissionHelper()
self._fsn_yield_helper = FissionYieldHelper(
self.chain.nuclides, len(self.burnable_mats))
self._fsn_yield_helper = ConstantFissionYieldHelper(
self.chain.nuclides, energy=fission_yield_energy)
def __call__(self, vec, power):
"""Runs a simulation.

42
tests/unit_tests/test_deplete_fission_yields.py Normal file
View file

@ -0,0 +1,42 @@
"""Test the FissionYieldHelpers"""
from collections import namedtuple
import pytest
from openmc.deplete.nuclide import Nuclide, FissionYieldDistribution
from openmc.deplete.helpers import ConstantFissionYieldHelper
@pytest.fixture(scope="module")
def nuclide_bundle():
u5yield_dict = {
0.0253: {"Xe135": 7.85e-4, "Gd155": 4.08e-12, "Sm149": 1.71e-12},
5.0e5: {"Xe135": 7.85e-4, "Sm149": 1.71e-12},
1.40e7: {"Xe135": 4.54e-3, "Gd155": 5.83e-8}}
u235 = Nuclide()
u235.name = "U235"
u235.yield_energies = (0.0253, 5.0e5, 1.40e7)
u235.yield_data = FissionYieldDistribution.from_dict(u5yield_dict)
u8yield_dict = {5.00e5: {"Xe135": 1.12e-3, "Gd155": 1.32e-12}}
u238 = Nuclide()
u238.name = "U238"
u238.yield_energies = (5.00e5, )
u238.yield_data = FissionYieldDistribution.from_dict(u8yield_dict)
xe135 = Nuclide()
xe135.name = "Xe135"
NuclideBundle = namedtuple("NuclideBundle", "u235 u238 xe135")
return NuclideBundle(u235, u238, xe135)
@pytest.mark.parametrize("input_energy, u5_yield_energy", (
(0.0253, 0.0253), (0.01, 0.0253), (4e5, 5e5)))
def test_constant_helper(nuclide_bundle, input_energy, u5_yield_energy):
helper = ConstantFissionYieldHelper(
nuclide_bundle, energy=input_energy)
assert helper.energy == input_energy
assert helper.constant_yields == {
"U235": nuclide_bundle.u235.yield_data[u5_yield_energy],
"U238": nuclide_bundle.u238.yield_data[5.00e5]} # only epithermal
assert helper.constant_yields == helper.weighted_yields(1)

96
tests/unit_tests/test_deplete_helpers.py
View file

@ -1,96 +0,0 @@
"""Test the Operator helpers"""
from unittest.mock import Mock
import pytest
import numpy
from numpy.testing import assert_array_equal
from openmc.deplete.nuclide import Nuclide, FissionYieldDistribution
from openmc.deplete.helpers import FissionYieldHelper
class FissionYieldHelperProxy(FissionYieldHelper):
def __init__(self, chain_nuclides, n_bmats):
super().__init__(chain_nuclides, n_bmats)
self._reaction_tally = Mock()
self.local_indexes = numpy.array([0])
def generate_tallies(self, *args, **kwargs):
# Avoid calls to the C-API
pass
def test_fission_yield_helper():
"""Test the collection of fission yield data using approximated tallies
"""
u5yield_dict = {
0.0253: {"Xe135": 7.85e-4, "Gd155": 4.08e-12, "Sm149": 1.71e-12},
1.40e7: {"Xe135": 4.54e-3, "Gd155": 5.83e-8}}
u235 = Nuclide()
u235.name = "U235"
u235.yield_energies = (0.0253, 1.40e7)
u235.yield_data = FissionYieldDistribution.from_dict(u5yield_dict)
u8yield_dict = {5.00e5: {"Xe135": 1.12e-3, "Gd155": 1.32e-12}}
u238 = Nuclide()
u238.name = "U238"
u238.yield_energies = (5.00e5, )
u238.yield_data = FissionYieldDistribution.from_dict(u8yield_dict)
xe135 = Nuclide()
xe135.name = "Xe135"
n_bmats = 2
helper = FissionYieldHelperProxy([u235, u238, xe135], n_bmats)
assert helper.energy_bounds == (0, 0.0253, 5.00e5, 1.40e7)
# test that tally must be created with nuclides with yields
with pytest.raises(ValueError, match="No overlap"):
helper.set_fissionable_nuclides(["Xe135", ])
act_nucs = helper.set_fissionable_nuclides(["U235", "U238", "Xe135"])
assert act_nucs == ("U235", "U238")
# Emulate getting tally data from transport run
# Test as if this Helper is responsible for one of two materials
# Tally results ordered [n_mat * n_ene, n_fiss_nuc, 3]
tally_res = numpy.zeros((3 * n_bmats, 2, 3))
u5_fiss_rates = numpy.array([1.0, 1.5, 2.0])
u8_fiss_rates = numpy.array([0.0, 1.0, 1.0])
tally_res[:3, 0, 1] = u5_fiss_rates
tally_res[:3, 1, 1] = u8_fiss_rates
helper._reaction_tally.results = tally_res
helper.unpack() # compute yield fractions
# Compare fraction fission rate from helper.reset
exp_results = numpy.empty((1, 3, 2))
# Fraction of fission events in each energy range
u5_vec = u5_fiss_rates / u5_fiss_rates.sum()
u8_vec = u8_fiss_rates / u8_fiss_rates.sum()
exp_results[..., 0] = u5_vec
exp_results[..., 1] = u8_vec
assert_array_equal(helper.results, exp_results)
# Compute and compare the library of fission yields
exp_lib = {
"U235": {
"Xe135": (u5yield_dict[0.0253]["Xe135"] * u5_vec[0]
+ u5yield_dict[1.40e7]["Xe135"] * u5_vec[2]),
"Gd155": (u5yield_dict[0.0253]["Gd155"] * u5_vec[0]
+ u5yield_dict[1.40e7]["Gd155"] * u5_vec[2]),
"Sm149": u5yield_dict[0.0253]["Sm149"] * u5_vec[0],
},
"U238": {
"Xe135": u8yield_dict[5.00e5]["Xe135"] * u8_vec[1],
"Gd155": u8yield_dict[5.00e5]["Gd155"] * u8_vec[1],
}
}
act_library = helper.compute_yields(0)
for parent, sub_yields in exp_lib.items():
assert act_library[parent] == pytest.approx(sub_yields)