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

added tracklength estimator capability for decay-rate tally and fixed other issues

Browse source
This commit is contained in:
samuelshaner 2016-12-15 20:23:07 +00:00
parent 9b94c4dbec
commit ff3d034275
9 changed files with 486 additions and 110 deletions
Download patch file Download diff file Expand all files Collapse all files

30
openmc/material.py
View file

@ -129,7 +129,7 @@ class Material(object):
string = 'Material\n'
string += '{0: <16}{1}{2}\n'.format('\tID', '=\t', self._id)
string += '{0: <16}{1}{2}\n'.format('\tName', '=\t', self._name)
string += '{0: <16}{1}{2}\n'.format('\Temperature', '=\t',
string += '{0: <16}{1}{2}\n'.format('\tTemperature', '=\t',
self._temperature)
string += '{0: <16}{1}{2}'.format('\tDensity', '=\t', self._density)
@ -251,7 +251,7 @@ class Material(object):
Parameters
----------
units : {'g/cm3', 'g/cc', 'km/cm3', 'atom/b-cm', 'atom/cm3', 'sum', 'macro'}
units : {'g/cm3', 'g/cc', 'kg/cm3', 'atom/b-cm', 'atom/cm3', 'sum', 'macro'}
Physical units of density.
density : float, optional
Value of the density. Must be specified unless units is given as
@ -641,37 +641,45 @@ class Material(object):
nucs = []
nuc_densities = []
nuc_density_types = []
for nuclide in nuclides.items():
nuc, nuc_density, nuc_density_type = nuclide[1]
nucs.append(nuc)
nuc_densities.append(nuc_density)
nuc_density_types.append(nuc_density_type)
nucs = np.array(nucs)
nuc_densities = np.array(nuc_densities)
nuc_density_types = np.array(nuc_density_types)
if sum_density:
density = np.sum(nuc_densities)
percent_in_atom = np.all(nuc_density_types == 'ao')
density_in_atom = density > 0.
sum_percent = 0.
awrs = []
for n, nuclide in enumerate(nuclides.items()):
for nuclide in nuclides.items():
awr = openmc.data.atomic_mass(nuclide[0])
if awr is not None:
awrs.append(awr / openmc.data.NEUTRON_MASS)
awrs.append(awr)
else:
raise ValueError(nuclide[0] + " is invalid")
# Convert the weight amounts to atomic amounts
if not percent_in_atom:
for n, nuc in enumerate(nucs):
nuc_densities[n] *= self.average_molar_mass / awrs[n]
# Now that we have the awr, lets finish calculating densities
sum_percent = np.sum(nuc_densities)
nuc_densities = nuc_densities / sum_percent
if not density_in_atom:
sum_percent = 0.
for n, nuc in enumerate(nucs):
x = nuc_densities[n]
sum_percent += x * awrs[n]
sum_percent = 1. / sum_percent
density = -density * sum_percent * \
openmc.data.AVOGADRO / openmc.data.NEUTRON_MASS * 1.E-24
density = -density / self.average_molar_mass * 1.E-24 \
* openmc.data.AVOGADRO
nuc_densities = density * nuc_densities
nuclides = OrderedDict()

9
openmc/mgxs/library.py
View file

@ -524,6 +524,15 @@ class Library(object):
for domain in self.domains:
for mgxs_type in self.mgxs_types:
mgxs = self.get_mgxs(domain, mgxs_type)
if mgxs_type in openmc.mgxs.MDGXS_TYPES:
if self.num_delayed_groups == 0:
mgxs.delayed_groups = None
else:
delayed_groups \
= list(range(1,self.num_delayed_groups+1))
mgxs.delayed_groups = delayed_groups
for tally in mgxs.tallies.values():
tallies_file.append(tally, merge=merge)

73
openmc/mgxs/mdgxs.py
View file

@ -966,6 +966,7 @@ class ChiDelayed(MDGXS):
super(ChiDelayed, self).__init__(domain, domain_type, energy_groups,
delayed_groups, by_nuclide, name)
self._rxn_type = 'chi-delayed'
self._estimator = 'analog'
@property
def scores(self):
@ -987,10 +988,6 @@ class ChiDelayed(MDGXS):
def tally_keys(self):
return ['delayed-nu-fission-in', 'delayed-nu-fission-out']
@property
def estimator(self):
return 'analog'
@property
def rxn_rate_tally(self):
if self._rxn_rate_tally is None:
@ -1018,6 +1015,29 @@ class ChiDelayed(MDGXS):
return self._xs_tally
def get_homogenized_mgxs(self, other_mgxs):
"""Construct a homogenized mgxs with other mgxs objects.
Parameters
----------
other_mgxs : openmc.mgxs.MGXS or Iterable of openmc.mgxs.MGXS
The MGXS to homogenize with this one.
Returns
-------
openmc.mgxs.MGXS
A new homogenized MGXS
Raises
------
ValueError
If the other_mgxs are of a different type.
"""
return self._get_homogenized_mgxs(other_mgxs, 'delayed-nu-fission-in')
def get_slice(self, nuclides=[], groups=[], delayed_groups=[]):
"""Build a sliced ChiDelayed for the specified nuclides and energy
groups.
@ -1586,6 +1606,28 @@ class Beta(MDGXS):
return self._xs_tally
def get_homogenized_mgxs(self, other_mgxs):
"""Construct a homogenized mgxs with other mgxs objects.
Parameters
----------
other_mgxs : openmc.mgxs.MGXS or Iterable of openmc.mgxs.MGXS
The MGXS to homogenize with this one.
Returns
-------
openmc.mgxs.MGXS
A new homogenized MGXS
Raises
------
ValueError
If the other_mgxs are of a different type.
"""
return self._get_homogenized_mgxs(other_mgxs, 'nu-fission')
class DecayRate(MDGXS):
r"""The decay rate for delayed neutron precursors.
@ -1703,7 +1745,6 @@ class DecayRate(MDGXS):
super(DecayRate, self).__init__(domain, domain_type, energy_groups,
delayed_groups, by_nuclide, name)
self._rxn_type = 'decay-rate'
self._estimator = 'analog'
@property
def scores(self):
@ -1738,6 +1779,28 @@ class DecayRate(MDGXS):
return self._xs_tally
def get_homogenized_mgxs(self, other_mgxs):
"""Construct a homogenized mgxs with other mgxs objects.
Parameters
----------
other_mgxs : openmc.mgxs.MGXS or Iterable of openmc.mgxs.MGXS
The MGXS to homogenize with this one.
Returns
-------
openmc.mgxs.MGXS
A new homogenized MGXS
Raises
------
ValueError
If the other_mgxs are of a different type.
"""
return self._get_homogenized_mgxs(other_mgxs, 'delayed-nu-fission')
@add_metaclass(ABCMeta)
class MatrixMDGXS(MDGXS):

118
openmc/mgxs/mgxs.py
View file

@ -996,6 +996,102 @@ class MGXS(object):
avg_xs.sparse = self.sparse
return avg_xs
def _get_homogenized_mgxs(self, other_mgxs, denom_score='flux'):
"""Construct a homogenized mgxs with other mgxs objects.
Parameters
----------
other_mgxs : openmc.mgxs.MGXS or Iterable of openmc.mgxs.MGXS
The MGXS to homogenize with this one.
denom_score : str
The denominator score in the denominator of computing the MGXS.
Returns
-------
openmc.mgxs.MGXS
A new homogenized MGXS
Raises
------
ValueError
If the other_mgxs are of a different type.
"""
# Check type of denom score
cv.check_type('denom_score', denom_score, str)
# Construct a collection of the subdomain filter bins to average across
if isinstance(other_mgxs, openmc.mgxs.MGXS):
other_mgxs = [other_mgxs]
cv.check_iterable_type('other_mgxs', other_mgxs, openmc.mgxs.MGXS)
for mgxs in other_mgxs:
if mgxs.rxn_type != self.rxn_type:
msg = 'Not able to homogenize two MGXS with different rxn types'
raise ValueError(msg)
# Clone this MGXS to initialize the homogenized version
homogenized_mgxs = copy.deepcopy(self)
homogenized_mgxs._derived = True
name = '{} + '.format(self.domain.name)
# Get the domain filter
filter_type = _DOMAIN_TO_FILTER[self.domain_type]
self_filter = self.rxn_rate_tally.find_filter(filter_type)
# Get the rxn rate and denom tallies
rxn_rate_tally = self.rxn_rate_tally
denom_tally = self.tallies[denom_score]
for mgxs in other_mgxs:
# Swap the domain filter bins for the other mgxs rxn rate tally
other_rxn_rate_tally = copy.deepcopy(mgxs.rxn_rate_tally)
other_filter = other_rxn_rate_tally.find_filter(filter_type)
other_filter._bins = self_filter._bins
# Swap the domain filter bins for the denom tally
other_denom_tally = copy.deepcopy(mgxs.tallies[denom_score])
other_filter = other_denom_tally.find_filter(filter_type)
other_filter._bins = self_filter._bins
# Add the rxn rate and denom tallies
rxn_rate_tally += other_rxn_rate_tally
denom_tally += other_denom_tally
# Update the name for the homogenzied MGXS
name += '{} + '.format(mgxs.domain.name)
# Set the properties of the homogenized MGXS
homogenized_mgxs._rxn_rate_tally = rxn_rate_tally
homogenized_mgxs.tallies[denom_score] = denom_tally
homogenized_mgxs._domain.name = name[:-3]
return homogenized_mgxs
def get_homogenized_mgxs(self, other_mgxs):
"""Construct a homogenized mgxs with other mgxs objects.
Parameters
----------
other_mgxs : openmc.mgxs.MGXS or Iterable of openmc.mgxs.MGXS
The MGXS to homogenize with this one.
Returns
-------
openmc.mgxs.MGXS
A new homogenized MGXS
Raises
------
ValueError
If the other_mgxs are of a different type.
"""
return self._get_homogenized_mgxs(other_mgxs, 'flux')
def get_slice(self, nuclides=[], groups=[]):
"""Build a sliced MGXS for the specified nuclides and energy groups.
@ -4561,6 +4657,28 @@ class Chi(MGXS):
return self._xs_tally
def get_homogenized_mgxs(self, other_mgxs):
"""Construct a homogenized mgxs with other mgxs objects.
Parameters
----------
other_mgxs : openmc.mgxs.MGXS or Iterable of openmc.mgxs.MGXS
The MGXS to homogenize with this one.
Returns
-------
openmc.mgxs.MGXS
A new homogenized MGXS
Raises
------
ValueError
If the other_mgxs are of a different type.
"""
return self._get_homogenized_mgxs(other_mgxs, 'nu-fission-in')
def get_slice(self, nuclides=[], groups=[]):
"""Build a sliced Chi for the specified nuclides and energy groups.

72
openmc/mgxs_library.py
View file

@ -15,8 +15,8 @@ from openmc.checkvalue import check_type, check_value, check_greater_than, \
# Supported incoming particle MGXS angular treatment representations
_REPRESENTATIONS = ['isotropic', 'angle']
_SCATTER_TYPES = ['tabular', 'legendre', 'histogram']
_XS_SHAPES = ["[G][G'][Order]", "[G]", "[G']", "[G][G']", "[DG]", "[G][DG]",
"[G'][DG]", "[G][G'][DG]"]
_XS_SHAPES = ["[G][G'][Order]", "[G]", "[G']", "[G][G']", "[DG]", "[DG][G]",
"[DG][G']", "[DG][G][G']"]
class XSdata(object):
@ -145,11 +145,11 @@ class XSdata(object):
[DG]: beta, decay_rate
[G][DG]: delayed_nu_fission, beta, decay_rate
[DG][G]: delayed_nu_fission, beta, decay_rate
[G'][DG]: chi_delayed
[DG][G']: chi_delayed
[G][G'][DG]: delayed_nu_fission
[DG][G][G']: delayed_nu_fission
"""
@ -296,11 +296,11 @@ class XSdata(object):
self._xs_shapes["[G][G']"] = (self.energy_groups.num_groups,
self.energy_groups.num_groups)
self._xs_shapes["[DG]"] = (self.num_delayed_groups,)
self._xs_shapes["[G][DG]"] = (self.energy_groups.num_groups,
self._xs_shapes["[DG][G]"] = (self.energy_groups.num_groups,
self.num_delayed_groups)
self._xs_shapes["[G'][DG]"] = (self.energy_groups.num_groups,
self._xs_shapes["[DG'][G']"] = (self.energy_groups.num_groups,
self.num_delayed_groups)
self._xs_shapes["[G][G'][DG]"] = (self.energy_groups.num_groups,
self._xs_shapes["[DG][G][G']"] = (self.energy_groups.num_groups,
self.energy_groups.num_groups,
self.num_delayed_groups)
@ -634,7 +634,7 @@ class XSdata(object):
"""
# Get the accepted shapes for this xs
shapes = [self.xs_shapes["[G']"], self.xs_shapes["[G'][DG]"]]
shapes = [self.xs_shapes["[G']"], self.xs_shapes["[DG][G']"]]
# Convert to a numpy array so we can easily get the shape for checking
chi_delayed = np.asarray(chi_delayed)
@ -664,7 +664,7 @@ class XSdata(object):
"""
# Get the accepted shapes for this xs
shapes = [self.xs_shapes["[DG]"], self.xs_shapes["[G][DG]"]]
shapes = [self.xs_shapes["[DG]"], self.xs_shapes["[DG][G]"]]
# Convert to a numpy array so we can easily get the shape for checking
beta = np.asarray(beta)
@ -694,7 +694,7 @@ class XSdata(object):
"""
# Get the accepted shapes for this xs
shapes = [self.xs_shapes["[DG]"], self.xs_shapes["[G][DG]"]]
shapes = [self.xs_shapes["[DG]"], self.xs_shapes["[DG][G]"]]
# Convert to a numpy array so we can easily get the shape for checking
decay_rate = np.asarray(decay_rate)
@ -856,7 +856,7 @@ class XSdata(object):
"""
# Get the accepted shapes for this xs
shapes = [self.xs_shapes["[G][DG]"], self.xs_shapes["[G][G'][DG]"]]
shapes = [self.xs_shapes["[G][DG]"], self.xs_shapes["[DG][G][G']"]]
# Convert to a numpy array so we can easily get the shape for checking
delayed_nu_fission = np.asarray(delayed_nu_fission)
@ -1644,6 +1644,54 @@ class XSdata(object):
self._multiplicity_matrix[i] = \
np.nan_to_num(self._multiplicity_matrix[i])
def set_inverse_velocity_mgxs(self, inverse_velocity, temperature=294.,
nuclide='total', xs_type='macro',
subdomain=None):
"""This method allows for an openmc.mgxs.InverseVelocity
to be used to set the absorption cross section for this XSdata object.
Parameters
----------
inverse_velocity : openmc.mgxs.InverseVelocity
MGXS Object containing the inverse velocity for the domain of
interest.
temperature : float
Temperature (in Kelvin) of the data. Defaults to room temperature
(294K).
nuclide : str
Individual nuclide (or 'total' if obtaining material-wise data)
to gather data for. Defaults to 'total'.
xs_type: {'macro', 'micro'}
Provide the macro or micro cross section in units of cm^-1 or
barns. Defaults to 'macro'.
subdomain : iterable of int
If the MGXS contains a mesh domain type, the subdomain parameter
specifies which mesh cell (i.e., [i, j, k] index) to use.
See also
--------
openmc.mgxs.Library.create_mg_library()
openmc.mgxs.Library.get_xsdata
"""
check_type('inverse_velocity', inverse_velocity, openmc.mgxs.InverseVelocity)
check_value('energy_groups', inverse_velocity.energy_groups,
[self.energy_groups])
check_value('domain_type', inverse_velocity.domain_type,
openmc.mgxs.DOMAIN_TYPES)
check_type('temperature', temperature, Real)
check_value('temperature', temperature, self.temperatures)
i = np.where(self.temperatures == temperature)[0][0]
if self.representation == 'isotropic':
self._inverse_velocity[i] = inverse_velocity.get_xs\
(nuclides=nuclide, xs_type=xs_type,
subdomains=subdomain)
elif self.representation == 'angle':
msg = 'Angular-Dependent MGXS have not yet been implemented'
raise ValueError(msg)
def to_hdf5(self, file):
"""Write XSdata to an HDF5 file

11
openmc/tallies.py
View file

@ -705,7 +705,7 @@ class Tally(object):
"""
# Two tallys must have the same number of filters
# Two tallies must have the same number of filters
if len(self.filters) != len(other.filters):
return False
@ -3486,9 +3486,14 @@ class Tallies(cv.CheckedList):
for d in derivs:
root_element.append(d.to_xml_element())
def export_to_xml(self):
def export_to_xml(self, path='tallies.xml'):
"""Create a tallies.xml file that can be used for a simulation.
Parameters
----------
path : str
Path to file to write. Defaults to 'tallies.xml'.
"""
root_element = ET.Element("tallies")
@ -3501,5 +3506,5 @@ class Tallies(cv.CheckedList):
# Write the XML Tree to the tallies.xml file
tree = ET.ElementTree(root_element)
tree.write("tallies.xml", xml_declaration=True,
tree.write(path, xml_declaration=True,
encoding='utf-8', method="xml")

7
src/input_xml.F90
View file

@ -3668,13 +3668,6 @@ contains
end if
case ('decay-rate')
t % score_bins(j) = SCORE_DECAY_RATE
! Set tally estimator to analog for CE mode
! (MG mode has all data available without a collision being
! necessary)
if (run_CE) then
t % estimator = ESTIMATOR_ANALOG
end if
case ('delayed-nu-fission')
t % score_bins(j) = SCORE_DELAYED_NU_FISSION
if (t % find_filter(FILTER_ENERGYOUT) > 0) then

3
src/nuclide_header.F90
View file

@ -621,7 +621,8 @@ module nuclide_header
case (EMISSION_DELAYED)
if (this % n_precursor > 0) then
if (present(group)) then
if (present(group) .and. group < &
size(this % reactions(this % index_fission(1)) % products)) then
! If delayed group specified, determine yield immediately
associate(p => this % reactions(this % index_fission(1)) % products(1 + group))
nu = p % yield % evaluate(E)

273
src/tally.F90
View file

@ -679,8 +679,8 @@ contains
end if
end if
case (SCORE_DECAY_RATE)
! make sure the correct energy is used
if (t % estimator == ESTIMATOR_TRACKLENGTH) then
E = p % E
@ -691,12 +691,134 @@ contains
! Set the delayedgroup filter index
dg_filter = t % find_filter(FILTER_DELAYEDGROUP)
if (survival_biasing) then
! No fission events occur if survival biasing is on -- need to
! calculate fraction of absorptions that would have resulted in
! delayed-nu-fission
if (micro_xs(p % event_nuclide) % absorption > ZERO .and. &
nuclides(p % event_nuclide) % fissionable) then
if (t % estimator == ESTIMATOR_ANALOG) then
if (survival_biasing) then
! No fission events occur if survival biasing is on -- need to
! calculate fraction of absorptions that would have resulted in
! delayed-nu-fission
if (micro_xs(p % event_nuclide) % absorption > ZERO .and. &
nuclides(p % event_nuclide) % fissionable) then
! Check if the delayed group filter is present
if (dg_filter > 0) then
select type(filt => t % filters(dg_filter) % obj)
type is (DelayedGroupFilter)
! Loop over all delayed group bins and tally to them
! individually
do d_bin = 1, filt % n_bins
! Get the delayed group for this bin
d = filt % groups(d_bin)
! Compute the yield for this delayed group
yield = nuclides(p % event_nuclide) &
% nu(E, EMISSION_DELAYED, d)
associate (rxn => nuclides(p % event_nuclide) % &
reactions(nuclides(p % event_nuclide) % index_fission(1)))
! Compute the score
score = p % absorb_wgt * yield * &
micro_xs(p % event_nuclide) % fission &
/ micro_xs(p % event_nuclide) % absorption &
* rxn % products(1 + d) % decay_rate
end associate
! Tally to bin
call score_fission_delayed_dg(t, d_bin, score, score_index)
end do
cycle SCORE_LOOP
end select
else
! If the delayed group filter is not present, compute the score
! by accumulating the absorbed weight times the decay rate times
! the fraction of the delayed-nu-fission xs to the absorption xs
! for all delayed groups.
score = ZERO
associate (rxn => nuclides(p % event_nuclide) % &
reactions(nuclides(p % event_nuclide) % index_fission(1)))
! We need to be careful not to overshoot the number of delayed
! groups since this could cause the range of the rxn % products
! array to be exceeded. Hence, we use the size of this array
! and not the MAX_DELAYED_GROUPS constant for this loop.
do d = 1, size(rxn % products) - 2
score = score + rxn % products(1 + d) % decay_rate * &
p % absorb_wgt * micro_xs(p % event_nuclide) % fission *&
nuclides(p % event_nuclide) % nu(E, EMISSION_DELAYED, d)&
/ micro_xs(p % event_nuclide) % absorption
end do
end associate
end if
end if
else
! Skip any non-fission events
if (.not. p % fission) cycle SCORE_LOOP
! If there is no outgoing energy filter, than we only need to
! score to one bin. For the score to be 'analog', we need to
! score the number of particles that were banked in the fission
! bank. Since this was weighted by 1/keff, we multiply by keff
! to get the proper score. Loop over the neutrons produced from
! fission and check which ones are delayed. If a delayed neutron is
! encountered, add its contribution to the fission bank to the
! score.
score = ZERO
! loop over number of particles banked
do k = 1, p % n_bank
! get the delayed group
g = fission_bank(n_bank - p % n_bank + k) % delayed_group
! Case for tallying delayed emissions
if (g /= 0) then
! Accumulate the decay rate times delayed nu fission score
associate (rxn => nuclides(p % event_nuclide) % &
reactions(nuclides(p % event_nuclide) % index_fission(1)))
! determine score based on bank site weight and keff.
score = score + keff * fission_bank(n_bank - p % n_bank + k) &
% wgt * rxn % products(1 + g) % decay_rate
end associate
! if the delayed group filter is present, tally to corresponding
! delayed group bin if it exists
if (dg_filter > 0) then
! declare the delayed group filter type
select type(filt => t % filters(dg_filter) % obj)
type is (DelayedGroupFilter)
! loop over delayed group bins until the corresponding bin is
! found
do d_bin = 1, filt % n_bins
d = filt % groups(d_bin)
! check whether the delayed group of the particle is equal to
! the delayed group of this bin
if (d == g) then
call score_fission_delayed_dg(t, d_bin, score, score_index)
end if
end do
end select
! Reset the score to zero
score = ZERO
end if
end if
end do
end if
else
! Check if tally is on a single nuclide
if (i_nuclide > 0) then
! Check if the delayed group filter is present
if (dg_filter > 0) then
@ -711,17 +833,14 @@ contains
d = filt % groups(d_bin)
! Compute the yield for this delayed group
yield = nuclides(p % event_nuclide) &
% nu(E, EMISSION_DELAYED, d)
yield = nuclides(i_nuclide) % nu(E, EMISSION_DELAYED, d)
associate (rxn => nuclides(p % event_nuclide) % &
reactions(nuclides(p % event_nuclide) % index_fission(1)))
associate (rxn => nuclides(i_nuclide) % &
reactions(nuclides(i_nuclide) % index_fission(1)))
! Compute the score
score = p % absorb_wgt * yield * &
micro_xs(p % event_nuclide) % fission &
/ micro_xs(p % event_nuclide) % absorption &
* rxn % products(1 + d) % decay_rate
! Compute the score and tally to bin
score = micro_xs(i_nuclide) % fission * yield * flux * &
atom_density * rxn % products(1 + d) % decay_rate
end associate
! Tally to bin
@ -746,78 +865,90 @@ contains
! and not the MAX_DELAYED_GROUPS constant for this loop.
do d = 1, size(rxn % products) - 2
score = score + rxn % products(1 + d) % decay_rate * &
p % absorb_wgt * micro_xs(p % event_nuclide) % fission *&
nuclides(p % event_nuclide) % nu(E, EMISSION_DELAYED, d)&
/ micro_xs(p % event_nuclide) % absorption
score = score + micro_xs(i_nuclide) % fission * flux * &
nuclides(i_nuclide) % nu(E, EMISSION_DELAYED) * &
atom_density * rxn % products(1 + d) % decay_rate
end do
end associate
end if
end if
else
! Skip any non-fission events
if (.not. p % fission) cycle SCORE_LOOP
! If there is no outgoing energy filter, than we only need to
! score to one bin. For the score to be 'analog', we need to
! score the number of particles that were banked in the fission
! bank. Since this was weighted by 1/keff, we multiply by keff
! to get the proper score. Loop over the neutrons produced from
! fission and check which ones are delayed. If a delayed neutron is
! encountered, add its contribution to the fission bank to the
! score.
! Tally is on total nuclides
else
score = ZERO
! Check if the delayed group filter is present
if (dg_filter > 0) then
select type(filt => t % filters(dg_filter) % obj)
type is (DelayedGroupFilter)
! loop over number of particles banked
do k = 1, p % n_bank
! Loop over all nuclides in the current material
do l = 1, materials(p % material) % n_nuclides
! get the delayed group
g = fission_bank(n_bank - p % n_bank + k) % delayed_group
! Get atom density
atom_density_ = materials(p % material) % atom_density(l)
! Case for tallying delayed emissions
if (g /= 0) then
! Get index in nuclides array
i_nuc = materials(p % material) % nuclide(l)
! Accumulate the decay rate times delayed nu fission score
associate (rxn => nuclides(p % event_nuclide) % &
reactions(nuclides(p % event_nuclide) % index_fission(1)))
if (nuclides(i_nuc) % fissionable) then
! determine score based on bank site weight and keff.
score = score + keff * fission_bank(n_bank - p % n_bank + k) &
% wgt * rxn % products(1 + g) % decay_rate
end associate
! Loop over all delayed group bins and tally to them
! individually
do d_bin = 1, filt % n_bins
! if the delayed group filter is present, tally to corresponding
! delayed group bin if it exists
if (dg_filter > 0) then
! Get the delayed group for this bin
d = filt % groups(d_bin)
! declare the delayed group filter type
select type(filt => t % filters(dg_filter) % obj)
type is (DelayedGroupFilter)
! Get the yield for the desired nuclide and delayed group
yield = nuclides(i_nuc) % nu(E, EMISSION_DELAYED, d)
! loop over delayed group bins until the corresponding bin is
! found
do d_bin = 1, filt % n_bins
d = filt % groups(d_bin)
associate (rxn => nuclides(i_nuc) % &
reactions(nuclides(i_nuc) % index_fission(1)))
! check whether the delayed group of the particle is equal to
! the delayed group of this bin
if (d == g) then
! Compute the score
score = micro_xs(i_nuc) % fission * yield * flux * &
atom_density_ * rxn % products(1 + d) % decay_rate
end associate
! Tally to bin
call score_fission_delayed_dg(t, d_bin, score, score_index)
end if
end do
end select
end do
end if
end do
cycle SCORE_LOOP
end select
else
! Reset the score to zero
score = ZERO
end if
score = ZERO
! Loop over all nuclides in the current material
do l = 1, materials(p % material) % n_nuclides
! Get atom density
atom_density_ = materials(p % material) % atom_density(l)
! Get index in nuclides array
i_nuc = materials(p % material) % nuclide(l)
if (nuclides(i_nuc) % fissionable) then
associate (rxn => nuclides(i_nuc) % &
reactions(nuclides(i_nuc) % index_fission(1)))
! We need to be careful not to overshoot the number of delayed
! groups since this could cause the range of the rxn % products
! array to be exceeded. Hence, we use the size of this array
! and not the MAX_DELAYED_GROUPS constant for this loop.
do d = 1, size(rxn % products) - 2
! Accumulate the contribution from each nuclide
score = score + micro_xs(i_nuc) % fission * nuclides(i_nuc) %&
nu(E, EMISSION_DELAYED) * atom_density_ * flux * &
rxn % products(1 + d) % decay_rate
end do
end associate
end if
end do
end if
end do
! If the delayed group filter is present, cycle because the
! score_fission_delayed_dg(...) has already tallied the score
if (dg_filter > 0) then
cycle SCORE_LOOP
end if
end if