crp/OpenMC
1
0
Fork 0
mirror of https://github.com/openmc-dev/openmc.git synced 2026-07-28 22:26:08 -04:00
Code Issues Projects Releases Packages Wiki Activity Actions 5

Merge remote-tracking branch 'shaner/tally-align' into sparse-tallies

Browse source
This commit is contained in:
wbinventor@gmail.com 2015-12-12 23:45:34 -05:00
commit e613658930
16 changed files with 637 additions and 215 deletions
Download patch file Download diff file Expand all files Collapse all files

8
CMakeLists.txt
View file

@ -235,6 +235,14 @@ if(NOT EXISTS ${CMAKE_CURRENT_SOURCE_DIR}/src/xml/fox/.git)
endif()
add_subdirectory(src/xml/fox)
#===============================================================================
# RPATH information
#===============================================================================
# add the automatically determined parts of the RPATH
# which point to directories outside the build tree to the install RPATH
set(CMAKE_INSTALL_RPATH_USE_LINK_PATH TRUE)
#===============================================================================
# Build OpenMC executable
#===============================================================================

2
docs/source/conf.py
View file

@ -55,7 +55,7 @@ copyright = u'2011-2015, Massachusetts Institute of Technology'
# The short X.Y version.
version = "0.7"
# The full version, including alpha/beta/rc tags.
release = "0.7.0"
release = "0.7.1"
# The language for content autogenerated by Sphinx. Refer to documentation
# for a list of supported languages.

63
docs/source/publications.rst
View file

@ -26,6 +26,10 @@ Overviews
Benchmarking
------------
- Khurrum S. Chaudri and Sikander M. Mirza, "Burnup dependent Monte Carlo
neutron physics calculations of IAEA MTR benchmark," *Prog. Nucl. Energy*,
**81**, 43-52 (2015). `<http://dx.doi.org/j.pnucene.2014.12.018>`_
- Daniel J. Kelly, Brian N. Aviles, Paul K. Romano, Bryan R. Herman,
Nicholas E. Horelik, and Benoit Forget, "Analysis of select BEAVRS PWR
benchmark cycle 1 results using MC21 and OpenMC," *Proc. PHYSOR*, Kyoto,
@ -57,13 +61,8 @@ Coupling and Multi-physics
- Bryan R. Herman, Benoit Forget, and Kord Smith, "Progress toward Monte
Carlo-thermal hydraulic coupling using low-order nonlinear diffusion
acceleration methods." In press, *Ann. Nucl. Energy*,
(2014). `<http://dx.doi.org/10.1016/j.anucene/2014.10.029>`_
- Adam G. Nelson and William R. Martin, "Improved Convergence of Monte Carlo
Generated Multi-Group Scattering Moments," *Proc. Int. Conf. Mathematics and
Computational Methods Applied to Nuclear Science and Engineering*, Sun Valley,
Idaho, May 5--9 (2013).
acceleration methods." *Ann. Nucl. Energy*, **84**, 63-72
(2015). `<http://dx.doi.org/10.1016/j.anucene.2014.10.029>`_
- Bryan R. Herman, Benoit Forget, and Kord Smith, "Utilizing CMFD in OpenMC to
Estimate Dominance Ratio and Adjoint," *Trans. Am. Nucl. Soc.*, **109**,
@ -81,19 +80,65 @@ Geometry
Miscellaneous
-------------
- William Boyd, Sterling Harper, and Paul K. Romano, "Equipping OpenMC for the
big data era," Accepted, *PHYSOR 2016*, Sun Valley, Idaho, May 1-5, 2016.
- Qicang Shen, William Boyd, Benoit Forget, and Kord Smith, "Tally precision
triggers for the OpenMC Monte Carlo code," *Trans. Am. Nucl. Soc.*, **112**,
637-640 (2015).
- Timothy P. Burke, Brian C. Kiedrowski, and William R. Martin, "Flux and
Reaction Rate Kernel Density Estimators in OpenMC," *Trans. Am. Nucl. Soc.*,
**109**, 683-686 (2013).
------------------------------------
Multi-group Cross Section Generation
------------------------------------
- Adam G. Nelson and William R. Martin, "Improved Monte Carlo tallying of
multi-group scattering moments using the NDPP code," *Trans. Am. Nucl. Soc.*,
**113**, 645-648 (2015)
- Adam G. Nelson and William R. Martin, "Improved Monte Carlo tallying of
multi-group scattering moment matrices," *Trans. Am. Nucl. Soc.*, **110**,
217-220 (2014).
- Adam G. Nelson and William R. Martin, "Improved Convergence of Monte Carlo
Generated Multi-Group Scattering Moments," *Proc. Int. Conf. Mathematics and
Computational Methods Applied to Nuclear Science and Engineering*, Sun Valley,
Idaho, May 5--9 (2013).
------------
Nuclear Data
------------
- Colin Josey, Pablo Ducru, Benoit Forget, and Kord Smith, "Windowed multipole
for cross section Doppler broadening," *J. Comput. Phys.*, In Press
(2016). `<http://dx.doi.org/10.1016/jcp.2015.08.013>`_
- Colin Josey, Benoit Forget, and Kord Smith, "Windowed multipole sensitivity to
target accuracy of the optimization procedure," *J. Nucl. Sci. Technol.*,
**52**, 987-992 (2015). `<http://dx.doi.org/10.1080/00223131.2015.1035353>`_
- Jonathan A. Walsh, Paul K. Romano, Benoit Forget, and Kord S. Smith,
"Optimizations of the energy grid search algorithm in continuous-energy Monte
Carlo particle transport codes", *Comput. Phys. Commun.*, **196**, 134-142
(2015). `<http://dx.doi.org/10.1016/j.cpc.2015.05.025>`_
- Jonathan A. Walsh, Benoit Forget, Kord S. Smith, Brian C. Kiedrowski, and
Forrest B. Brown, "Direct, on-the-fly calculation of unresolved resonance
region cross sections in Monte Carlo simulations," *Proc. Joint
Int. Conf. M&C+SNA+MC*, Nashville, Tennessee, Apr. 19--23 (2015).
- Amanda L. Lund, Andrew R. Siegel, Benoit Forget, Colin Josey, and
Paul K. Romano, "Using fractional cascading to accelerate cross section
lookups in Monte Carlo particle transport calculations," *Proc. Joint
Int. Conf. M&C+SNA+MC*, Nashville, Tennessee, Apr. 19--23 (2015).
- Ronald O. Rahaman, Andrew R. Siegel, and Paul K. Romano, "Monte Carlo
performance analysis for varying cross section parameter regimes,"
*Proc. Joint Int. Conf. M&C+SNA+MC*, Nashville, Tennessee, Apr. 19--23 (2015).
- Paul K. Romano and Timothy H. Trumbull, "Comparison of algorithms for Doppler
broadening pointwise tabulated cross sections," *Ann. Nucl. Energy*, **75**,
358--364 (2015). `<http://dx.doi.org/10.1016/j.anucene.2014.08.046>`_
@ -114,6 +159,10 @@ Nuclear Data
Parallelism
-----------
- Paul K. Romano, John R. Tramm, and Andrew R. Siegel, "Efficacy of hardware
threading for Monte Carlo particle transport calculations on multi- and
many-core systems," Accepted, *PHYSOR 2016*, Sun Valley, Idaho, May 1-5, 2016.
- David Ozog, Allen D. Malony, and Andrew R. Siegel, "A performance analysis of
SIMD algorithms for Monte Carlo simulations of nuclear reactor cores,"
*Proc. IEEE Int. Parallel and Distributed Processing Symposium*, Hyderabad,

83
docs/source/releasenotes.rst
View file

@ -1,9 +1,30 @@
.. _releasenotes:
==============================
Release Notes for OpenMC 0.7.0
Release Notes for OpenMC 0.7.1
==============================
This release of OpenMC provides some substantial improvements over version
0.7.0. Non-simple cell regions can now be defined through the ``|`` (union) and
``~`` (complement) operators. Similar changes in the Python API also allow
complex cell regions to be defined. A true secondary particle bank now exists;
this is crucial for photon transport (to be added in the next minor release). A
rich API for multi-group cross section generation has been added via the
``openmc.mgxs`` Python module.
Various improvements to tallies have also been made. It is now possible to
explicitly specify that a collision estimator be used in a tally. A new
``delayedgroup`` filter and ``delayed-nu-fission`` score allow a user to obtain
delayed fission neutron production rates filtered by delayed group. Finally, the
new ``inverse-velocity`` score may be useful for calculating kinetics
parameters.
.. caution:: In previous versions, depending on how OpenMC was compiled binary
output was either given in HDF5 or a flat binary format. With this
version, all binary output is now HDF5 which means you **must**
have HDF5 in order to install OpenMC. Please consult the user's
guide for instructions on how to compile with HDF5.
-------------------
System Requirements
-------------------
@ -17,36 +38,41 @@ the problem at hand (mostly on the number of nuclides in the problem).
New Features
------------
- Complete Python API
- Python 3 compatability for all scripts
- All scripts consistently named openmc-* and installed together
- New 'distribcell' tally filter for repeated cells
- Ability to specify outer lattice universe
- XML input validation utility (openmc-validate-xml)
- Support for hexagonal lattices
- Material union energy grid method
- Tally triggers
- Remove dependence on PETSc
- Significant OpenMP performance improvements
- Support for Fortran 2008 MPI interface
- Use of Travis CI for continuous integration
- Simplifications and improvements to test suite
- Support for complex cell regions (union and complement operators)
- Generic quadric surface type
- Improved handling of secondary particles
- Binary output is now solely HDF5
- ``openmc.mgxs`` Python module enabling multi-group cross section generation
- Collision estimator for tallies
- Delayed fission neutron production tallies with ability to filter by delayed
group
- Inverse velocity tally score
- Performance improvements for binary search
- Performance improvements for reaction rate tallies
---------
Bug Fixes
---------
- b5f712_: Fix bug in spherical harmonics tallies
- e6675b_: Ensure all constants are double precision
- 04e2c1_: Fix potential bug in sample_nuclide routine
- 6121d9_: Fix bugs related to particle track files
- 2f0e89_: Fixes for nuclide specification in tallies
- 299322_: Bug with material filter when void material present
- d74840_: Fix triggers on tallies with multiple filters
- c29a81_: Correctly handle maximum transport energy
- 3edc23_: Fixes in the nu-scatter score
- 629e3b_: Assume unspecified surface coefficients are zero in Python API
- 5dbe8b_: Fix energy filters for openmc-plot-mesh-tally
- ff66f4_: Fixes in the openmc-plot-mesh-tally script
- 441fd4_: Fix bug in kappa-fission score
- 7e5974_: Allow fixed source simulations from Python API
.. _b5f712: https://github.com/mit-crpg/openmc/commit/b5f712
.. _e6675b: https://github.com/mit-crpg/openmc/commit/e6675b
.. _04e2c1: https://github.com/mit-crpg/openmc/commit/04e2c1
.. _6121d9: https://github.com/mit-crpg/openmc/commit/6121d9
.. _2f0e89: https://github.com/mit-crpg/openmc/commit/2f0e89
.. _299322: https://github.com/mit-crpg/openmc/commit/299322
.. _d74840: https://github.com/mit-crpg/openmc/commit/d74840
.. _c29a81: https://github.com/mit-crpg/openmc/commit/c29a81
.. _3edc23: https://github.com/mit-crpg/openmc/commit/3edc23
.. _629e3b: https://github.com/mit-crpg/openmc/commit/629e3b
.. _5dbe8b: https://github.com/mit-crpg/openmc/commit/5dbe8b
.. _ff66f4: https://github.com/mit-crpg/openmc/commit/ff66f4
.. _441fd4: https://github.com/mit-crpg/openmc/commit/441fd4
.. _7e5974: https://github.com/mit-crpg/openmc/commit/7e5974
------------
Contributors
@ -55,13 +81,10 @@ Contributors
This release contains new contributions from the following people:
- `Will Boyd <wbinventor@gmail.com>`_
- `Matt Ellis <mellis13@mit.edu>`_
- `Sterling Harper <sterlingmharper@mit.edu>`_
- `Bryan Herman <bherman@mit.edu>`_
- `Nicholas Horelik <nicholas.horelik@gmail.com>`_
- `Colin Josey <cjosey@mit.edu>`_
- `William Lyu <PaleNeutron@users.noreply.github.com>`_
- `Adam Nelson <nelsonag@umich.edu>`_
- `Paul Romano <paul.k.romano@gmail.com>`_
- `Anthony Scopatz <scopatz@gmail.com>`_
- `Kelly Rowland <kellylynnerowland@gmail.com>`_
- `Sam Shaner <samuelshaner@gmail.com>`_
- `Jon Walsh <walshjon@mit.edu>`_

11
docs/source/usersguide/install.rst
View file

@ -8,7 +8,7 @@ Installation and Configuration
Installing on Ubuntu with PPA
-----------------------------
For users with Ubuntu 11.10 or later, a binary package for OpenMC is available
For users with Ubuntu 15.04 or later, a binary package for OpenMC is available
through a Personal Package Archive (PPA) and can be installed through the APT
package manager. First, add the following PPA to the repository sources:
@ -28,6 +28,9 @@ Now OpenMC should be recognized within the repository and can be installed:
sudo apt-get install openmc
Binary packages from this PPA may exist for earlier versions of Ubuntu, but they
are no longer supported.
--------------------
Building from Source
--------------------
@ -74,6 +77,12 @@ Prerequisites
You may omit ``--enable-parallel`` if you want to compile HDF5_ in serial.
.. important::
OpenMC uses various parts of the HDF5 Fortran 2003 API; as such you
must include ``--enable-fortran2003`` or else OpenMC will not be able
to compile.
On Debian derivatives, HDF5 and/or parallel HDF5 can be installed through
the APT package manager:

2
openmc/clean_xml.py
View file

@ -1,7 +1,7 @@
def sort_xml_elements(tree):
# Retrieve all children of the root XML node in the tree
elements = tree.getchildren()
elements = list(tree)
# Initialize empty lists for the sorted and comment elements
sorted_elements = []

12
openmc/mgxs/mgxs.py
View file

@ -47,7 +47,7 @@ _DOMAINS = [openmc.Cell,
class MGXS(object):
"""An abstract multi-group cross section for some energy group structure
"""An abstract multi-group cross section for some energy group structure
within some spatial domain.
This class can be used for both OpenMC input generation and tally data
@ -95,7 +95,7 @@ class MGXS(object):
num_sumbdomains : Integral
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
for 'distribcell' domain types (it is equal to unity prior to loading
tally data from a statepoint file).
num_nuclides : Integral
The number of nuclides for which the multi-group cross section is
@ -1420,8 +1420,8 @@ class AbsorptionXS(MGXS):
class CaptureXS(MGXS):
"""A capture multi-group cross section.
The Neutron capture reaction rate is defined as the difference between
The Neutron capture reaction rate is defined as the difference between
OpenMC's 'absorption' and 'fission' reaction rate score types. This includes
not only radiative capture, but all forms of neutron disappearance aside
from fission (e.g., MT > 100).
@ -1723,8 +1723,8 @@ class ScatterMatrixXS(MGXS):
if self.correction == 'P0':
scatter_p1 = self.tallies['scatter-P1']
scatter_p1 = scatter_p1.get_slice(scores=['scatter-P1'])
energy_filter = openmc.Filter(type='energy')
energy_filter.bins = self.energy_groups.group_edges
energy_filter = copy.deepcopy(self.tallies['scatter'].
find_filter('energy'))
scatter_p1 = scatter_p1.diagonalize_filter(energy_filter)
rxn_tally = self.tallies['scatter'] - scatter_p1
else:

503
openmc/tallies.py
View file

@ -25,6 +25,7 @@ if sys.version_info[0] >= 3:
# "Static" variable for auto-generated Tally IDs
AUTO_TALLY_ID = 10000
_PRODUCT_TYPES = ['tensor', 'entrywise']
def reset_auto_tally_id():
global AUTO_TALLY_ID
@ -1255,12 +1256,12 @@ class Tally(object):
# Tile the nuclide bins into a DataFrame column
nuclides = np.repeat(nuclides, len(self.scores))
tile_factor = data_size / len(nuclides)
df['nuclide'] = np.tile(nuclides, tile_factor)
df['nuclide'] = np.tile(nuclides, int(tile_factor))
# Include column for scores if user requested it
if scores:
tile_factor = data_size / len(self.scores)
df['score'] = np.tile(self.scores, tile_factor)
df['score'] = np.tile(self.scores, int(tile_factor))
# Append columns with mean, std. dev. for each tally bin
df['mean'] = self.mean.ravel()
@ -1480,25 +1481,47 @@ class Tally(object):
# Pickle the Tally results to a file
pickle.dump(tally_results, open(filename, 'wb'))
def _outer_product(self, other, binary_op):
def hybrid_product(self, other, binary_op, filter_product='None',
nuclide_product='None', score_product='None'):
"""Combines filters, scores and nuclides with another tally.
This is a helper method for the tally arithmetic methods. The filters,
scores and nuclides from both tallies are enumerated into all possible
combinations and expressed as CrossFilter, CrossScore and
CrossNuclide objects in the new derived tally.
This is a helper method for the tally arithmetic methods. It is called a
"hybrid product" because it performs a combination of tensor
(or Kronecker) and entrywise (or Hadamard) products. The filters,
nuclides, and scores from both tallies are combined using an entrywise
(or Hadamard) product on matching filters. By default, if all nuclides
are identical in the two tallies, the entrywise product is performed
across nuclides; else the tensor product is performed. By default, if all
scores are identical in the two tallies, the entrywise product is
performed across scores; else the tensor product is performed. Users can
also call the method explicitly and specify the desired product.
Parameters
----------
other : Tally
The tally on the right hand side of the outer product
The tally on the right hand side of the hybrid product
binary_op : {'+', '-', '*', '/', '^'}
The binary operation in the outer product
The binary operation in the hybrid product
filter_product : str, optional
The type of product (tensor or entrywise) to be performed between
filter data. The default is the entrywise product. Currently only
the entrywise product is supported since a tally cannot contain
two of the same tallies.
nuclide_product : str, optional
The type of product (tensor or entrywise) to be performed between
nuclide data. The default is the entrywise product if all nuclides
between the two tallies are the same; otherwise the default is
the tensor product.
score_product : str, optional
The type of product (tensor or entrywise) to be performed between
score data. The default is the entrywise product if all scores
between the two tallies are the same; otherwise the default is
the tensor product.
Returns
-------
Tally
A new Tally that is the outer product with this one.
A new Tally that is the hybrid product with this one.
Raises
------
@ -1508,6 +1531,33 @@ class Tally(object):
"""
# Set default value for filter product if it was not set
if filter_product == 'None':
filter_product = 'entrywise'
elif filter_product == 'tensor':
msg = 'Unable to perform Tally arithmetic with a tensor product' \
'for the filter data as this not currently supported.'
raise ValueError(msg)
# Set default value for nuclide product if it was not set
if nuclide_product == 'None':
if self.nuclides == other.nuclides:
nuclide_product = 'entrywise'
else:
nuclide_product = 'tensor'
# Set default value for score product if it was not set
if score_product == 'None':
if self.scores == other.scores:
score_product = 'entrywise'
else:
score_product = 'tensor'
# Check product types
cv.check_value('filter product', filter_product, _PRODUCT_TYPES)
cv.check_value('nuclide product', nuclide_product, _PRODUCT_TYPES)
cv.check_value('score product', score_product, _PRODUCT_TYPES)
# Check that results have been read
if not other.derived and other.sum is None:
msg = 'Unable to use tally arithmetic with Tally ID="{0}" ' \
@ -1523,32 +1573,23 @@ class Tally(object):
new_name = '({0} {1} {2})'.format(self.name, binary_op, other.name)
new_tally.name = new_name
# Query the mean and std dev so the tally data is read in from file
# if it has not already been read in.
self.mean
self.std_dev
other.mean
other.std_dev
# Create copies of self and other tallies to rearrange for tally
# arithmetic
self_copy = copy.deepcopy(self)
other_copy = copy.deepcopy(other)
# Find any shared filters between the two tallies
filter_intersect = []
for filter in self_copy.filters:
if filter in other_copy.filters:
filter_intersect.append(filter)
# Align the shared filters in successive order
for i, filter in enumerate(filter_intersect):
self_index = self_copy.filters.index(filter)
other_index = other_copy.filters.index(filter)
# If necessary, swap self filter
if self_index != i:
self_copy.swap_filters(filter, self_copy.filters[i], inplace=True)
# If necessary, swap other filter
if other_index != i:
other_copy.swap_filters(filter, other_copy.filters[i], inplace=True)
data = self_copy._align_tally_data(other_copy)
# Align the tally data based on desired hybrid product
data = self_copy._align_tally_data(other_copy, filter_product,
nuclide_product, score_product)
# Perform tally arithmetic operation
if binary_op == '+':
new_tally._mean = data['self']['mean'] + data['other']['mean']
new_tally._std_dev = np.sqrt(data['self']['std. dev.']**2 +
@ -1578,6 +1619,13 @@ class Tally(object):
new_tally._std_dev = np.abs(new_tally.mean) * \
np.sqrt(first_term**2 + second_term**2)
# Convert any infs and nans to zero
new_tally._mean[np.isinf(new_tally._mean)] = 0
new_tally._mean = np.nan_to_num(new_tally._mean)
new_tally._std_dev[np.isinf(new_tally._std_dev)] = 0
new_tally._std_dev = np.nan_to_num(new_tally._std_dev)
# Set tally attributes
if self_copy.estimator == other_copy.estimator:
new_tally.estimator = self_copy.estimator
if self_copy.with_summary and other_copy.with_summary:
@ -1585,43 +1633,28 @@ class Tally(object):
if self_copy.num_realizations == other_copy.num_realizations:
new_tally.num_realizations = self_copy.num_realizations
# If filters are identical, simply reuse them in derived tally
if self_copy.filters == other_copy.filters:
# Add filters to the new tally
if filter_product == 'entrywise':
for self_filter in self_copy.filters:
new_tally.add_filter(self_filter)
# Generate filter "outer products" for non-identical filters
new_tally.filters.append(self_filter)
else:
all_filters = [self_copy.filters, other_copy.filters]
for self_filter, other_filter in itertools.product(*all_filters):
new_filter = CrossFilter(self_filter, other_filter, binary_op)
new_tally.add_filter(new_filter)
# Find the common longest sequence of shared filters
match = 0
for self_filter, other_filter in zip(self_copy.filters, other_copy.filters):
if self_filter == other_filter:
match += 1
else:
break
# Add nuclides to the new tally
if nuclide_product == 'entrywise':
for self_nuclide in self_copy.nuclides:
new_tally.nuclides.append(self_nuclide)
else:
all_nuclides = [self_copy.nuclides, other_copy.nuclides]
for self_nuclide, other_nuclide in itertools.product(*all_nuclides):
new_nuclide = CrossNuclide(self_nuclide, other_nuclide, binary_op)
new_tally.add_nuclide(new_nuclide)
match_filters = self_copy.filters[:match]
cross_filters = [self_copy.filters[match:], other_copy.filters[match:]]
# Simply reuse shared filters in derived tally
for filter in match_filters:
new_tally.add_filter(filter)
# Use cross filters to combine non-shared filters in derived tally
if len(self_copy.filters) != match and len(other_copy.filters) == match:
for filter in cross_filters[0]:
new_tally.add_filter(filter)
elif len(self_copy.filters) == match and len(other_copy.filters) != match:
for filter in cross_filters[1]:
new_tally.add_filter(filter)
else:
for self_filter, other_filter in itertools.product(*cross_filters):
new_filter = CrossFilter(self_filter, other_filter, binary_op)
new_tally.add_filter(new_filter)
# Generate score "outer products"
if self_copy.scores == other_copy.scores:
# Add scores to the new tally
if score_product == 'entrywise':
new_tally.num_score_bins = self_copy.num_score_bins
for self_score in self_copy.scores:
new_tally.add_score(self_score)
@ -1632,16 +1665,6 @@ class Tally(object):
new_score = CrossScore(self_score, other_score, binary_op)
new_tally.add_score(new_score)
# Generate nuclide "outer products"
if self_copy.nuclides == other_copy.nuclides:
for self_nuclide in self_copy.nuclides:
new_tally.nuclides.append(self_nuclide)
else:
all_nuclides = [self_copy.nuclides, other_copy.nuclides]
for self_nuclide, other_nuclide in itertools.product(*all_nuclides):
new_nuclide = CrossNuclide(self_nuclide, other_nuclide, binary_op)
new_tally.add_nuclide(new_nuclide)
# Correct each Filter's stride
stride = new_tally.num_nuclides * new_tally.num_score_bins
for filter in reversed(new_tally.filters):
@ -1650,12 +1673,13 @@ class Tally(object):
return new_tally
def _align_tally_data(self, other):
def _align_tally_data(self, other, filter_product, nuclide_product,
score_product):
"""Aligns data from two tallies for tally arithmetic.
This is a helper method to construct a dict of dicts of the "aligned"
data arrays from each tally for tally arithmetic. The method analyzes
the filters, scores and nuclides in both tally's and determines how to
the filters, scores and nuclides in both tallies and determines how to
appropriately align the data for vectorized arithmetic. For example,
if the two tallies have different filters, this method will use NumPy
'tile' and 'repeat' operations to the new data arrays such that all
@ -1666,6 +1690,15 @@ class Tally(object):
----------
other : Tally
The tally to outer product with this tally
filter_product : str
The type of product (tensor or entrywise) to be performed between
filter data.
nuclide_product : str
The type of product (tensor or entrywise) to be performed between
nuclide data.
score_product : str
The type of product (tensor or entrywise) to be performed between
score data.
Returns
-------
@ -1675,98 +1708,124 @@ class Tally(object):
"""
# Get the set of filters that each tally is missing
other_missing_filters = set(self.filters).difference(set(other.filters))
self_missing_filters = set(other.filters).difference(set(self.filters))
# Add filters present in self but not in other to other
for filter in other_missing_filters:
filter = copy.deepcopy(filter)
other._mean = np.repeat(other.mean, filter.num_bins, axis=0)
other._std_dev = np.repeat(other.std_dev, filter.num_bins, axis=0)
other.add_filter(filter)
# Add filters present in other but not in self to self
for filter in self_missing_filters:
filter = copy.deepcopy(filter)
self._mean = np.repeat(self.mean, filter.num_bins, axis=0)
self._std_dev = np.repeat(self.std_dev, filter.num_bins, axis=0)
self.add_filter(filter)
# Align other filters with self filters
for i, filter in enumerate(self.filters):
other_index = other.filters.index(filter)
# If necessary, swap other filter
if other_index != i:
other.swap_filters(filter, other.filters[i], inplace=True)
# Repeat and tile the data by nuclide in preparation for performing
# the tensor product across nuclides.
if nuclide_product == 'tensor':
self._mean = np.repeat(self.mean, other.num_nuclides, axis=1)
self._std_dev = np.repeat(self.std_dev, other.num_nuclides, axis=1)
other._mean = np.tile(other.mean, (1, self.num_nuclides, 1))
other._std_dev = np.tile(other.std_dev, (1, self.num_nuclides, 1))
# Add nuclides to each tally such that each tally contains the complete
# set of nuclides necessary to perform an entrywise product. New nuclides
# added to a tally will have all their scores set to zero.
else:
# Get the set of nuclides that each tally is missing
other_missing_nuclides = set(self.nuclides).difference(set(other.nuclides))
self_missing_nuclides = set(other.nuclides).difference(set(self.nuclides))
# Add nuclides present in self but not in other to other
for nuclide in other_missing_nuclides:
other._mean = np.insert(other.mean, other.num_nuclides, 0, axis=1)
other._std_dev = np.insert(other.std_dev, other.num_nuclides, 0, axis=1)
other.add_nuclide(nuclide)
# Add nuclides present in other but not in self to self
for nuclide in self_missing_nuclides:
self._mean = np.insert(self.mean, self.num_nuclides, 0, axis=1)
self._std_dev = np.insert(self.std_dev, self.num_nuclides, 0, axis=1)
self.add_nuclide(nuclide)
# Align other nuclides with self nuclides
for i, nuclide in enumerate(self.nuclides):
other_index = other.nuclides.index(nuclide)
# If necessary, swap other nuclide
if other_index != i:
other.swap_nuclides(nuclide, other.nuclides[i])
# Repeat and tile the data by score in preparation for performing
# the tensor product across scores.
if score_product == 'tensor':
self._mean = np.repeat(self.mean, other.num_score_bins, axis=2)
self._std_dev = np.repeat(self.std_dev, other.num_score_bins, axis=2)
other._mean = np.tile(other.mean, (1, 1, self.num_score_bins))
other._std_dev = np.tile(other.std_dev, (1, 1, self.num_score_bins))
# Add scores to each tally such that each tally contains the complete set
# of scores necessary to perform an entrywise product. New scores added
# to a tally will be set to zero.
else:
# Get the set of scores that each tally is missing
other_missing_scores = set(self.scores).difference(set(other.scores))
self_missing_scores = set(other.scores).difference(set(self.scores))
# Add scores present in self but not in other to other
for score in other_missing_scores:
other._mean = np.insert(other.mean, other.num_score_bins, 0, axis=2)
other._std_dev = np.insert(other.std_dev, other.num_score_bins, 0, axis=2)
other.add_score(score)
# Add scores present in other but not in self to self
for score in self_missing_scores:
self._mean = np.insert(self.mean, self.num_score_bins, 0, axis=2)
self._std_dev = np.insert(self.std_dev, self.num_score_bins, 0, axis=2)
self.add_score(score)
# Align other scores with self scores
for i, score in enumerate(self.scores):
other_index = other.scores.index(score)
# If necessary, swap other score
if other_index != i:
other.swap_scores(score, other.scores[i])
# Correct the stride for other filters
stride = other.num_nuclides * other.num_score_bins
for filter in reversed(other.filters):
filter.stride = stride
stride *= filter.num_bins
# Correct the stride for self filters
stride = self.num_nuclides * self.num_score_bins
for filter in reversed(self.filters):
filter.stride = stride
stride *= filter.num_bins
# Deep copy the mean and std dev data
self_mean = copy.deepcopy(self.mean)
self_std_dev = copy.deepcopy(self.std_dev)
other_mean = copy.deepcopy(other.mean)
other_std_dev = copy.deepcopy(other.std_dev)
if self.filters != other.filters:
# Determine the number of paired combinations of filter bins
# between the two tallies and repeat arrays along filter axes
diff1 = list(set(self.filters).difference(set(other.filters)))
diff2 = list(set(other.filters).difference(set(self.filters)))
# Determine the factors by which each tally operands' data arrays
# must be tiled or repeated for the tally outer product
other_tile_factor = 1
self_repeat_factor = 1
for filter in diff1:
other_tile_factor *= filter.num_bins
for filter in diff2:
self_repeat_factor *= filter.num_bins
# Tile / repeat the tally data for the tally outer product
self_shape = list(self_mean.shape)
other_shape = list(other_mean.shape)
self_shape[0] *= self_repeat_factor
self_mean = np.repeat(self_mean, self_repeat_factor)
self_std_dev = np.repeat(self_std_dev, self_repeat_factor)
if self_repeat_factor == 1:
other_shape[0] *= other_tile_factor
other_mean = np.repeat(other_mean, other_tile_factor, axis=0)
other_std_dev = np.repeat(other_std_dev, other_tile_factor,
axis=0)
else:
other_mean = np.tile(other_mean, (other_tile_factor, 1, 1))
other_std_dev = np.tile(other_std_dev, (other_tile_factor, 1, 1))
# NumPy repeat and tile routines return 1D flattened arrays
# Reshape arrays as 3D with filters, nuclides and scores axes
self_mean.shape = tuple(self_shape)
self_std_dev.shape = tuple(self_shape)
other_mean.shape = tuple(other_shape)
other_std_dev.shape = tuple(other_shape)
if self.nuclides != other.nuclides:
# Determine the number of paired combinations of nuclides
# between the two tallies and repeat arrays along nuclide axes
self_repeat_factor = other.num_nuclides
other_tile_factor = self.num_nuclides
# Tile / repeat the tally data for the tally outer product
self_shape = list(self_mean.shape)
other_shape = list(other_mean.shape)
self_shape[1] *= self_repeat_factor
other_shape[1] *= other_tile_factor
self_mean = np.repeat(self_mean, self_repeat_factor, axis=1)
other_mean = np.tile(other_mean, (1, other_tile_factor, 1))
self_std_dev = np.repeat(self_std_dev, self_repeat_factor, axis=1)
other_std_dev = np.tile(other_std_dev, (1, other_tile_factor, 1))
# NumPy repeat and tile routines return 1D flattened arrays
# Reshape arrays as 3D with filters, nuclides and scores axes
self_mean.shape = tuple(self_shape)
self_std_dev.shape = tuple(self_shape)
other_mean.shape = tuple(other_shape)
other_std_dev.shape = tuple(other_shape)
if self.scores != other.scores:
# Determine the number of paired combinations of score bins
# between the two tallies and repeat arrays along score axes
self_repeat_factor = other.num_score_bins
other_tile_factor = self.num_score_bins
# Tile / repeat the tally data for the tally outer product
self_shape = list(self_mean.shape)
other_shape = list(other_mean.shape)
self_shape[2] *= self_repeat_factor
other_shape[2] *= other_tile_factor
self_mean = np.repeat(self_mean, self_repeat_factor, axis=2)
other_mean = np.tile(other_mean, (1, 1, other_tile_factor))
self_std_dev = np.repeat(self_std_dev, self_repeat_factor, axis=2)
other_std_dev = np.tile(other_std_dev, (1, 1, other_tile_factor))
# NumPy repeat and tile routines return 1D flattened arrays
# Reshape arrays as 3D with filters, nuclides and scores axes
self_mean.shape = tuple(self_shape)
self_std_dev.shape = tuple(self_shape)
other_mean.shape = tuple(other_shape)
other_std_dev.shape = tuple(other_shape)
data = {}
data['self'] = {}
data['other'] = {}
@ -1903,6 +1962,116 @@ class Tally(object):
if not inplace:
return swap_tally
def swap_nuclides(self, nuclide1, nuclide2):
"""Reverse the ordering of two nuclides in this tally
This is a helper method for tally arithmetic which helps align the data
in two tallies with shared nuclides. This method reverses the order of
the two nuclides in place.
Parameters
----------
nuclide1 : Nuclide
The nuclide to swap with nuclide2
nuclide2 : Nuclide
The nuclide to swap with nuclide1
Raises
------
ValueError
If this is a derived tally or this method is called before the tally
is populated with data by the StatePoint.read_results() method.
"""
# Check that results have been read
if not self.derived and self.sum is None:
msg = 'Unable to use tally arithmetic with Tally ID="{0}" ' \
'since it does not contain any results.'.format(self.id)
raise ValueError(msg)
# Swap the nuclides in the Tally
nuclide1_index = self.nuclides.index(nuclide1)
nuclide2_index = self.nuclides.index(nuclide2)
self.nuclides[nuclide1_index] = nuclide2
self.nuclides[nuclide2_index] = nuclide1
# Copy the tally data
self_mean = copy.deepcopy(self.mean)
self_std_dev = copy.deepcopy(self.std_dev)
# Swap nuclide 1 in place of nuclide 2
self._mean = np.delete(self.mean, nuclide2_index, axis=1)
self._mean = np.insert(self.mean, nuclide2_index,
self_mean[:,nuclide1_index,:], axis=1)
self._std_dev = np.delete(self.std_dev, nuclide2_index, axis=1)
self._std_dev = np.insert(self.std_dev, nuclide2_index,
self_mean[:,nuclide1_index,:], axis=1)
# Swap nuclide 2 in place of nuclide 1
self._mean = np.delete(self.mean, nuclide1_index, axis=1)
self._mean = np.insert(self.mean, nuclide1_index,
self_mean[:,nuclide2_index,:], axis=1)
self._std_dev = np.delete(self.std_dev, nuclide1_index, axis=1)
self._std_dev = np.insert(self.std_dev, nuclide1_index,
self_mean[:,nuclide2_index,:], axis=1)
def swap_scores(self, score1, score2):
"""Reverse the ordering of two scores in this tally
This is a helper method for tally arithmetic which helps align the data
in two tallies with shared scores. This method copies reverses the order
of the two scores in place.
Parameters
----------
score1 : Score
The score to swap with score2
score2 : Score
The score to swap with score1
Raises
------
ValueError
If this is a derived tally or this method is called before the tally
is populated with data by the StatePoint.read_results() method.
"""
# Check that results have been read
if not self.derived and self.sum is None:
msg = 'Unable to use tally arithmetic with Tally ID="{0}" ' \
'since it does not contain any results.'.format(self.id)
raise ValueError(msg)
# Swap the scores in the Tally
score1_index = self.scores.index(score1)
score2_index = self.scores.index(score2)
self.scores[score1_index] = score2
self.scores[score2_index] = score1
# Copy the tally data
self_mean = copy.deepcopy(self.mean)
self_std_dev = copy.deepcopy(self.std_dev)
# Swap score 1 in place of score 2
self._mean = np.delete(self.mean, score2_index, axis=2)
self._mean = np.insert(self.mean, score2_index,
self_mean[:,:,score1_index], axis=2)
self._std_dev = np.delete(self.std_dev, score2_index, axis=2)
self._std_dev = np.insert(self.std_dev, score2_index,
self_mean[:,:,score1_index], axis=2)
# Swap score 2 in place of score 1
self._mean = np.delete(self.mean, score1_index, axis=2)
self._mean = np.insert(self.mean, score1_index,
self_mean[:,:,score2_index], axis=2)
self._std_dev = np.delete(self.std_dev, score1_index, axis=2)
self._std_dev = np.insert(self.std_dev, score1_index,
self_mean[:,:,score2_index], axis=2)
def __add__(self, other):
"""Adds this tally to another tally or scalar value.
@ -1946,7 +2115,7 @@ class Tally(object):
raise ValueError(msg)
if isinstance(other, Tally):
new_tally = self._outer_product(other, binary_op='+')
new_tally = self.hybrid_product(other, binary_op='+')
elif isinstance(other, Real):
new_tally = Tally(name='derived')
@ -2016,7 +2185,7 @@ class Tally(object):
raise ValueError(msg)
if isinstance(other, Tally):
new_tally = self._outer_product(other, binary_op='-')
new_tally = self.hybrid_product(other, binary_op='-')
elif isinstance(other, Real):
new_tally = Tally(name='derived')
@ -2086,7 +2255,7 @@ class Tally(object):
raise ValueError(msg)
if isinstance(other, Tally):
new_tally = self._outer_product(other, binary_op='*')
new_tally = self.hybrid_product(other, binary_op='*')
elif isinstance(other, Real):
new_tally = Tally(name='derived')
@ -2156,7 +2325,7 @@ class Tally(object):
raise ValueError(msg)
if isinstance(other, Tally):
new_tally = self._outer_product(other, binary_op='/')
new_tally = self.hybrid_product(other, binary_op='/')
elif isinstance(other, Real):
new_tally = Tally(name='derived')
@ -2229,7 +2398,7 @@ class Tally(object):
raise ValueError(msg)
if isinstance(power, Tally):
new_tally = self._outer_product(power, binary_op='^')
new_tally = self.hybrid_product(power, binary_op='^')
elif isinstance(power, Real):
new_tally = Tally(name='derived')

2
setup.py
View file

@ -10,7 +10,7 @@ except ImportError:
have_setuptools = False
kwargs = {'name': 'openmc',
'version': '0.7.0',
'version': '0.7.1',
'packages': ['openmc', 'openmc.mgxs'],
'scripts': glob.glob('scripts/openmc-*'),

2
src/constants.F90
View file

@ -8,7 +8,7 @@ module constants
! OpenMC major, minor, and release numbers
integer, parameter :: VERSION_MAJOR = 0
integer, parameter :: VERSION_MINOR = 7
integer, parameter :: VERSION_RELEASE = 0
integer, parameter :: VERSION_RELEASE = 1
! Revision numbers for binary files
integer, parameter :: REVISION_STATEPOINT = 14

8
tests/test_tally_arithmetic/geometry.xml Normal file
View file

@ -0,0 +1,8 @@
<?xml version="1.0"?>
<geometry>
<!-- Sphere with radius 10 -->
<surface id="1" type="sphere" coeffs="0 0 0 10" boundary="vacuum"/>
<cell id="1" material="1" region="-1" />
</geometry>

11
tests/test_tally_arithmetic/materials.xml Normal file
View file

@ -0,0 +1,11 @@
<?xml version="1.0"?>
<materials>
<material id="1">
<density value="4.5" units="g/cc" />
<nuclide name="U-238" xs="71c" ao="1.0" />
<nuclide name="U-235" xs="71c" ao="1.0" />
<nuclide name="Pu-239" xs="71c" ao="1.0" />
</material>
</materials>

53
tests/test_tally_arithmetic/results_true.dat Normal file
View file

@ -0,0 +1,53 @@
Tally
ID = 10000
Name = (tally 1 + tally 2)
Filters =
cell [1]
energy [ 0. 20.]
material [1]
Nuclides = (U-235 + U-235) (U-235 + Pu-239) (U-238 + U-235) (U-238 + Pu-239)
Scores = [(fission + fission), (fission + absorption), (nu-fission + fission), (nu-fission + absorption)]
Estimator = tracklength
[[[ 0.07510122 0.07839105 0.13713377 0.1404236 ]
[ 0.0916553 0.09321683 0.15368785 0.15524938]
[ 0.04596277 0.0492526 0.06126275 0.06455259]
[ 0.06251685 0.06407838 0.07781683 0.07937836]]]Tally
ID = 10001
Name = (tally 1 - tally 2)
Filters =
cell [1]
energy [ 0. 20.]
material [1]
Nuclides = (U-235 - U-235) (U-235 - Pu-239) (U-238 - U-235) (U-238 - Pu-239)
Scores = [(fission - fission), (fission - absorption), (nu-fission - fission), (nu-fission - absorption)]
Estimator = tracklength
[[[ 0. -0.00328983 0.06203255 0.05874271]
[-0.01655408 -0.01811561 0.04547847 0.04391694]
[-0.02913845 -0.03242829 -0.01383847 -0.0171283 ]
[-0.04569253 -0.04725406 -0.03039255 -0.03195408]]]Tally
ID = 10002
Name = (tally 1 * tally 2)
Filters =
cell [1]
energy [ 0. 20.]
material [1]
Nuclides = (U-235 * U-235) (U-235 * Pu-239) (U-238 * U-235) (U-238 * Pu-239)
Scores = [(fission * fission), (fission * absorption), (nu-fission * fission), (nu-fission * absorption)]
Estimator = tracklength
[[[ 0.00141005 0.00153358 0.00373941 0.00406702]
[ 0.00203166 0.0020903 0.00538792 0.00554342]
[ 0.00031588 0.00034356 0.00089041 0.00096841]
[ 0.00045514 0.00046827 0.00128294 0.00131997]]]Tally
ID = 10003
Name = (tally 1 / tally 2)
Filters =
cell [1]
energy [ 0. 20.]
material [1]
Nuclides = (U-235 / U-235) (U-235 / Pu-239) (U-238 / U-235) (U-238 / Pu-239)
Scores = [(fission / fission), (fission / absorption), (nu-fission / fission), (nu-fission / absorption)]
Estimator = tracklength
[[[ 1. 0.91944666 2.65197177 2.4383466 ]
[ 0.69403611 0.67456727 1.84056418 1.78893335]
[ 0.22402189 0.20597618 0.63147153 0.58060439]
[ 0.15547928 0.15111784 0.43826405 0.42597002]]]

18
tests/test_tally_arithmetic/settings.xml Normal file
View file

@ -0,0 +1,18 @@
<?xml version="1.0"?>
<settings>
<eigenvalue>
<batches>10</batches>
<inactive>5</inactive>
<particles>1000</particles>
</eigenvalue>
<source>
<space type="box">
<parameters>-4 -4 -4 4 4 4</parameters>
</space>
</source>
<output summary="true"/>
</settings>

15
tests/test_tally_arithmetic/tallies.xml Normal file
View file

@ -0,0 +1,15 @@
<?xml version='1.0' encoding='utf-8'?>
<tallies>
<tally id="10000" name="tally 1">
<filter bins="1" type="cell" />
<filter bins="0.0 20.0" type="energy" />
<nuclides>U-235 U-238</nuclides>
<scores>fission nu-fission</scores>
</tally>
<tally id="10001" name="tally 2">
<filter bins="1" type="material" />
<filter bins="0.0 20.0" type="energy" />
<nuclides>U-235 Pu-239</nuclides>
<scores>fission absorption</scores>
</tally>
</tallies>

59
tests/test_tally_arithmetic/test_tally_arithmetic.py Normal file
View file

@ -0,0 +1,59 @@
#!/usr/bin/env python
import os
import sys
import glob
import hashlib
sys.path.insert(0, os.pardir)
from testing_harness import TestHarness
import openmc
class TallyArithmeticTestHarness(TestHarness):
def _get_results(self, hash_output=False):
"""Digest info in the statepoint and return as a string."""
# Read the statepoint file.
statepoint = glob.glob(os.path.join(os.getcwd(), self._sp_name))[0]
sp = openmc.StatePoint(statepoint)
# Read the summary file.
summary = glob.glob(os.path.join(os.getcwd(), 'summary.h5'))[0]
su = openmc.Summary(summary)
sp.link_with_summary(su)
# Load the tallies
tally_1 = sp.get_tally(name='tally 1')
tally_2 = sp.get_tally(name='tally 2')
# Perform all the tally arithmetic operations and output results
outstr = ''
tally_3 = tally_1 + tally_2
outstr += repr(tally_3)
outstr += str(tally_3.mean)
tally_3 = tally_1 - tally_2
outstr += repr(tally_3)
outstr += str(tally_3.mean)
tally_3 = tally_1 * tally_2
outstr += repr(tally_3)
outstr += str(tally_3.mean)
tally_3 = tally_1 / tally_2
outstr += repr(tally_3)
outstr += str(tally_3.mean)
print(outstr)
# Hash the results if necessary
if hash_output:
sha512 = hashlib.sha512()
sha512.update(outstr.encode('utf-8'))
outstr = sha512.hexdigest()
return outstr
if __name__ == '__main__':
harness = TallyArithmeticTestHarness('statepoint.10.h5', True)
harness.main()