Implement FluxCollapseHelper that is used by Operator

This commit is contained in:
Paul Romano 2020-08-04 14:35:12 -05:00
parent b7df2d95e1
commit 76fc6fb3e7
3 changed files with 150 additions and 15 deletions

View file

@ -225,11 +225,11 @@ class TransportOperator(ABC):
class ReactionRateHelper(ABC):
"""Abstract class for generating reaction rates for operators
Responsible for generating reaction rate tallies for burnable
materials, given nuclides and scores from the operator.
Responsible for generating reaction rate tallies for burnable materials,
given nuclides and scores from the operator.
Reaction rates are passed back to the operator for be used in
an :class:`openmc.deplete.OperatorResult` instance
Reaction rates are passed back to the operator for be used in an
:class:`openmc.deplete.OperatorResult` instance.
Parameters
----------
@ -246,7 +246,6 @@ class ReactionRateHelper(ABC):
def __init__(self, n_nucs, n_react):
self._nuclides = None
self._rate_tally = None
self._results_cache = empty((n_nucs, n_react))
@abstractmethod
@ -262,7 +261,6 @@ class ReactionRateHelper(ABC):
def nuclides(self, nuclides):
check_type("nuclides", nuclides, list, str)
self._nuclides = nuclides
self._rate_tally.nuclides = nuclides
@abstractmethod
def get_material_rates(self, mat_id, nuc_index, react_index):
@ -281,8 +279,7 @@ class ReactionRateHelper(ABC):
def divide_by_adens(self, number):
"""Normalize reaction rates by number of nuclides
Acts on the current material examined by
:meth:`get_material_rates`
Acts on the current material examined by :meth:`get_material_rates`
Parameters
----------
@ -559,7 +556,6 @@ class TalliedFissionYieldHelper(FissionYieldHelper):
self._fission_rate_tally = Tally()
self._fission_rate_tally.writable = False
self._fission_rate_tally.scores = ['fission']
self._fission_rate_tally.filters = [MaterialFilter(materials)]
def update_tally_nuclides(self, nuclides):

View file

@ -7,13 +7,14 @@ from numbers import Real
import bisect
from collections import defaultdict
from numpy import dot, zeros, newaxis
from numpy import dot, zeros, newaxis, asarray
from . import comm
from openmc.checkvalue import check_type, check_greater_than
from openmc.data import JOULE_PER_EV
from openmc.data import JOULE_PER_EV, REACTION_NAME
from openmc.lib import (
Tally, MaterialFilter, EnergyFilter, EnergyFunctionFilter)
import openmc.lib
from .abc import (
ReactionRateHelper, NormalizationHelper, FissionYieldHelper,
TalliedFissionYieldHelper)
@ -21,7 +22,7 @@ from .abc import (
__all__ = (
"DirectReactionRateHelper", "ChainFissionHelper", "EnergyScoreHelper"
"SourceRateHelper", "ConstantFissionYieldHelper", "FissionYieldCutoffHelper",
"AveragedFissionYieldHelper")
"AveragedFissionYieldHelper", "FluxCollapseHelper")
# -------------------------------------
# Helpers for generating reaction rates
@ -43,6 +44,14 @@ class DirectReactionRateHelper(ReactionRateHelper):
nuclides : list of str
All nuclides with desired reaction rates.
"""
def __init__(self, n_nuc, n_react):
super().__init__(n_nuc, n_react)
self._rate_tally = None
@ReactionRateHelper.nuclides.setter
def nuclides(self, nuclides):
ReactionRateHelper.nuclides.fset(self, nuclides)
self._rate_tally.nuclides = nuclides
def generate_tallies(self, materials, scores):
"""Produce one-group reaction rate tally
@ -92,6 +101,107 @@ class DirectReactionRateHelper(ReactionRateHelper):
return self._results_cache
class FluxCollapseHelper(ReactionRateHelper):
"""Class that generates tallies for one-group rates
.. versionadded:: 0.12.1
Parameters
----------
n_nucs : int
Number of burnable nuclides tracked by :class:`openmc.deplete.Operator`
n_react : int
Number of reactions tracked by :class:`openmc.deplete.Operator`
energies : iterable of float
Energy group boundaries for flux spectrum in [eV]
Attributes
----------
nuclides : list of str
All nuclides with desired reaction rates.
"""
def __init__(self, n_nucs, n_reacts, energies):
super().__init__(n_nucs, n_reacts)
self._energies = asarray(energies)
def generate_tallies(self, materials, scores):
"""Produce multigroup flux spectrum tally
Uses the :mod:`openmc.lib` module to generate a multigroup flux tally
for each burnable material.
Parameters
----------
materials : iterable of :class:`openmc.Material`
Burnable materials in the problem. Used to construct a
:class:`openmc.MaterialFilter`
scores : iterable of str
Reaction identifiers, e.g. ``"(n, fission)"``, ``"(n, gamma)"``,
needed for the reaction rate tally.
"""
self._materials = materials
# Convert reactions to MT values (needed when collapsing)
mt_values = {v: k for k, v in REACTION_NAME.items()}
mt_values['fission'] = 18
self._mts = [mt_values[x] for x in scores]
# Create flux tally with material and energy filters
self._flux_tally = Tally()
self._flux_tally.writable = False
self._flux_tally.filters = [
MaterialFilter(materials),
EnergyFilter(self._energies)
]
self._flux_tally.scores = ['flux']
def get_material_rates(self, mat_index, nuc_index, react_index):
"""Return an array of reaction rates for a material
Parameters
----------
mat_index : int
Index for material
nuc_index : iterable of int
Index for each nuclide in :attr:`nuclides` in the
desired reaction rate matrix
react_index : iterable of int
Index for each reaction scored in the tally
Returns
-------
rates : numpy.ndarray
Array with shape ``(n_nuclides, n_rxns)`` with the reaction rates in
this material
"""
self._results_cache.fill(0.0)
# Get flux for specified material
shape = (len(self._materials), len(self._energies) - 1)
mean_value = self._flux_tally.mean.reshape(shape)
flux = mean_value[mat_index]
mat = self._materials[mat_index]
for name, i_nuc in zip(self.nuclides, nuc_index):
# Determine density of nuclide
# TODO: O(N) search for nuclide is not ideal
for nucname, density in zip(mat.nuclides, mat.densities):
if nucname == name:
break
else:
raise ValueError("Couldn't find {} in material {}".format(name, mat_id))
for mt, i_react in zip(self._mts, react_index):
# Use flux to collapse reaction rate (per N)
nuc = openmc.lib.nuclides[name]
rate_per_nuc = nuc.collapse_rate(mt, self._energies, flux)
# Multiply by density to get absolute reaction rate
self._results_cache[i_nuc, i_react] = rate_per_nuc * density
return self._results_cache
# ------------------------------------------
# Helpers for obtaining normalization factor
# ------------------------------------------

View file

@ -27,7 +27,7 @@ from .results_list import ResultsList
from .helpers import (
DirectReactionRateHelper, ChainFissionHelper, ConstantFissionYieldHelper,
FissionYieldCutoffHelper, AveragedFissionYieldHelper, EnergyScoreHelper,
SourceRateHelper)
SourceRateHelper, FluxCollapseHelper)
__all__ = ["Operator", "OperatorResult"]
@ -113,6 +113,18 @@ class Operator(TransportOperator):
``fission_yield_mode``. Will be passed directly on to the
helper. Passing a value of None will use the defaults for
the associated helper.
reaction_rate_mode : {"direct", "flux"}
Indicate how one-group reaction rates should be calculated. The "direct"
method tallies transmutation reaction rates directly. The "flux" method
tallies a multigroup flux spectrum and then collapses one-group reaction
rates after a transport solve.
.. versionadded:: 0.12.1
reaction_rate_energies : iterable of float
Energy group boundaries that are to be used for calculating a multigroup
flux spectrum when the "flux" based ``reaction_rate_mode`` is being used.
.. versionadded:: 0.12.1
reduce_chain : bool, optional
If True, use :meth:`openmc.deplete.Chain.reduce` to reduce the
depletion chain up to ``reduce_chain_level``. Default is False.
@ -170,6 +182,7 @@ class Operator(TransportOperator):
diff_burnable_mats=False, normalization_mode="fission-q",
fission_q=None, dilute_initial=1.0e3,
fission_yield_mode="constant", fission_yield_opts=None,
reaction_rate_mode="direct", reaction_rate_energies=None,
reduce_chain=False, reduce_chain_level=None):
check_value('fission yield mode', fission_yield_mode,
self._fission_helpers.keys())
@ -239,8 +252,24 @@ class Operator(TransportOperator):
self.local_mats, self._burnable_nucs, self.chain.reactions)
# Get classes to assist working with tallies
self._rate_helper = DirectReactionRateHelper(
self.reaction_rates.n_nuc, self.reaction_rates.n_react)
if reaction_rate_mode == "direct":
self._rate_helper = DirectReactionRateHelper(
self.reaction_rates.n_nuc, self.reaction_rates.n_react)
elif reaction_rate_mode == "flux":
# Ensure energy group boundaries were specified
if reaction_rate_energies is None:
raise ValueError(
"Energy group boundaries must be specified in the "
"reaction_rate_energies argument when reaction_rate_mode is"
"set to 'flux'.")
self._rate_helper = FluxCollapseHelper(
self.reaction_rates.n_nuc,
self.reaction_rates.n_react,
reaction_rate_energies
)
else:
raise ValueError("Invalid reaction rate mode.")
if normalization_mode == "fission-q":
self._normalization_helper = ChainFissionHelper()