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

Merging with upstream

Browse source
This commit is contained in:
Adam Nelson 2016-11-19 05:25:32 -05:00
commit 4459c93932
54 changed files with 3468 additions and 557 deletions
Download patch file Download diff file Expand all files Collapse all files

26
docs/source/io_formats/statepoint.rst
View file

@ -162,6 +162,18 @@ if run_mode == 'k-eigenvalue':
Width of each mesh cell in each dimension.
**/tallies/derivatives/derivative <id>/independent variable** (*char[]*)
Independent variable of tally derivative
**/tallies/derivatives/derivative <id>/material** (*int*)
ID of the perturbed material
**/tallies/derivatives/derivative <id>/nuclide** (*char[]*)
Alias of the perturbed nuclide
**/tallies/n_tallies** (*int*)
Number of user-defined tallies.
@ -204,6 +216,10 @@ if run_mode == 'k-eigenvalue':
Array of nuclides to tally. Note that if no nuclide is specified in the user
input, a single 'total' nuclide appears here.
**/tallies/tally <uid>/derivative** (*int*)
ID of the derivative applied to the tally.
**/tallies/tally <uid>/n_score_bins** (*int*)
Number of scoring bins for a single nuclide. In general, this can be greater
@ -224,12 +240,12 @@ if run_mode == 'k-eigenvalue':
Tallying moment orders for Legendre and spherical harmonic tally expansions
(*e.g.*, 'P2', 'Y1,2', etc.).
**/tallies/tally <uid>/results** (Compound type)
**/tallies/tally <uid>/results** (*double[][][2]*)
Accumulated sum and sum-of-squares for each bin of the i-th tally. This is a
two-dimensional array, the first dimension of which represents combinations
of filter bins and the second dimensions of which represents scoring
bins. Each element of the array has fields 'sum' and 'sum_sq'.
Accumulated sum and sum-of-squares for each bin of the i-th tally. The first
dimension represents combinations of filter bins, the second dimensions
represents scoring bins, and the third dimension has two entries for the sum
and the sum-of-squares.
**/source_present** (*int*)

2
docs/source/pythonapi/index.rst
View file

@ -293,6 +293,7 @@ Multi-group Cross Sections
openmc.mgxs.NuScatterXS
openmc.mgxs.NuScatterMatrixXS
openmc.mgxs.PromptNuFissionXS
openmc.mgxs.PromptNuFissionMatrixXS
openmc.mgxs.ScatterXS
openmc.mgxs.ScatterMatrixXS
openmc.mgxs.TotalXS
@ -309,6 +310,7 @@ Multi-delayed-group Cross Sections
openmc.mgxs.MDGXS
openmc.mgxs.ChiDelayed
openmc.mgxs.DelayedNuFissionXS
openmc.mgxs.DelayedNuFissionMatrixXS
openmc.mgxs.Beta
openmc.mgxs.DecayRate

40
docs/source/usersguide/input.rst
View file

@ -1916,6 +1916,13 @@ The ``<tally>`` element accepts the following sub-elements:
*Default*: "all"
:derivative:
The id of a ``derivative`` element. This derivative will be applied to all
scores in the tally. Differential tallies are currently only implemented
for collision and analog estimators.
*Default*: None
``<mesh>`` Element
------------------
@ -1944,6 +1951,39 @@ attributes/sub-elements:
One of ``<upper_right>`` or ``<width>`` must be specified, but not both
(even if they are consistent with one another).
``<derivative>`` Element
------------------------
OpenMC can take the first-order derivative of many tallies with respect to
material perturbations. It works by propagating a derivative through the
transport equation. Essentially, OpenMC keeps track of how each particle's
weight would change as materials are perturbed, and then accounts for that
weight change in the tallies. Note that this assumes material perturbations are
small enough not to change the distribution of fission sites. This element has
the following attributes/sub-elements:
:id:
A unique integer that can be used to identify the derivative.
:variable:
The independent variable of the derivative. Accepted options are "density",
"nuclide_density", and "temperature". A "density" derivative will give the
derivative with respect to the density of the material in [g / cm^3]. A
"nuclide_density" derivative will give the derivative with respect to the
density of a particular nuclide in units of [atom / b / cm]. A
"temperature" derivative is with respect to a material temperature in units
of [K]. The temperature derivative requires windowed multipole to be
turned on. Note also that the temperature derivative only accounts for
resolved resonance Doppler broadening. It does not account for thermal
expansion, S(a, b) scattering, resonance scattering, or unresolved Doppler
broadening.
:material:
The perturbed material. (Necessary for all derivative types)
:nuclide:
The perturbed nuclide. (Necessary only for "nuclide_density")
``<assume_separate>`` Element
-----------------------------

1
openmc/__init__.py
View file

@ -16,6 +16,7 @@ from openmc.universe import *
from openmc.mesh import *
from openmc.filter import *
from openmc.trigger import *
from openmc.tally_derivative import *
from openmc.tallies import *
from openmc.mgxs_library import *
from openmc.cmfd import *

2
openmc/element.py
View file

@ -212,7 +212,7 @@ class Element(object):
# its natural nuclides
else:
for nuclide in natural_nuclides:
abundances[nuclide] = NATURAL_ABUNDNACE[nuclide]
abundances[nuclide] = NATURAL_ABUNDANCE[nuclide]
# Modify mole fractions if enrichment provided
if enrichment is not None:

14
openmc/mgxs/library.py
View file

@ -582,7 +582,7 @@ class Library(object):
----------
domain : Material or Cell or Universe or Integral
The material, cell, or universe object of interest (or its ID)
mgxs_type : {'total', 'transport', 'nu-transport', 'absorption', 'capture', 'fission', 'nu-fission', 'kappa-fission', 'scatter', 'nu-scatter', 'scatter matrix', 'nu-scatter matrix', 'multiplicity matrix', 'nu-fission matrix', chi', 'chi-prompt', 'inverse-velocity', 'prompt-nu-fission', 'delayed-nu-fission', 'chi-delayed', 'beta'}
mgxs_type : {'total', 'transport', 'nu-transport', 'absorption', 'capture', 'fission', 'nu-fission', 'kappa-fission', 'scatter', 'nu-scatter', 'scatter matrix', 'nu-scatter matrix', 'multiplicity matrix', 'nu-fission matrix', chi', 'chi-prompt', 'inverse-velocity', 'prompt-nu-fission', 'prompt-nu-fission matrix', 'delayed-nu-fission', 'delayed-nu-fission matrix', 'chi-delayed', 'beta'}
The type of multi-group cross section object to return
Returns
@ -1026,12 +1026,24 @@ class Library(object):
nuclide=[nuclide],
subdomain=subdomain)
if 'prompt-nu-fission matrix' in self.mgxs_types:
mymgxs = self.get_mgxs(domain, 'prompt-nu-fission matrix')
xsdata.set_prompt_nu_fission_mgxs(mymgxs, xs_type=xs_type,
nuclide=[nuclide],
subdomain=subdomain)
if 'delayed-nu-fission' in self.mgxs_types:
mymgxs = self.get_mgxs(domain, 'delayed-nu-fission')
xsdata.set_delayed_nu_fission_mgxs(mymgxs, xs_type=xs_type,
nuclide=[nuclide],
subdomain=subdomain)
if 'delayed-nu-fission matrix' in self.mgxs_types:
mymgxs = self.get_mgxs(domain, 'delayed-nu-fission matrix')
xsdata.set_delayed_nu_fission_mgxs(mymgxs, xs_type=xs_type,
nuclide=[nuclide],
subdomain=subdomain)
if 'beta' in self.mgxs_types:
mymgxs = self.get_mgxs(domain, 'nu-fission')
xsdata.set_beta_mgxs(mymgxs, xs_type=xs_type, nuclide=[nuclide],

614
openmc/mgxs/mdgxs.py
View file

@ -21,7 +21,8 @@ import openmc.checkvalue as cv
MDGXS_TYPES = ['delayed-nu-fission',
'chi-delayed',
'beta',
'decay-rate']
'decay-rate',
'delayed-nu-fission matrix']
# Maximum number of delayed groups, from src/constants.F90
MAX_DELAYED_GROUPS = 8
@ -211,7 +212,7 @@ class MDGXS(MGXS):
Parameters
----------
mdgxs_type : {'delayed-nu-fission', 'chi-delayed', 'beta', 'decay-rate'}
mdgxs_type : {'delayed-nu-fission', 'chi-delayed', 'beta', 'decay-rate', 'delayed-nu-fission matrix'}
The type of multi-delayed-group cross section object to return
domain : openmc.Material or openmc.Cell or openmc.Universe or
openmc.Mesh
@ -249,6 +250,9 @@ class MDGXS(MGXS):
mdgxs = Beta(domain, domain_type, energy_groups, delayed_groups)
elif mdgxs_type == 'decay-rate':
mdgxs = DecayRate(domain, domain_type, energy_groups, delayed_groups)
elif mdgxs_type == 'delayed-nu-fission matrix':
mdgxs = DelayedNuFissionMatrixXS(domain, domain_type, energy_groups,
delayed_groups)
mdgxs.by_nuclide = by_nuclide
mdgxs.name = name
@ -1733,3 +1737,609 @@ class DecayRate(MDGXS):
super(DecayRate, self)._compute_xs()
return self._xs_tally
@add_metaclass(ABCMeta)
class MatrixMDGXS(MDGXS):
"""An abstract multi-delayed-group cross section for some energy group and
delayed group structure within some spatial domain. This class is
specifically intended for cross sections which depend on both the incoming
and outgoing energy groups and are therefore represented by matrices.
An example of this is the delayed-nu-fission matrix.
This class can be used for both OpenMC input generation and tally data
post-processing to compute spatially-homogenized and energy-integrated
multi-group and multi-delayed-group cross sections for downstream neutronics
calculations.
NOTE: Users should instantiate the subclasses of this abstract class.
Parameters
----------
domain : openmc.Material or openmc.Cell or openmc.Universe or openmc.Mesh
The domain for spatial homogenization
domain_type : {'material', 'cell', 'distribcell', 'universe', 'mesh'}
The domain type for spatial homogenization
energy_groups : openmc.mgxs.EnergyGroups
The energy group structure for energy condensation
by_nuclide : bool
If true, computes cross sections for each nuclide in domain
name : str, optional
Name of the multi-group cross section. Used as a label to identify
tallies in OpenMC 'tallies.xml' file.
delayed_groups : list of int
Delayed groups to filter out the xs
Attributes
----------
name : str, optional
Name of the multi-group cross section
rxn_type : str
Reaction type (e.g., 'total', 'nu-fission', etc.)
by_nuclide : bool
If true, computes cross sections for each nuclide in domain
domain : Material or Cell or Universe or Mesh
Domain for spatial homogenization
domain_type : {'material', 'cell', 'distribcell', 'universe', 'mesh'}
Domain type for spatial homogenization
energy_groups : openmc.mgxs.EnergyGroups
Energy group structure for energy condensation
delayed_groups : list of int
Delayed groups to filter out the xs
tally_trigger : openmc.Trigger
An (optional) tally precision trigger given to each tally used to
compute the cross section
scores : list of str
The scores in each tally used to compute the multi-group cross section
filters : list of openmc.Filter
The filters in each tally used to compute the multi-group cross section
tally_keys : list of str
The keys into the tallies dictionary for each tally used to compute
the multi-group cross section
estimator : {'tracklength', 'collision', 'analog'}
The tally estimator used to compute the multi-group cross section
tallies : collections.OrderedDict
OpenMC tallies needed to compute the multi-group cross section
rxn_rate_tally : openmc.Tally
Derived tally for the reaction rate tally used in the numerator to
compute the multi-group cross section. This attribute is None
unless the multi-group cross section has been computed.
xs_tally : openmc.Tally
Derived tally for the multi-group cross section. This attribute
is None unless the multi-group cross section has been computed.
num_subdomains : int
The number of subdomains is unity for 'material', 'cell' and 'universe'
domain types. This is equal to the number of cell instances
for 'distribcell' domain types (it is equal to unity prior to loading
tally data from a statepoint file) and the number of mesh cells for
'mesh' domain types.
num_nuclides : int
The number of nuclides for which the multi-group cross section is
being tracked. This is unity if the by_nuclide attribute is False.
nuclides : Iterable of str or 'sum'
The optional user-specified nuclides for which to compute cross
sections (e.g., 'U238', 'O16'). If by_nuclide is True but nuclides
are not specified by the user, all nuclides in the spatial domain
are included. This attribute is 'sum' if by_nuclide is false.
sparse : bool
Whether or not the MGXS' tallies use SciPy's LIL sparse matrix format
for compressed data storage
loaded_sp : bool
Whether or not a statepoint file has been loaded with tally data
derived : bool
Whether or not the MGXS is merged from one or more other MGXS
hdf5_key : str
The key used to index multi-group cross sections in an HDF5 data store
"""
@property
def filters(self):
# Create the non-domain specific Filters for the Tallies
group_edges = self.energy_groups.group_edges
energy = openmc.EnergyFilter(group_edges)
energyout = openmc.EnergyoutFilter(group_edges)
if self.delayed_groups is not None:
delayed = openmc.DelayedGroupFilter(self.delayed_groups)
return [[energy], [delayed, energy, energyout]]
else:
return [[energy], [energy, energyout]]
def get_xs(self, in_groups='all', out_groups='all',
subdomains='all', nuclides='all',
xs_type='macro', order_groups='increasing',
row_column='inout', value='mean', delayed_groups='all',
squeeze=True, **kwargs):
"""Returns an array of multi-group cross sections.
This method constructs a 4D NumPy array for the requested
multi-group cross section data for one or more subdomains
(1st dimension), delayed groups (2nd dimension), energy groups in
(3rd dimension), energy groups out (4th dimension), and nuclides
(5th dimension).
Parameters
----------
in_groups : Iterable of Integral or 'all'
Incoming energy groups of interest. Defaults to 'all'.
out_groups : Iterable of Integral or 'all'
Outgoing energy groups of interest. Defaults to 'all'.
subdomains : Iterable of Integral or 'all'
Subdomain IDs of interest. Defaults to 'all'.
nuclides : Iterable of str or 'all' or 'sum'
A list of nuclide name strings (e.g., ['U235', 'U238']). The
special string 'all' will return the cross sections for all
nuclides in the spatial domain. The special string 'sum' will
return the cross section summed over all nuclides. Defaults to
'all'.
xs_type: {'macro', 'micro'}
Return the macro or micro cross section in units of cm^-1 or barns.
Defaults to 'macro'.
order_groups: {'increasing', 'decreasing'}
Return the cross section indexed according to increasing or
decreasing energy groups (decreasing or increasing energies).
Defaults to 'increasing'.
row_column: {'inout', 'outin'}
Return the cross section indexed first by incoming group and
second by outgoing group ('inout'), or vice versa ('outin').
Defaults to 'inout'.
value : {'mean', 'std_dev', 'rel_err'}
A string for the type of value to return. Defaults to 'mean'.
delayed_groups : list of int or 'all'
Delayed groups of interest. Defaults to 'all'.
squeeze : bool
A boolean representing whether to eliminate the extra dimensions
of the multi-dimensional array to be returned. Defaults to True.
Returns
-------
numpy.ndarray
A NumPy array of the multi-group cross section indexed in the order
each group and subdomain is listed in the parameters.
Raises
------
ValueError
When this method is called before the multi-group cross section is
computed from tally data.
"""
cv.check_value('value', value, ['mean', 'std_dev', 'rel_err'])
cv.check_value('xs_type', xs_type, ['macro', 'micro'])
# FIXME: Unable to get microscopic xs for mesh domain because the mesh
# cells do not know the nuclide densities in each mesh cell.
if self.domain_type == 'mesh' and xs_type == 'micro':
msg = 'Unable to get micro xs for mesh domain since the mesh ' \
'cells do not know the nuclide densities in each mesh cell.'
raise ValueError(msg)
filters = []
filter_bins = []
# Construct a collection of the domain filter bins
if not isinstance(subdomains, string_types):
cv.check_iterable_type('subdomains', subdomains, Integral,
max_depth=3)
for subdomain in subdomains:
filters.append(_DOMAIN_TO_FILTER[self.domain_type])
filter_bins.append((subdomain,))
# Construct list of energy group bounds tuples for all requested groups
if not isinstance(in_groups, string_types):
cv.check_iterable_type('groups', in_groups, Integral)
for group in in_groups:
filters.append(openmc.EnergyFilter)
filter_bins.append((
self.energy_groups.get_group_bounds(group),))
# Construct list of energy group bounds tuples for all requested groups
if not isinstance(out_groups, string_types):
cv.check_iterable_type('groups', out_groups, Integral)
for group in out_groups:
filters.append(openmc.EnergyoutFilter)
filter_bins.append((
self.energy_groups.get_group_bounds(group),))
# Construct list of delayed group tuples for all requested groups
if not isinstance(delayed_groups, string_types):
cv.check_type('delayed groups', delayed_groups, list, int)
for delayed_group in delayed_groups:
filters.append(openmc.DelayedGroupFilter)
filter_bins.append((delayed_group,))
# Construct a collection of the nuclides to retrieve from the xs tally
if self.by_nuclide:
if nuclides == 'all' or nuclides == 'sum' or nuclides == ['sum']:
query_nuclides = self.get_nuclides()
else:
query_nuclides = nuclides
else:
query_nuclides = ['total']
# Use tally summation if user requested the sum for all nuclides
if nuclides == 'sum' or nuclides == ['sum']:
xs_tally = self.xs_tally.summation(nuclides=query_nuclides)
xs = xs_tally.get_values(filters=filters, filter_bins=filter_bins,
value=value)
else:
xs = self.xs_tally.get_values(filters=filters,
filter_bins=filter_bins,
nuclides=query_nuclides, value=value)
# Divide by atom number densities for microscopic cross sections
if xs_type == 'micro':
if self.by_nuclide:
densities = self.get_nuclide_densities(nuclides)
else:
densities = self.get_nuclide_densities('sum')
if value == 'mean' or value == 'std_dev':
xs /= densities[np.newaxis, :, np.newaxis]
# Eliminate the trivial score dimension
xs = np.squeeze(xs, axis=len(xs.shape) - 1)
xs = np.nan_to_num(xs)
if in_groups == 'all':
num_in_groups = self.num_groups
else:
num_in_groups = len(in_groups)
if out_groups == 'all':
num_out_groups = self.num_groups
else:
num_out_groups = len(out_groups)
if delayed_groups == 'all':
num_delayed_groups = self.num_delayed_groups
else:
num_delayed_groups = len(delayed_groups)
# Reshape tally data array with separate axes for domain and energy
num_subdomains = int(xs.shape[0] / (num_in_groups * num_out_groups *
num_delayed_groups))
new_shape = (num_subdomains, num_delayed_groups, num_in_groups,
num_out_groups)
new_shape += xs.shape[1:]
xs = np.reshape(xs, new_shape)
# Transpose the matrix if requested by user
if row_column == 'outin':
xs = np.swapaxes(xs, 2, 3)
# Reverse data if user requested increasing energy groups since
# tally data is stored in order of increasing energies
if order_groups == 'increasing':
xs = xs[:, :, ::-1, ::-1, :]
if squeeze:
xs = np.squeeze(xs)
xs = np.atleast_2d(xs)
return xs
def get_slice(self, nuclides=[], in_groups=[], out_groups=[],
delayed_groups=[]):
"""Build a sliced MatrixMDGXS object for the specified nuclides and
energy groups.
This method constructs a new MdGXS to encapsulate a subset of the data
represented by this MdGXS. The subset of data to include in the tally
slice is determined by the nuclides, energy groups, and delayed groups
specified in the input parameters.
Parameters
----------
nuclides : list of str
A list of nuclide name strings
(e.g., ['U235', 'U238']; default is [])
in_groups : list of int
A list of incoming energy group indices starting at 1 for the high
energies (e.g., [1, 2, 3]; default is [])
out_groups : list of int
A list of outgoing energy group indices starting at 1 for the high
energies (e.g., [1, 2, 3]; default is [])
delayed_groups : list of int
A list of delayed group indices
(e.g., [1, 2, 3]; default is [])
Returns
-------
openmc.mgxs.MatrixMDGXS
A new MatrixMDGXS object which encapsulates the subset of data
requested for the nuclide(s) and/or energy group(s) requested in
the parameters.
"""
# Call super class method and null out derived tallies
slice_xs = super(MatrixMDGXS, self).get_slice(nuclides, in_groups,
delayed_groups)
slice_xs._rxn_rate_tally = None
slice_xs._xs_tally = None
# Slice outgoing energy groups if needed
if len(out_groups) != 0:
filter_bins = []
for group in out_groups:
group_bounds = self.energy_groups.get_group_bounds(group)
filter_bins.append(group_bounds)
filter_bins = [tuple(filter_bins)]
# Slice each of the tallies across energyout groups
for tally_type, tally in slice_xs.tallies.items():
if tally.contains_filter(openmc.EnergyoutFilter):
tally_slice = tally.get_slice(
filters=[openmc.EnergyoutFilter],
filter_bins=filter_bins)
slice_xs.tallies[tally_type] = tally_slice
slice_xs.sparse = self.sparse
return slice_xs
def print_xs(self, subdomains='all', nuclides='all', xs_type='macro'):
"""Prints a string representation for the multi-group cross section.
Parameters
----------
subdomains : Iterable of Integral or 'all'
The subdomain IDs of the cross sections to include in the report.
Defaults to 'all'.
nuclides : Iterable of str or 'all' or 'sum'
The nuclides of the cross-sections to include in the report. This
may be a list of nuclide name strings (e.g., ['U235', 'U238']).
The special string 'all' will report the cross sections for all
nuclides in the spatial domain. The special string 'sum' will
report the cross sections summed over all nuclides. Defaults to
'all'.
xs_type: {'macro', 'micro'}
Return the macro or micro cross section in units of cm^-1 or barns.
Defaults to 'macro'.
"""
# Construct a collection of the subdomains to report
if not isinstance(subdomains, string_types):
cv.check_iterable_type('subdomains', subdomains, Integral)
elif self.domain_type == 'distribcell':
subdomains = np.arange(self.num_subdomains, dtype=np.int)
elif self.domain_type == 'mesh':
xyz = [range(1, x+1) for x in self.domain.dimension]
subdomains = list(itertools.product(*xyz))
else:
subdomains = [self.domain.id]
# Construct a collection of the nuclides to report
if self.by_nuclide:
if nuclides == 'all':
nuclides = self.get_nuclides()
if nuclides == 'sum':
nuclides = ['sum']
else:
cv.check_iterable_type('nuclides', nuclides, string_types)
else:
nuclides = ['sum']
cv.check_value('xs_type', xs_type, ['macro', 'micro'])
# Build header for string with type and domain info
string = 'Multi-Delayed-Group XS\n'
string += '{0: <16}=\t{1}\n'.format('\tReaction Type', self.rxn_type)
string += '{0: <16}=\t{1}\n'.format('\tDomain Type', self.domain_type)
string += '{0: <16}=\t{1}\n'.format('\tDomain ID', self.domain.id)
# Generate the header for an individual XS
xs_header = '\tCross Sections [{0}]:'.format(self.get_units(xs_type))
# If cross section data has not been computed, only print string header
if self.tallies is None:
print(string)
return
string += '{0: <16}\n'.format('\tEnergy Groups:')
template = '{0: <12}Group {1} [{2: <10} - {3: <10}MeV]\n'
# Loop over energy groups ranges
for group in range(1, self.num_groups + 1):
bounds = self.energy_groups.get_group_bounds(group)
string += template.format('', group, bounds[0], bounds[1])
# Loop over all subdomains
for subdomain in subdomains:
if self.domain_type == 'distribcell':
string += '{: <16}=\t{}\n'.format('\tSubdomain', subdomain)
# Loop over all Nuclides
for nuclide in nuclides:
# Build header for nuclide type
if xs_type != 'sum':
string += '{: <16}=\t{}\n'.format('\tNuclide', nuclide)
# Build header for cross section type
string += '{: <16}\n'.format(xs_header)
if self.delayed_groups is not None:
for delayed_group in self.delayed_groups:
template = '{0: <12}Delayed Group {1}:\t'
string += template.format('', delayed_group)
string += '\n'
template = '{0: <12}Group {1} -> Group {2}:\t\t'
# Loop over incoming/outgoing energy groups ranges
for in_group in range(1, self.num_groups + 1):
for out_group in range(1, self.num_groups + 1):
string += template.format('', in_group, out_group)
average = self.get_xs([in_group], [out_group],
[subdomain], [nuclide],
xs_type=xs_type,
value='mean',
delayed_groups=[delayed_group])
rel_err = self.get_xs([in_group], [out_group],
[subdomain], [nuclide],
xs_type=xs_type,
value='rel_err',
delayed_groups=[delayed_group])
average = average.flatten()[0]
rel_err = rel_err.flatten()[0] * 100.
string += '{:.2e} +/- {:.2e}%'.format(average,
rel_err)
string += '\n'
string += '\n'
string += '\n'
else:
template = '{0: <12}Group {1} -> Group {2}:\t\t'
# Loop over incoming/outgoing energy groups ranges
for in_group in range(1, self.num_groups + 1):
for out_group in range(1, self.num_groups + 1):
string += template.format('', in_group, out_group)
average = self.get_xs([in_group], [out_group],
[subdomain], [nuclide],
xs_type=xs_type, value='mean')
rel_err = self.get_xs([in_group], [out_group],
[subdomain], [nuclide],
xs_type=xs_type, value='rel_err')
average = average.flatten()[0]
rel_err = rel_err.flatten()[0] * 100.
string += '{:.2e} +/- {:.2e}%'.format(average,
rel_err)
string += '\n'
string += '\n'
string += '\n'
string += '\n'
print(string)
class DelayedNuFissionMatrixXS(MatrixMDGXS):
r"""A fission delayed neutron production matrix multi-group cross section.
This class can be used for both OpenMC input generation and tally data
post-processing to compute spatially-homogenized and energy-integrated
multi-group fission neutron production cross sections for multi-group
neutronics calculations. At a minimum, one needs to set the
:attr:`DelayedNuFissionMatrixXS.energy_groups` and
:attr:`DelayedNuFissionMatrixXS.domain` properties. Tallies for the flux and
appropriate reaction rates over the specified domain are generated
automatically via the :attr:`DelayedNuFissionMatrixXS.tallies` property,
which can then be appended to a :class:`openmc.Tallies` instance.
For post-processing, the :meth:`MGXS.load_from_statepoint` will pull in the
necessary data to compute multi-group cross sections from a
:class:`openmc.StatePoint` instance. The derived multi-group cross section
can then be obtained from the :attr:`DelayedNuFissionMatrixXS.xs_tally`
property.
For a spatial domain :math:`V`, energy group :math:`[E_g,E_{g-1}]`, and
delayed group :math:`d`, the fission delayed neutron production cross
section is calculated as:
.. math::
\langle \nu\sigma_{f,g'\rightarrow g} \phi \rangle &= \int_{r \in V} dr
\int_{4\pi} d\Omega' \int_{E_{g'}}^{E_{g'-1}} dE' \int_{E_g}^{E_{g-1}} dE
\; \chi(E) \nu\sigma_f^d (r, E') \psi(r, E', \Omega')\\
\langle \phi \rangle &= \int_{r \in V} dr \int_{4\pi} d\Omega
\int_{E_g}^{E_{g-1}} dE \; \psi (r, E, \Omega) \\
\nu\sigma_{f,g'\rightarrow g} &= \frac{\langle \nu\sigma_{f,g'\rightarrow
g}^d \phi \rangle}{\langle \phi \rangle}
Parameters
----------
domain : openmc.Material or openmc.Cell or openmc.Universe or openmc.Mesh
The domain for spatial homogenization
domain_type : {'material', 'cell', 'distribcell', 'universe', 'mesh'}
The domain type for spatial homogenization
groups : openmc.mgxs.EnergyGroups
The energy group structure for energy condensation
by_nuclide : bool
If true, computes cross sections for each nuclide in domain
name : str, optional
Name of the multi-group cross section. Used as a label to identify
tallies in OpenMC 'tallies.xml' file.
delayed_groups : list of int
Delayed groups to filter out the xs
Attributes
----------
name : str, optional
Name of the multi-group cross section
rxn_type : str
Reaction type (e.g., 'total', 'nu-fission', etc.)
by_nuclide : bool
If true, computes cross sections for each nuclide in domain
domain : Material or Cell or Universe or Mesh
Domain for spatial homogenization
domain_type : {'material', 'cell', 'distribcell', 'universe', 'mesh'}
Domain type for spatial homogenization
energy_groups : openmc.mgxs.EnergyGroups
Energy group structure for energy condensation
delayed_groups : list of int
Delayed groups to filter out the xs
tally_trigger : openmc.Trigger
An (optional) tally precision trigger given to each tally used to
compute the cross section
scores : list of str
The scores in each tally used to compute the multi-group cross section
filters : list of openmc.Filter
The filters in each tally used to compute the multi-group cross section
tally_keys : list of str
The keys into the tallies dictionary for each tally used to compute
the multi-group cross section
estimator : {'tracklength', 'analog'}
The tally estimator used to compute the multi-group cross section
tallies : collections.OrderedDict
OpenMC tallies needed to compute the multi-group cross section. The keys
are strings listed in the :attr:`DelayedNuFissionXS.tally_keys` property
and values are instances of :class:`openmc.Tally`.
rxn_rate_tally : openmc.Tally
Derived tally for the reaction rate tally used in the numerator to
compute the multi-group cross section. This attribute is None
unless the multi-group cross section has been computed.
xs_tally : openmc.Tally
Derived tally for the multi-group cross section. This attribute
is None unless the multi-group cross section has been computed.
num_subdomains : int
The number of subdomains is unity for 'material', 'cell' and 'universe'
domain types. When the This is equal to the number of cell instances
for 'distribcell' domain types (it is equal to unity prior to loading
tally data from a statepoint file).
num_nuclides : int
The number of nuclides for which the multi-group cross section is
being tracked. This is unity if the by_nuclide attribute is False.
nuclides : Iterable of str or 'sum'
The optional user-specified nuclides for which to compute cross
sections (e.g., 'U-238', 'O-16'). If by_nuclide is True but nuclides
are not specified by the user, all nuclides in the spatial domain
are included. This attribute is 'sum' if by_nuclide is false.
sparse : bool
Whether or not the MGXS' tallies use SciPy's LIL sparse matrix format
for compressed data storage
loaded_sp : bool
Whether or not a statepoint file has been loaded with tally data
derived : bool
Whether or not the MGXS is merged from one or more other MGXS
hdf5_key : str
The key used to index multi-group cross sections in an HDF5 data store
"""
def __init__(self, domain=None, domain_type=None, energy_groups=None,
delayed_groups=None, by_nuclide=False, name=''):
super(DelayedNuFissionMatrixXS, self).__init__(domain, domain_type,
energy_groups,
delayed_groups,
by_nuclide, name)
self._rxn_type = 'delayed-nu-fission'
self._hdf5_key = 'delayed-nu-fission matrix'
self._estimator = 'analog'
self._valid_estimators = ['analog']

128
openmc/mgxs/mgxs.py
View file

@ -36,7 +36,8 @@ MGXS_TYPES = ['total',
'chi',
'chi-prompt',
'inverse-velocity',
'prompt-nu-fission']
'prompt-nu-fission',
'prompt-nu-fission matrix']
# Supported domain types
DOMAIN_TYPES = ['cell',
@ -451,7 +452,7 @@ class MGXS(object):
Parameters
----------
mgxs_type : {'total', 'transport', 'nu-transport', 'absorption', 'capture', 'fission', 'nu-fission', 'kappa-fission', 'scatter', 'nu-scatter', 'scatter matrix', 'nu-scatter matrix', 'multiplicity matrix', 'nu-fission matrix', 'chi', 'chi-prompt', 'inverse-velocity', 'prompt-nu-fission'}
mgxs_type : {'total', 'transport', 'nu-transport', 'absorption', 'capture', 'fission', 'nu-fission', 'kappa-fission', 'scatter', 'nu-scatter', 'scatter matrix', 'nu-scatter matrix', 'multiplicity matrix', 'nu-fission matrix', 'chi', 'chi-prompt', 'inverse-velocity', 'prompt-nu-fission', 'prompt-nu-fission matrix'}
The type of multi-group cross section object to return
domain : openmc.Material or openmc.Cell or openmc.Universe or openmc.Mesh
The domain for spatial homogenization
@ -512,6 +513,8 @@ class MGXS(object):
mgxs = InverseVelocity(domain, domain_type, energy_groups)
elif mgxs_type == 'prompt-nu-fission':
mgxs = PromptNuFissionXS(domain, domain_type, energy_groups)
elif mgxs_type == 'prompt-nu-fission matrix':
mgxs = PromptNuFissionMatrixXS(domain, domain_type, energy_groups)
mgxs.by_nuclide = by_nuclide
mgxs.name = name
@ -1755,7 +1758,7 @@ class MatrixMGXS(MGXS):
Returns
-------
ndarray
numpy.ndarray
A NumPy array of the multi-group cross section indexed in the order
each group and subdomain is listed in the parameters.
@ -3654,7 +3657,7 @@ class ScatterMatrixXS(MatrixMGXS):
Returns
-------
ndarray
numpy.ndarray
A NumPy array of the multi-group cross section indexed in the order
each group and subdomain is listed in the parameters.
@ -5232,3 +5235,120 @@ class PromptNuFissionXS(MGXS):
super(PromptNuFissionXS, self).__init__(domain, domain_type, groups,
by_nuclide, name)
self._rxn_type = 'prompt-nu-fission'
class PromptNuFissionMatrixXS(MatrixMGXS):
r"""A prompt fission neutron production matrix multi-group cross section.
This class can be used for both OpenMC input generation and tally data
post-processing to compute spatially-homogenized and energy-integrated
multi-group cross sections for multi-group neutronics calculations. At a
minimum, one needs to set the :attr:`PromptNuFissionMatrixXS.energy_groups`
and :attr:`PromptNuFissionMatrixXS.domain` properties. Tallies for the flux
and appropriate reaction rates over the specified domain are generated
automatically via the :attr:`PromptNuFissionMatrixXS.tallies` property,
which can then be appended to a :class:`openmc.Tallies` instance.
For post-processing, the :meth:`MGXS.load_from_statepoint` will pull in the
necessary data to compute multi-group cross sections from a
:class:`openmc.StatePoint` instance. The derived multi-group cross section
can then be obtained from the :attr:`PromptNuFissionMatrixXS.xs_tally`
property.
For a spatial domain :math:`V` and energy group :math:`[E_g,E_{g-1}]`, the
fission spectrum is calculated as:
.. math::
\langle \nu\sigma_{f,g'\rightarrow g} \phi \rangle &= \int_{r \in V} dr
\int_{4\pi} d\Omega' \int_{E_{g'}}^{E_{g'-1}} dE' \int_{E_g}^{E_{g-1}} dE
\; \chi(E) \nu\sigma_f^p (r, E') \psi(r, E', \Omega')\\
\langle \phi \rangle &= \int_{r \in V} dr \int_{4\pi} d\Omega
\int_{E_g}^{E_{g-1}} dE \; \psi (r, E, \Omega) \\
\nu\sigma_{f,g'\rightarrow g} &= \frac{\langle \nu\sigma_{f,g'\rightarrow
g}^p \phi \rangle}{\langle \phi \rangle}
Parameters
----------
domain : openmc.Material or openmc.Cell or openmc.Universe or openmc.Mesh
The domain for spatial homogenization
domain_type : {'material', 'cell', 'distribcell', 'universe', 'mesh'}
The domain type for spatial homogenization
groups : openmc.mgxs.EnergyGroups
The energy group structure for energy condensation
by_nuclide : bool
If true, computes cross sections for each nuclide in domain
name : str, optional
Name of the multi-group cross section. Used as a label to identify
tallies in OpenMC 'tallies.xml' file.
Attributes
----------
name : str, optional
Name of the multi-group cross section
rxn_type : str
Reaction type (e.g., 'total', 'nu-fission', etc.)
by_nuclide : bool
If true, computes cross sections for each nuclide in domain
domain : Material or Cell or Universe or Mesh
Domain for spatial homogenization
domain_type : {'material', 'cell', 'distribcell', 'universe', 'mesh'}
Domain type for spatial homogenization
energy_groups : openmc.mgxs.EnergyGroups
Energy group structure for energy condensation
tally_trigger : openmc.Trigger
An (optional) tally precision trigger given to each tally used to
compute the cross section
scores : list of str
The scores in each tally used to compute the multi-group cross section
filters : list of openmc.Filter
The filters in each tally used to compute the multi-group cross section
tally_keys : list of str
The keys into the tallies dictionary for each tally used to compute
the multi-group cross section
estimator : {'tracklength', 'collision', 'analog'}
The tally estimator used to compute the multi-group cross section
tallies : collections.OrderedDict
OpenMC tallies needed to compute the multi-group cross section. The keys
are strings listed in the :attr:`PromptNuFissionXS.tally_keys` property
and values are instances of :class:`openmc.Tally`.
rxn_rate_tally : openmc.Tally
Derived tally for the reaction rate tally used in the numerator to
compute the multi-group cross section. This attribute is None
unless the multi-group cross section has been computed.
xs_tally : openmc.Tally
Derived tally for the multi-group cross section. This attribute
is None unless the multi-group cross section has been computed.
num_subdomains : int
The number of subdomains is unity for 'material', 'cell' and 'universe'
domain types. This is equal to the number of cell instances
for 'distribcell' domain types (it is equal to unity prior to loading
tally data from a statepoint file).
num_nuclides : int
The number of nuclides for which the multi-group cross section is
being tracked. This is unity if the by_nuclide attribute is False.
nuclides : Iterable of str or 'sum'
The optional user-specified nuclides for which to compute cross
sections (e.g., 'U-238', 'O-16'). If by_nuclide is True but nuclides
are not specified by the user, all nuclides in the spatial domain
are included. This attribute is 'sum' if by_nuclide is false.
sparse : bool
Whether or not the MGXS' tallies use SciPy's LIL sparse matrix format
for compressed data storage
loaded_sp : bool
Whether or not a statepoint file has been loaded with tally data
derived : bool
Whether or not the MGXS is merged from one or more other MGXS
hdf5_key : str
The key used to index multi-group cross sections in an HDF5 data store
"""
def __init__(self, domain=None, domain_type=None,
groups=None, by_nuclide=False, name=''):
super(PromptNuFissionMatrixXS, self).__init__(domain, domain_type,
groups, by_nuclide, name)
self._rxn_type = 'prompt-nu-fission'
self._hdf5_key = 'prompt-nu-fission matrix'
self._estimator = 'analog'
self._valid_estimators = ['analog']

35
openmc/mgxs_library.py
View file

@ -137,11 +137,11 @@ class XSdata(object):
[G][G'][Order]: scatter_matrix
[G]: total, absorption, fission, kappa_fission, nu_fission,
prompt_nu_fission, inverse_velocity
prompt_nu_fission, delayed_nu_fission, inverse_velocity
[G']: chi, chi_prompt, chi_delayed
[G][G']: multiplicity_matrix, nu_fission
[G][G']: multiplicity_matrix, nu_fission, prompt_nu_fission
[DG]: beta, decay_rate
@ -149,6 +149,8 @@ class XSdata(object):
[G'][DG]: chi_delayed
[G][G'][DG]: delayed_nu_fission
"""
def __init__(self, name, energy_groups, temperatures=[294.],
@ -298,6 +300,9 @@ class XSdata(object):
self.num_delayed_groups)
self._xs_shapes["[G'][DG]"] = (self.energy_groups.num_groups,
self.num_delayed_groups)
self._xs_shapes["[G][G'][DG]"] = (self.energy_groups.num_groups,
self.energy_groups.num_groups,
self.num_delayed_groups)
self._xs_shapes["[G][G'][Order]"] \
= (self.energy_groups.num_groups,
self.energy_groups.num_groups, self.num_orders)
@ -816,7 +821,7 @@ class XSdata(object):
"""
# Get the accepted shapes for this xs
shapes = [self.xs_shapes["[G]"]]
shapes = [self.xs_shapes["[G]"], self.xs_shapes["[G][G']"]]
# Convert to a numpy array so we can easily get the shape for checking
prompt_nu_fission = np.asarray(prompt_nu_fission)
@ -850,7 +855,7 @@ class XSdata(object):
"""
# Get the accepted shapes for this xs
shapes = [self.xs_shapes["[G][DG]"]]
shapes = [self.xs_shapes["[G][DG]"], self.xs_shapes["[G][G'][DG]"]]
# Convert to a numpy array so we can easily get the shape for checking
delayed_nu_fission = np.asarray(delayed_nu_fission)
@ -1086,12 +1091,15 @@ class XSdata(object):
def set_prompt_nu_fission_mgxs(self, prompt_nu_fission, temperature=294.,
nuclide='total', xs_type='macro',
subdomain=None):
"""This method allows for an openmc.mgxs.PromptNuFissionXS to be used to
set the prompt-nu-fission cross section for this XSdata object.
"""Sets the prompt-nu-fission cross section.
This method allows for an openmc.mgxs.PromptNuFissionXS or
openmc.mgxs.PromptNuFissionMatrixXS to be used to set the
prompt-nu-fission cross section for this XSdata object.
Parameters
----------
prompt_nu_fission: openmc.mgxs.PromptNuFissionXS
prompt_nu_fission: openmc.mgxs.PromptNuFissionXS or openmc.mgxs.PromptNuFissionMatrixXS
MGXS Object containing the prompt-nu-fission cross section
for the domain of interest.
temperature : float
@ -1115,7 +1123,8 @@ class XSdata(object):
"""
check_type('prompt_nu_fission', prompt_nu_fission,
(openmc.mgxs.PromptNuFissionXS,))
(openmc.mgxs.PromptNuFissionXS,
openmc.mgxs.PromptNuFissionMatrixXS))
check_value('energy_groups', prompt_nu_fission.energy_groups,
[self.energy_groups])
check_value('domain_type', prompt_nu_fission.domain_type,
@ -1139,12 +1148,13 @@ class XSdata(object):
def set_delayed_nu_fission_mgxs(self, delayed_nu_fission, temperature=294.,
nuclide='total', xs_type='macro',
subdomain=None):
"""This method allows for an openmc.mgxs.DelayedNuFissionXS to be used
to set the delayed-nu-fission cross section for this XSdata object.
"""This method allows for an openmc.mgxs.DelayedNuFissionXS or
openmc.mgxs.DelayedNuFissionMatrixXS to be used to set the
delayed-nu-fission cross section for this XSdata object.
Parameters
----------
delayed_nu_fission: openmc.mgxs.DelayedNuFissionXS
delayed_nu_fission: openmc.mgxs.DelayedNuFissionXS or openmc.mgxs.DelayedNuFissionMatrixXS
MGXS Object containing the delayed-nu-fission cross section
for the domain of interest.
temperature : float
@ -1168,7 +1178,8 @@ class XSdata(object):
"""
check_type('delayed_nu_fission', delayed_nu_fission,
(openmc.mgxs.DelayedNuFissionXS,))
(openmc.mgxs.DelayedNuFissionXS,
openmc.mgxs.DelayedNuFissionMatrixXS))
check_value('energy_groups', delayed_nu_fission.energy_groups,
[self.energy_groups])
check_value('num_delayed_groups', delayed_nu_fission.num_delayed_groups,

8
openmc/particle_restart.py
View file

@ -1,3 +1,5 @@
import h5py
class Particle(object):
"""Information used to restart a specific particle that caused a simulation to
fail.
@ -33,12 +35,6 @@ class Particle(object):
"""
def __init__(self, filename):
import h5py
if h5py.__version__ == '2.6.0':
raise ImportError("h5py 2.6.0 has a known bug which makes it "
"incompatible with OpenMC's HDF5 files. "
"Please switch to a different version.")
self._f = h5py.File(filename, 'r')
# Ensure filetype and revision are correct

60
openmc/statepoint.py
View file

@ -5,6 +5,7 @@ import warnings
import glob
import numpy as np
import h5py
import openmc
import openmc.checkvalue as cv
@ -96,6 +97,9 @@ class StatePoint(object):
Dictionary whose keys are tally IDs and whose values are Tally objects
tallies_present : bool
Indicate whether user-defined tallies are present
tally_derivatives : dict
Dictionary whose keys are tally derivative IDs and whose values are
TallyDerivative objects
version: tuple of Integral
Version of OpenMC
summary : None or openmc.Summary
@ -104,12 +108,6 @@ class StatePoint(object):
"""
def __init__(self, filename, autolink=True):
import h5py
if h5py.__version__ == '2.6.0':
raise ImportError("h5py 2.6.0 has a known bug which makes it "
"incompatible with OpenMC's HDF5 files. "
"Please switch to a different version.")
self._f = h5py.File(filename, 'r')
# Ensure filetype and revision are correct
@ -134,6 +132,7 @@ class StatePoint(object):
self._summary = None
self._global_tallies = None
self._sparse = False
self._derivs_read = False
# Automatically link in a summary file if one exists
if autolink:
@ -209,13 +208,13 @@ class StatePoint(object):
def global_tallies(self):
if self._global_tallies is None:
data = self._f['global_tallies'].value
gt = np.zeros_like(data, dtype=[
gt = np.zeros(data.shape[0], dtype=[
('name', 'a14'), ('sum', 'f8'), ('sum_sq', 'f8'),
('mean', 'f8'), ('std_dev', 'f8')])
gt['name'] = ['k-collision', 'k-absorption', 'k-tracklength',
'leakage']
gt['sum'] = data['sum']
gt['sum_sq'] = data['sum_sq']
gt['sum'] = data[:,1]
gt['sum_sq'] = data[:,2]
# Calculate mean and sample standard deviation of mean
n = self.n_realizations
@ -396,6 +395,12 @@ class StatePoint(object):
base, tally_key)].value.decode()
tally.num_realizations = n_realizations
# Read derivative information.
if 'derivative' in self._f['{0}{1}'.format(base, tally_key)]:
deriv_id = self._f['{0}{1}/derivative'.format(
base, tally_key)].value
tally.derivative = self.tally_derivatives[deriv_id]
# Read the number of Filters
n_filters = \
self._f['{0}{1}/n_filters'.format(base, tally_key)].value
@ -457,6 +462,43 @@ class StatePoint(object):
def tallies_present(self):
return self._f['tallies/tallies_present'].value
@property
def tally_derivatives(self):
if not self._derivs_read:
# Initialize dictionaries for the Meshes
# Keys - Derivative IDs
# Values - TallyDerivative objects
self._derivs = {}
# Populate the dictionary if any derivatives are present.
if 'derivatives' in self._f['tallies']:
# Read the derivative ids.
base = 'tallies/derivatives'
deriv_ids = [int(k.split(' ')[1]) for k in self._f[base]]
# Create each derivative object and add it to the dictionary.
for d_id in deriv_ids:
base = 'tallies/derivatives/derivative {:d}'.format(d_id)
deriv = openmc.TallyDerivative(derivative_id=d_id)
deriv.variable = \
self._f[base + '/independent variable'].value.decode()
if deriv.variable == 'density':
deriv.material = self._f[base + '/material'].value
elif deriv.variable == 'nuclide_density':
deriv.material = self._f[base + '/material'].value
deriv.nuclide = \
self._f[base + '/nuclide'].value.decode()
elif deriv.variable == 'temperature':
deriv.material = self._f[base + '/material'].value
else:
raise RuntimeError('Unrecognized tally differential '
'variable')
self._derivs[d_id] = deriv
self._derivs_read = True
return self._derivs
@property
def version(self):
return (self._f['version_major'].value,

10
openmc/summary.py
View file

@ -2,6 +2,7 @@ from collections import Iterable
import re
import numpy as np
import h5py
import openmc
from openmc.region import Region
@ -23,15 +24,6 @@ class Summary(object):
"""
def __init__(self, filename):
# A user may not have h5py, but they can still use the rest of the
# Python API so we'll only try to import h5py if the user actually inits
# a Summary object.
import h5py
if h5py.__version__ == '2.6.0':
raise ImportError("h5py 2.6.0 has a known bug which makes it "
"incompatible with OpenMC's HDF5 files. "
"Please switch to a different version.")
openmc.reset_auto_ids()
if not filename.endswith(('.h5', '.hdf5')):

83
openmc/tallies.py
View file

@ -13,6 +13,7 @@ from xml.etree import ElementTree as ET
from six import string_types
import numpy as np
import h5py
import openmc
import openmc.checkvalue as cv
@ -103,6 +104,8 @@ class Tally(object):
sparse : bool
Whether or not the tally uses SciPy's LIL sparse matrix format for
compressed data storage
derivative : openmc.TallyDerivative
A material perturbation derivative to apply to all scores in the tally.
"""
@ -115,6 +118,7 @@ class Tally(object):
self._scores = cv.CheckedList(_SCORE_CLASSES, 'tally scores')
self._estimator = None
self._triggers = cv.CheckedList(openmc.Trigger, 'tally triggers')
self._derivative = None
self._num_realizations = 0
self._with_summary = False
@ -150,6 +154,10 @@ class Tally(object):
if nuclide not in other.nuclides:
return False
# Check derivatives
if self.derivative != other.derivative:
return False
# Check all scores
if len(self.scores) != len(other.scores):
return False
@ -171,27 +179,31 @@ class Tally(object):
def __repr__(self):
string = 'Tally\n'
string += '{0: <16}{1}{2}\n'.format('\tID', '=\t', self.id)
string += '{0: <16}{1}{2}\n'.format('\tName', '=\t', self.name)
string += '{: <16}=\t{}\n'.format('\tID', self.id)
string += '{: <16}=\t{}\n'.format('\tName', self.name)
string += '{0: <16}{1}\n'.format('\tFilters', '=\t')
if self.derivative is not None:
string += '{: <16}=\t{}\n'.format('\tDerivative ID',
str(self.derivative.id))
string += '{: <16}=\n'.format('\tFilters')
for self_filter in self.filters:
string += '{0: <16}\t\t{1}\t{2}\n'.format('',
string += '{: <16}\t\t{}\t{}\n'.format('',
type(self_filter).__name__, self_filter.bins)
string += '{0: <16}{1}'.format('\tNuclides', '=\t')
string += '{: <16}=\t'.format('\tNuclides')
for nuclide in self.nuclides:
if isinstance(nuclide, openmc.Nuclide):
string += '{0} '.format(nuclide.name)
string += nuclide.name + ' '
else:
string += '{0} '.format(nuclide)
string += nuclide + ' '
string += '\n'
string += '{0: <16}{1}{2}\n'.format('\tScores', '=\t', self.scores)
string += '{0: <16}{1}{2}\n'.format('\tEstimator', '=\t', self.estimator)
string += '{: <16}=\t{}\n'.format('\tScores', self.scores)
string += '{: <16}=\t{}\n'.format('\tEstimator', self.estimator)
return string
@ -269,20 +281,14 @@ class Tally(object):
return None
if not self._results_read:
import h5py
if h5py.__version__ == '2.6.0':
raise ImportError("h5py 2.6.0 has a known bug which makes it "
"incompatible with OpenMC's HDF5 files. "
"Please switch to a different version.")
# Open the HDF5 statepoint file
f = h5py.File(self._sp_filename, 'r')
# Extract Tally data from the file
data = f['tallies/tally {0}/results'.format(
self.id)].value
sum = data['sum']
sum_sq = data['sum_sq']
sum = data[:,:,0]
sum_sq = data[:,:,1]
# Reshape the results arrays
sum = np.reshape(sum, self.shape)
@ -380,6 +386,10 @@ class Tally(object):
def derived(self):
return self._derived
@property
def derivative(self):
return self._derivative
@property
def sparse(self):
return self._sparse
@ -434,6 +444,12 @@ class Tally(object):
else:
self._name = ''
@derivative.setter
def derivative(self, deriv):
if deriv is not None:
cv.check_type('tally derivative', deriv, openmc.TallyDerivative)
self._derivative = deriv
@filters.setter
def filters(self, filters):
cv.check_type('tally filters', filters, MutableSequence)
@ -1060,6 +1076,11 @@ class Tally(object):
for trigger in self.triggers:
trigger.get_trigger_xml(element)
# Optional derivatives
if self.derivative is not None:
subelement = ET.SubElement(element, "derivative")
subelement.text = str(self.derivative.id)
return element
def contains_filter(self, filter_type):
@ -1495,7 +1516,8 @@ class Tally(object):
return data
def get_pandas_dataframe(self, filters=True, nuclides=True, scores=True,
distribcell_paths=True, float_format='{:.2e}'):
derivative=True, distribcell_paths=True,
float_format='{:.2e}'):
"""Build a Pandas DataFrame for the Tally data.
This method constructs a Pandas DataFrame object for the Tally data
@ -1513,6 +1535,8 @@ class Tally(object):
Include columns with nuclide bin information (default is True).
scores : bool
Include columns with score bin information (default is True).
derivative : bool
Include columns with differential tally info (default is True).
distribcell_paths : bool, optional
Construct columns for distribcell tally filters (default is True).
The geometric information in the Summary object is embedded into a
@ -1593,6 +1617,14 @@ class Tally(object):
tile_factor = data_size / len(self.scores)
df[column_name] = np.tile(scores, int(tile_factor))
# Include columns for derivatives if user requested it
if derivative and (self.derivative is not None):
df['d_variable'] = self.derivative.variable
if self.derivative.material is not None:
df['d_material'] = self.derivative.material
if self.derivative.nuclide is not None:
df['d_nuclide'] = self.derivative.nuclide
# Append columns with mean, std. dev. for each tally bin
df['mean'] = self.mean.ravel()
df['std. dev.'] = self.std_dev.ravel()
@ -1727,8 +1759,6 @@ class Tally(object):
# HDF5 binary file
if format == 'hdf5':
import h5py
filename = directory + '/' + filename + '.h5'
if append:
@ -3562,6 +3592,18 @@ class Tallies(cv.CheckedList):
self._tallies_file.append(xml_element)
already_written.add(f.mesh)
def _create_derivative_subelements(self):
# Get a list of all derivatives referenced in a tally.
derivs = []
for tally in self:
deriv = tally.derivative
if deriv is not None and deriv not in derivs:
derivs.append(deriv)
# Add the derivatives to the XML tree.
for d in derivs:
self._tallies_file.append(d.to_xml_element())
def export_to_xml(self):
"""Create a tallies.xml file that can be used for a simulation.
@ -3572,6 +3614,7 @@ class Tallies(cv.CheckedList):
self._create_mesh_subelements()
self._create_tally_subelements()
self._create_derivative_subelements()
# Clean the indentation in the file to be user-readable
clean_xml_indentation(self._tallies_file)

141
openmc/tally_derivative.py Normal file
View file

@ -0,0 +1,141 @@
from __future__ import division
import sys
from numbers import Integral
from xml.etree import ElementTree as ET
from six import string_types
import openmc.checkvalue as cv
from openmc.mixin import EqualityMixin
# "Static" variable for auto-generated TallyDerivative IDs
AUTO_TALLY_DERIV_ID = 10000
def reset_auto_tally_deriv_id():
global AUTO_TALLY_ID
AUTO_TALLY_DERIV_ID = 10000
class TallyDerivative(EqualityMixin):
"""A material perturbation derivative to apply to a tally.
Parameters
----------
derivative_id : Integral, optional
Unique identifier for the tally derivative. If none is specified, an
identifier will automatically be assigned
variable : str, optional
Accepted values are 'density', 'nuclide_density', and 'temperature'
material : Integral, optional
The perturubed material ID
nuclide : str, optional
The perturbed nuclide. Only needed for 'nuclide_density' derivatives.
Ex: 'Xe135'
Attributes
----------
id : Integral
Unique identifier for the tally derivative
variable : str
Accepted values are 'density', 'nuclide_density', and 'temperature'
material : Integral
The perturubed material ID
nuclide : str
The perturbed nuclide. Only needed for 'nuclide_density' derivatives.
Ex: 'Xe135'
"""
def __init__(self, derivative_id=None, variable=None, material=None,
nuclide=None):
# Initialize Tally class attributes
self.id = derivative_id
self.variable = variable
self.material = material
self.nuclide = nuclide
def __hash__(self):
return hash(repr(self))
def __repr__(self):
string = 'Tally Derivative\n'
string += '{: <16}=\t{}\n'.format('\tID', self.id)
string += '{: <16}=\t{}\n'.format('\tVariable', self.variable)
if self.variable == 'density':
string += '{: <16}=\t{}\n'.format('\tMaterial', self.material)
elif self.variable == 'nuclide_density':
string += '{: <16}=\t{}\n'.format('\tMaterial', self.material)
string += '{: <16}=\t{}\n'.format('\tNuclide', self.nuclide)
elif self.variable == 'temperature':
string += '{: <16}=\t{}\n'.format('\tMaterial', self.material)
return string
@property
def id(self):
return self._id
@property
def variable(self):
return self._variable
@property
def material(self):
return self._material
@property
def nuclide(self):
return self._nuclide
@id.setter
def id(self, deriv_id):
if deriv_id is None:
global AUTO_TALLY_DERIV_ID
self._id = AUTO_TALLY_DERIV_ID
AUTO_TALLY_DERIV_ID += 1
else:
cv.check_type('tally derivative ID', deriv_id, Integral)
cv.check_greater_than('tally derivative ID', deriv_id, 0,
equality=True)
self._id = deriv_id
@variable.setter
def variable(self, var):
if var is not None:
cv.check_type('derivative variable', var, string_types)
cv.check_value('derivative variable', var,
('density', 'nuclide_density', 'temperature'))
self._variable = var
@material.setter
def material(self, mat):
if mat is not None:
cv.check_type('derivative material', mat, Integral)
self._material = mat
@nuclide.setter
def nuclide(self, nuc):
if nuc is not None:
cv.check_type('derivative nuclide', nuc, string_types)
self._nuclide = nuc
def to_xml_element(self):
"""Return XML representation of the tally derivative
Returns
-------
element : xml.etree.ElementTree.Element
XML element containing derivative data
"""
element = ET.Element("derivative")
element.set("id", str(self.id))
element.set("variable", self.variable)
element.set("material", str(self.material))
if self.variable == 'nuclide_density':
element.set("nuclide", self.nuclide)
return element

3
openmc/volume.py
View file

@ -5,6 +5,7 @@ from warnings import warn
import numpy as np
import pandas as pd
import h5py
import openmc
import openmc.checkvalue as cv
@ -187,8 +188,6 @@ class VolumeCalculation(object):
Results of the stochastic volume calculation
"""
import h5py
with h5py.File(filename, 'r') as f:
domain_type = f.attrs['domain_type'].decode()
samples = f.attrs['samples']

39
src/cmfd_data.F90
View file

@ -164,7 +164,7 @@ contains
* t%stride) + 1
! Get flux
flux = t % results(1,score_index) % sum
flux = t % results(RESULT_SUM,1,score_index)
cmfd % flux(h,i,j,k) = flux
! Detect zero flux, abort if located
@ -175,10 +175,10 @@ contains
end if
! Get total rr and convert to total xs
cmfd % totalxs(h,i,j,k) = t % results(2,score_index) % sum / flux
cmfd % totalxs(h,i,j,k) = t % results(RESULT_SUM,2,score_index) / flux
! Get p1 scatter rr and convert to p1 scatter xs
cmfd % p1scattxs(h,i,j,k) = t % results(3,score_index) % sum / flux
cmfd % p1scattxs(h,i,j,k) = t % results(RESULT_SUM,3,score_index) / flux
! Calculate diffusion coefficient
cmfd % diffcof(h,i,j,k) = ONE/(3.0_8*(cmfd % totalxs(h,i,j,k) - &
@ -211,19 +211,18 @@ contains
* t%stride) + 1
! Get scattering
cmfd % scattxs(h,g,i,j,k) = t % results(1,score_index) % sum /&
cmfd % scattxs(h,g,i,j,k) = t % results(RESULT_SUM,1,score_index) /&
cmfd % flux(h,i,j,k)
! Get nu-fission
cmfd % nfissxs(h,g,i,j,k) = t % results(2,score_index) % sum /&
cmfd % nfissxs(h,g,i,j,k) = t % results(RESULT_SUM,2,score_index) /&
cmfd % flux(h,i,j,k)
! Bank source
cmfd % openmc_src(g,i,j,k) = cmfd % openmc_src(g,i,j,k) + &
t % results(2,score_index) % sum
t % results(RESULT_SUM,2,score_index)
cmfd % keff_bal = cmfd % keff_bal + &
t % results(2,score_index) % sum / &
dble(t % n_realizations)
t % results(RESULT_SUM,2,score_index) / t % n_realizations
end do INGROUP
@ -243,67 +242,67 @@ contains
matching_bins(i_filter_surf) = OUT_LEFT
score_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
cmfd % current(1,h,i,j,k) = t % results(1,score_index) % sum
cmfd % current(1,h,i,j,k) = t % results(RESULT_SUM,1,score_index)
matching_bins(i_filter_surf) = IN_LEFT
score_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
cmfd % current(2,h,i,j,k) = t % results(1,score_index) % sum
cmfd % current(2,h,i,j,k) = t % results(RESULT_SUM,1,score_index)
! Right surface
matching_bins(i_filter_surf) = IN_RIGHT
score_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
cmfd % current(3,h,i,j,k) = t % results(1,score_index) % sum
cmfd % current(3,h,i,j,k) = t % results(RESULT_SUM,1,score_index)
matching_bins(i_filter_surf) = OUT_RIGHT
score_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
cmfd % current(4,h,i,j,k) = t % results(1,score_index) % sum
cmfd % current(4,h,i,j,k) = t % results(RESULT_SUM,1,score_index)
! Back surface
matching_bins(i_filter_surf) = OUT_BACK
score_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
cmfd % current(5,h,i,j,k) = t % results(1,score_index) % sum
cmfd % current(5,h,i,j,k) = t % results(RESULT_SUM,1,score_index)
matching_bins(i_filter_surf) = IN_BACK
score_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
cmfd % current(6,h,i,j,k) = t % results(1,score_index) % sum
cmfd % current(6,h,i,j,k) = t % results(RESULT_SUM,1,score_index)
! Front surface
matching_bins(i_filter_surf) = IN_FRONT
score_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
cmfd % current(7,h,i,j,k) = t % results(1,score_index) % sum
cmfd % current(7,h,i,j,k) = t % results(RESULT_SUM,1,score_index)
matching_bins(i_filter_surf) = OUT_FRONT
score_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
cmfd % current(8,h,i,j,k) = t % results(1,score_index) % sum
cmfd % current(8,h,i,j,k) = t % results(RESULT_SUM,1,score_index)
! Bottom surface
matching_bins(i_filter_surf) = OUT_BOTTOM
score_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
cmfd % current(9,h,i,j,k) = t % results(1,score_index) % sum
cmfd % current(9,h,i,j,k) = t % results(RESULT_SUM,1,score_index)
matching_bins(i_filter_surf) = IN_BOTTOM
score_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
cmfd % current(10,h,i,j,k) = t % results(1,score_index) % sum
cmfd % current(10,h,i,j,k) = t % results(RESULT_SUM,1,score_index)
! Top surface
matching_bins(i_filter_surf) = IN_TOP
score_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
cmfd % current(11,h,i,j,k) = t % results(1,score_index) % sum
cmfd % current(11,h,i,j,k) = t % results(RESULT_SUM,1,score_index)
matching_bins(i_filter_surf) = OUT_TOP
score_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
cmfd % current(12,h,i,j,k) = t % results(1,score_index) % sum
cmfd % current(12,h,i,j,k) = t % results(RESULT_SUM,1,score_index)
end if TALLY

12
src/cmfd_execute.F90
View file

@ -365,22 +365,18 @@ contains
subroutine cmfd_tally_reset()
use global, only: n_cmfd_tallies, cmfd_tallies
use global, only: cmfd_tallies
use output, only: write_message
use tally, only: reset_result
integer :: i ! loop counter
! Print message
call write_message("CMFD tallies reset", 7)
! Begin loop around CMFD tallies
do i = 1, n_cmfd_tallies
! Reset that tally
! Reset CMFD tallies
do i = 1, size(cmfd_tallies)
cmfd_tallies(i) % n_realizations = 0
call reset_result(cmfd_tallies(i) % results)
cmfd_tallies(i) % results(:,:,:) = ZERO
end do
end subroutine cmfd_tally_reset

12
src/constants.F90
View file

@ -269,6 +269,12 @@ module constants
! ============================================================================
! TALLY-RELATED CONSTANTS
! Tally result entries
integer, parameter :: &
RESULT_VALUE = 1, &
RESULT_SUM = 2, &
RESULT_SUM_SQ = 3
! Tally type
integer, parameter :: &
TALLY_VOLUME = 1, &
@ -385,6 +391,12 @@ module constants
K_TRACKLENGTH = 3, &
LEAKAGE = 4
! Differential tally independent variables
integer, parameter :: &
DIFF_DENSITY = 1, &
DIFF_NUCLIDE_DENSITY = 2, &
DIFF_TEMPERATURE = 3
! ============================================================================
! RANDOM NUMBER STREAM CONSTANTS

90
src/cross_section.F90
View file

@ -7,7 +7,7 @@ module cross_section
use global
use list_header, only: ListElemInt
use material_header, only: Material
use math, only: faddeeva, broaden_wmp_polynomials
use math, only: faddeeva, w_derivative, broaden_wmp_polynomials
use multipole_header, only: FORM_RM, FORM_MLBW, MP_EA, RM_RT, RM_RA, RM_RF, &
MLBW_RT, MLBW_RX, MLBW_RA, MLBW_RF, FIT_T, FIT_A,&
FIT_F, MultipoleArray
@ -706,6 +706,94 @@ contains
end if
end subroutine multipole_eval
!===============================================================================
! MULTIPOLE_DERIV_EVAL evaluates the windowed multipole equations for the
! derivative of cross sections in the resolved resonance regions with respect to
! temperature.
!===============================================================================
subroutine multipole_deriv_eval(multipole, E, sqrtkT, sigT, sigA, sigF)
type(MultipoleArray), intent(in) :: multipole ! The windowed multipole
! object to process.
real(8), intent(in) :: E ! The energy at which to
! evaluate the cross section
real(8), intent(in) :: sqrtkT ! The temperature in the form
! sqrt(kT), at which to
! evaluate the XS.
real(8), intent(out) :: sigT ! Total cross section
real(8), intent(out) :: sigA ! Absorption cross section
real(8), intent(out) :: sigF ! Fission cross section
complex(8) :: w_val ! The faddeeva function evaluated at Z
complex(8) :: Z ! sqrt(atomic weight ratio / kT) * (sqrt(E) - pole)
complex(8) :: sigT_factor(multipole % num_l)
real(8) :: sqrtE ! sqrt(E), eV
real(8) :: invE ! 1/E, eV
real(8) :: dopp ! sqrt(atomic weight ratio / kT)
integer :: i_pole ! index of pole
integer :: i_window ! index of window
integer :: startw ! window start pointer (for poles)
integer :: endw ! window end pointer
real(8) :: T
! ==========================================================================
! Bookkeeping
! Define some frequently used variables.
sqrtE = sqrt(E)
invE = ONE / E
dopp = multipole % sqrtAWR / sqrtkT
T = sqrtkT**2 / K_BOLTZMANN
if (sqrtkT == ZERO) call fatal_error("Windowed multipole temperature &
&derivatives are not implemented for 0 Kelvin cross sections.")
! Locate us
i_window = floor((sqrtE - sqrt(multipole % start_E)) / multipole % spacing &
+ ONE)
startw = multipole % w_start(i_window)
endw = multipole % w_end(i_window)
! Fill in factors.
if (startw <= endw) then
call compute_sigT_factor(multipole, sqrtE, sigT_factor)
end if
! Initialize the ouptut cross sections.
sigT = ZERO
sigA = ZERO
sigF = ZERO
! TODO Polynomials: Some of the curvefit polynomials Doppler broaden so
! rigorously we should be computing the derivative of those. But in
! practice, those derivatives are only large at very low energy and they
! have no effect on reactor calculations.
! ==========================================================================
! Add the contribution from the poles in this window.
if (endw >= startw) then
do i_pole = startw, endw
Z = (sqrtE - multipole % data(MP_EA, i_pole)) * dopp
w_val = -invE * SQRT_PI * HALF * w_derivative(Z, 2)
if (multipole % formalism == FORM_MLBW) then
sigT = sigT + real((multipole % data(MLBW_RT, i_pole) * &
sigT_factor(multipole%l_value(i_pole)) + &
multipole % data(MLBW_RX, i_pole)) * w_val)
sigA = sigA + real(multipole % data(MLBW_RA, i_pole) * w_val)
sigF = sigF + real(multipole % data(MLBW_RF, i_pole) * w_val)
else if (multipole % formalism == FORM_RM) then
sigT = sigT + real(multipole % data(RM_RT, i_pole) * w_val * &
sigT_factor(multipole % l_value(i_pole)))
sigA = sigA + real(multipole % data(RM_RA, i_pole) * w_val)
sigF = sigF + real(multipole % data(RM_RF, i_pole) * w_val)
end if
end do
sigT = -HALF*multipole % sqrtAWR / sqrt(K_BOLTZMANN) * T**(-1.5) * sigT
sigA = -HALF*multipole % sqrtAWR / sqrt(K_BOLTZMANN) * T**(-1.5) * sigA
sigF = -HALF*multipole % sqrtAWR / sqrt(K_BOLTZMANN) * T**(-1.5) * sigF
end if
end subroutine multipole_deriv_eval
!===============================================================================
! COMPUTE_SIGT_FACTOR calculates the sigT_factor, a factor inside of the sigT
! equation not present in the sigA and sigF equations.

14
src/eigenvalue.F90
View file

@ -374,7 +374,7 @@ contains
subroutine calculate_generation_keff()
! Get keff for this generation by subtracting off the starting value
keff_generation = global_tallies(K_TRACKLENGTH) % value - keff_generation
keff_generation = global_tallies(RESULT_VALUE, K_TRACKLENGTH) - keff_generation
#ifdef MPI
! Combine values across all processors
@ -466,14 +466,14 @@ contains
k_combined = ZERO
! Copy estimates of k-effective and its variance (not variance of the mean)
kv(1) = global_tallies(K_COLLISION) % sum / n
kv(2) = global_tallies(K_ABSORPTION) % sum / n
kv(3) = global_tallies(K_TRACKLENGTH) % sum / n
cov(1,1) = (global_tallies(K_COLLISION) % sum_sq - &
kv(1) = global_tallies(RESULT_SUM, K_COLLISION) / n
kv(2) = global_tallies(RESULT_SUM, K_ABSORPTION) / n
kv(3) = global_tallies(RESULT_SUM, K_TRACKLENGTH) / n
cov(1,1) = (global_tallies(RESULT_SUM_SQ, K_COLLISION) - &
n * kv(1) * kv(1)) / (n - 1)
cov(2,2) = (global_tallies(K_ABSORPTION) % sum_sq - &
cov(2,2) = (global_tallies(RESULT_SUM_SQ, K_ABSORPTION) - &
n * kv(2) * kv(2)) / (n - 1)
cov(3,3) = (global_tallies(K_TRACKLENGTH) % sum_sq - &
cov(3,3) = (global_tallies(RESULT_SUM_SQ, K_TRACKLENGTH) - &
n * kv(3) * kv(3)) / (n - 1)
! Calculate covariances based on sums with Bessel's correction

4
src/finalize.F90
View file

@ -9,7 +9,7 @@ module finalize
use message_passing
#endif
use hdf5_interface, only: hdf5_bank_t, hdf5_tallyresult_t
use hdf5_interface, only: hdf5_bank_t
use hdf5, only: h5tclose_f, h5close_f
implicit none
@ -53,7 +53,6 @@ contains
call free_memory()
! Release compound datatypes
call h5tclose_f(hdf5_tallyresult_t, hdf5_err)
call h5tclose_f(hdf5_bank_t, hdf5_err)
! Close FORTRAN interface.
@ -62,7 +61,6 @@ contains
#ifdef MPI
! Free all MPI types
call MPI_TYPE_FREE(MPI_BANK, mpi_err)
call MPI_TYPE_FREE(MPI_TALLYRESULT, mpi_err)
! If MPI is in use and enabled, terminate it
call MPI_FINALIZE(mpi_err)

20
src/global.F90
View file

@ -1,5 +1,11 @@
module global
use, intrinsic :: ISO_C_BINDING
#ifdef MPIF08
use mpi_f08
#endif
use bank_header, only: Bank
use cmfd_header
use constants
@ -14,15 +20,11 @@ module global
use set_header, only: SetInt
use surface_header, only: SurfaceContainer
use source_header, only: SourceDistribution
use tally_header, only: TallyObject, TallyResult
use tally_header, only: TallyObject, TallyDerivative
use trigger_header, only: KTrigger
use timer_header, only: Timer
use volume_header, only: VolumeCalculation
#ifdef MPIF08
use mpi_f08
#endif
implicit none
! ============================================================================
@ -164,7 +166,7 @@ module global
! 3) track-length estimate of k-eff
! 4) leakage fraction
type(TallyResult), allocatable, target :: global_tallies(:)
real(C_DOUBLE), allocatable, target :: global_tallies(:,:)
! It is possible to protect accumulate operations on global tallies by using
! an atomic update. However, when multiple threads accumulate to the same
@ -183,6 +185,10 @@ module global
integer :: n_tallies = 0 ! # of tallies
integer :: n_user_tallies = 0 ! # of user tallies
! Tally derivatives
type(TallyDerivative), allocatable :: tally_derivs(:)
!$omp threadprivate(tally_derivs)
! Normalization for statistics
integer :: n_realizations = 0 ! # of independent realizations
real(8) :: total_weight ! total starting particle weight in realization
@ -272,10 +278,8 @@ module global
integer :: mpi_err ! MPI error code
#ifdef MPIF08
type(MPI_Datatype) :: MPI_BANK
type(MPI_Datatype) :: MPI_TALLYRESULT
#else
integer :: MPI_BANK ! MPI datatype for fission bank
integer :: MPI_TALLYRESULT ! MPI datatype for TallyResult
#endif
#ifdef _OPENMP

130
src/hdf5_interface.F90
View file

@ -16,7 +16,6 @@ module hdf5_interface
use h5lt
use error, only: fatal_error
use tally_header, only: TallyResult
#ifdef PHDF5
use message_passing, only: MPI_COMM_WORLD, MPI_INFO_NULL
#endif
@ -24,7 +23,6 @@ module hdf5_interface
implicit none
private
integer(HID_T), public :: hdf5_tallyresult_t ! Compound type for TallyResult
integer(HID_T), public :: hdf5_bank_t ! Compound type for Bank
integer(HID_T), public :: hdf5_integer8_t ! type for integer(8)
@ -42,8 +40,6 @@ module hdf5_interface
module procedure write_long
module procedure write_string
module procedure write_string_1D
module procedure write_tally_result_1D
module procedure write_tally_result_2D
end interface write_dataset
interface read_dataset
@ -60,8 +56,6 @@ module hdf5_interface
module procedure read_long
module procedure read_string
module procedure read_string_1D
module procedure read_tally_result_1D
module procedure read_tally_result_2D
module procedure read_complex_2D
end interface read_dataset
@ -2063,130 +2057,6 @@ contains
call h5ltset_attribute_string_f(group_id, var, attr_type, attr_str, hdf5_err)
end subroutine write_attribute_string
!===============================================================================
! WRITE_TALLY_RESULT writes an OpenMC TallyResult type
!===============================================================================
subroutine write_tally_result_1D(group_id, name, buffer)
integer(HID_T), intent(in) :: group_id
character(*), intent(in) :: name ! name of data
type(TallyResult), intent(in), target :: buffer(:) ! data to write
integer(HSIZE_T) :: dims(1)
dims(:) = shape(buffer)
call write_tally_result_1D_explicit(group_id, dims, name, buffer)
end subroutine write_tally_result_1D
subroutine write_tally_result_1D_explicit(group_id, dims, name, buffer)
integer(HID_T), intent(in) :: group_id
integer(HSIZE_T), intent(in) :: dims(1)
character(*), intent(in) :: name ! name of data
type(TallyResult), intent(in), target :: buffer(dims(1))
integer :: hdf5_err
integer(HID_T) :: dset ! data set handle
integer(HID_T) :: dspace ! data or file space handle
type(c_ptr) :: f_ptr
call h5screate_simple_f(1, dims, dspace, hdf5_err)
call h5dcreate_f(group_id, trim(name), hdf5_tallyresult_t, &
dspace, dset, hdf5_err)
f_ptr = c_loc(buffer)
call h5dwrite_f(dset, hdf5_tallyresult_t, f_ptr, hdf5_err)
call h5dclose_f(dset, hdf5_err)
call h5sclose_f(dspace, hdf5_err)
end subroutine write_tally_result_1D_explicit
subroutine write_tally_result_2D(group_id, name, buffer)
integer(HID_T), intent(in) :: group_id
character(*), intent(in) :: name ! name of data
type(TallyResult), intent(in), target :: buffer(:,:) ! data to write
integer(HSIZE_T) :: dims(2)
dims(:) = shape(buffer)
call write_tally_result_2D_explicit(group_id, dims, name, buffer)
end subroutine write_tally_result_2D
subroutine write_tally_result_2D_explicit(group_id, dims, name, buffer)
integer(HID_T), intent(in) :: group_id
integer(HSIZE_T), intent(in) :: dims(2)
character(*), intent(in) :: name ! name of data
type(TallyResult), intent(in), target :: buffer(dims(1),dims(2))
integer :: hdf5_err
integer(HID_T) :: dset ! data set handle
integer(HID_T) :: dspace ! data or file space handle
type(c_ptr) :: f_ptr
call h5screate_simple_f(2, dims, dspace, hdf5_err)
call h5dcreate_f(group_id, trim(name), hdf5_tallyresult_t, &
dspace, dset, hdf5_err)
f_ptr = c_loc(buffer)
call h5dwrite_f(dset, hdf5_tallyresult_t, f_ptr, hdf5_err)
call h5dclose_f(dset, hdf5_err)
call h5sclose_f(dspace, hdf5_err)
end subroutine write_tally_result_2D_explicit
!===============================================================================
! READ_TALLY_RESULT reads OpenMC TallyResult data
!===============================================================================
subroutine read_tally_result_1D(group_id, name, buffer)
integer(HID_T), intent(in) :: group_id
character(*), intent(in) :: name ! name of data
type(TallyResult), intent(inout), target :: buffer(:) ! read data here
integer(HSIZE_T) :: dims(1)
dims(:) = shape(buffer)
call read_tally_result_1D_explicit(group_id, dims, name, buffer)
end subroutine read_tally_result_1D
subroutine read_tally_result_1D_explicit(group_id, dims, name, buffer)
integer(HID_T), intent(in) :: group_id
integer(HSIZE_T), intent(in) :: dims(1)
character(*), intent(in) :: name ! name of data
type(TallyResult), intent(inout), target :: buffer(dims(1))
integer :: hdf5_err
integer(HID_T) :: dset ! data set handle
type(c_ptr) :: f_ptr
call h5dopen_f(group_id, trim(name), dset, hdf5_err)
f_ptr = c_loc(buffer)
call h5dread_f(dset, hdf5_tallyresult_t, f_ptr, hdf5_err)
call h5dclose_f(dset, hdf5_err)
end subroutine read_tally_result_1D_explicit
subroutine read_tally_result_2D(group_id, name, buffer)
integer(HID_T), intent(in) :: group_id
character(*), intent(in) :: name ! name of data
type(TallyResult), intent(inout), target :: buffer(:,:)
integer(HSIZE_T) :: dims(2)
dims(:) = shape(buffer)
call read_tally_result_2D_explicit(group_id, dims, name, buffer)
end subroutine read_tally_result_2D
subroutine read_tally_result_2D_explicit(group_id, dims, name, buffer)
integer(HID_T), intent(in) :: group_id
integer(HSIZE_T), intent(in) :: dims(2)
character(*), intent(in) :: name ! name of data
type(TallyResult), intent(inout), target :: buffer(dims(1),dims(2))
integer :: hdf5_err
integer(HID_T) :: dset ! data set handle
type(c_ptr) :: f_ptr
call h5dopen_f(group_id, trim(name), dset, hdf5_err)
f_ptr = c_loc(buffer)
call h5dread_f(dset, hdf5_tallyresult_t, f_ptr, hdf5_err)
call h5dclose_f(dset, hdf5_err)
end subroutine read_tally_result_2D_explicit
subroutine read_attribute_double(buffer, obj_id, name)
real(8), intent(inout), target :: buffer
integer(HID_T), intent(in) :: obj_id

51
src/initialize.F90
View file

@ -12,7 +12,7 @@ module initialize
&BASE_UNIVERSE
use global
use hdf5_interface, only: file_open, read_dataset, file_close, hdf5_bank_t,&
hdf5_tallyresult_t, hdf5_integer8_t
hdf5_integer8_t
use input_xml, only: read_input_xml, cells_in_univ_dict, read_plots_xml
use material_header, only: Material
use mgxs_data, only: read_mgxs, create_macro_xs
@ -22,7 +22,7 @@ module initialize
use state_point, only: load_state_point
use string, only: to_str, starts_with, ends_with, str_to_int
use summary, only: write_summary
use tally_header, only: TallyObject, TallyResult
use tally_header, only: TallyObject
use tally_initialize,only: configure_tallies
use tally_filter
use tally, only: init_tally_routines
@ -169,21 +169,11 @@ contains
integer :: bank_blocks(5) ! Count for each datatype
#ifdef MPIF08
type(MPI_Datatype) :: bank_types(5)
type(MPI_Datatype) :: result_types(1)
type(MPI_Datatype) :: temp_type
#else
integer :: bank_types(5) ! Datatypes
integer :: result_types(1) ! Datatypes
integer :: temp_type ! temporary derived type
#endif
integer(MPI_ADDRESS_KIND) :: bank_disp(5) ! Displacements
integer :: result_blocks(1) ! Count for each datatype
integer(MPI_ADDRESS_KIND) :: result_disp(1) ! Displacements
integer(MPI_ADDRESS_KIND) :: result_base_disp ! Base displacement
integer(MPI_ADDRESS_KIND) :: lower_bound ! Lower bound for TallyResult
integer(MPI_ADDRESS_KIND) :: extent ! Extent for TallyResult
type(Bank) :: b
type(TallyResult) :: tr
! Indicate that MPI is turned on
mpi_enabled = .true.
@ -222,34 +212,6 @@ contains
bank_types, MPI_BANK, mpi_err)
call MPI_TYPE_COMMIT(MPI_BANK, mpi_err)
! ==========================================================================
! CREATE MPI_TALLYRESULT TYPE
! Determine displacements for MPI_BANK type
call MPI_GET_ADDRESS(tr%value, result_base_disp, mpi_err)
call MPI_GET_ADDRESS(tr%sum, result_disp(1), mpi_err)
! Adjust displacements
result_disp = result_disp - result_base_disp
! Define temporary type for TallyResult
result_blocks = (/ 2 /)
result_types = (/ MPI_REAL8 /)
call MPI_TYPE_CREATE_STRUCT(1, result_blocks, result_disp, result_types, &
temp_type, mpi_err)
! Adjust lower-bound and extent of type for tally score
lower_bound = 0
extent = result_disp(1) + 16
call MPI_TYPE_CREATE_RESIZED(temp_type, lower_bound, extent, &
MPI_TALLYRESULT, mpi_err)
! Commit derived type for tally scores
call MPI_TYPE_COMMIT(MPI_TALLYRESULT, mpi_err)
! Free temporary MPI type
call MPI_TYPE_FREE(temp_type, mpi_err)
end subroutine initialize_mpi
#endif
@ -259,7 +221,6 @@ contains
subroutine hdf5_initialize()
type(TallyResult), target :: tmp(2) ! temporary TallyResult
type(Bank), target :: tmpb(2) ! temporary Bank
integer :: hdf5_err
integer(HID_T) :: coordinates_t ! HDF5 type for 3 reals
@ -268,14 +229,6 @@ contains
! Initialize FORTRAN interface.
call h5open_f(hdf5_err)
! Create the compound datatype for TallyResult
call h5tcreate_f(H5T_COMPOUND_F, h5offsetof(c_loc(tmp(1)), &
c_loc(tmp(2))), hdf5_tallyresult_t, hdf5_err)
call h5tinsert_f(hdf5_tallyresult_t, "sum", h5offsetof(c_loc(tmp(1)), &
c_loc(tmp(1)%sum)), H5T_NATIVE_DOUBLE, hdf5_err)
call h5tinsert_f(hdf5_tallyresult_t, "sum_sq", h5offsetof(c_loc(tmp(1)), &
c_loc(tmp(1)%sum_sq)), H5T_NATIVE_DOUBLE, hdf5_err)
! Create compound type for xyz and uvw
call h5tarray_create_f(H5T_NATIVE_DOUBLE, 1, dims, coordinates_t, hdf5_err)

160
src/input_xml.F90
View file

@ -2039,7 +2039,6 @@ contains
type(Library), allocatable :: libraries(:)
type(VectorReal), allocatable :: nuc_temps(:) ! List of T to read for each nuclide
type(VectorReal), allocatable :: sab_temps(:) ! List of T to read for each S(a,b)
character(MAX_LINE_LEN) :: temp_str
real(8), allocatable :: material_temps(:)
logical :: file_exists
character(MAX_FILE_LEN) :: env_variable
@ -2629,11 +2628,13 @@ contains
type(Node), pointer :: node_mesh => null()
type(Node), pointer :: node_tal => null()
type(Node), pointer :: node_filt => null()
type(Node), pointer :: node_trigger=>null()
type(Node), pointer :: node_trigger => null()
type(Node), pointer :: node_deriv => null()
type(NodeList), pointer :: node_mesh_list => null()
type(NodeList), pointer :: node_tal_list => null()
type(NodeList), pointer :: node_filt_list => null()
type(NodeList), pointer :: node_trigger_list => null()
type(NodeList), pointer :: node_deriv_list => null()
type(ElemKeyValueCI), pointer :: scores
type(ElemKeyValueCI), pointer :: next
@ -2641,6 +2642,12 @@ contains
filename = trim(path_input) // "tallies.xml"
inquire(FILE=filename, EXIST=file_exists)
if (.not. file_exists) then
! We need to allocate tally_derivs to avoid segfaults. Also needs to be
! done in parallel because tally derivs are threadprivate.
!$omp parallel
allocate(tally_derivs(0))
!$omp end parallel
! Since a tallies.xml file is optional, no error is issued here
return
end if
@ -2819,6 +2826,94 @@ contains
! We only need the mesh info for plotting
if (run_mode == MODE_PLOTTING) return
! ==========================================================================
! READ DATA FOR DERIVATIVES
! Get pointer list to XML <derivative> nodes and allocate global array.
! The array is threadprivate so it must be allocated in parallel.
call get_node_list(doc, "derivative", node_deriv_list)
!$omp parallel
allocate(tally_derivs(get_list_size(node_deriv_list)))
!$omp end parallel
! Make sure this is not an MG run.
if (.not. run_CE .and. get_list_size(node_deriv_list) > 0) then
call fatal_error("Differential tallies not supported in multi-group mode")
end if
! Read derivative attributes.
do i = 1, get_list_size(node_deriv_list)
associate(deriv => tally_derivs(i))
! Get pointer to derivative node.
call get_list_item(node_deriv_list, i, node_deriv)
! Copy the derivative id.
if (check_for_node(node_deriv, "id")) then
call get_node_value(node_deriv, "id", deriv % id)
else
call fatal_error("Must specify an ID for <derivative> elements in the&
& tally XML file")
end if
! Make sure the id is > 0.
if (deriv % id <= 0) then
call fatal_error("<derivative> IDs must be an integer greater than &
&zero")
end if
! Make sure this id has not already been used.
do j = 1, i-1
if (tally_derivs(j) % id == deriv % id) then
call fatal_error("Two or more <derivative>'s use the same unique &
&ID: " // trim(to_str(deriv % id)))
end if
end do
! Read the independent variable name.
temp_str = ""
call get_node_value(node_deriv, "variable", temp_str)
temp_str = to_lower(temp_str)
select case(temp_str)
case("density")
deriv % variable = DIFF_DENSITY
call get_node_value(node_deriv, "material", deriv % diff_material)
case("nuclide_density")
deriv % variable = DIFF_NUCLIDE_DENSITY
call get_node_value(node_deriv, "material", deriv % diff_material)
call get_node_value(node_deriv, "nuclide", word)
word = trim(to_lower(word))
pair_list => nuclide_dict % keys()
do while (associated(pair_list))
if (starts_with(pair_list % key, word)) then
word = pair_list % key(1:150)
exit
end if
! Advance to next
pair_list => pair_list % next
end do
! Check if no nuclide was found
if (.not. associated(pair_list)) then
call fatal_error("Could not find the nuclide " &
// trim(word) // " specified in derivative " &
// trim(to_str(deriv % id)) // " in any material.")
end if
deallocate(pair_list)
deriv % diff_nuclide = nuclide_dict % get_key(word)
case("temperature")
deriv % variable = DIFF_TEMPERATURE
call get_node_value(node_deriv, "material", deriv % diff_material)
end select
end associate
end do
! ==========================================================================
! READ TALLY DATA
@ -2840,7 +2935,7 @@ contains
t % estimator = ESTIMATOR_TRACKLENGTH
! Copy material id
! Copy tally id
if (check_for_node(node_tal, "id")) then
call get_node_value(node_tal, "id", t % id)
else
@ -3430,6 +3525,7 @@ contains
call fatal_error("Cannot tally flux with an outgoing energy &
&filter.")
end if
case ('flux-yn')
! Prohibit user from tallying flux for an individual nuclide
if (.not. (t % n_nuclide_bins == 1 .and. &
@ -3767,6 +3863,48 @@ contains
// trim(to_str(t % id)) // ".")
end if
! Check for a tally derivative.
if (check_for_node(node_tal, "derivative")) then
! Temporarily store the derivative id.
call get_node_value(node_tal, "derivative", t % deriv)
! Find the derivative with the given id, and store it's index.
do j = 1, size(tally_derivs)
if (tally_derivs(j) % id == t % deriv) then
t % deriv = j
! Only analog or collision estimators are supported for differential
! tallies.
if (t % estimator == ESTIMATOR_TRACKLENGTH) then
t % estimator = ESTIMATOR_COLLISION
end if
! We found the derivative we were looking for; exit the do loop.
exit
end if
if (j == size(tally_derivs)) then
call fatal_error("Could not find derivative " &
// trim(to_str(t % deriv)) // " specified on tally " &
// trim(to_str(t % id)))
end if
end do
if (tally_derivs(t % deriv) % variable == DIFF_NUCLIDE_DENSITY &
.or. tally_derivs(t % deriv) % variable == DIFF_TEMPERATURE) then
if (any(t % nuclide_bins == -1)) then
if (t % find_filter(FILTER_ENERGYOUT) > 0) then
call fatal_error("Error on tally " // trim(to_str(t % id)) &
// ": Cannot use a 'nuclide_density' or 'temperature' &
&derivative on a tally with an outgoing energy filter and &
&'total' nuclide rate. Instead, tally each nuclide in the &
&material individually.")
! Note that diff tallies with these characteristics would work
! correctly if no tally events occur in the perturbed material
! (e.g. pertrubing moderator but only tallying fuel), but this
! case would be hard to check for by only reading inputs.
end if
end if
end if
end if
! If settings.xml trigger is turned on, create tally triggers
if (trigger_on) then
@ -4788,8 +4926,8 @@ contains
sum_percent = sum(mat % atom_density)
mat % atom_density = mat % atom_density / sum_percent
! Change density in g/cm^3 to atom/b-cm. Since all values are now in atom
! percent, the sum needs to be re-evaluated as 1/sum(x*awr)
! Change density in g/cm^3 to atom/b-cm. Since all values are now in
! atom percent, the sum needs to be re-evaluated as 1/sum(x*awr)
if (.not. density_in_atom) then
sum_percent = ZERO
do j = 1, mat % n_nuclides
@ -4808,6 +4946,18 @@ contains
! Calculate nuclide atom densities
mat % atom_density = mat % density * mat % atom_density
! Calculate density in g/cm^3.
mat % density_gpcc = ZERO
do j = 1, mat % n_nuclides
if (run_CE) then
awr = nuclides(mat % nuclide(j)) % awr
else
awr = ONE
end if
mat % density_gpcc = mat % density_gpcc &
+ mat % atom_density(j) * awr * MASS_NEUTRON / N_AVOGADRO
end do
end associate
end do

1
src/material_header.F90
View file

@ -13,6 +13,7 @@ module material_header
integer, allocatable :: nuclide(:) ! index in nuclides array
real(8) :: density ! total atom density in atom/b-cm
real(8), allocatable :: atom_density(:) ! nuclide atom density in atom/b-cm
real(8) :: density_gpcc ! total density in g/cm^3
! Energy grid information
integer :: n_grid ! # of union material grid points

19
src/math.F90
View file

@ -1,8 +1,9 @@
module math
use, intrinsic :: ISO_C_BINDING
use constants
use random_lcg, only: prn
use ISO_C_BINDING
implicit none
@ -752,6 +753,22 @@ contains
end function faddeeva
recursive function w_derivative(z, order) result(wv)
complex(C_DOUBLE_COMPLEX), intent(in) :: z ! The point to evaluate Z at
integer, intent(in) :: order
complex(8) :: wv ! The resulting w(z) value
select case(order)
case (0)
wv = faddeeva(z)
case (1)
wv = -TWO * z * faddeeva(z) + TWO * ONEI / SQRT_PI
case default
wv = -TWO * z * w_derivative(z, order-1) &
- TWO * (order-1) * w_derivative(z, order-2)
end select
end function w_derivative
!===============================================================================
! BROADEN_WMP_POLYNOMIALS Doppler broadens the windowed multipole curvefit. The
! curvefit is a polynomial of the form

231
src/mgxs_header.F90
View file

@ -434,7 +434,7 @@ module mgxs_header
integer :: ndims
integer(HSIZE_T) :: dims(2)
real(8), allocatable :: temp_arr(:), temp_2d(:, :)
real(8), allocatable :: temp_beta(:, :)
real(8), allocatable :: temp_beta(:, :), temp_3d(:, :, :)
real(8) :: dmu, mu, norm, chi_sum
integer :: order, order_dim, gin, gout, l, imu, length
type(VectorInt) :: temps_to_read
@ -767,9 +767,57 @@ module mgxs_header
! If prompt-nu-fission present, set prompt-nu-fission
if (object_exists(xsdata_grp, "prompt-nu-fission")) then
! Set prompt-nu-fission
call read_dataset(xs % prompt_nu_fission, xsdata_grp, &
"prompt-nu-fission")
! Get the dimensions of the prompt-nu-fission dataset
xsdata = open_dataset(xsdata_grp, "prompt-nu-fission")
call get_ndims(xsdata, ndims)
! If prompt-nu-fission is a vector
if (ndims == 1) then
! Set prompt_nu_fission
call read_dataset(xs % prompt_nu_fission, xsdata_grp, &
"prompt-nu-fission")
! If prompt-nu-fission is a matrix, set prompt_nu_fission and
! chi_prompt.
else if (ndims == 2) then
! chi_prompt is embedded in prompt_nu_fission -> extract
! chi_prompt
allocate(temp_arr(energy_groups * energy_groups))
call read_dataset(temp_arr, xsdata_grp, "prompt-nu-fission")
allocate(temp_2d(energy_groups, energy_groups))
temp_2d = reshape(temp_arr, (/energy_groups, energy_groups/))
! Deallocate temporary 1D array for prompt_nu_fission matrix
deallocate(temp_arr)
! Set the vector prompt-nu-fission from the matrix
! prompt-nu-fission
do gin = 1, energy_groups
xs % prompt_nu_fission(gin) = sum(temp_2d(:, gin))
end do
! Now pull out information needed for chi
xs % chi_prompt(:, :) = temp_2d
! Deallocate temporary 2D array for nu_fission matrix
deallocate(temp_2d)
! Normalize chi so its CDF goes to 1
do gin = 1, energy_groups
chi_sum = sum(xs % chi_prompt(:, gin))
if (chi_sum == ZERO) then
call fatal_error("Encountered chi prompt for a group &
&that sums to zero")
else
xs % chi_prompt(:, gin) = xs % chi_prompt(:, gin) / chi_sum
end if
end do
else
call fatal_error("prompt-nu-fission must be provided as a 1D &
&or 2D array")
end if
end if
! If delayed-nu-fission provided, set delayed-nu-fission. If
@ -848,9 +896,52 @@ module mgxs_header
! Deallocate temporary array for delayed-nu-fission matrix
deallocate(temp_arr)
! If delayed nu-fission is a 3D matrix, set delayed_nu_fission
! and chi_delayed.
else if (ndims == 3) then
! chi_delayed is embedded in delayed_nu_fission -> extract
! chi_delayed
allocate(temp_arr(delayed_groups * energy_groups * &
energy_groups))
call read_dataset(temp_arr, xsdata_grp, "delayed-nu-fission")
allocate(temp_3d(delayed_groups, energy_groups, energy_groups))
temp_3d = reshape(temp_arr, (/delayed_groups, energy_groups, &
energy_groups/))
! Deallocate temporary 1D array for delayed_nu_fission matrix
deallocate(temp_arr)
! Set the 2D delayed-nu-fission matrix and 3D chi_dealyed matrix
! from the 3D delayed-nu-fission matrix
do dg = 1, delayed_groups
do gin = 1, energy_groups
xs % delayed_nu_fission(dg, gin) = sum(temp_3d(dg, :, gin))
do gout = 1, energy_groups
xs % chi_delayed(dg, gout, gin) = temp_3d(dg, gout, gin)
end do
end do
end do
! Normalize chi_delayed so its CDF goes to 1
do dg = 1, delayed_groups
do gin = 1, energy_groups
chi_sum = sum(xs % chi_delayed(dg, :, gin))
if (chi_sum == ZERO) then
call fatal_error("Encountered chi delayed for a group &
&that sums to zero")
else
xs % chi_delayed(dg, :, gin) = &
xs % chi_delayed(dg, :, gin) / chi_sum
end if
end do
end do
! Deallocate temporary 3D matrix for delayed_nu_fission
deallocate(temp_3d)
else
call fatal_error("delayed-nu-fission must be provided as a &
&1D or 2D array")
&1D, 2D, or 3D array")
end if
end if
@ -1119,7 +1210,8 @@ module mgxs_header
integer(HSIZE_T) :: dims(4)
integer, allocatable :: int_arr(:)
real(8), allocatable :: temp_1d(:), temp_3d(:, :, :)
real(8), allocatable :: temp_4d(:, :, :, :), temp_beta(:, :, :, :)
real(8), allocatable :: temp_4d(:, :, :, :), temp_5d(:, :, :, :, :)
real(8), allocatable :: temp_beta(:, :, :, :)
real(8) :: dmu, mu, norm, chi_sum
integer :: order, order_dim, gin, gout, l, imu, dg
type(VectorInt) :: temps_to_read
@ -1540,16 +1632,70 @@ module mgxs_header
! If prompt-nu-fission present, set prompt-nu-fission
if (object_exists(xsdata_grp, "prompt-nu-fission")) then
! Allocate temporary array for prompt-nu-fission
allocate(temp_1d(energy_groups * this % n_azi * this % n_pol))
! Get the dimensions of the prompt-nu-fission dataset
xsdata = open_dataset(xsdata_grp, "prompt-nu-fission")
call get_ndims(xsdata, ndims)
! Read prompt-nu-fission
call read_dataset(temp_1d, xsdata_grp, "prompt-nu-fission")
xs % prompt_nu_fission = reshape(temp_1d, (/energy_groups, &
this % n_azi, this % n_pol/))
! If prompt-nu-fission is a vector for each azi and pol
if (ndims == 3) then
! Deallocate temporary array for prompt-nu-fission
deallocate(temp_1d)
! Set prompt_nu_fission
call read_dataset(xs % prompt_nu_fission, xsdata_grp, &
"prompt-nu-fission")
! If prompt-nu-fission is a matrix for each azi and pol,
! set prompt_nu_fission and chi_prompt.
else if (ndims == 4) then
! chi_prompt is embedded in prompt_nu_fission -> extract
! chi_prompt
allocate(temp_1d(energy_groups * energy_groups &
* this % n_azi * this % n_pol))
allocate(temp_4d(energy_groups, energy_groups, this % n_azi, &
this % n_pol))
call read_dataset(temp_1d, xsdata_grp, "prompt-nu-fission")
temp_4d = reshape(temp_1d, (/energy_groups, energy_groups, &
this % n_azi, this % n_pol/))
! Deallocate temporary 1D array for prompt_nu_fission matrix
deallocate(temp_1d)
! Set the vector prompt-nu-fission from the matrix
! prompt-nu-fission
do ipol = 1, this % n_pol
do iazi = 1, this % n_azi
do gin = 1, energy_groups
xs % prompt_nu_fission(gin, iazi, ipol) = &
sum(temp_4d(:, gin, iazi, ipol))
end do
end do
end do
! Now pull out information needed for chi
xs % chi_prompt(:, :, :, :) = temp_4d
! Deallocate temporary 4D array for nu_fission matrix
deallocate(temp_4d)
! Normalize chi so its CDF goes to 1
do ipol = 1, this % n_pol
do iazi = 1, this % n_azi
do gin = 1, energy_groups
chi_sum = sum(xs % chi_prompt(:, gin, iazi, ipol))
if (chi_sum == ZERO) then
call fatal_error("Encountered chi prompt for a group &
&that sums to zero")
else
xs % chi_prompt(:, gin, iazi, ipol) = &
xs % chi_prompt(:, gin, iazi, ipol) / chi_sum
end if
end do
end do
end do
else
call fatal_error("prompt-nu-fission must be provided as a 3D &
&or 4D array")
end if
end if
! If delayed-nu-fission provided, set delayed-nu-fission. If
@ -1639,9 +1785,64 @@ module mgxs_header
! Deallocate temporary array for delayed-nu-fission matrix
deallocate(temp_1d)
! If delayed nu-fission is a 5D matrix, set delayed_nu_fission
! and chi_delayed.
else if (ndims == 5) then
! chi_delayed is embedded in delayed_nu_fission -> extract
! chi_delayed
allocate(temp_1d(delayed_groups * energy_groups * &
energy_groups * this % n_azi * this % n_pol))
allocate(temp_5d(delayed_groups, energy_groups, energy_groups, &
this % n_azi, this % n_pol))
call read_dataset(temp_1d, xsdata_grp, "delayed-nu-fission")
temp_5d = reshape(temp_1d, (/delayed_groups, energy_groups, &
energy_groups, this % n_azi, this % n_pol/))
! Deallocate temporary 1D array for delayed_nu_fission matrix
deallocate(temp_1d)
! Set the 4D delayed-nu-fission matrix and 5D chi_delayed matrix
! from the 5D delayed-nu-fission matrix
do ipol = 1, this % n_pol
do iazi = 1, this % n_azi
do dg = 1, delayed_groups
do gin = 1, energy_groups
xs % delayed_nu_fission(dg, gin, iazi, ipol) = &
sum(temp_5d(dg, :, gin, iazi, ipol))
do gout = 1, energy_groups
xs % chi_delayed(dg, gout, gin, iazi, ipol) = &
temp_5d(dg, gout, gin, iazi, ipol)
end do
end do
end do
end do
end do
! Normalize chi_delayed so its CDF goes to 1
do ipol = 1, this % n_pol
do iazi = 1, this % n_azi
do dg = 1, delayed_groups
do gin = 1, energy_groups
chi_sum = sum(xs % chi_delayed(dg, :, gin, iazi, ipol))
if (chi_sum == ZERO) then
call fatal_error("Encountered chi delayed for a group&
& that sums to zero")
else
xs % chi_delayed(dg, :, gin, iazi, ipol) = &
xs % chi_delayed(dg, :, gin, iazi, ipol) / &
chi_sum
end if
end do
end do
end do
end do
! Deallocate temporary 5D matrix for delayed_nu_fission
deallocate(temp_5d)
else
call fatal_error("delayed-nu-fission must be provided as a &
&1D or 2D array")
&3D, 4D, or 5D array")
end if
end if

117
src/output.F90
View file

@ -611,7 +611,7 @@ contains
t_value = t_percentile(ONE - alpha/TWO, n_realizations - 1)
! Adjust sum_sq
global_tallies(:) % sum_sq = t_value * global_tallies(:) % sum_sq
global_tallies(RESULT_SUM_SQ,:) = t_value * global_tallies(RESULT_SUM_SQ,:)
! Adjust combined estimator
if (n_realizations > 3) then
@ -623,26 +623,26 @@ contains
! write global tallies
if (n_realizations > 1) then
if (run_mode == MODE_EIGENVALUE) then
write(ou,102) "k-effective (Collision)", global_tallies(K_COLLISION) &
% sum, global_tallies(K_COLLISION) % sum_sq
write(ou,102) "k-effective (Track-length)", global_tallies(K_TRACKLENGTH) &
% sum, global_tallies(K_TRACKLENGTH) % sum_sq
write(ou,102) "k-effective (Absorption)", global_tallies(K_ABSORPTION) &
% sum, global_tallies(K_ABSORPTION) % sum_sq
write(ou,102) "k-effective (Collision)", global_tallies(RESULT_SUM, &
K_COLLISION), global_tallies(RESULT_SUM_SQ, K_COLLISION)
write(ou,102) "k-effective (Track-length)", global_tallies(RESULT_SUM, &
K_TRACKLENGTH), global_tallies(RESULT_SUM_SQ, K_TRACKLENGTH)
write(ou,102) "k-effective (Absorption)", global_tallies(RESULT_SUM, &
K_ABSORPTION), global_tallies(RESULT_SUM_SQ, K_ABSORPTION)
if (n_realizations > 3) write(ou,102) "Combined k-effective", k_combined
end if
write(ou,102) "Leakage Fraction", global_tallies(LEAKAGE) % sum, &
global_tallies(LEAKAGE) % sum_sq
write(ou,102) "Leakage Fraction", global_tallies(RESULT_SUM, LEAKAGE), &
global_tallies(RESULT_SUM_SQ, LEAKAGE)
else
if (master) call warning("Could not compute uncertainties -- only one &
&active batch simulated!")
if (run_mode == MODE_EIGENVALUE) then
write(ou,103) "k-effective (Collision)", global_tallies(K_COLLISION) % sum
write(ou,103) "k-effective (Track-length)", global_tallies(K_TRACKLENGTH) % sum
write(ou,103) "k-effective (Absorption)", global_tallies(K_ABSORPTION) % sum
write(ou,103) "k-effective (Collision)", global_tallies(RESULT_SUM, K_COLLISION)
write(ou,103) "k-effective (Track-length)", global_tallies(RESULT_SUM, K_TRACKLENGTH)
write(ou,103) "k-effective (Absorption)", global_tallies(RESULT_SUM, K_ABSORPTION)
end if
write(ou,103) "Leakage Fraction", global_tallies(LEAKAGE) % sum
write(ou,103) "Leakage Fraction", global_tallies(RESULT_SUM, LEAKAGE)
end if
write(ou,*)
@ -765,7 +765,7 @@ contains
end if
! Multiply uncertainty by t-value
t % results % sum_sq = t_value * t % results % sum_sq
t % results(RESULT_SUM_SQ,:,:) = t_value * t % results(RESULT_SUM_SQ,:,:)
end if
! Write header block
@ -777,6 +777,28 @@ contains
// trim(t % name), unit=unit_tally, level=3)
endif
! Write derivative information.
if (t % deriv /= NONE) then
associate(deriv => tally_derivs(t % deriv))
select case (deriv % variable)
case (DIFF_DENSITY)
write(unit=unit_tally, fmt="(' Density derivative Material ',A)") &
to_str(deriv % diff_material)
case (DIFF_NUCLIDE_DENSITY)
write(unit=unit_tally, fmt="(' Nuclide density derivative &
&Material ',A,' Nuclide ',A)") &
trim(to_str(deriv % diff_material)), &
trim(nuclides(deriv % diff_nuclide) % name)
case (DIFF_TEMPERATURE)
write(unit=unit_tally, fmt="(' Temperature derivative Material ',&
&A)") to_str(deriv % diff_material)
case default
call fatal_error("Differential tally dependent variable for tally "&
// trim(to_str(t % id)) // " not defined in output.F90.")
end select
end associate
end if
! Handle surface current tallies separately
if (t % type == TALLY_SURFACE_CURRENT) then
call write_surface_current(t, unit_tally)
@ -876,8 +898,8 @@ contains
score_names(abs(t % score_bins(k)))
write(UNIT=unit_tally, FMT='(1X,2A,1X,A,"+/- ",A)') &
repeat(" ", indent), score_name, &
to_str(t % results(score_index,filter_index) % sum), &
trim(to_str(t % results(score_index,filter_index) % sum_sq))
to_str(t % results(RESULT_SUM,score_index,filter_index)), &
trim(to_str(t % results(RESULT_SUM_SQ,score_index,filter_index)))
case (SCORE_SCATTER_PN, SCORE_NU_SCATTER_PN)
score_index = score_index - 1
do n_order = 0, t % moment_order(k)
@ -886,9 +908,8 @@ contains
score_names(abs(t % score_bins(k)))
write(UNIT=unit_tally, FMT='(1X,2A,1X,A,"+/- ",A)') &
repeat(" ", indent), score_name, &
to_str(t % results(score_index,filter_index) % sum), &
trim(to_str(t % results(score_index,filter_index) &
% sum_sq))
to_str(t % results(RESULT_SUM,score_index,filter_index)), &
trim(to_str(t % results(RESULT_SUM_SQ,score_index,filter_index)))
end do
k = k + t % moment_order(k)
case (SCORE_SCATTER_YN, SCORE_NU_SCATTER_YN, SCORE_FLUX_YN, &
@ -902,9 +923,9 @@ contains
// score_names(abs(t % score_bins(k)))
write(UNIT=unit_tally, FMT='(1X,2A,1X,A,"+/- ",A)') &
repeat(" ", indent), score_name, &
to_str(t % results(score_index,filter_index) % sum), &
trim(to_str(t % results(score_index,filter_index)&
% sum_sq))
to_str(t % results(RESULT_SUM,score_index,filter_index)), &
trim(to_str(t % results(RESULT_SUM_SQ,score_index,&
filter_index)))
end do
end do
k = k + (t % moment_order(k) + 1)**2 - 1
@ -916,8 +937,8 @@ contains
end if
write(UNIT=unit_tally, FMT='(1X,2A,1X,A,"+/- ",A)') &
repeat(" ", indent), score_name, &
to_str(t % results(score_index,filter_index) % sum), &
trim(to_str(t % results(score_index,filter_index) % sum_sq))
to_str(t % results(RESULT_SUM,score_index,filter_index)), &
trim(to_str(t % results(RESULT_SUM_SQ,score_index,filter_index)))
end select
end do
indent = indent - 2
@ -1020,16 +1041,16 @@ contains
* t % stride) + 1
write(UNIT=unit_tally, FMT='(5X,A,T35,A,"+/- ",A)') &
"Outgoing Current on Left", &
to_str(t % results(1,filter_index) % sum), &
trim(to_str(t % results(1,filter_index) % sum_sq))
to_str(t % results(RESULT_SUM,1,filter_index)), &
trim(to_str(t % results(RESULT_SUM_SQ,1,filter_index)))
matching_bins(i_filter_surf) = IN_LEFT
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
write(UNIT=unit_tally, FMT='(5X,A,T35,A,"+/- ",A)') &
"Incoming Current on Left", &
to_str(t % results(1,filter_index) % sum), &
trim(to_str(t % results(1,filter_index) % sum_sq))
to_str(t % results(RESULT_SUM,1,filter_index)), &
trim(to_str(t % results(RESULT_SUM_SQ,1,filter_index)))
! Right Surface
matching_bins(i_filter_surf) = OUT_RIGHT
@ -1037,16 +1058,16 @@ contains
* t % stride) + 1
write(UNIT=unit_tally, FMT='(5X,A,T35,A,"+/- ",A)') &
"Outgoing Current on Right", &
to_str(t % results(1,filter_index) % sum), &
trim(to_str(t % results(1,filter_index) % sum_sq))
to_str(t % results(RESULT_SUM,1,filter_index)), &
trim(to_str(t % results(RESULT_SUM_SQ,1,filter_index)))
matching_bins(i_filter_surf) = IN_RIGHT
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
write(UNIT=unit_tally, FMT='(5X,A,T35,A,"+/- ",A)') &
"Incoming Current on Right", &
to_str(t % results(1,filter_index) % sum), &
trim(to_str(t % results(1,filter_index) % sum_sq))
to_str(t % results(RESULT_SUM,1,filter_index)), &
trim(to_str(t % results(RESULT_SUM_SQ,1,filter_index)))
if (n_dim >= 2) then
@ -1056,16 +1077,16 @@ contains
* t % stride) + 1
write(UNIT=unit_tally, FMT='(5X,A,T35,A,"+/- ",A)') &
"Outgoing Current on Back", &
to_str(t % results(1,filter_index) % sum), &
trim(to_str(t % results(1,filter_index) % sum_sq))
to_str(t % results(RESULT_SUM,1,filter_index)), &
trim(to_str(t % results(RESULT_SUM_SQ,1,filter_index)))
matching_bins(i_filter_surf) = IN_BACK
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
write(UNIT=unit_tally, FMT='(5X,A,T35,A,"+/- ",A)') &
"Incoming Current on Back", &
to_str(t % results(1,filter_index) % sum), &
trim(to_str(t % results(1,filter_index) % sum_sq))
to_str(t % results(RESULT_SUM,1,filter_index)), &
trim(to_str(t % results(RESULT_SUM_SQ,1,filter_index)))
! Front Surface
matching_bins(i_filter_surf) = OUT_FRONT
@ -1073,16 +1094,16 @@ contains
* t % stride) + 1
write(UNIT=unit_tally, FMT='(5X,A,T35,A,"+/- ",A)') &
"Net Current on Front", &
to_str(t % results(1,filter_index) % sum), &
trim(to_str(t % results(1,filter_index) % sum_sq))
to_str(t % results(RESULT_SUM,1,filter_index)), &
trim(to_str(t % results(RESULT_SUM_SQ,1,filter_index)))
matching_bins(i_filter_surf) = IN_FRONT
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
write(UNIT=unit_tally, FMT='(5X,A,T35,A,"+/- ",A)') &
"Net Current on Front", &
to_str(t % results(1,filter_index) % sum), &
trim(to_str(t % results(1,filter_index) % sum_sq))
to_str(t % results(RESULT_SUM,1,filter_index)), &
trim(to_str(t % results(RESULT_SUM_SQ,1,filter_index)))
end if
if (n_dim == 3) then
@ -1092,16 +1113,16 @@ contains
* t % stride) + 1
write(UNIT=unit_tally, FMT='(5X,A,T35,A,"+/- ",A)') &
"Outgoing Current on Bottom", &
to_str(t % results(1,filter_index) % sum), &
trim(to_str(t % results(1,filter_index) % sum_sq))
to_str(t % results(RESULT_SUM,1,filter_index)), &
trim(to_str(t % results(RESULT_SUM_SQ,1,filter_index)))
matching_bins(i_filter_surf) = IN_BOTTOM
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
write(UNIT=unit_tally, FMT='(5X,A,T35,A,"+/- ",A)') &
"Incoming Current on Bottom", &
to_str(t % results(1,filter_index) % sum), &
trim(to_str(t % results(1,filter_index) % sum_sq))
to_str(t % results(RESULT_SUM,1,filter_index)), &
trim(to_str(t % results(RESULT_SUM_SQ,1,filter_index)))
! Top Surface
matching_bins(i_filter_surf) = OUT_TOP
@ -1109,16 +1130,16 @@ contains
* t % stride) + 1
write(UNIT=unit_tally, FMT='(5X,A,T35,A,"+/- ",A)') &
"Outgoing Current on Top", &
to_str(t % results(1,filter_index) % sum), &
trim(to_str(t % results(1,filter_index) % sum_sq))
to_str(t % results(RESULT_SUM,1,filter_index)), &
trim(to_str(t % results(RESULT_SUM_SQ,1,filter_index)))
matching_bins(i_filter_surf) = IN_TOP
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
write(UNIT=unit_tally, FMT='(5X,A,T35,A,"+/- ",A)') &
"Incoming Current on Top", &
to_str(t % results(1,filter_index) % sum), &
trim(to_str(t % results(1,filter_index) % sum_sq))
to_str(t % results(RESULT_SUM,1,filter_index)), &
trim(to_str(t % results(RESULT_SUM_SQ,1,filter_index)))
end if
end do
end do

8
src/physics_mg.F90
View file

@ -119,16 +119,16 @@ contains
! Score implicit absorption estimate of keff
!$omp atomic
global_tallies(K_ABSORPTION) % value = &
global_tallies(K_ABSORPTION) % value + p % absorb_wgt * &
global_tallies(RESULT_VALUE, K_ABSORPTION) = &
global_tallies(RESULT_VALUE, K_ABSORPTION) + p % absorb_wgt * &
material_xs % nu_fission / material_xs % absorption
else
! See if disappearance reaction happens
if (material_xs % absorption > prn() * material_xs % total) then
! Score absorption estimate of keff
!$omp atomic
global_tallies(K_ABSORPTION) % value = &
global_tallies(K_ABSORPTION) % value + p % wgt * &
global_tallies(RESULT_VALUE, K_ABSORPTION) = &
global_tallies(RESULT_VALUE, K_ABSORPTION) + p % wgt * &
material_xs % nu_fission / material_xs % absorption
p % alive = .false.

20
src/relaxng/tallies.rnc
View file

@ -15,6 +15,23 @@ element tallies {
)
}* &
element derivative {
(element id { xsd:int } | attribute id { xsd:int }) &
(element material { xsd:int } | attribute material { xsd:int }) &
( (element variable { ( "density") }
| attribute variable { ( "density" ) } ) |
(
(element variable { ( "nuclide_density" ) }
| attribute variable { ( "nuclide_density" ) } )
&
(element nuclide { xsd:string { maxLength = "12" } }
| attribute nuclide { xsd:string { maxLength = "12" } } ) |
)
(element variable { ( "temperature") }
| attribute variable { ( "temperature" ) } )
)
}* &
element tally {
(element id { xsd:int } | attribute id { xsd:int }) &
(element name { xsd:string { maxLength="52" } } |
@ -41,7 +58,8 @@ element tallies {
(element type { xsd:string } | attribute type { xsd:string }) &
(element threshold { xsd:double} | attribute threshold { xsd:double }) &
(element scores { list { xsd:string { maxLength = "20" }+ } } | attribute scores { list { xsd:string { maxLength = "20"}+ } } )?
}*
}* &
(element derivative { xsd:int } | attribute derivative { xsd:int } )?
}* &
element assume_separate { xsd:boolean }?

72
src/relaxng/tallies.rng
View file

@ -89,6 +89,68 @@
</interleave>
</element>
</zeroOrMore>
<zeroOrMore>
<element name="derivative">
<interleave>
<choice>
<element name="id">
<data type="int"/>
</element>
<attribute name="id">
<data type="int"/>
</attribute>
</choice>
<choice>
<element name="material">
<data type="int"/>
</element>
<attribute name="material">
<data type="int"/>
</attribute>
</choice>
<choice>
<choice>
<element name="variable">
<value>density</value>
</element>
<attribute name="variable">
<value>density</value>
</attribute>
</choice>
<interleave>
<choice>
<element name="variable">
<value>nuclide_density</value>
</element>
<attribute name="variable">
<value>nuclide_density</value>
</attribute>
</choice>
<choice>
<element name="nuclide">
<data type="string">
<param name="maxLength">12</param>
</data>
</element>
<attribute name="nuclide">
<data type="string">
<param name="maxLength">12</param>
</data>
</attribute>
</choice>
</interleave>
<choice>
<element name="variable">
<value>temperature</value>
</element>
<attribute name="variable">
<value>temperature</value>
</attribute>
</choice>
</choice>
</interleave>
</element>
</zeroOrMore>
<zeroOrMore>
<element name="tally">
<interleave>
@ -254,6 +316,16 @@
</interleave>
</element>
</zeroOrMore>
<optional>
<choice>
<element name="derivative">
<data type="int"/>
</element>
<attribute name="derivative">
<data type="int"/>
</attribute>
</choice>
</optional>
</interleave>
</element>
</zeroOrMore>

33
src/simulation.F90
View file

@ -20,8 +20,7 @@ module simulation
use source, only: initialize_source, sample_external_source
use state_point, only: write_state_point, write_source_point
use string, only: to_str
use tally, only: synchronize_tallies, setup_active_usertallies, &
reset_result
use tally, only: synchronize_tallies, setup_active_usertallies
use trigger, only: check_triggers
use tracking, only: transport
use volume_calc, only: run_volume_calculations
@ -84,7 +83,7 @@ contains
! ====================================================================
! LOOP OVER PARTICLES
!$omp parallel do schedule(static) firstprivate(p)
!$omp parallel do schedule(static) firstprivate(p) copyin(tally_derivs)
PARTICLE_LOOP: do i_work = 1, work
current_work = i_work
@ -220,7 +219,7 @@ contains
if (ufs) call count_source_for_ufs()
! Store current value of tracklength k
keff_generation = global_tallies(K_TRACKLENGTH) % value
keff_generation = global_tallies(RESULT_VALUE, K_TRACKLENGTH)
end if
end subroutine initialize_generation
@ -237,24 +236,24 @@ contains
!$omp parallel
!$omp critical
if (run_mode == MODE_EIGENVALUE) then
global_tallies(K_COLLISION) % value = &
global_tallies(K_COLLISION) % value + global_tally_collision
global_tallies(K_ABSORPTION) % value = &
global_tallies(K_ABSORPTION) % value + global_tally_absorption
global_tallies(K_TRACKLENGTH) % value = &
global_tallies(K_TRACKLENGTH) % value + global_tally_tracklength
global_tallies(RESULT_VALUE, K_COLLISION) = &
global_tallies(RESULT_VALUE, K_COLLISION) + global_tally_collision
global_tallies(RESULT_VALUE, K_ABSORPTION) = &
global_tallies(RESULT_VALUE, K_ABSORPTION) + global_tally_absorption
global_tallies(RESULT_VALUE, K_TRACKLENGTH) = &
global_tallies(RESULT_VALUE, K_TRACKLENGTH) + global_tally_tracklength
end if
global_tallies(LEAKAGE) % value = &
global_tallies(LEAKAGE) % value + global_tally_leakage
global_tallies(RESULT_VALUE, LEAKAGE) = &
global_tallies(RESULT_VALUE, LEAKAGE) + global_tally_leakage
!$omp end critical
! reset private tallies
if (run_mode == MODE_EIGENVALUE) then
global_tally_collision = 0
global_tally_absorption = 0
global_tally_tracklength = 0
global_tally_collision = ZERO
global_tally_absorption = ZERO
global_tally_tracklength = ZERO
end if
global_tally_leakage = 0
global_tally_leakage = ZERO
!$omp end parallel
if (run_mode == MODE_EIGENVALUE) then
@ -302,7 +301,7 @@ contains
! Reset global tally results
if (.not. active_batches) then
call reset_result(global_tallies)
global_tallies(:,:) = ZERO
n_realizations = 0
end if

92
src/state_point.F90
View file

@ -50,7 +50,8 @@ contains
integer, allocatable :: key_array(:)
integer(HID_T) :: file_id
integer(HID_T) :: cmfd_group, tallies_group, tally_group, meshes_group, &
mesh_group, filter_group, runtime_group
mesh_group, filter_group, derivs_group, deriv_group, &
runtime_group
character(MAX_WORD_LEN), allocatable :: str_array(:)
character(MAX_FILE_LEN) :: filename
type(RegularMesh), pointer :: meshp
@ -198,6 +199,38 @@ contains
call close_group(meshes_group)
! Write information for derivatives.
if (size(tally_derivs) > 0) then
derivs_group = create_group(tallies_group, "derivatives")
do i = 1, size(tally_derivs)
associate(deriv => tally_derivs(i))
deriv_group = create_group(derivs_group, "derivative " &
// trim(to_str(deriv % id)))
select case (deriv % variable)
case (DIFF_DENSITY)
call write_dataset(deriv_group, "independent variable", "density")
call write_dataset(deriv_group, "material", deriv % diff_material)
case (DIFF_NUCLIDE_DENSITY)
call write_dataset(deriv_group, "independent variable", &
"nuclide_density")
call write_dataset(deriv_group, "material", deriv % diff_material)
call write_dataset(deriv_group, "nuclide", &
nuclides(deriv % diff_nuclide) % name)
case (DIFF_TEMPERATURE)
call write_dataset(deriv_group, "independent variable", &
"temperature")
call write_dataset(deriv_group, "material", deriv % diff_material)
case default
call fatal_error("Independent variable for derivative " &
// trim(to_str(deriv % id)) // " not defined in &
&state_point.F90.")
end select
call close_group(deriv_group)
end associate
end do
call close_group(derivs_group)
end if
! Write number of tallies
call write_dataset(tallies_group, "n_tallies", n_tallies)
@ -273,6 +306,13 @@ contains
call write_dataset(tally_group, "nuclides", str_array)
deallocate(str_array)
! Write derivative information.
if (tally % deriv /= NONE) then
call write_dataset(tally_group, "derivative", &
tally_derivs(tally % deriv) % id)
end if
! Write scores.
call write_dataset(tally_group, "n_score_bins", tally % n_score_bins)
allocate(str_array(size(tally % score_bins)))
do j = 1, size(tally % score_bins)
@ -353,7 +393,7 @@ contains
! Write sum and sum_sq for each bin
tally_group = open_group(tallies_group, "tally " &
// to_str(tally % id))
call write_dataset(tally_group, "results", tally % results)
call tally % write_results_hdf5(tally_group)
call close_group(tally_group)
end do TALLY_RESULTS
@ -481,7 +521,7 @@ contains
integer :: n_bins ! total number of bins
integer(HID_T) :: tallies_group, tally_group
real(8), allocatable :: tally_temp(:,:,:) ! contiguous array of results
real(8), target :: global_temp(2,N_GLOBAL_TALLIES)
real(8), target :: global_temp(3,N_GLOBAL_TALLIES)
#ifdef MPI
real(8) :: dummy ! temporary receive buffer for non-root reduces
#endif
@ -489,7 +529,7 @@ contains
type(ElemKeyValueII), pointer :: current
type(ElemKeyValueII), pointer :: next
type(TallyObject), pointer :: tally
type(TallyResult), allocatable :: tallyresult_temp(:,:)
type(TallyObject) :: dummy_tally
! ==========================================================================
! COLLECT AND WRITE GLOBAL TALLIES
@ -505,9 +545,8 @@ contains
end if
! Copy global tallies into temporary array for reducing
n_bins = 2 * N_GLOBAL_TALLIES
global_temp(1,:) = global_tallies(:)%sum
global_temp(2,:) = global_tallies(:)%sum_sq
n_bins = 3 * N_GLOBAL_TALLIES
global_temp(:,:) = global_tallies(:,:)
if (master) then
! The MPI_IN_PLACE specifier allows the master to copy values into a
@ -519,20 +558,11 @@ contains
! Transfer values to value on master
if (current_batch == n_max_batches .or. satisfy_triggers) then
global_tallies(:)%sum = global_temp(1,:)
global_tallies(:)%sum_sq = global_temp(2,:)
global_tallies(:,:) = global_temp(:,:)
end if
! Put reduced value in temporary tally result
allocate(tallyresult_temp(N_GLOBAL_TALLIES, 1))
tallyresult_temp(:,1)%sum = global_temp(1,:)
tallyresult_temp(:,1)%sum_sq = global_temp(2,:)
! Write out global tallies sum and sum_sq
call write_dataset(file_id, "global_tallies", tallyresult_temp)
! Deallocate temporary tally result
deallocate(tallyresult_temp)
call write_dataset(file_id, "global_tallies", global_temp)
else
! Receive buffer not significant at other processors
#ifdef MPI
@ -568,15 +598,15 @@ contains
tally => tallies(i)
! Determine size of tally results array
m = size(tally%results, 1)
n = size(tally%results, 2)
m = size(tally%results, 2)
n = size(tally%results, 3)
n_bins = m*n*2
! Allocate array for storing sums and sums of squares, but
! contiguously in memory for each
allocate(tally_temp(2,m,n))
tally_temp(1,:,:) = tally%results(:,:)%sum
tally_temp(2,:,:) = tally%results(:,:)%sum_sq
tally_temp(1,:,:) = tally%results(RESULT_SUM,:,:)
tally_temp(2,:,:) = tally%results(RESULT_SUM_SQ,:,:)
if (master) then
tally_group = open_group(tallies_group, "tally " // &
@ -592,20 +622,20 @@ contains
! At the end of the simulation, store the results back in the
! regular TallyResults array
if (current_batch == n_max_batches .or. satisfy_triggers) then
tally%results(:,:)%sum = tally_temp(1,:,:)
tally%results(:,:)%sum_sq = tally_temp(2,:,:)
tally%results(RESULT_SUM,:,:) = tally_temp(1,:,:)
tally%results(RESULT_SUM_SQ,:,:) = tally_temp(2,:,:)
end if
! Put in temporary tally result
allocate(tallyresult_temp(m,n))
tallyresult_temp(:,:)%sum = tally_temp(1,:,:)
tallyresult_temp(:,:)%sum_sq = tally_temp(2,:,:)
allocate(dummy_tally % results(3,m,n))
dummy_tally % results(RESULT_SUM,:,:) = tally_temp(1,:,:)
dummy_tally % results(RESULT_SUM_SQ,:,:) = tally_temp(2,:,:)
! Write reduced tally results to file
call write_dataset(tally_group, "results", tally%results)
call dummy_tally % write_results_hdf5(tally_group)
! Deallocate temporary tally result
deallocate(tallyresult_temp)
deallocate(dummy_tally % results)
else
! Receive buffer not significant at other processors
#ifdef MPI
@ -771,7 +801,7 @@ contains
call read_dataset(n_realizations, file_id, "n_realizations", indep=.true.)
! Read global tally data
call read_dataset(file_id, "global_tallies", global_tallies)
call read_dataset(global_tallies, file_id, "global_tallies")
! Check if tally results are present
tallies_group = open_group(file_id, "tallies")
@ -787,7 +817,7 @@ contains
! Read sum, sum_sq, and N for each bin
tally_group = open_group(tallies_group, "tally " // &
trim(to_str(tally % id)))
call read_dataset(tally_group, "results", tally % results)
call tally % read_results_hdf5(tally_group)
call read_dataset(tally % n_realizations, tally_group, &
"n_realizations")
call close_group(tally_group)

989
src/tally.F90
View file

File diff suppressed because it is too large Load diff

99
src/tally_header.F90
View file

@ -1,22 +1,27 @@
module tally_header
use, intrinsic :: ISO_C_BINDING
use hdf5
use constants, only: NONE, N_FILTER_TYPES
use tally_filter_header, only: TallyFilterContainer
use trigger_header, only: TriggerObject
use, intrinsic :: ISO_C_BINDING
implicit none
!===============================================================================
! TALLYRESULT provides accumulation of results in a particular tally bin
! TALLYDERIVATIVE describes a first-order derivative that can be applied to
! tallies.
!===============================================================================
type, bind(C) :: TallyResult
real(C_DOUBLE) :: value = 0.
real(C_DOUBLE) :: sum = 0.
real(C_DOUBLE) :: sum_sq = 0.
end type TallyResult
type TallyDerivative
integer :: id
integer :: variable
integer :: diff_material
integer :: diff_nuclide
real(8) :: flux_deriv
end type TallyDerivative
!===============================================================================
! TALLYOBJECT describes a user-specified tally. The region of phase space to
@ -68,7 +73,7 @@ module tally_header
integer :: total_filter_bins
integer :: total_score_bins
type(TallyResult), allocatable :: results(:,:)
real(C_DOUBLE), allocatable :: results(:,:,:)
! reset property - allows a tally to be reset after every batch
logical :: reset = .false.
@ -79,6 +84,82 @@ module tally_header
! Tally precision triggers
integer :: n_triggers = 0 ! # of triggers
type(TriggerObject), allocatable :: triggers(:) ! Array of triggers
! Index for the TallyDerivative for differential tallies.
integer :: deriv = NONE
contains
procedure :: write_results_hdf5
procedure :: read_results_hdf5
end type TallyObject
contains
subroutine write_results_hdf5(this, group_id)
class(TallyObject), intent(in) :: this
integer(HID_T), intent(in) :: group_id
integer :: hdf5_err
integer(HID_T) :: dset, dspace
integer(HID_T) :: memspace
integer(HSIZE_T) :: dims(3)
integer(HSIZE_T) :: dims_slab(3)
integer(HSIZE_T) :: offset(3) = [1,0,0]
! Create file dataspace
dims_slab(:) = shape(this % results)
dims_slab(1) = 2
call h5screate_simple_f(3, dims_slab, dspace, hdf5_err)
! Create memory dataspace that contains only SUM and SUM_SQ values
dims(:) = shape(this % results)
call h5screate_simple_f(3, dims, memspace, hdf5_err)
call h5sselect_hyperslab_f(memspace, H5S_SELECT_SET_F, offset, dims_slab, &
hdf5_err)
! Create and write to dataset
call h5dcreate_f(group_id, "results", H5T_NATIVE_DOUBLE, dspace, dset, &
hdf5_err)
call h5dwrite_f(dset, H5T_NATIVE_DOUBLE, this % results, dims_slab, &
hdf5_err, mem_space_id=memspace)
! Close identifiers
call h5dclose_f(dset, hdf5_err)
call h5sclose_f(memspace, hdf5_err)
call h5sclose_f(dspace, hdf5_err)
end subroutine write_results_hdf5
subroutine read_results_hdf5(this, group_id)
class(TallyObject), intent(inout) :: this
integer(HID_T), intent(in) :: group_id
integer :: hdf5_err
integer(HID_T) :: dset, dspace
integer(HID_T) :: memspace
integer(HSIZE_T) :: dims(3)
integer(HSIZE_T) :: dims_slab(3)
integer(HSIZE_T) :: offset(3) = [1,0,0]
! Create file dataspace
dims_slab(:) = shape(this % results)
dims_slab(1) = 2
call h5screate_simple_f(3, dims_slab, dspace, hdf5_err)
! Create memory dataspace that contains only SUM and SUM_SQ values
dims(:) = shape(this % results)
call h5screate_simple_f(3, dims, memspace, hdf5_err)
call h5sselect_hyperslab_f(memspace, H5S_SELECT_SET_F, offset, dims_slab, &
hdf5_err)
! Create and write to dataset
call h5dopen_f(group_id, "results", dset, hdf5_err)
call h5dread_f(dset, H5T_NATIVE_DOUBLE, this % results, dims_slab, &
hdf5_err, mem_space_id=memspace)
! Close identifiers
call h5dclose_f(dset, hdf5_err)
call h5sclose_f(memspace, hdf5_err)
call h5sclose_f(dspace, hdf5_err)
end subroutine read_results_hdf5
end module tally_header

6
src/tally_initialize.F90
View file

@ -20,7 +20,8 @@ contains
subroutine configure_tallies()
! Allocate global tallies
allocate(global_tallies(N_GLOBAL_TALLIES))
allocate(global_tallies(3, N_GLOBAL_TALLIES))
global_tallies(:,:) = ZERO
call setup_tally_arrays()
@ -62,7 +63,8 @@ contains
t % total_score_bins = t % n_score_bins * t % n_nuclide_bins
! Allocate results array
allocate(t % results(t % total_score_bins, t % total_filter_bins))
allocate(t % results(3, t % total_score_bins, t % total_filter_bins))
t % results(:,:,:) = ZERO
end do TALLY_LOOP

13
src/tracking.F90
View file

@ -14,7 +14,9 @@ module tracking
use random_lcg, only: prn
use string, only: to_str
use tally, only: score_analog_tally, score_tracklength_tally, &
score_collision_tally, score_surface_current
score_collision_tally, score_surface_current, &
score_track_derivative, &
score_collision_derivative, zero_flux_derivs
use track_output, only: initialize_particle_track, write_particle_track, &
add_particle_track, finalize_particle_track
@ -63,6 +65,9 @@ contains
call initialize_particle_track()
endif
! Every particle starts with no accumulated flux derivative.
if (active_tallies % size() > 0) call zero_flux_derivs()
EVENT_LOOP: do
! If the cell hasn't been determined based on the particle's location,
! initiate a search for the current cell. This generally happens at the
@ -137,6 +142,9 @@ contains
distance * material_xs % nu_fission
end if
! Score flux derivative accumulators for differential tallies.
if (active_tallies % size() > 0) call score_track_derivative(p, distance)
if (d_collision > d_boundary) then
! ====================================================================
! PARTICLE CROSSES SURFACE
@ -213,6 +221,9 @@ contains
p % coord(j + 1) % uvw = p % coord(j) % uvw
end if
end do
! Score flux derivative accumulators for differential tallies.
if (active_tallies % size() > 0) call score_collision_derivative(p)
end if
! Save coordinates for tallying purposes

16
src/trigger.F90
View file

@ -432,17 +432,19 @@ contains
real(8), intent(inout) :: rel_err ! tally relative error
integer, intent(in) :: score_index ! tally results score index
integer, intent(in) :: filter_index ! tally results filter index
integer :: n ! number of realizations
real(8) :: mean ! tally mean
type(TallyResult) :: tally_result ! pointer to TallyResult
type(TallyObject), pointer :: t ! tally pointer
type(TallyObject), intent(in) :: t ! tally
integer :: n ! number of realizations
real(8) :: mean ! tally mean
real(8) :: tally_sum, tally_sum_sq ! results for a single tally bin
n = t % n_realizations
tally_result = t % results(score_index, filter_index)
tally_sum = t % results(RESULT_SUM, score_index, filter_index)
tally_sum_sq = t % results(RESULT_SUM_SQ, score_index, filter_index)
! Compute the tally mean and standard deviation
mean = tally_result % sum / n
std_dev = sqrt((tally_result % sum_sq / n - mean * mean) / (n - 1))
mean = tally_sum / n
std_dev = sqrt((tally_sum_sq / n - mean * mean) / (n - 1))
! Compute the relative error if the mean is non-zero
if (mean == ZERO) then

1
tests/test_diff_tally/inputs_true.dat Normal file
View file

@ -0,0 +1 @@
df5a0ee7d353fe25344496058cea42e3244a7ab32c13c34c03b4b2f6adf88579f08a71ec0d22d8d264e50e3663717068ff54308e2d16aed8e085c222a8c2fdbd

177
tests/test_diff_tally/results_true.dat Normal file
View file

@ -0,0 +1,177 @@
d_material,d_nuclide,d_variable,score,mean,std. dev.
1,,density,nu-fission,0.0000000e+00,0.0000000e+00
1,,density,scatter,3.9902949e-02,9.1258428e-03
1,,density,nu-fission,0.0000000e+00,0.0000000e+00
1,,density,scatter,9.8213745e-04,9.8213745e-04
1,,density,nu-fission,2.7945389e-02,2.0999368e-02
1,,density,scatter,4.0626155e-01,1.7017674e-02
1,,density,nu-fission,2.4595279e-02,2.1224443e-02
1,,density,scatter,2.6505962e-03,2.3698970e-03
1,,density,nu-fission,0.0000000e+00,0.0000000e+00
1,,density,scatter,-1.2721728e-01,1.8810986e-01
1,,density,nu-fission,0.0000000e+00,0.0000000e+00
1,,density,scatter,0.0000000e+00,0.0000000e+00
1,,density,nu-fission,0.0000000e+00,0.0000000e+00
1,,density,scatter,4.1950592e-02,7.3680059e-02
1,,density,nu-fission,0.0000000e+00,0.0000000e+00
1,,density,scatter,0.0000000e+00,0.0000000e+00
1,O16,nuclide_density,nu-fission,0.0000000e+00,0.0000000e+00
1,O16,nuclide_density,scatter,-1.9798576e-02,1.9798576e-02
1,O16,nuclide_density,nu-fission,9.3632921e-01,2.0322286e+00
1,O16,nuclide_density,scatter,-2.2042881e-01,1.6781767e-01
1,O16,nuclide_density,nu-fission,0.0000000e+00,0.0000000e+00
1,O16,nuclide_density,scatter,0.0000000e+00,0.0000000e+00
1,O16,nuclide_density,nu-fission,0.0000000e+00,0.0000000e+00
1,O16,nuclide_density,scatter,0.0000000e+00,0.0000000e+00
1,U235,nuclide_density,nu-fission,0.0000000e+00,0.0000000e+00
1,U235,nuclide_density,scatter,2.3761276e+00,2.3761276e+00
1,U235,nuclide_density,nu-fission,7.7441690e+02,8.6460933e+01
1,U235,nuclide_density,scatter,9.6917103e+01,1.0750582e+01
1,U235,nuclide_density,nu-fission,0.0000000e+00,0.0000000e+00
1,U235,nuclide_density,scatter,0.0000000e+00,0.0000000e+00
1,U235,nuclide_density,nu-fission,0.0000000e+00,0.0000000e+00
1,U235,nuclide_density,scatter,0.0000000e+00,0.0000000e+00
1,,temperature,nu-fission,0.0000000e+00,0.0000000e+00
1,,temperature,scatter,-2.0607998e-06,2.0607998e-06
1,,temperature,nu-fission,3.8845881e-06,3.5405083e-05
1,,temperature,scatter,-6.6772431e-07,2.0750826e-07
1,,temperature,nu-fission,0.0000000e+00,0.0000000e+00
1,,temperature,scatter,0.0000000e+00,0.0000000e+00
1,,temperature,nu-fission,0.0000000e+00,0.0000000e+00
1,,temperature,scatter,0.0000000e+00,0.0000000e+00
3,,density,flux,-7.6179774e+00,5.0326447e+00
3,,density,flux,-1.6242019e+01,6.6231774e+00
1,,density,flux,-2.2394746e-01,5.4534297e-02
1,,density,flux,-5.7132854e-02,1.6898134e-01
1,O16,nuclide_density,flux,-1.3403658e+01,1.3017127e+01
1,O16,nuclide_density,flux,-1.0499715e+01,2.1684953e+01
1,U235,nuclide_density,flux,-2.4741168e+03,5.7986812e+01
1,U235,nuclide_density,flux,-1.9955533e+03,4.3254820e+02
1,,temperature,flux,1.0601815e-04,3.3076550e-04
1,,temperature,flux,1.6664450e-04,2.8761112e-04
3,,density,total,-3.3161585e+00,2.5615932e+00
3,,density,absorption,-4.3491248e-01,5.1559021e-01
3,,density,scatter,-2.8812460e+00,2.0540305e+00
3,,density,fission,-2.3060366e-01,3.1191109e-01
3,,density,nu-fission,-5.6817321e-01,7.6143389e-01
3,,density,total,-2.8908007e-01,3.9383184e-01
3,,density,absorption,-2.5237485e-01,3.6565053e-01
3,,density,scatter,-3.6705221e-02,2.9275371e-02
3,,density,fission,-2.1271120e-01,3.0873698e-01
3,,density,nu-fission,-5.1866339e-01,7.5235533e-01
3,,density,total,2.7320569e+00,8.7709465e+00
3,,density,absorption,8.4322466e-02,1.2807898e-01
3,,density,scatter,2.6477344e+00,8.6432567e+00
3,,density,fission,0.0000000e+00,0.0000000e+00
3,,density,nu-fission,0.0000000e+00,0.0000000e+00
3,,density,total,0.0000000e+00,0.0000000e+00
3,,density,absorption,0.0000000e+00,0.0000000e+00
3,,density,scatter,0.0000000e+00,0.0000000e+00
3,,density,fission,0.0000000e+00,0.0000000e+00
3,,density,nu-fission,0.0000000e+00,0.0000000e+00
1,,density,total,4.7523439e-01,2.5652004e-02
1,,density,absorption,3.1301234e-02,7.6793142e-03
1,,density,scatter,4.4393315e-01,1.8023252e-02
1,,density,fission,1.5559355e-02,3.2310294e-03
1,,density,nu-fission,3.8658109e-02,7.8012233e-03
1,,density,total,2.2062939e-02,4.3689061e-03
1,,density,absorption,1.6845626e-02,4.0750601e-03
1,,density,scatter,5.2173132e-03,2.9663452e-04
1,,density,fission,1.3503649e-02,3.2642875e-03
1,,density,nu-fission,3.2945474e-02,7.9507685e-03
1,,density,total,-9.4735691e-02,2.5584665e-01
1,,density,absorption,-4.0246397e-03,5.2111770e-03
1,,density,scatter,-9.0711052e-02,2.5073762e-01
1,,density,fission,0.0000000e+00,0.0000000e+00
1,,density,nu-fission,0.0000000e+00,0.0000000e+00
1,,density,total,0.0000000e+00,0.0000000e+00
1,,density,absorption,0.0000000e+00,0.0000000e+00
1,,density,scatter,0.0000000e+00,0.0000000e+00
1,,density,fission,0.0000000e+00,0.0000000e+00
1,,density,nu-fission,0.0000000e+00,0.0000000e+00
1,O16,nuclide_density,total,4.4488215e+01,5.3981603e+00
1,O16,nuclide_density,absorption,-2.9372611e-01,9.3282813e-01
1,O16,nuclide_density,scatter,4.4781941e+01,4.5257555e+00
1,O16,nuclide_density,fission,9.9827872e-02,4.0698927e-01
1,O16,nuclide_density,nu-fission,2.3342880e-01,9.8749212e-01
1,O16,nuclide_density,total,4.0914051e-02,6.2301812e-01
1,O16,nuclide_density,absorption,8.8289354e-02,5.5570922e-01
1,O16,nuclide_density,scatter,-4.7375303e-02,6.7714071e-02
1,O16,nuclide_density,fission,1.3678733e-01,4.4006413e-01
1,O16,nuclide_density,nu-fission,3.3257520e-01,1.0719313e+00
1,O16,nuclide_density,total,-1.6055918e+01,2.3439335e+01
1,O16,nuclide_density,absorption,-3.9601231e-01,3.0131622e-01
1,O16,nuclide_density,scatter,-1.5659906e+01,2.3140407e+01
1,O16,nuclide_density,fission,0.0000000e+00,0.0000000e+00
1,O16,nuclide_density,nu-fission,0.0000000e+00,0.0000000e+00
1,O16,nuclide_density,total,0.0000000e+00,0.0000000e+00
1,O16,nuclide_density,absorption,0.0000000e+00,0.0000000e+00
1,O16,nuclide_density,scatter,0.0000000e+00,0.0000000e+00
1,O16,nuclide_density,fission,0.0000000e+00,0.0000000e+00
1,O16,nuclide_density,nu-fission,0.0000000e+00,0.0000000e+00
1,U235,nuclide_density,total,-6.1444902e+02,6.0111692e+01
1,U235,nuclide_density,absorption,2.3977134e+02,4.5868405e+01
1,U235,nuclide_density,scatter,-8.5422036e+02,2.4877920e+01
1,U235,nuclide_density,fission,3.1319132e+02,3.1465126e+01
1,U235,nuclide_density,nu-fission,7.6415465e+02,7.6501604e+01
1,U235,nuclide_density,total,5.0916179e+02,3.5499319e+01
1,U235,nuclide_density,absorption,3.9477360e+02,3.4326165e+01
1,U235,nuclide_density,scatter,1.1438818e+02,2.7378217e+00
1,U235,nuclide_density,fission,3.1360739e+02,3.2243605e+01
1,U235,nuclide_density,nu-fission,7.6527228e+02,7.8683697e+01
1,U235,nuclide_density,total,-4.1815711e+03,8.7796611e+02
1,U235,nuclide_density,absorption,-1.0563226e+02,2.5308627e+01
1,U235,nuclide_density,scatter,-4.0759388e+03,8.5332185e+02
1,U235,nuclide_density,fission,0.0000000e+00,0.0000000e+00
1,U235,nuclide_density,nu-fission,0.0000000e+00,0.0000000e+00
1,U235,nuclide_density,total,0.0000000e+00,0.0000000e+00
1,U235,nuclide_density,absorption,0.0000000e+00,0.0000000e+00
1,U235,nuclide_density,scatter,0.0000000e+00,0.0000000e+00
1,U235,nuclide_density,fission,0.0000000e+00,0.0000000e+00
1,U235,nuclide_density,nu-fission,0.0000000e+00,0.0000000e+00
1,,temperature,total,2.1367130e-04,1.8632816e-04
1,,temperature,absorption,6.9456249e-05,3.2199390e-05
1,,temperature,scatter,1.4421505e-04,1.5754180e-04
1,,temperature,fission,3.0674981e-06,1.8048377e-05
1,,temperature,nu-fission,7.4768877e-06,4.3976766e-05
1,,temperature,total,5.8223116e-06,2.3786474e-05
1,,temperature,absorption,5.2142698e-06,2.1901911e-05
1,,temperature,scatter,6.0804185e-07,1.9834890e-06
1,,temperature,fission,3.0674200e-06,1.8048542e-05
1,,temperature,nu-fission,7.4767028e-06,4.3977153e-05
1,,temperature,total,2.1510156e-04,4.6388707e-04
1,,temperature,absorption,2.2409527e-06,8.2901060e-06
1,,temperature,scatter,2.1286061e-04,4.5560034e-04
1,,temperature,fission,0.0000000e+00,0.0000000e+00
1,,temperature,nu-fission,0.0000000e+00,0.0000000e+00
1,,temperature,total,0.0000000e+00,0.0000000e+00
1,,temperature,absorption,0.0000000e+00,0.0000000e+00
1,,temperature,scatter,0.0000000e+00,0.0000000e+00
1,,temperature,fission,0.0000000e+00,0.0000000e+00
1,,temperature,nu-fission,0.0000000e+00,0.0000000e+00
3,,density,absorption,-1.6517201e-01,3.0984738e-01
3,,density,absorption,8.0402344e-03,1.4689335e-01
1,,density,absorption,2.9069882e-02,2.2139609e-03
1,,density,absorption,-9.4690065e-03,8.6605392e-03
1,O16,nuclide_density,absorption,7.6911962e-01,4.1945687e-01
1,O16,nuclide_density,absorption,-6.3724795e-01,6.6341488e-01
1,U235,nuclide_density,absorption,1.4109543e+02,1.8518890e+01
1,U235,nuclide_density,absorption,-1.2052822e+02,3.2084002e+01
1,,temperature,absorption,3.9995404e-05,2.1703384e-05
1,,temperature,absorption,2.7274361e-06,8.9339257e-06
3,,density,nu-fission,0.0000000e+00,0.0000000e+00
3,,density,scatter,-8.7934551e-01,1.0210652e+00
3,,density,nu-fission,0.0000000e+00,0.0000000e+00
3,,density,scatter,-3.5236516e-03,3.5236516e-03
3,,density,nu-fission,7.2953012e-02,2.9725863e-01
3,,density,scatter,-2.2716410e+00,1.2546347e+00
3,,density,nu-fission,8.9171481e-02,3.0634856e-01
3,,density,scatter,-1.0151747e-03,9.6296775e-03
3,,density,nu-fission,0.0000000e+00,0.0000000e+00
3,,density,scatter,2.3163923e+00,4.8766690e+00
3,,density,nu-fission,0.0000000e+00,0.0000000e+00
3,,density,scatter,0.0000000e+00,0.0000000e+00
3,,density,nu-fission,0.0000000e+00,0.0000000e+00
3,,density,scatter,4.0762434e-01,3.8504141e+00
3,,density,nu-fission,0.0000000e+00,0.0000000e+00
3,,density,scatter,0.0000000e+00,0.0000000e+00

140
tests/test_diff_tally/test_diff_tally.py Normal file
View file

@ -0,0 +1,140 @@
#!/usr/bin/env python
import glob
import os
import sys
import pandas as pd
sys.path.insert(0, os.pardir)
from testing_harness import PyAPITestHarness
import openmc
class DiffTallyTestHarness(PyAPITestHarness):
def _build_inputs(self):
# Build default materials/geometry
self._input_set.build_default_materials_and_geometry()
# Set settings explicitly
self._input_set.settings.batches = 3
self._input_set.settings.inactive = 0
self._input_set.settings.particles = 100
self._input_set.settings.source = openmc.Source(space=openmc.stats.Box(
[-160, -160, -183], [160, 160, 183]))
self._input_set.settings.temperature['multipole'] = True
self._input_set.tallies = openmc.Tallies()
filt_mats = openmc.MaterialFilter((1, 3))
filt_eout = openmc.EnergyoutFilter((0.0, 0.625, 20.0e6))
# We want density derivatives for both water and fuel to get coverage
# for both fissile and non-fissile materials.
d1 = openmc.TallyDerivative(derivative_id=1)
d1.variable = 'density'
d1.material = 3
d2 = openmc.TallyDerivative(derivative_id=2)
d2.variable = 'density'
d2.material = 1
# O-16 is a good nuclide to test against because it is present in both
# water and fuel. Some routines need to recognize that they have the
# perturbed nuclide but not the perturbed material.
d3 = openmc.TallyDerivative(derivative_id=3)
d3.variable = 'nuclide_density'
d3.material = 1
d3.nuclide = 'O16'
# A fissile nuclide, just for good measure.
d4 = openmc.TallyDerivative(derivative_id=4)
d4.variable = 'nuclide_density'
d4.material = 1
d4.nuclide = 'U235'
# Temperature derivatives.
d5 = openmc.TallyDerivative(derivative_id=5)
d5.variable = 'temperature'
d5.material = 1
derivs = [d1, d2, d3, d4, d5]
# Cover the flux score.
for i in range(5):
t = openmc.Tally()
t.add_score('flux')
t.add_filter(filt_mats)
t.derivative = derivs[i]
self._input_set.tallies.append(t)
# Cover supported scores with a collision estimator.
for i in range(5):
t = openmc.Tally()
t.add_score('total')
t.add_score('absorption')
t.add_score('scatter')
t.add_score('fission')
t.add_score('nu-fission')
t.add_filter(filt_mats)
t.add_nuclide('total')
t.add_nuclide('U235')
t.derivative = derivs[i]
self._input_set.tallies.append(t)
# Cover an analog estimator.
for i in range(5):
t = openmc.Tally()
t.add_score('absorption')
t.add_filter(filt_mats)
t.estimator = 'analog'
t.derivative = derivs[i]
self._input_set.tallies.append(t)
# Energyout filter and total nuclide for the density derivatives.
for i in range(2):
t = openmc.Tally()
t.add_score('nu-fission')
t.add_score('scatter')
t.add_filter(filt_mats)
t.add_filter(filt_eout)
t.add_nuclide('total')
t.add_nuclide('U235')
t.derivative = derivs[i]
self._input_set.tallies.append(t)
# Energyout filter without total nuclide for other derivatives.
for i in range(2, 5):
t = openmc.Tally()
t.add_score('nu-fission')
t.add_score('scatter')
t.add_filter(filt_mats)
t.add_filter(filt_eout)
t.add_nuclide('U235')
t.derivative = derivs[i]
self._input_set.tallies.append(t)
self._input_set.export()
def _get_results(self):
# Read the statepoint and summary files.
statepoint = glob.glob(os.path.join(os.getcwd(), self._sp_name))[0]
sp = openmc.StatePoint(statepoint)
# Extract the tally data as a Pandas DataFrame.
df = pd.DataFrame()
for t in sp.tallies.values():
df = df.append(t.get_pandas_dataframe(), ignore_index=True)
# Extract the relevant data as a CSV string.
cols = ('d_material', 'd_nuclide', 'd_variable', 'score', 'mean',
'std. dev.')
return df.to_csv(None, columns=cols, index=False, float_format='%.7e')
def _cleanup(self):
super(DiffTallyTestHarness, self)._cleanup()
f = os.path.join(os.getcwd(), 'tallies.xml')
if os.path.exists(f): os.remove(f)
if __name__ == '__main__':
harness = DiffTallyTestHarness('statepoint.3.h5', True)
harness.main()

2
tests/test_mgxs_library_condense/inputs_true.dat
View file

@ -1 +1 @@
e86f24e20f37096c7898f459fc5bf336f3e0670fb30f321443f2bb3a01a154ef3d234bb48d4cee8e0d3730fd182b6a2051ec1b4c09d771420886a56ba928fdd9
4291b40470e7d59383c9e51c4178ca923b698cb1aaea16c1982fe3789ca980df10a65b84fb5021dacd4d290ebc232c2579f99b6c990fb6b8f28f67eef2aabcb2

27
tests/test_mgxs_library_condense/results_true.dat
View file

@ -40,6 +40,8 @@
0 10000 1 total 4.996730e-07 3.650635e-08
material group in nuclide mean std. dev.
0 10000 1 total 0.090004 0.006367
material group in group out nuclide mean std. dev.
0 10000 1 1 total 0.084542 0.005716
material delayedgroup group in nuclide mean std. dev.
0 10000 1 1 total 0.000021 0.000001
1 10000 2 1 total 0.000110 0.000008
@ -68,6 +70,13 @@
3 10000 4 1 total 0.302780 0.109110
4 10000 5 1 total 0.000000 0.000000
5 10000 6 1 total 0.000000 0.000000
material delayedgroup group in group out nuclide mean std. dev.
0 10000 1 1 1 total 0.000000 0.000000
1 10000 2 1 1 total 0.000384 0.000236
2 10000 3 1 1 total 0.000179 0.000180
3 10000 4 1 1 total 0.000730 0.000188
4 10000 5 1 1 total 0.000000 0.000000
5 10000 6 1 1 total 0.000000 0.000000
material group in nuclide mean std. dev.
0 10001 1 total 0.311594 0.013793
material group in nuclide mean std. dev.
@ -110,6 +119,8 @@
0 10001 1 total 5.454760e-07 4.949800e-08
material group in nuclide mean std. dev.
0 10001 1 total 0.0 0.0
material group in group out nuclide mean std. dev.
0 10001 1 1 total 0.0 0.0
material delayedgroup group in nuclide mean std. dev.
0 10001 1 1 total 0.0 0.0
1 10001 2 1 total 0.0 0.0
@ -138,6 +149,13 @@
3 10001 4 1 total 0.0 0.0
4 10001 5 1 total 0.0 0.0
5 10001 6 1 total 0.0 0.0
material delayedgroup group in group out nuclide mean std. dev.
0 10001 1 1 1 total 0.0 0.0
1 10001 2 1 1 total 0.0 0.0
2 10001 3 1 1 total 0.0 0.0
3 10001 4 1 1 total 0.0 0.0
4 10001 5 1 1 total 0.0 0.0
5 10001 6 1 1 total 0.0 0.0
material group in nuclide mean std. dev.
0 10002 1 total 0.904999 0.043964
material group in nuclide mean std. dev.
@ -180,6 +198,8 @@
0 10002 1 total 5.773006e-07 5.322132e-08
material group in nuclide mean std. dev.
0 10002 1 total 0.0 0.0
material group in group out nuclide mean std. dev.
0 10002 1 1 total 0.0 0.0
material delayedgroup group in nuclide mean std. dev.
0 10002 1 1 total 0.0 0.0
1 10002 2 1 total 0.0 0.0
@ -208,3 +228,10 @@
3 10002 4 1 total 0.0 0.0
4 10002 5 1 total 0.0 0.0
5 10002 6 1 total 0.0 0.0
material delayedgroup group in group out nuclide mean std. dev.
0 10002 1 1 1 total 0.0 0.0
1 10002 2 1 1 total 0.0 0.0
2 10002 3 1 1 total 0.0 0.0
3 10002 4 1 1 total 0.0 0.0
4 10002 5 1 1 total 0.0 0.0
5 10002 6 1 1 total 0.0 0.0

2
tests/test_mgxs_library_distribcell/inputs_true.dat
View file

@ -1 +1 @@
67ee414eed54f596464831e734ac365e43b88bf03297e2b08adfb97510e1d8e631ff45061d5b99aa0ad85ec08d5b51e341866e60f35cdbfae078030bcde6c794
df187239f7481867cc09138709da90bc28eeb647b4bcbb08b894e7fdf6ff45510341b77e34754c2358f0e0665f8e552bb1eafed8c42cd2c50595b60366320ce4

9
tests/test_mgxs_library_distribcell/results_true.dat
View file

@ -40,6 +40,8 @@
0 (0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13,... 1 total 5.139437e-07 2.133314e-08
avg(distribcell) group in nuclide mean std. dev.
0 (0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13,... 1 total 0.091725 0.003604
avg(distribcell) group in group out nuclide mean std. dev.
0 (0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13,... 1 1 total 0.093985 0.005872
avg(distribcell) delayedgroup group in nuclide mean std. dev.
0 (0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13,... 1 1 total 0.000021 8.253907e-07
1 (0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13,... 2 1 total 0.000112 4.284000e-06
@ -68,3 +70,10 @@
3 (0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13,... 4 1 total 0.000000 0.000000
4 (0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13,... 5 1 total 0.000000 0.000000
5 (0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13,... 6 1 total 2.853000 4.034751
avg(distribcell) delayedgroup group in group out nuclide mean std. dev.
0 (0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13,... 1 1 1 total 0.000000 0.000000
1 (0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13,... 2 1 1 total 0.000175 0.000175
2 (0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13,... 3 1 1 total 0.000178 0.000178
3 (0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13,... 4 1 1 total 0.000000 0.000000
4 (0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13,... 5 1 1 total 0.000000 0.000000
5 (0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13,... 6 1 1 total 0.000178 0.000178

2
tests/test_mgxs_library_hdf5/inputs_true.dat
View file

@ -1 +1 @@
e86f24e20f37096c7898f459fc5bf336f3e0670fb30f321443f2bb3a01a154ef3d234bb48d4cee8e0d3730fd182b6a2051ec1b4c09d771420886a56ba928fdd9
4291b40470e7d59383c9e51c4178ca923b698cb1aaea16c1982fe3789ca980df10a65b84fb5021dacd4d290ebc232c2579f99b6c990fb6b8f28f67eef2aabcb2

120
tests/test_mgxs_library_hdf5/results_true.dat
View file

@ -72,6 +72,11 @@ domain=10000 type=inverse-velocity
domain=10000 type=prompt-nu-fission
[1.92392215e-02 4.66719027e-01]
[1.30950595e-03 4.14108704e-02]
domain=10000 type=prompt-nu-fission matrix
[[2.01424282e-02 0.00000000e+00]
[4.45819177e-01 0.00000000e+00]]
[[3.14909168e-03 0.00000000e+00]
[2.86750787e-02 0.00000000e+00]]
domain=10000 type=delayed-nu-fission
[[2.29808234e-05 1.06974158e-04]
[1.43606337e-04 5.52167907e-04]
@ -124,6 +129,41 @@ domain=10000 type=decay-rate
[0.00000000e+00 1.09109511e-01]
[0.00000000e+00 0.00000000e+00]
[0.00000000e+00 0.00000000e+00]]
domain=10000 type=delayed-nu-fission matrix
[[[0.00000000e+00 0.00000000e+00]
[0.00000000e+00 0.00000000e+00]]
[[0.00000000e+00 0.00000000e+00]
[2.53814542e-03 0.00000000e+00]]
[[0.00000000e+00 0.00000000e+00]
[1.18579166e-03 0.00000000e+00]]
[[0.00000000e+00 0.00000000e+00]
[8.59787018e-04 0.00000000e+00]]
[[0.00000000e+00 0.00000000e+00]
[0.00000000e+00 0.00000000e+00]]
[[0.00000000e+00 0.00000000e+00]
[0.00000000e+00 0.00000000e+00]]]
[[[0.00000000e+00 0.00000000e+00]
[0.00000000e+00 0.00000000e+00]]
[[0.00000000e+00 0.00000000e+00]
[1.56094584e-03 0.00000000e+00]]
[[0.00000000e+00 0.00000000e+00]
[1.18610401e-03 0.00000000e+00]]
[[0.00000000e+00 0.00000000e+00]
[2.22194634e-04 0.00000000e+00]]
[[0.00000000e+00 0.00000000e+00]
[0.00000000e+00 0.00000000e+00]]
[[0.00000000e+00 0.00000000e+00]
[0.00000000e+00 0.00000000e+00]]]
domain=10001 type=total
[3.13737671e-01 3.00821402e-01]
[1.55819024e-02 2.80524484e-02]
@ -198,6 +238,11 @@ domain=10001 type=inverse-velocity
domain=10001 type=prompt-nu-fission
[0.00000000e+00 0.00000000e+00]
[0.00000000e+00 0.00000000e+00]
domain=10001 type=prompt-nu-fission matrix
[[0.00000000e+00 0.00000000e+00]
[0.00000000e+00 0.00000000e+00]]
[[0.00000000e+00 0.00000000e+00]
[0.00000000e+00 0.00000000e+00]]
domain=10001 type=delayed-nu-fission
[[0.00000000e+00 0.00000000e+00]
[0.00000000e+00 0.00000000e+00]
@ -250,6 +295,41 @@ domain=10001 type=decay-rate
[0.00000000e+00 0.00000000e+00]
[0.00000000e+00 0.00000000e+00]
[0.00000000e+00 0.00000000e+00]]
domain=10001 type=delayed-nu-fission matrix
[[[0.00000000e+00 0.00000000e+00]
[0.00000000e+00 0.00000000e+00]]
[[0.00000000e+00 0.00000000e+00]
[0.00000000e+00 0.00000000e+00]]
[[0.00000000e+00 0.00000000e+00]
[0.00000000e+00 0.00000000e+00]]
[[0.00000000e+00 0.00000000e+00]
[0.00000000e+00 0.00000000e+00]]
[[0.00000000e+00 0.00000000e+00]
[0.00000000e+00 0.00000000e+00]]
[[0.00000000e+00 0.00000000e+00]
[0.00000000e+00 0.00000000e+00]]]
[[[0.00000000e+00 0.00000000e+00]
[0.00000000e+00 0.00000000e+00]]
[[0.00000000e+00 0.00000000e+00]
[0.00000000e+00 0.00000000e+00]]
[[0.00000000e+00 0.00000000e+00]
[0.00000000e+00 0.00000000e+00]]
[[0.00000000e+00 0.00000000e+00]
[0.00000000e+00 0.00000000e+00]]
[[0.00000000e+00 0.00000000e+00]
[0.00000000e+00 0.00000000e+00]]
[[0.00000000e+00 0.00000000e+00]
[0.00000000e+00 0.00000000e+00]]]
domain=10002 type=total
[6.64572261e-01 2.05238401e+00]
[3.12147519e-02 2.24342907e-01]
@ -324,6 +404,11 @@ domain=10002 type=inverse-velocity
domain=10002 type=prompt-nu-fission
[0.00000000e+00 0.00000000e+00]
[0.00000000e+00 0.00000000e+00]
domain=10002 type=prompt-nu-fission matrix
[[0.00000000e+00 0.00000000e+00]
[0.00000000e+00 0.00000000e+00]]
[[0.00000000e+00 0.00000000e+00]
[0.00000000e+00 0.00000000e+00]]
domain=10002 type=delayed-nu-fission
[[0.00000000e+00 0.00000000e+00]
[0.00000000e+00 0.00000000e+00]
@ -376,3 +461,38 @@ domain=10002 type=decay-rate
[0.00000000e+00 0.00000000e+00]
[0.00000000e+00 0.00000000e+00]
[0.00000000e+00 0.00000000e+00]]
domain=10002 type=delayed-nu-fission matrix
[[[0.00000000e+00 0.00000000e+00]
[0.00000000e+00 0.00000000e+00]]
[[0.00000000e+00 0.00000000e+00]
[0.00000000e+00 0.00000000e+00]]
[[0.00000000e+00 0.00000000e+00]
[0.00000000e+00 0.00000000e+00]]
[[0.00000000e+00 0.00000000e+00]
[0.00000000e+00 0.00000000e+00]]
[[0.00000000e+00 0.00000000e+00]
[0.00000000e+00 0.00000000e+00]]
[[0.00000000e+00 0.00000000e+00]
[0.00000000e+00 0.00000000e+00]]]
[[[0.00000000e+00 0.00000000e+00]
[0.00000000e+00 0.00000000e+00]]
[[0.00000000e+00 0.00000000e+00]
[0.00000000e+00 0.00000000e+00]]
[[0.00000000e+00 0.00000000e+00]
[0.00000000e+00 0.00000000e+00]]
[[0.00000000e+00 0.00000000e+00]
[0.00000000e+00 0.00000000e+00]]
[[0.00000000e+00 0.00000000e+00]
[0.00000000e+00 0.00000000e+00]]
[[0.00000000e+00 0.00000000e+00]
[0.00000000e+00 0.00000000e+00]]]

2
tests/test_mgxs_library_mesh/inputs_true.dat
View file

@ -1 +1 @@
a0d62fc011ae33756cd87202432f5c49f847793a40ec50efc6e4366755a70196f9ae970b24244766c011e2991961ea92317b53f6839b89363262b4e67ed4b289
03d894a7995ac40f7971b17349b460f4b563cb1c94abaa7e7241e6f7edad7f0cb80d3f80829db14ea6b7d70d18197cb15e308047cd10f699d834178eeb6396be

32
tests/test_mgxs_library_mesh/results_true.dat
View file

@ -130,6 +130,12 @@
1 1 2 1 1 total 0.020397 0.008086
2 2 1 1 1 total 0.025824 0.003192
3 2 2 1 1 total 0.020865 0.004879
mesh 1 group in group out nuclide mean std. dev.
x y z
0 1 1 1 1 1 total 0.020874 0.002977
1 1 2 1 1 1 total 0.017348 0.008786
2 2 1 1 1 1 total 0.020409 0.003354
3 2 2 1 1 1 total 0.011105 0.003806
mesh 1 delayedgroup group in nuclide mean std. dev.
x y z
0 1 1 1 1 1 total 0.000005 1.004627e-06
@ -234,3 +240,29 @@
21 2 2 1 4 1 total 0.00000 0.00000
22 2 2 1 5 1 total 0.00000 0.00000
23 2 2 1 6 1 total 0.00000 0.00000
mesh 1 delayedgroup group in group out nuclide mean std. dev.
x y z
0 1 1 1 1 1 1 total 0.000000 0.000000
1 1 1 1 2 1 1 total 0.000000 0.000000
2 1 1 1 3 1 1 total 0.000000 0.000000
3 1 1 1 4 1 1 total 0.000000 0.000000
4 1 1 1 5 1 1 total 0.000185 0.000186
5 1 1 1 6 1 1 total 0.000000 0.000000
6 1 2 1 1 1 1 total 0.000000 0.000000
7 1 2 1 2 1 1 total 0.000000 0.000000
8 1 2 1 3 1 1 total 0.000000 0.000000
9 1 2 1 4 1 1 total 0.000000 0.000000
10 1 2 1 5 1 1 total 0.000000 0.000000
11 1 2 1 6 1 1 total 0.000000 0.000000
12 2 1 1 1 1 1 total 0.000000 0.000000
13 2 1 1 2 1 1 total 0.000000 0.000000
14 2 1 1 3 1 1 total 0.000000 0.000000
15 2 1 1 4 1 1 total 0.000000 0.000000
16 2 1 1 5 1 1 total 0.000000 0.000000
17 2 1 1 6 1 1 total 0.000000 0.000000
18 2 2 1 1 1 1 total 0.000000 0.000000
19 2 2 1 2 1 1 total 0.000000 0.000000
20 2 2 1 3 1 1 total 0.000000 0.000000
21 2 2 1 4 1 1 total 0.000000 0.000000
22 2 2 1 5 1 1 total 0.000000 0.000000
23 2 2 1 6 1 1 total 0.000000 0.000000

2
tests/test_mgxs_library_no_nuclides/inputs_true.dat
View file

@ -1 +1 @@
e86f24e20f37096c7898f459fc5bf336f3e0670fb30f321443f2bb3a01a154ef3d234bb48d4cee8e0d3730fd182b6a2051ec1b4c09d771420886a56ba928fdd9
4291b40470e7d59383c9e51c4178ca923b698cb1aaea16c1982fe3789ca980df10a65b84fb5021dacd4d290ebc232c2579f99b6c990fb6b8f28f67eef2aabcb2

90
tests/test_mgxs_library_no_nuclides/results_true.dat
View file

@ -84,6 +84,11 @@
material group in nuclide mean std. dev.
1 10000 1 total 0.019239 0.001310
0 10000 2 total 0.466719 0.041411
material group in group out nuclide mean std. dev.
3 10000 1 1 total 0.020142 0.003149
2 10000 1 2 total 0.000000 0.000000
1 10000 2 1 total 0.445819 0.028675
0 10000 2 2 total 0.000000 0.000000
material delayedgroup group in nuclide mean std. dev.
1 10000 1 1 total 0.000023 0.000002
3 10000 2 1 total 0.000144 0.000011
@ -136,6 +141,31 @@
6 10000 4 2 total 0.302780 0.109110
8 10000 5 2 total 0.000000 0.000000
10 10000 6 2 total 0.000000 0.000000
material delayedgroup group in group out nuclide mean std. dev.
3 10000 1 1 1 total 0.000000 0.000000
7 10000 2 1 1 total 0.000000 0.000000
11 10000 3 1 1 total 0.000000 0.000000
15 10000 4 1 1 total 0.000000 0.000000
19 10000 5 1 1 total 0.000000 0.000000
23 10000 6 1 1 total 0.000000 0.000000
2 10000 1 1 2 total 0.000000 0.000000
6 10000 2 1 2 total 0.000000 0.000000
10 10000 3 1 2 total 0.000000 0.000000
14 10000 4 1 2 total 0.000000 0.000000
18 10000 5 1 2 total 0.000000 0.000000
22 10000 6 1 2 total 0.000000 0.000000
1 10000 1 2 1 total 0.000000 0.000000
5 10000 2 2 1 total 0.002538 0.001561
9 10000 3 2 1 total 0.001186 0.001186
13 10000 4 2 1 total 0.000860 0.000222
17 10000 5 2 1 total 0.000000 0.000000
21 10000 6 2 1 total 0.000000 0.000000
0 10000 1 2 2 total 0.000000 0.000000
4 10000 2 2 2 total 0.000000 0.000000
8 10000 3 2 2 total 0.000000 0.000000
12 10000 4 2 2 total 0.000000 0.000000
16 10000 5 2 2 total 0.000000 0.000000
20 10000 6 2 2 total 0.000000 0.000000
material group in nuclide mean std. dev.
1 10001 1 total 0.313738 0.015582
0 10001 2 total 0.300821 0.028052
@ -222,6 +252,11 @@
material group in nuclide mean std. dev.
1 10001 1 total 0.0 0.0
0 10001 2 total 0.0 0.0
material group in group out nuclide mean std. dev.
3 10001 1 1 total 0.0 0.0
2 10001 1 2 total 0.0 0.0
1 10001 2 1 total 0.0 0.0
0 10001 2 2 total 0.0 0.0
material delayedgroup group in nuclide mean std. dev.
1 10001 1 1 total 0.0 0.0
3 10001 2 1 total 0.0 0.0
@ -274,6 +309,31 @@
6 10001 4 2 total 0.0 0.0
8 10001 5 2 total 0.0 0.0
10 10001 6 2 total 0.0 0.0
material delayedgroup group in group out nuclide mean std. dev.
3 10001 1 1 1 total 0.0 0.0
7 10001 2 1 1 total 0.0 0.0
11 10001 3 1 1 total 0.0 0.0
15 10001 4 1 1 total 0.0 0.0
19 10001 5 1 1 total 0.0 0.0
23 10001 6 1 1 total 0.0 0.0
2 10001 1 1 2 total 0.0 0.0
6 10001 2 1 2 total 0.0 0.0
10 10001 3 1 2 total 0.0 0.0
14 10001 4 1 2 total 0.0 0.0
18 10001 5 1 2 total 0.0 0.0
22 10001 6 1 2 total 0.0 0.0
1 10001 1 2 1 total 0.0 0.0
5 10001 2 2 1 total 0.0 0.0
9 10001 3 2 1 total 0.0 0.0
13 10001 4 2 1 total 0.0 0.0
17 10001 5 2 1 total 0.0 0.0
21 10001 6 2 1 total 0.0 0.0
0 10001 1 2 2 total 0.0 0.0
4 10001 2 2 2 total 0.0 0.0
8 10001 3 2 2 total 0.0 0.0
12 10001 4 2 2 total 0.0 0.0
16 10001 5 2 2 total 0.0 0.0
20 10001 6 2 2 total 0.0 0.0
material group in nuclide mean std. dev.
1 10002 1 total 0.664572 0.031215
0 10002 2 total 2.052384 0.224343
@ -360,6 +420,11 @@
material group in nuclide mean std. dev.
1 10002 1 total 0.0 0.0
0 10002 2 total 0.0 0.0
material group in group out nuclide mean std. dev.
3 10002 1 1 total 0.0 0.0
2 10002 1 2 total 0.0 0.0
1 10002 2 1 total 0.0 0.0
0 10002 2 2 total 0.0 0.0
material delayedgroup group in nuclide mean std. dev.
1 10002 1 1 total 0.0 0.0
3 10002 2 1 total 0.0 0.0
@ -412,3 +477,28 @@
6 10002 4 2 total 0.0 0.0
8 10002 5 2 total 0.0 0.0
10 10002 6 2 total 0.0 0.0
material delayedgroup group in group out nuclide mean std. dev.
3 10002 1 1 1 total 0.0 0.0
7 10002 2 1 1 total 0.0 0.0
11 10002 3 1 1 total 0.0 0.0
15 10002 4 1 1 total 0.0 0.0
19 10002 5 1 1 total 0.0 0.0
23 10002 6 1 1 total 0.0 0.0
2 10002 1 1 2 total 0.0 0.0
6 10002 2 1 2 total 0.0 0.0
10 10002 3 1 2 total 0.0 0.0
14 10002 4 1 2 total 0.0 0.0
18 10002 5 1 2 total 0.0 0.0
22 10002 6 1 2 total 0.0 0.0
1 10002 1 2 1 total 0.0 0.0
5 10002 2 2 1 total 0.0 0.0
9 10002 3 2 1 total 0.0 0.0
13 10002 4 2 1 total 0.0 0.0
17 10002 5 2 1 total 0.0 0.0
21 10002 6 2 1 total 0.0 0.0
0 10002 1 2 2 total 0.0 0.0
4 10002 2 2 2 total 0.0 0.0
8 10002 3 2 2 total 0.0 0.0
12 10002 4 2 2 total 0.0 0.0
16 10002 5 2 2 total 0.0 0.0
20 10002 6 2 2 total 0.0 0.0

2
tests/test_mgxs_library_nuclides/inputs_true.dat
View file

@ -1 +1 @@
1e8d0f408ba9d47fc9278e750dbf68d3104d248f7783ca1185c8b7666c9ee2ba3ec6fcf4a49bb614956bb8eb3a57dc1c104725b4160171dfbe6dc972a268df56
eed0190893105747f0146472dfea86ca58b9fc2e8d039961d9fef8f235f79da32215f9b59b840419e6b69c5a04544b4659f2634b63e3a76b6b0f8c75d05e416c

2
tests/test_mgxs_library_nuclides/results_true.dat
View file

@ -1 +1 @@
ce682f577fd65dd9b6e5da58763cf2ae4cf8e5041b38b5d1ed64d3ba7a3a4df7338aadfa8830a0764fc9ffb737f05fc989a494158c714047186dd18605a7d2b8
8bc6694ee99cc05ec59143dac065f3f51026713044c48a36ac486f46014289c1f18d4d8e14ef0bb468050aea0bac5628425fd76154a01c7438c7d85139aa2f06