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 'upstream/develop' into multipole

Browse source
This commit is contained in:
Sterling Harper 2016-06-09 13:32:26 -05:00
commit 3280d5465e
87 changed files with 3983 additions and 1619 deletions
Download patch file Download diff file Expand all files Collapse all files

8
docs/source/_templates/myclassinherit.rst Normal file
View file

@ -0,0 +1,8 @@
{{ fullname }}
{{ underline }}
.. currentmodule:: {{ module }}
.. autoclass:: {{ objname }}
:members:
:inherited-members:

4
docs/source/io_formats/mgxs_library.rst
View file

@ -172,7 +172,7 @@ attributes/sub-elements required to describe the meta-data:
during the scattering process. Specifically, the options are to either
convert the Legendre expansion to a tabular representation or leave it as
a set of Legendre coefficients. Converting to a tabular representation will
cost memory but is likely to decrease runtime compared to leaving as a
cost memory but can allow for a decrease in runtime compared to leaving as a
set of Legendre coefficients. This element has the following
attributes/sub-elements:
@ -181,7 +181,7 @@ attributes/sub-elements required to describe the meta-data:
tabular format should be performed or not. A value of "true" means
the conversion should be performed, "false" means it should not.
*Default*: "true"
*Default*: "false"
:num_points:
If the conversion is to take place the number of tabular points is

2
docs/source/methods/geometry.rst
View file

@ -437,6 +437,8 @@ where :math:`(x_0, y_0, z_0)` are the coordinates to the lower-left-bottom
corner of the lattice, and :math:`p_0, p_1, p_2` are the pitches along the
:math:`x`, :math:`y`, and :math:`z` axes, respectively.
.. _hexagonal_indexing:
Hexagonal Lattice Indexing
--------------------------

14
docs/source/pythonapi/examples/mgxs-part-iii.ipynb
View file

@ -260,7 +260,6 @@
"source": [
"# Create fuel assembly Lattice\n",
"assembly = openmc.RectLattice(name='1.6% Fuel Assembly')\n",
"assembly.dimension = (17, 17)\n",
"assembly.pitch = (1.26, 1.26)\n",
"assembly.lower_left = [-1.26 * 17. / 2.0] * 2"
]
@ -939,7 +938,8 @@
"cell_type": "markdown",
"metadata": {},
"source": [
"The `NuFissionXS` object supports all of the methods described previously the `openmc.mgxs` tutorials, such as [Pandas](http://pandas.pydata.org/) `DataFrames`:"
"The `NuFissionXS` object supports all of the methods described previously in the `openmc.mgxs` tutorials, such as [Pandas](http://pandas.pydata.org/) `DataFrames`:\n",
"Note that since so few histories were simulated, we should expect a few division-by-error errors as some tallies have not yet scored any results."
]
},
{
@ -1596,21 +1596,21 @@
],
"metadata": {
"kernelspec": {
"display_name": "Python 2",
"display_name": "Python 3",
"language": "python",
"name": "python2"
"name": "python3"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 2
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython2",
"version": "2.7.11"
"pygments_lexer": "ipython3",
"version": "3.5.1"
}
},
"nbformat": 4,

144
docs/source/pythonapi/examples/mgxs-part-iv.ipynb
View file

File diff suppressed because one or more lines are too long

11
docs/source/pythonapi/examples/pandas-dataframes.ipynb
View file

@ -199,7 +199,6 @@
"source": [
"# Create fuel assembly Lattice\n",
"assembly = openmc.RectLattice(name='1.6% Fuel - 0BA')\n",
"assembly.dimension = (17, 17)\n",
"assembly.pitch = (1.26, 1.26)\n",
"assembly.lower_left = [-1.26 * 17. / 2.0] * 2\n",
"assembly.universes = [[pin_cell_universe] * 17] * 17"
@ -2194,21 +2193,21 @@
],
"metadata": {
"kernelspec": {
"display_name": "Python 2",
"display_name": "Python 3",
"language": "python",
"name": "python2"
"name": "python3"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 2
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython2",
"version": "2.7.6"
"pygments_lexer": "ipython3",
"version": "3.5.1"
}
},
"nbformat": 4,

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

@ -271,14 +271,17 @@ Multi-group Cross Sections
.. autosummary::
:toctree: generated
:nosignatures:
:template: myclass.rst
:template: myclassinherit.rst
openmc.mgxs.MGXS
openmc.mgxs.AbsorptionXS
openmc.mgxs.CaptureXS
openmc.mgxs.Chi
openmc.mgxs.FissionXS
openmc.mgxs.KappaFissionXS
openmc.mgxs.MultiplicityMatrixXS
openmc.mgxs.NuFissionXS
openmc.mgxs.NuFissionMatrixXS
openmc.mgxs.NuScatterXS
openmc.mgxs.NuScatterMatrixXS
openmc.mgxs.ScatterXS
@ -296,6 +299,33 @@ Multi-group Cross Section Libraries
openmc.mgxs.Library
-------------------------------------
:mod:`openmc.model` -- Model Building
-------------------------------------
TRISO Fuel Modeling
-------------------
Classes
+++++++
.. autosummary::
:toctree: generated
:nosignatures:
:template: myclass.rst
openmc.model.TRISO
Functions
+++++++++
.. autosummary::
:toctree: generated
:nosignatures:
openmc.model.create_triso_lattice
.. _Jupyter: https://jupyter.org/
.. _NumPy: http://www.numpy.org/
.. _Codecademy: https://www.codecademy.com/tracks/python

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

@ -921,11 +921,19 @@ Each ``<surface>`` element can have the following attributes or sub-elements:
*Default*: None
:boundary:
The boundary condition for the surface. This can be "transmission",
"vacuum", or "reflective".
The boundary condition for the surface. This can be "transmission",
"vacuum", "reflective", or "periodic". Periodic boundary conditions can
only be applied to x-, y-, and z-planes. Only axis-aligned periodicity is
supported, i.e., x-planes can only be paired with x-planes. Specify which
planes are periodic and the code will automatically identify which planes
are paired together.
*Default*: "transmission"
:periodic_surface_id:
If a periodic boundary condition is applied, this attribute identifies the
``id`` of the corresponding periodic sufrace.
The following quadratic surfaces can be modeled:
:x-plane:
@ -1623,7 +1631,8 @@ The ``<tally>`` element accepts the following sub-elements:
|Score | Description |
+======================+===================================================+
|absorption |Total absorption rate. This accounts for all |
| |reactions which do not produce secondary neutrons. |
| |reactions which do not produce secondary neutrons |
| |as well as fission. |
+----------------------+---------------------------------------------------+
|elastic |Elastic scattering reaction rate. |
+----------------------+---------------------------------------------------+

3
examples/python/basic/build-xml.py
View file

@ -74,8 +74,7 @@ universe1.add_cells([cell2, cell3])
root.add_cells([cell1, cell4])
# Instantiate a Geometry, register the root Universe, and export to XML
geometry = openmc.Geometry()
geometry.root_universe = root
geometry = openmc.Geometry(root)
geometry.export_to_xml()

3
examples/python/boxes/build-xml.py
View file

@ -97,8 +97,7 @@ root = openmc.Universe(universe_id=0, name='root universe')
root.add_cells([inner_box, middle_box, outer_box])
# Instantiate a Geometry, register the root Universe, and export to XML
geometry = openmc.Geometry()
geometry.root_universe = root
geometry = openmc.Geometry(root)
geometry.export_to_xml()

3
examples/python/lattice/hexagonal/build-xml.py
View file

@ -105,8 +105,7 @@ lattice.outer = univ2
cell1.fill = lattice
# Instantiate a Geometry, register the root Universe, and export to XML
geometry = openmc.Geometry()
geometry.root_universe = root
geometry = openmc.Geometry(root)
geometry.export_to_xml()

5
examples/python/lattice/nested/build-xml.py
View file

@ -98,14 +98,12 @@ univ4.add_cell(cell2)
# Instantiate nested Lattices
lattice1 = openmc.RectLattice(lattice_id=4, name='4x4 assembly')
lattice1.dimension = [2, 2]
lattice1.lower_left = [-1., -1.]
lattice1.pitch = [1., 1.]
lattice1.universes = [[univ1, univ2],
[univ2, univ3]]
lattice2 = openmc.RectLattice(lattice_id=6, name='4x4 core')
lattice2.dimension = [2, 2]
lattice2.lower_left = [-2., -2.]
lattice2.pitch = [2., 2.]
lattice2.universes = [[univ4, univ4],
@ -116,8 +114,7 @@ cell1.fill = lattice2
cell2.fill = lattice1
# Instantiate a Geometry, register the root Universe, and export to XML
geometry = openmc.Geometry()
geometry.root_universe = root
geometry = openmc.Geometry(root)
geometry.export_to_xml()

4
examples/python/lattice/simple/build-xml.py
View file

@ -94,7 +94,6 @@ root.add_cell(cell1)
# Instantiate a Lattice
lattice = openmc.RectLattice(lattice_id=5)
lattice.dimension = [4, 4]
lattice.lower_left = [-2., -2.]
lattice.pitch = [1., 1.]
lattice.universes = [[univ1, univ2, univ1, univ2],
@ -106,8 +105,7 @@ lattice.universes = [[univ1, univ2, univ1, univ2],
cell1.fill = lattice
# Instantiate a Geometry, register the root Universe, and export to XML
geometry = openmc.Geometry()
geometry.root_universe = root
geometry = openmc.Geometry(root)
geometry.export_to_xml()

3
examples/python/pincell/build-xml.py
View file

@ -149,8 +149,7 @@ root = openmc.Universe(universe_id=0, name='root universe')
root.add_cells([fuel, gap, clad, water])
# Instantiate a Geometry, register the root Universe, and export to XML
geometry = openmc.Geometry()
geometry.root_universe = root
geometry = openmc.Geometry(root)
geometry.export_to_xml()

3
examples/python/pincell_multigroup/build-xml.py
View file

@ -119,8 +119,7 @@ root = openmc.Universe(universe_id=0, name='root universe')
root.add_cells([fuel, moderator])
# Instantiate a Geometry, register the root Universe, and export to XML
geometry = openmc.Geometry()
geometry.root_universe = root
geometry = openmc.Geometry(root)
geometry.export_to_xml()

3
examples/python/reflective/build-xml.py
View file

@ -64,8 +64,7 @@ root = openmc.Universe(universe_id=0, name='root universe')
root.add_cell(cell)
# Instantiate a Geometry, register the root Universe, and export to XML
geometry = openmc.Geometry()
geometry.root_universe = root
geometry = openmc.Geometry(root)
geometry.export_to_xml()

112
openmc/arithmetic.py
View file

@ -67,24 +67,6 @@ class CrossScore(object):
def __ne__(self, other):
return not self == other
def __deepcopy__(self, memo):
existing = memo.get(id(self))
# If this is the first time we have tried to copy this object, create a copy
if existing is None:
clone = type(self).__new__(type(self))
clone._left_score = self.left_score
clone._right_score = self.right_score
clone._binary_op = self.binary_op
memo[id(self)] = clone
return clone
# If this object has been copied before, return the first copy made
else:
return existing
def __repr__(self):
string = '({0} {1} {2})'.format(self.left_score,
self.binary_op, self.right_score)
@ -169,28 +151,9 @@ class CrossNuclide(object):
def __ne__(self, other):
return not self == other
def __deepcopy__(self, memo):
existing = memo.get(id(self))
# If this is the first time we have tried to copy this object, create a copy
if existing is None:
clone = type(self).__new__(type(self))
clone._left_nuclide = self.left_nuclide
clone._right_nuclide = self.right_nuclide
clone._binary_op = self.binary_op
memo[id(self)] = clone
return clone
# If this object has been copied before, return the first copy made
else:
return existing
def __repr__(self):
return self.name
@property
def left_nuclide(self):
return self._left_nuclide
@ -325,27 +288,6 @@ class CrossFilter(object):
string += '{0: <16}{1}{2}\n'.format('\tBins', '=\t', filter_bins)
return string
def __deepcopy__(self, memo):
existing = memo.get(id(self))
# If this is the first time we have tried to copy this object, create a copy
if existing is None:
clone = type(self).__new__(type(self))
clone._left_filter = self.left_filter
clone._right_filter = self.right_filter
clone._binary_op = self.binary_op
clone._type = self.type
clone._bins = self._bins
clone._stride = self.stride
memo[id(self)] = clone
return clone
# If this object has been copied before, return the first copy made
else:
return existing
@property
def left_filter(self):
return self._left_filter
@ -532,23 +474,6 @@ class AggregateScore(object):
def __ne__(self, other):
return not self == other
def __deepcopy__(self, memo):
existing = memo.get(id(self))
# If this is the first time we have tried to copy this object, create a copy
if existing is None:
clone = type(self).__new__(type(self))
clone._scores = self.scores
clone._aggregate_op = self.aggregate_op
memo[id(self)] = clone
return clone
# If this object has been copied before, return the first copy made
else:
return existing
def __repr__(self):
string = ', '.join(map(str, self.scores))
string = '{0}({1})'.format(self.aggregate_op, string)
@ -622,23 +547,6 @@ class AggregateNuclide(object):
def __ne__(self, other):
return not self == other
def __deepcopy__(self, memo):
existing = memo.get(id(self))
# If this is the first time we have tried to copy this object, create a copy
if existing is None:
clone = type(self).__new__(type(self))
clone._nuclides = self.nuclides
clone._aggregate_op = self._aggregate_op
memo[id(self)] = clone
return clone
# If this object has been copied before, return the first copy made
else:
return existing
def __repr__(self):
# Append each nuclide in the aggregate to the string
@ -757,26 +665,6 @@ class AggregateFilter(object):
string += '{0: <16}{1}{2}\n'.format('\tBins', '=\t', self.bins)
return string
def __deepcopy__(self, memo):
existing = memo.get(id(self))
# If this is the first time we have tried to copy this object, create a copy
if existing is None:
clone = type(self).__new__(type(self))
clone._type = self.type
clone._aggregate_filter = self.aggregate_filter
clone._aggregate_op = self.aggregate_op
clone._bins = self._bins
clone._stride = self.stride
memo[id(self)] = clone
return clone
# If this object has been copied before, return the first copy made
else:
return existing
@property
def aggregate_filter(self):
return self._aggregate_filter

175
openmc/cell.py
View file

@ -1,9 +1,12 @@
from collections import OrderedDict, Iterable
from math import cos, sin, pi
from numbers import Real, Integral
from xml.etree import ElementTree as ET
import sys
import warnings
import numpy as np
import openmc
import openmc.checkvalue as cv
from openmc.surface import Halfspace
@ -33,7 +36,7 @@ class Cell(object):
automatically be assigned.
name : str, optional
Name of the cell. If not specified, the name is the empty string.
fill : openmc.Material or openmc.Universe or openmc.Lattice or 'void' or iterable of openmc.Material, optional
fill : openmc.Material or openmc.Universe or openmc.Lattice or None or iterable of openmc.Material, optional
Indicates what the region of space is filled with
region : openmc.Region, optional
Region of space that is assigned to the cell.
@ -44,9 +47,13 @@ class Cell(object):
Unique identifier for the cell
name : str
Name of the cell
fill : openmc.Material or openmc.Universe or openmc.Lattice or 'void' or iterable of openmc.Material
Indicates what the region of space is filled with
region : openmc.Region
fill : openmc.Material or openmc.Universe or openmc.Lattice or None or iterable of openmc.Material
Indicates what the region of space is filled with. If None, the cell is
treated as a void. An iterable of materials is used to fill repeated
instances of a cell with different materials.
fill_type : {'material', 'universe', 'lattice', 'distribmat', 'void'}
Indicates what the cell is filled with.
region : openmc.Region or None
Region of space that is assigned to the cell.
rotation : Iterable of float
If the cell is filled with a universe, this array specifies the angles
@ -63,6 +70,9 @@ class Cell(object):
\sin\phi \sin\theta \sin\psi & -\sin\phi \cos\psi + \cos\phi
\sin\theta \sin\psi \\ -\sin\theta & \sin\phi \cos\theta & \cos\phi
\cos\theta \end{array} \right ]
rotation_matrix : numpy.ndarray
The rotation matrix defined by the angles specified in the
:attr:`Cell.rotation` property.
temperature : float or iterable of float
Temperature of the cell in Kelvin. Multiple temperatures can be given
to give each distributed cell instance a unique temperature.
@ -80,19 +90,20 @@ class Cell(object):
# Initialize Cell class attributes
self.id = cell_id
self.name = name
self._fill = None
self._type = None
self._region = None
self._temperature = None
self.fill = fill
self.region = region
self._rotation = None
self._rotation_matrix = None
self._temperature = None
self._translation = None
self._offsets = None
self._distribcell_index = None
if fill is not None:
self.fill = fill
if region is not None:
self.region = region
def __contains__(self, point):
if self.region is None:
return True
else:
return point in self.region
def __eq__(self, other):
if not isinstance(other, Cell):
@ -122,36 +133,27 @@ class Cell(object):
def __repr__(self):
string = 'Cell\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)
if isinstance(self._fill, openmc.Material):
string += '{0: <16}{1}{2}\n'.format('\tMaterial', '=\t',
self._fill._id)
elif isinstance(self._fill, Iterable):
string += '{0: <16}{1}'.format('\tMaterial', '=\t')
string += '['
string += ', '.join(['void' if m == 'void' else str(m.id)
for m in self.fill])
string += ']\n'
elif isinstance(self._fill, (openmc.Universe, openmc.Lattice)):
string += '{0: <16}{1}{2}\n'.format('\tFill', '=\t',
self._fill._id)
if self.fill_type == 'material':
string += '{: <16}=\tMaterial {}\n'.format('\tFill', self.fill.id)
elif self.fill_type == 'void':
string += '{: <16}=\tNone\n'.format('\tFill')
elif self.fill_type == 'distribmat':
string += '{: <16}=\t{}\n'.format('\tFill', list(map(
lambda m: m if m is None else m.id, self.fill)))
else:
string += '{0: <16}{1}{2}\n'.format('\tFill', '=\t', self._fill)
string += '{: <16}=\t{}\n'.format('\tFill', self.fill.id)
string += '{0: <16}{1}{2}\n'.format('\tRegion', '=\t', self._region)
string += '{0: <16}{1}{2}\n'.format('\tRotation', '=\t',
self._rotation)
string += '{: <16}=\t{}\n'.format('\tRegion', self.region)
string += '{: <16}=\t{}\n'.format('\tRotation', self.rotation)
if self.fill_type == 'material':
string += '\t{0: <15}=\t{1}\n'.format('Temperature',
self.temperature)
string += '{0: <16}{1}{2}\n'.format('\tTranslation', '=\t',
self._translation)
string += '{0: <16}{1}{2}\n'.format('\tOffset', '=\t', self._offsets)
string += '{0: <16}{1}{2}\n'.format('\tDistribcell index', '=\t',
self._distribcell_index)
string += '{: <16}=\t{}\n'.format('\tTranslation', self.translation)
string += '{: <16}=\t{}\n'.format('\tOffset', self.offsets)
string += '{: <16}=\t{}\n'.format('\tDistribcell index', self.distribcell_index)
return string
@ -175,8 +177,10 @@ class Cell(object):
return 'universe'
elif isinstance(self.fill, openmc.Lattice):
return 'lattice'
elif isinstance(self.fill, Iterable):
return 'distribmat'
else:
return None
return 'void'
@property
def region(self):
@ -186,6 +190,10 @@ class Cell(object):
def rotation(self):
return self._rotation
@property
def rotation_matrix(self):
return self._rotation_matrix
@property
def temperature(self):
return self._temperature
@ -223,33 +231,25 @@ class Cell(object):
@fill.setter
def fill(self, fill):
if isinstance(fill, basestring):
if fill.strip().lower() == 'void':
self._type = 'void'
else:
if fill is not None:
if isinstance(fill, basestring):
if fill.strip().lower() != 'void':
msg = 'Unable to set Cell ID="{0}" to use a non-Material ' \
'or Universe fill "{1}"'.format(self._id, fill)
raise ValueError(msg)
fill = None
elif isinstance(fill, Iterable):
for i, f in enumerate(fill):
if f is not None:
cv.check_type('cell.fill[i]', f, openmc.Material)
elif not isinstance(fill, (openmc.Material, openmc.Lattice,
openmc.Universe)):
msg = 'Unable to set Cell ID="{0}" to use a non-Material or ' \
'Universe fill "{1}"'.format(self._id, fill)
'Universe fill "{1}"'.format(self._id, fill)
raise ValueError(msg)
elif isinstance(fill, openmc.Material):
self._type = 'normal'
elif isinstance(fill, Iterable):
cv.check_type('cell.fill', fill, Iterable,
(openmc.Material, basestring))
self._type = 'normal'
elif isinstance(fill, openmc.Universe):
self._type = 'fill'
elif isinstance(fill, openmc.Lattice):
self._type = 'lattice'
else:
msg = 'Unable to set Cell ID="{0}" to use a non-Material or ' \
'Universe fill "{1}"'.format(self._id, fill)
raise ValueError(msg)
self._fill = fill
@rotation.setter
@ -260,13 +260,23 @@ class Cell(object):
cv.check_type('cell rotation', rotation, Iterable, Real)
cv.check_length('cell rotation', rotation, 3)
self._rotation = rotation
self._rotation = np.asarray(rotation)
# Save rotation matrix
phi, theta, psi = self.rotation*(-pi/180.)
c3, s3 = cos(phi), sin(phi)
c2, s2 = cos(theta), sin(theta)
c1, s1 = cos(psi), sin(psi)
self._rotation_matrix = np.array([
[c1*c2, c1*s2*s3 - c3*s1, s1*s3 + c1*c3*s2],
[c2*s1, c1*c3 + s1*s2*s3, c3*s1*s2 - c1*s3],
[-s2, c2*s3, c2*c3]])
@translation.setter
def translation(self, translation):
cv.check_type('cell translation', translation, Iterable, Real)
cv.check_length('cell translation', translation, 3)
self._translation = translation
self._translation = np.asarray(translation)
@temperature.setter
def temperature(self, temperature):
@ -286,7 +296,8 @@ class Cell(object):
@region.setter
def region(self, region):
cv.check_type('cell region', region, Region)
if region is not None:
cv.check_type('cell region', region, Region)
self._region = region
@distribcell_index.setter
@ -341,11 +352,11 @@ class Cell(object):
def get_cell_instance(self, path, distribcell_index):
# If the Cell is filled by a Material
if self._type == 'normal' or self._type == 'void':
if self.fill_type in ('material', 'distribmat', 'void'):
offset = 0
# If the Cell is filled by a Universe
elif self._type == 'fill':
elif self.fill_type == 'universe':
offset = self.offsets[distribcell_index-1]
offset += self.fill.get_cell_instance(path, distribcell_index)
@ -368,8 +379,8 @@ class Cell(object):
nuclides = OrderedDict()
if self._type != 'void':
nuclides.update(self._fill.get_all_nuclides())
if self.fill_type != 'void':
nuclides.update(self.fill.get_all_nuclides())
return nuclides
@ -387,8 +398,8 @@ class Cell(object):
cells = OrderedDict()
if self._type == 'fill' or self._type == 'lattice':
cells.update(self._fill.get_all_cells())
if self.fill_type in ('universe', 'lattice'):
cells.update(self.fill.get_all_cells())
return cells
@ -428,11 +439,11 @@ class Cell(object):
universes = OrderedDict()
if self._type == 'fill':
universes[self._fill._id] = self._fill
universes.update(self._fill.get_all_universes())
elif self._type == 'lattice':
universes.update(self._fill.get_all_universes())
if self.fill_type == 'universe':
universes[self.fill.id] = self.fill
universes.update(self.fill.get_all_universes())
elif self.fill_type == 'lattice':
universes.update(self.fill.get_all_universes())
return universes
@ -443,24 +454,20 @@ class Cell(object):
if len(self._name) > 0:
element.set("name", str(self.name))
if isinstance(self.fill, basestring):
if self.fill_type == 'void':
element.set("material", "void")
elif isinstance(self.fill, openmc.Material):
elif self.fill_type == 'material':
element.set("material", str(self.fill.id))
elif isinstance(self.fill, Iterable):
element.set("material", ' '.join([m if m == 'void' else str(m.id)
elif self.fill_type == 'distribmat':
element.set("material", ' '.join(['void' if m is None else str(m.id)
for m in self.fill]))
elif isinstance(self.fill, (openmc.Universe, openmc.Lattice)):
elif self.fill_type in ('universe', 'lattice'):
element.set("fill", str(self.fill.id))
self.fill.create_xml_subelement(xml_element)
else:
element.set("fill", str(self.fill))
self.fill.create_xml_subelement(xml_element)
if self.region is not None:
# Set the region attribute with the region specification
element.set("region", str(self.region))

41
openmc/checkvalue.py
View file

@ -4,33 +4,6 @@ from numbers import Integral, Real
import numpy as np
def _isinstance(value, expected_type):
"""A Numpy-aware replacement for isinstance
This function will be obsolete when Numpy v. >= 1.9 is established.
"""
# Declare numpy numeric types.
np_ints = (np.int_, np.intc, np.intp, np.int8, np.int16, np.int32, np.int64,
np.uint8, np.uint16, np.uint32, np.uint64)
np_floats = (np.float_, np.float16, np.float32, np.float64)
# Include numpy integers, if necessary.
if type(expected_type) is tuple:
if Integral in expected_type:
expected_type = expected_type + np_ints
elif expected_type is Integral:
expected_type = (Integral, ) + np_ints
# Include numpy floats, if necessary.
if type(expected_type) is tuple:
if Real in expected_type:
expected_type = expected_type + np_floats
elif expected_type is Real:
expected_type = (Real, ) + np_floats
# Now, make the instance check.
return isinstance(value, expected_type)
def check_type(name, value, expected_type, expected_iter_type=None):
"""Ensure that an object is of an expected type. Optionally, if the object is
@ -50,7 +23,7 @@ def check_type(name, value, expected_type, expected_iter_type=None):
"""
if not _isinstance(value, expected_type):
if not isinstance(value, expected_type):
if isinstance(expected_type, Iterable):
msg = 'Unable to set "{0}" to "{1}" which is not one of the ' \
'following types: "{2}"'.format(name, value, ', '.join(
@ -61,8 +34,16 @@ def check_type(name, value, expected_type, expected_iter_type=None):
raise TypeError(msg)
if expected_iter_type:
if isinstance(value, np.ndarray):
if not issubclass(value.dtype.type, expected_iter_type):
msg = 'Unable to set "{0}" to "{1}" since each item must be ' \
'of type "{2}"'.format(name, value,
expected_iter_type.__name__)
else:
return
for item in value:
if not _isinstance(item, expected_iter_type):
if not isinstance(item, expected_iter_type):
if isinstance(expected_iter_type, Iterable):
msg = 'Unable to set "{0}" to "{1}" since each item must be ' \
'one of the following types: "{2}"'.format(
@ -118,7 +99,7 @@ def check_iterable_type(name, value, expected_type, min_depth=1, max_depth=1):
# If this item is of the expected type, then we've reached the bottom
# level of this branch.
if _isinstance(current_item, expected_type):
if isinstance(current_item, expected_type):
# Is this deep enough?
if len(tree) < min_depth:
msg = 'Error setting "{0}": The item at {1} does not meet the '\

4
openmc/element.py
View file

@ -126,8 +126,8 @@ class Element(object):
"""
isotopes = []
for isotope, abundance in natural_abundance.items():
if isotope.startswith(self.name):
for isotope, abundance in sorted(natural_abundance.items()):
if isotope.startswith(self.name + '-'):
nuc = openmc.Nuclide(isotope, self.xs)
isotopes.append((nuc, abundance))
return isotopes

25
openmc/filter.py
View file

@ -104,27 +104,6 @@ class Filter(object):
def __hash__(self):
return hash(repr(self))
def __deepcopy__(self, memo):
existing = memo.get(id(self))
# If this is the first time we have tried to copy this object, create a copy
if existing is None:
clone = type(self).__new__(type(self))
clone._type = self.type
clone._bins = copy.deepcopy(self.bins, memo)
clone._num_bins = self.num_bins
clone._mesh = copy.deepcopy(self.mesh, memo)
clone._stride = self.stride
clone._distribcell_paths = copy.deepcopy(self.distribcell_paths)
memo[id(self)] = clone
return clone
# If this object has been copied before, return the first copy made
else:
return existing
def __repr__(self):
string = 'Filter\n'
string += '{0: <16}{1}{2}\n'.format('\tType', '=\t', self.type)
@ -196,7 +175,7 @@ class Filter(object):
elif self.type in ['energy', 'energyout']:
for edge in bins:
if not cv._isinstance(edge, Real):
if not isinstance(edge, Real):
msg = 'Unable to add bin edge "{0}" to a "{1}" Filter ' \
'since it is a non-integer or floating point ' \
'value'.format(edge, self.type)
@ -220,7 +199,7 @@ class Filter(object):
msg = 'Unable to add bins "{0}" to a mesh Filter since ' \
'only a single mesh can be used per tally'.format(bins)
raise ValueError(msg)
elif not cv._isinstance(bins[0], Integral):
elif not isinstance(bins[0], Integral):
msg = 'Unable to add bin "{0}" to mesh Filter since it ' \
'is a non-integer'.format(bins[0])
raise ValueError(msg)

81
openmc/geometry.py
View file

@ -15,6 +15,11 @@ def reset_auto_ids():
class Geometry(object):
"""Geometry representing a collection of surfaces, cells, and universes.
Parameters
----------
root_universe : openmc.Universe, optional
Root universe which contains all others
Attributes
----------
root_universe : openmc.Universe
@ -22,9 +27,11 @@ class Geometry(object):
"""
def __init__(self):
def __init__(self, root_universe=None):
self._root_universe = None
self._offsets = {}
if root_universe is not None:
self.root_universe = root_universe
@property
def root_universe(self):
@ -62,6 +69,23 @@ class Geometry(object):
tree.write("geometry.xml", xml_declaration=True, encoding='utf-8',
method="xml")
def find(self, point):
"""Find cells/universes/lattices which contain a given point
Parameters
----------
point : 3-tuple of float
Cartesian coordinates of the point
Returns
-------
list
Sequence of universes, cells, and lattices which are traversed to
find the given point
"""
return self.root_universe.find(point)
def get_cell_instance(self, path):
"""Return the instance number for the final cell in a geometry path.
@ -120,14 +144,8 @@ class Geometry(object):
"""
all_cells = self._root_universe.get_all_cells()
cells = set()
for cell in all_cells.values():
if cell._type == 'normal':
cells.add(cell)
cells = list(cells)
all_cells = self.root_universe.get_all_cells()
cells = list(set(all_cells.values()))
cells.sort(key=lambda x: x.id)
return cells
@ -142,12 +160,7 @@ class Geometry(object):
"""
all_universes = self._root_universe.get_all_universes()
universes = set()
for universe in all_universes.values():
universes.add(universe)
universes = list(universes)
universes = list(set(all_universes.values()))
universes.sort(key=lambda x: x.id)
return universes
@ -180,15 +193,17 @@ class Geometry(object):
"""
material_cells = self.get_all_material_cells()
materials = set()
materials = []
for cell in material_cells:
if isinstance(cell.fill, Iterable):
for m in cell.fill: materials.add(m)
else:
materials.add(cell.fill)
if cell.fill_type == 'distribmat':
for m in cell.fill:
if m is not None and m not in materials:
materials.append(m)
elif cell.fill_type == 'material':
if cell.fill not in materials:
materials.append(cell.fill)
materials = list(materials)
materials.sort(key=lambda x: x.id)
return materials
@ -203,13 +218,13 @@ class Geometry(object):
"""
all_cells = self.get_all_cells()
material_cells = set()
material_cells = []
for cell in all_cells:
if cell._type == 'normal':
material_cells.add(cell)
if cell.fill_type in ('material', 'distribmat'):
if cell not in material_cells:
material_cells.append(cell)
material_cells = list(material_cells)
material_cells.sort(key=lambda x: x.id)
return material_cells
@ -224,15 +239,15 @@ class Geometry(object):
"""
all_universes = self.get_all_universes()
material_universes = set()
material_universes = []
for universe in all_universes:
cells = universe.cells
for cell in cells:
if cell._type == 'normal':
material_universes.add(universe)
if cell.fill_type in ('material', 'distribmat', 'void'):
if universe not in material_universes:
material_universes.append(universe)
material_universes = list(material_universes)
material_universes.sort(key=lambda x: x.id)
return material_universes
@ -247,13 +262,13 @@ class Geometry(object):
"""
cells = self.get_all_cells()
lattices = set()
lattices = []
for cell in cells:
if isinstance(cell.fill, openmc.Lattice):
lattices.add(cell.fill)
if cell.fill_type == 'lattice':
if cell.fill not in lattices:
lattices.append(cell.fill)
lattices = list(lattices)
lattices.sort(key=lambda x: x.id)
return lattices

461
openmc/lattice.py
View file

@ -1,8 +1,12 @@
from __future__ import division
import abc
from collections import OrderedDict, Iterable
from math import sqrt, floor
from numbers import Real, Integral
from xml.etree import ElementTree as ET
import sys
import warnings
import numpy as np
@ -113,12 +117,6 @@ class Lattice(object):
cv.check_type('outer universe', outer, openmc.Universe)
self._outer = outer
@universes.setter
def universes(self, universes):
cv.check_iterable_type('lattice universes', universes, openmc.Universe,
min_depth=2, max_depth=3)
self._universes = np.asarray(universes)
def get_unique_universes(self):
"""Determine all unique universes in the lattice
@ -239,6 +237,19 @@ class Lattice(object):
class RectLattice(Lattice):
"""A lattice consisting of rectangular prisms.
To completely define a rectangular lattice, the
:attr:`RectLattice.lower_left` :attr:`RectLattice.pitch`,
:attr:`RectLattice.outer`, and :attr:`RectLattice.universes` properties need
to be set.
Most methods for this class use a natural indexing scheme wherein elements
are assigned an index corresponding to their position relative to the
(x,y,z) axes in a Cartesian coordinate system, i.e., an index of (0,0,0) in
the lattice gives the element whose x, y, and z coordinates are the
smallest. However, note that when universes are assigned to lattice elements
using the :attr:`RectLattice.universes` property, the array indices do not
correspond to natural indices.
Parameters
----------
lattice_id : int, optional
@ -253,12 +264,6 @@ class RectLattice(Lattice):
Unique identifier for the lattice
name : str
Name of the lattice
dimension : Iterable of int
An array of two or three integers representing the number of lattice
cells in the x- and y- (and z-) directions, respectively.
lower_left : Iterable of float
The coordinates of the lower-left corner of the lattice. If the lattice
is two-dimensional, only the x- and y-coordinates are specified.
pitch : Iterable of float
Pitch of the lattice in the x, y, and (if applicable) z directions in
cm.
@ -266,7 +271,25 @@ class RectLattice(Lattice):
A universe to fill all space outside the lattice
universes : Iterable of Iterable of openmc.Universe
A two- or three-dimensional list/array of universes filling each element
of the lattice
of the lattice. The first dimension corresponds to the z-direction (if
applicable), the second dimension corresponds to the y-direction, and
the third dimension corresponds to the x-direction. Note that for the
y-direction, a higher index corresponds to a lower physical
y-value. Each z-slice in the array can be thought of as a top-down view
of the lattice.
lower_left : Iterable of float
The Cartesian coordinates of the lower-left corner of the lattice. If
the lattice is two-dimensional, only the x- and y-coordinates are
specified.
indices : list of tuple
A list of all possible (z,y,x) or (y,x) lattice element indices. These
indices correspond to indices in the :attr:`RectLattice.universes`
property.
ndim : int
The number of dimensions of the lattice
shape : Iterable of int
An array of two or three integers representing the number of lattice
cells in the x- and y- (and z-) directions, respectively.
"""
@ -274,7 +297,6 @@ class RectLattice(Lattice):
super(RectLattice, self).__init__(lattice_id, name)
# Initialize Lattice class attributes
self._dimension = None
self._lower_left = None
self._offsets = None
@ -283,7 +305,7 @@ class RectLattice(Lattice):
return False
elif not super(RectLattice, self).__eq__(other):
return False
elif self.dimension != other.dimension:
elif self.shape != other.shape:
return False
elif self.lower_left != other.lower_left:
return False
@ -300,8 +322,8 @@ class RectLattice(Lattice):
string = 'RectLattice\n'
string += '{0: <16}{1}{2}\n'.format('\tID', '=\t', self._id)
string += '{0: <16}{1}{2}\n'.format('\tName', '=\t', self._name)
string += '{0: <16}{1}{2}\n'.format('\tDimension', '=\t',
self._dimension)
string += '{0: <16}{1}{2}\n'.format('\tShape', '=\t',
self.shape)
string += '{0: <16}{1}{2}\n'.format('\tLower Left', '=\t',
self._lower_left)
string += '{0: <16}{1}{2}\n'.format('\tPitch', '=\t', self._pitch)
@ -320,7 +342,7 @@ class RectLattice(Lattice):
string += '{0} '.format(universe._id)
# Add a newline character every time we reach end of row of cells
if (i+1) % self._dimension[-1] == 0:
if (i+1) % self.shape[0] == 0:
string += '\n'
string = string.rstrip('\n')
@ -333,7 +355,7 @@ class RectLattice(Lattice):
string += '{0} '.format(offset)
# Add a newline character when we reach end of row of cells
if (i+1) % self._dimension[-1] == 0:
if (i+1) % self.shape[0] == 0:
string += '\n'
string = string.rstrip('\n')
@ -341,24 +363,29 @@ class RectLattice(Lattice):
return string
@property
def dimension(self):
return self._dimension
def indices(self):
if self.ndim == 2:
return list(np.broadcast(*np.ogrid[
:self.shape[1], :self.shape[0]]))
else:
return list(np.broadcast(*np.ogrid[
:self.shape[2], :self.shape[1], :self.shape[0]]))
@property
def lower_left(self):
return self._lower_left
@property
def ndim(self):
return len(self.pitch)
@property
def offsets(self):
return self._offsets
@dimension.setter
def dimension(self, dimension):
cv.check_type('lattice dimension', dimension, Iterable, Integral)
cv.check_length('lattice dimension', dimension, 2, 3)
for dim in dimension:
cv.check_greater_than('lattice dimension', dim, 0)
self._dimension = dimension
@property
def shape(self):
return self._universes.shape[::-1]
@lower_left.setter
def lower_left(self, lower_left):
@ -379,8 +406,13 @@ class RectLattice(Lattice):
cv.check_greater_than('lattice pitch', dim, 0.0)
self._pitch = pitch
def get_cell_instance(self, path, distribcell_index):
@Lattice.universes.setter
def universes(self, universes):
cv.check_iterable_type('lattice universes', universes, openmc.Universe,
min_depth=2, max_depth=3)
self._universes = np.asarray(universes)
def get_cell_instance(self, path, distribcell_index):
# Extract the lattice element from the path
next_index = path.index('-')
lat_id_indices = path[:next_index]
@ -395,7 +427,7 @@ class RectLattice(Lattice):
lat_z = int(i.split(',')[2]) - 1
# For 2D Lattices
if len(self._dimension) == 2:
if self.ndim == 2:
offset = self._offsets[lat_z, lat_y, lat_x, distribcell_index-1]
offset += self._universes[lat_x][lat_y].get_cell_instance(path,
distribcell_index)
@ -408,6 +440,128 @@ class RectLattice(Lattice):
return offset
def find_element(self, point):
"""Determine index of lattice element and local coordinates for a point
Parameters
----------
point : Iterable of float
Cartesian coordinates of point
Returns
-------
2- or 3-tuple of int
A tuple of the corresponding (x,y,z) lattice element indices
3-tuple of float
Carestian coordinates of the point in the corresponding lattice
element coordinate system
"""
ix = floor((point[0] - self.lower_left[0])/self.pitch[0])
iy = floor((point[1] - self.lower_left[1])/self.pitch[1])
if self.ndim == 2:
idx = (ix, iy)
else:
iz = floor((point[2] - self.lower_left[2])/self.pitch[2])
idx = (ix, iy, iz)
return idx, self.get_local_coordinates(point, idx)
def get_local_coordinates(self, point, idx):
"""Determine local coordinates of a point within a lattice element
Parameters
----------
point : Iterable of float
Cartesian coordinates of point
idx : Iterable of int
(x,y,z) indices of lattice element. If the lattice is 2D, the z
index can be omitted.
Returns
-------
3-tuple of float
Cartesian coordinates of point in the lattice element coordinate
system
"""
x = point[0] - (self.lower_left[0] + (idx[0] + 0.5)*self.pitch[0])
y = point[1] - (self.lower_left[1] + (idx[1] + 0.5)*self.pitch[1])
if self.ndim == 2:
z = point[2]
else:
z = point[2] - (self.lower_left[2] + (idx[2] + 0.5)*self.pitch[2])
return (x, y, z)
def get_universe_index(self, idx):
"""Return index in the universes array corresponding to a lattice element index
Parameters
----------
idx : Iterable of int
Lattice element indices in the :math:`(x,y,z)` coordinate system
Returns
-------
2- or 3-tuple of int
Indices used when setting the :attr:`RectLattice.universes` property
"""
max_y = self.shape[1] - 1
if self.ndim == 2:
x, y = idx
return (max_y - y, x)
else:
x, y, z = idx
return (z, max_y - y, x)
def is_valid_index(self, idx):
"""Determine whether lattice element index is within defined range
Parameters
----------
idx : Iterable of int
Lattice element indices in the :math:`(x,y,z)` coordinate system
Returns
-------
bool
Whether index is valid
"""
if self.ndim == 2:
return (0 <= idx[0] < self.shape[0] and
0 <= idx[1] < self.shape[1])
else:
return (0 <= idx[0] < self.shape[0] and
0 <= idx[1] < self.shape[1] and
0 <= idx[2] < self.shape[2])
def find(self, point):
"""Find cells/universes/lattices which contain a given point
Parameters
----------
point : 3-tuple of float
Cartesian coordinatesof the point
Returns
-------
list
Sequence of universes, cells, and lattices which are traversed to
find the given point
"""
idx, p = self.find_element(point)
if self.is_valid_index(idx):
idx_u = self.get_universe_index(idx)
u = self.universes[idx_u]
else:
if self.outer is not None:
u = self.outer
else:
return []
return [(self, idx)] + u.find(p)
def create_xml_subelement(self, xml_element):
# Determine if XML element already contains subelement for this Lattice
@ -436,7 +590,7 @@ class RectLattice(Lattice):
# Export Lattice cell dimensions
dimension = ET.SubElement(lattice_subelement, "dimension")
dimension.text = ' '.join(map(str, self._dimension))
dimension.text = ' '.join(map(str, self.shape))
# Export Lattice lower left
lower_left = ET.SubElement(lattice_subelement, "lower_left")
@ -446,10 +600,10 @@ class RectLattice(Lattice):
universe_ids = '\n'
# 3D Lattices
if len(self._dimension) == 3:
for z in range(self._dimension[2]):
for y in range(self._dimension[1]):
for x in range(self._dimension[0]):
if self.ndim == 3:
for z in range(self.shape[2]):
for y in range(self.shape[1]):
for x in range(self.shape[0]):
universe = self._universes[z][y][x]
# Append Universe ID to the Lattice XML subelement
@ -466,8 +620,8 @@ class RectLattice(Lattice):
# 2D Lattices
else:
for y in range(self._dimension[1]):
for x in range(self._dimension[0]):
for y in range(self.shape[1]):
for x in range(self.shape[0]):
universe = self._universes[y][x]
# Append Universe ID to Lattice XML subelement
@ -490,7 +644,18 @@ class RectLattice(Lattice):
class HexLattice(Lattice):
"""A lattice consisting of hexagonal prisms.
r"""A lattice consisting of hexagonal prisms.
To completely define a hexagonal lattice, the :attr:`HexLattice.center`,
:attr:`HexLattice.pitch`, :attr:`HexLattice.universes`, and
:attr:`HexLattice.outer` properties need to be set.
Most methods for this class use a natural indexing scheme wherein elements
are assigned an index corresponding to their position relative to skewed
:math:`(x,\alpha,z)` axes as described fully in
:ref:`hexagonal_indexing`. However, note that when universes are assigned to
lattice elements using the :attr:`HexLattice.universes` property, the array
indices do not correspond to natural indices.
Parameters
----------
@ -506,26 +671,31 @@ class HexLattice(Lattice):
Unique identifier for the lattice
name : str
Name of the lattice
num_rings : int
Number of radial ring positions in the xy-plane
num_axial : int
Number of positions along the z-axis.
center : Iterable of float
Coordinates of the center of the lattice. If the lattice does not have
axial sections then only the x- and y-coordinates are specified
pitch : Iterable of float
Pitch of the lattice in cm. The first item in the iterable specifies the
pitch in the radial direction and, if the lattice is 3D, the second item
in the iterable specifies the pitch in the axial direction.
outer : openmc.Universe
A universe to fill all space outside the lattice
universes : Iterable of Iterable of openmc.Universe
universes : Nested Iterable of openmc.Universe
A two- or three-dimensional list/array of universes filling each element
of the lattice. Each sub-list corresponds to one ring of universes and
should be ordered from outermost ring to innermost ring. The universes
within each sub-list are ordered from the "top" and proceed in a
clockwise fashion. The :meth:`HexLattice.show_indices` method can be
used to help figure out indices for this property.
center : Iterable of float
Coordinates of the center of the lattice. If the lattice does not have
axial sections then only the x- and y-coordinates are specified
indices : list of tuple
A list of all possible (z,r,i) or (r,i) lattice element indices that are
possible, where z is the axial index, r is in the ring index (starting
from the outermost ring), and i is the index with a ring starting from
the top and proceeding clockwise.
num_rings : int
Number of radial ring positions in the xy-plane
num_axial : int
Number of positions along the z-axis.
"""
@ -597,17 +767,15 @@ class HexLattice(Lattice):
def center(self):
return self._center
@num_rings.setter
def num_rings(self, num_rings):
cv.check_type('number of rings', num_rings, Integral)
cv.check_greater_than('number of rings', num_rings, 0)
self._num_rings = num_rings
@num_axial.setter
def num_axial(self, num_axial):
cv.check_type('number of axial', num_axial, Integral)
cv.check_greater_than('number of axial', num_axial, 0)
self._num_axial = num_axial
@property
def indices(self):
if self.num_axial is None:
return [(r, i) for r in range(self.num_rings)
for i in range(max(6*(self.num_rings - 1 - r), 1))]
else:
return [(z, r, i) for z in range(self.num_axial)
for r in range(self.num_rings)
for i in range(max(6*(self.num_rings - 1 - r), 1))]
@center.setter
def center(self, center):
@ -625,8 +793,9 @@ class HexLattice(Lattice):
@Lattice.universes.setter
def universes(self, universes):
# Call Lattice.universes parent class setter property
Lattice.universes.fset(self, universes)
cv.check_iterable_type('lattice universes', universes, openmc.Universe,
min_depth=2, max_depth=3)
self._universes = universes
# NOTE: This routine assumes that the user creates a "ragged" list of
# lists, where each sub-list corresponds to one ring of Universes.
@ -649,14 +818,14 @@ class HexLattice(Lattice):
# Set the number of axial positions.
if n_dims == 3:
self.num_axial = len(self._universes)
self._num_axial = len(self._universes)
else:
self._num_axial = None
# Set the number of rings and make sure this number is consistent for
# all axial positions.
if n_dims == 3:
self.num_rings = len(self._universes[0])
self._num_rings = len(self._universes[0])
for rings in self._universes:
if len(rings) != self._num_rings:
msg = 'HexLattice ID={0:d} has an inconsistent number of ' \
@ -664,7 +833,7 @@ class HexLattice(Lattice):
raise ValueError(msg)
else:
self.num_rings = len(self._universes)
self._num_rings = len(self._universes)
# Make sure there are the correct number of elements in each ring.
if n_dims == 3:
@ -705,6 +874,170 @@ class HexLattice(Lattice):
6*(self._num_rings - 1 - r))
raise ValueError(msg)
def find_element(self, point):
r"""Determine index of lattice element and local coordinates for a point
Parameters
----------
point : Iterable of float
Cartesian coordinates of point
Returns
-------
3-tuple of int
Indices of corresponding lattice element in :math:`(x,\alpha,z)`
bases
numpy.ndarray
Carestian coordinates of the point in the corresponding lattice
element coordinate system
"""
# Convert coordinates to skewed bases
x = point[0] - self.center[0]
y = point[1] - self.center[1]
if self._num_axial is None:
iz = 1
else:
z = point[2] - self.center[2]
iz = floor(z/self.pitch[1] + 0.5*self.num_axial)
alpha = y - x/sqrt(3.)
ix = floor(x/(sqrt(0.75) * self.pitch[0]))
ia = floor(alpha/self.pitch[0])
# Check four lattice elements to see which one is closest based on local
# coordinates
d_min = np.inf
for idx in [(ix, ia, iz), (ix + 1, ia, iz), (ix, ia + 1, iz),
(ix + 1, ia + 1, iz)]:
p = self.get_local_coordinates(point, idx)
d = p[0]**2 + p[1]**2
if d < d_min:
d_min = d
idx_min = idx
p_min = p
return idx_min, p_min
def get_local_coordinates(self, point, idx):
r"""Determine local coordinates of a point within a lattice element
Parameters
----------
point : Iterable of float
Cartesian coordinates of point
idx : Iterable of int
Indices of lattice element in :math:`(x,\alpha,z)` bases
Returns
-------
3-tuple of float
Cartesian coordinates of point in the lattice element coordinate
system
"""
x = point[0] - (self.center[0] + sqrt(0.75)*self.pitch[0]*idx[0])
y = point[1] - (self.center[1] + (0.5*idx[0] + idx[1])*self.pitch[0])
if self._num_axial is None:
z = point[2]
else:
z = point[2] - (self.center[2] + (idx[2] + 0.5 - 0.5*self.num_axial)*
self.pitch[1])
return (x, y, z)
def get_universe_index(self, idx):
r"""Return index in the universes array corresponding to a lattice element index
Parameters
----------
idx : Iterable of int
Lattice element indices in the :math:`(x,\alpha,z)` coordinate
system
Returns
-------
2- or 3-tuple of int
Indices used when setting the :attr:`HexLattice.universes` property
"""
# First we determine which ring the index corresponds to.
x = idx[0]
a = idx[1]
z = -a - x
g = max(abs(x), abs(a), abs(z))
# Next we use a clever method to figure out where along the ring we are.
i_ring = self._num_rings - 1 - g
if x >= 0:
if a >= 0:
i_within = x
else:
i_within = 2*g + z
else:
if a <= 0:
i_within = 3*g - x
else:
i_within = 5*g - z
if self.num_axial is None:
return (i_ring, i_within)
else:
return (idx[2], i_ring, i_within)
def is_valid_index(self, idx):
r"""Determine whether lattice element index is within defined range
Parameters
----------
idx : Iterable of int
Lattice element indices in the :math:`(x,\alpha,z)` coordinate
system
Returns
-------
bool
Whether index is valid
"""
x = idx[0]
y = idx[1]
z = 0 - y - x
g = max(abs(x), abs(y), abs(z))
if self.num_axial is None:
return g < self.num_rings
else:
return g < self.num_rings and 0 <= idx[2] < self.num_axial
def find(self, point):
"""Find cells/universes/lattices which contain a given point
Parameters
----------
point : 3-tuple of float
Cartesian coordinatesof the point
Returns
-------
list
Sequence of universes, cells, and lattices which are traversed to
find the given point
"""
idx, p = self.find_element(point)
if self.is_valid_index(idx):
idx_u = self.get_universe_index(idx)
if self.num_axial is None:
u = self.universes[idx_u[0]][idx_u[1]]
else:
u = self.universes[idx_u[0]][idx_u[1]][idx_u[2]]
else:
if self.outer is not None:
u = self.outer
else:
return []
return [(self, idx)] + u.find(p)
def create_xml_subelement(self, xml_element):
# Determine if XML element already contains subelement for this Lattice
path = './hex_lattice[@id=\'{0}\']'.format(self._id)
@ -736,8 +1069,8 @@ class HexLattice(Lattice):
lattice_subelement.set("n_axial", str(self._num_axial))
# Export Lattice cell center
dimension = ET.SubElement(lattice_subelement, "center")
dimension.text = ' '.join(map(str, self._center))
center = ET.SubElement(lattice_subelement, "center")
center.text = ' '.join(map(str, self._center))
# Export the Lattice nested Universe IDs.

120
openmc/material.py
View file

@ -50,14 +50,14 @@ class Material(object):
Units used for `density`. Can be one of 'g/cm3', 'g/cc', 'kg/cm3',
'atom/b-cm', 'atom/cm3', 'sum', or 'macro'. The 'macro' unit only
applies in the case of a multi-group calculation.
elements : collections.OrderedDict
Dictionary whose keys are element names and values are 3-tuples
consisting of an :class:`openmc.Element` instance, the percent density,
and the percent type (atom or weight fraction).
nuclides : collections.OrderedDict
Dictionary whose keys are nuclide names and values are 3-tuples
consisting of an :class:`openmc.Nuclide` instance, the percent density,
and the percent type (atom or weight fraction).
elements : list of tuple
List in which each item is a 3-tuple consisting of an
:class:`openmc.Element` instance, the percent density, and the percent
type ('ao' or 'wo').
nuclides : list of tuple
List in which each item is a 3-tuple consisting of an
:class:`openmc.Nuclide` instance, the percent density, and the percent
type ('ao' or 'wo').
"""
@ -68,19 +68,15 @@ class Material(object):
self._density = None
self._density_units = ''
# An ordered dictionary of Nuclides (order affects OpenMC results)
# Keys - Nuclide names
# Values - tuple (nuclide, percent, percent type)
self._nuclides = OrderedDict()
# A list of tuples (nuclide, percent, percent type)
self._nuclides = []
# The single instance of Macroscopic data present in this material
# (only one is allowed, hence this is different than _nuclides, etc)
self._macroscopic = None
# An ordered dictionary of Elements (order affects OpenMC results)
# Keys - Element names
# Values - tuple (element, percent, percent type)
self._elements = OrderedDict()
# A list of tuples (element, percent, percent type)
self._elements = []
# If specified, a list of tuples of (table name, xs identifier)
self._sab = []
@ -134,10 +130,8 @@ class Material(object):
string += '{0: <16}\n'.format('\tNuclides')
for nuclide in self._nuclides:
percent = self._nuclides[nuclide][1]
percent_type = self._nuclides[nuclide][2]
string += '{0: <16}'.format('\t{0}'.format(nuclide))
for nuclide, percent, percent_type in self._nuclides:
string += '{0: <16}'.format('\t{0.name}.{0.xs}'.format(nuclide))
string += '=\t{0: <12} [{1}]\n'.format(percent, percent_type)
if self._macroscopic is not None:
@ -146,39 +140,12 @@ class Material(object):
string += '{0: <16}\n'.format('\tElements')
for element in self._elements:
percent = self._elements[element][1]
percent_type = self._elements[element][2]
string += '{0: <16}'.format('\t{0}'.format(element))
for element, percent, percent_type in self._elements:
string += '{0: <16}'.format('\t{0.name}.{0.xs}'.format(element))
string += '=\t{0: <12} [{1}]\n'.format(percent, percent_type)
return string
def __deepcopy__(self, memo):
existing = memo.get(id(self))
if existing is None:
# If this is the first time we have tried to copy this object, create a copy
clone = type(self).__new__(type(self))
clone._id = self._id
clone._name = self._name
clone._density = self._density
clone._density_units = self._density_units
clone._nuclides = deepcopy(self._nuclides, memo)
clone._macroscopic = self._macroscopic
clone._elements = deepcopy(self._elements, memo)
clone._sab = deepcopy(self._sab, memo)
clone._convert_to_distrib_comps = self._convert_to_distrib_comps
clone._distrib_otf_file = self._distrib_otf_file
memo[id(self)] = clone
return clone
else:
# If this object has been copied before, return the first copy made
return existing
@property
def id(self):
return self._id
@ -326,7 +293,7 @@ class Material(object):
else:
nuclide = openmc.Nuclide(nuclide)
self._nuclides[nuclide._name] = (nuclide, percent, percent_type)
self._nuclides.append((nuclide, percent, percent_type))
def remove_nuclide(self, nuclide):
"""Remove a nuclide from the material
@ -344,8 +311,9 @@ class Material(object):
raise ValueError(msg)
# If the Material contains the Nuclide, delete it
if nuclide._name in self._nuclides:
del self._nuclides[nuclide._name]
for nuc in self._nuclides:
if nuclide == nuc:
self._nuclides.remove(nuc)
def add_macroscopic(self, macroscopic):
"""Add a macroscopic to the material. This will also set the
@ -464,10 +432,9 @@ class Material(object):
raise NotImplementedError('Expanding natural element based on '
'weight percent is not yet supported.')
for isotope, abundance in element.expand():
self._nuclides[isotope.name] = (
isotope, percent*abundance, percent_type)
self._nuclides.append((isotope, percent*abundance, percent_type))
else:
self._elements[element.name] = (element, percent, percent_type)
self._elements.append((element, percent, percent_type))
def remove_element(self, element):
"""Remove a natural element from the material
@ -479,9 +446,15 @@ class Material(object):
"""
# If the Material contains the Element, delete it
if element._name in self._elements:
del self._elements[element._name]
if not isinstance(nuclide, openmc.Element):
msg = 'Unable to remove "{0}" in Material ID="{1}" ' \
'since it is not an Element'.format(self.id, element)
raise ValueError(msg)
# If the Material contains the Nuclide, delete it
for elm in self._elements:
if element == elm:
self._nuclides.remove(elm)
def add_s_alpha_beta(self, name, xs):
r"""Add an :math:`S(\alpha,\beta)` table to the material
@ -513,10 +486,10 @@ class Material(object):
self._sab.append((name, xs))
def make_isotropic_in_lab(self):
for nuclide_name in self._nuclides:
self._nuclides[nuclide_name][0].scattering = 'iso-in-lab'
for element_name in self._elements:
self._elements[element_name][0].scattering = 'iso-in-lab'
for nuclide, percent, percent_type in self._nuclides:
nuclide.scattering = 'iso-in-lab'
for element, percent, percent_type in self._elements:
element.scattering = 'iso-in-lab'
def get_all_nuclides(self):
"""Returns all nuclides in the material
@ -531,15 +504,10 @@ class Material(object):
nuclides = OrderedDict()
for nuclide_name, nuclide_tuple in self._nuclides.items():
nuclide = nuclide_tuple[0]
density = nuclide_tuple[1]
nuclides[nuclide._name] = (nuclide, density)
for element_name, element_tuple in self._elements.items():
element = element_tuple[0]
density = element_tuple[1]
for nuclide, density, density_type in self._nuclides:
nuclides[nuclide.name] = (nuclide, density)
for element, density, density_type in self._elements:
# Expand natural element into isotopes
for isotope, abundance in element.expand():
nuclides[isotope.name] = (isotope, density*abundance)
@ -548,7 +516,7 @@ class Material(object):
def _get_nuclide_xml(self, nuclide, distrib=False):
xml_element = ET.Element("nuclide")
xml_element.set("name", nuclide[0]._name)
xml_element.set("name", nuclide[0].name)
if not distrib:
if nuclide[2] == 'ao':
@ -575,7 +543,7 @@ class Material(object):
def _get_element_xml(self, element, distrib=False):
xml_element = ET.Element("element")
xml_element.set("name", str(element[0]._name))
xml_element.set("name", str(element[0].name))
if not distrib:
if element[2] == 'ao':
@ -594,7 +562,7 @@ class Material(object):
def _get_nuclides_xml(self, nuclides, distrib=False):
xml_elements = []
for nuclide in nuclides.values():
for nuclide in nuclides:
xml_elements.append(self._get_nuclide_xml(nuclide, distrib))
return xml_elements
@ -602,7 +570,7 @@ class Material(object):
def _get_elements_xml(self, elements, distrib=False):
xml_elements = []
for element in elements.values():
for element in elements:
xml_elements.append(self._get_element_xml(element, distrib))
return xml_elements
@ -650,7 +618,7 @@ class Material(object):
subelement = ET.SubElement(element, "compositions")
comps = []
allnucs = self._nuclides.values() + self._elements.values()
allnucs = self._nuclides + self._elements
dist_per_type = allnucs[0][2]
for nuc, per, typ in allnucs:
if not typ == dist_per_type:
@ -845,4 +813,4 @@ class Materials(cv.CheckedList):
# Write the XML Tree to the materials.xml file
tree = ET.ElementTree(self._materials_file)
tree.write("materials.xml", xml_declaration=True,
encoding='utf-8', method="xml")
encoding='utf-8', method="xml")

23
openmc/mesh.py
View file

@ -1,5 +1,4 @@
from collections import Iterable
import copy
from numbers import Real, Integral
from xml.etree import ElementTree as ET
import sys
@ -83,28 +82,6 @@ class Mesh(object):
else:
return True
def __deepcopy__(self, memo):
existing = memo.get(id(self))
# If this is the first time we have tried to copy this object, create a copy
if existing is None:
clone = type(self).__new__(type(self))
clone._id = self._id
clone._name = self._name
clone._type = self._type
clone._dimension = copy.deepcopy(self._dimension, memo)
clone._lower_left = copy.deepcopy(self._lower_left, memo)
clone._upper_right = copy.deepcopy(self._upper_right, memo)
clone._width = copy.deepcopy(self._width, memo)
memo[id(self)] = clone
return clone
# If this object has been copied before, return the first copy made
else:
return existing
@property
def id(self):
return self._id

116
openmc/mgxs/library.py
View file

@ -34,7 +34,7 @@ class Library(object):
Parameters
----------
openmc_geometry : openmc.Geometry
An geometry which has been initialized with a root universe
A geometry which has been initialized with a root universe
by_nuclide : bool
If true, computes cross sections for each nuclide in each domain
mgxs_types : Iterable of str
@ -460,7 +460,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', 'chi'}
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'}
The type of multi-group cross section object to return
Returns
@ -484,7 +484,7 @@ class Library(object):
cv.check_type('domain', domain, (openmc.Universe, Integral))
# Check that requested domain is included in library
if cv._isinstance(domain, Integral):
if isinstance(domain, Integral):
domain_id = domain
for domain in self.domains:
if domain_id == domain.id:
@ -853,16 +853,29 @@ class Library(object):
mymgxs = self.get_mgxs(domain, 'kappa-fission')
xsdata.set_kappa_fission_mgxs(mymgxs, xs_type=xs_type,
nuclide=[nuclide])
if 'chi' in self.mgxs_types:
mymgxs = self.get_mgxs(domain, 'chi')
xsdata.set_chi_mgxs(mymgxs, xs_type=xs_type, nuclide=[nuclide])
if 'nu-fission' in self.mgxs_types:
mymgxs = self.get_mgxs(domain, 'nu-fission')
# For chi and nu-fission we can either have only a nu-fission matrix
# provided, or vectors of chi and nu-fission provided
if 'nu-fission matrix' in self.mgxs_types:
mymgxs = self.get_mgxs(domain, 'nu-fission matrix')
xsdata.set_nu_fission_mgxs(mymgxs, xs_type=xs_type,
nuclide=[nuclide])
# multiplicity requires scatter and nu-scatter
if ((('scatter matrix' in self.mgxs_types) and
('nu-scatter matrix' in self.mgxs_types))):
else:
if 'chi' in self.mgxs_types:
mymgxs = self.get_mgxs(domain, 'chi')
xsdata.set_chi_mgxs(mymgxs, xs_type=xs_type, nuclide=[nuclide])
if 'nu-fission' in self.mgxs_types:
mymgxs = self.get_mgxs(domain, 'nu-fission')
xsdata.set_nu_fission_mgxs(mymgxs, xs_type=xs_type,
nuclide=[nuclide])
# If multiplicity matrix is available, prefer that
if 'multiplicity matrix' in self.mgxs_types:
mymgxs = self.get_mgxs(domain, 'multiplicity matrix')
xsdata.set_multiplicity_mgxs(mymgxs, xs_type=xs_type,
nuclide=[nuclide])
using_multiplicity = True
# multiplicity wil fall back to using scatter and nu-scatter
elif ((('scatter matrix' in self.mgxs_types) and
('nu-scatter matrix' in self.mgxs_types))):
scatt_mgxs = self.get_mgxs(domain, 'scatter matrix')
nuscatt_mgxs = self.get_mgxs(domain, 'nu-scatter matrix')
xsdata.set_multiplicity_mgxs(nuscatt_mgxs, scatt_mgxs,
@ -926,6 +939,7 @@ class Library(object):
See also
--------
Library.dump_to_file()
Library.create_mg_mode()
"""
@ -975,14 +989,14 @@ class Library(object):
return mgxs_file
def create_mg_library_and_materials(self, xsdata_names=None, xs_ids=None,
material_ids=None):
def create_mg_mode(self, xsdata_names=None, xs_ids=None):
"""Creates an openmc.MGXSLibrary object to contain the MGXS data for the
Multi-Group mode of OpenMC as well as the associated openmc.Materials
objects. This method cannot be used for Library objects with
`Library.by_nuclide == True` since the materials to output would be
problem dependent and thus any Materials object produced by this method
would not be useful.
and openmc.Geometry objects. The created Geometry is the same as that
used to generate the MGXS data, with the only differences being
modifications to point to newly-created Materials which point to the
multi-group data. This method only creates a macroscopic
MGXS Library even if nuclidic tallies are specified in the Library.
Parameters
----------
@ -993,18 +1007,16 @@ class Library(object):
Cross section set identifier (i.e., '71c') for all
data sets (if only str) or for each individual one
(if iterable of str). Defaults to '1m'.
material_ids : None or Iterable of Integral
An optional list of material IDs to pass to the materials in
materials_file. Defaults to `None` implying the materials will be
given an ID number which matches the index of the domain in
`self.domains`
Returns
-------
mgxs_file : openmc.MGXSLibrary
Multi-Group Cross Section File that is ready to be printed to the
file of choice by the user.
materials_file : openmc.Materials
materials : openmc.Materials
Materials file ready to be printed with all the macroscopic data
present within this Library.
geometry : openmc.Geometry
Materials file ready to be printed with all the macroscopic data
present within this Library.
@ -1036,42 +1048,51 @@ class Library(object):
cv.check_iterable_type('xs_ids', xs_ids, basestring)
else:
xs_ids = ['1m' for i in range(len(self.domains))]
if material_ids is not None:
cv.check_iterable_type('material_ids', material_ids, Integral)
xs_type = 'macro'
# Initialize files
# Initialize MGXS File
mgxs_file = openmc.MGXSLibrary(self.energy_groups)
materials = []
macroscopics = []
nuclide = 'total'
# Create a copy of the Geometry to differentiate for these Macroscopics
geometry = copy.deepcopy(self.openmc_geometry)
materials = openmc.Materials()
# Get all Cells from the Geometry for differentiation
all_cells = geometry.get_all_material_cells()
# Create the xsdata object and add it to the mgxs_file
for i, domain in enumerate(self.domains):
# Build & add metadata to XSdata object
if xsdata_names is None:
xsdata_name = 'set' + str(i + 1)
else:
xsdata_name = xsdata_names[i]
xsdata = self.get_xsdata(domain, xsdata_name, nuclide=nuclide,
# Create XSdata and Macroscopic for this domain
xsdata = self.get_xsdata(domain, xsdata_name, nuclide='total',
xs_type=xs_type, xs_id=xs_ids[i])
mgxs_file.add_xsdata(xsdata)
macroscopic = openmc.Macroscopic(name=xsdata_name, xs=xs_ids[i])
macroscopics.append(openmc.Macroscopic(name=xsdata_name,
xs=xs_ids[i]))
if material_ids is not None:
mat_id = material_ids[i]
else:
mat_id = i
materials.append(openmc.Material(name=xsdata_name + '.' +
xs_ids[i], material_id=mat_id))
materials[-1].add_macroscopic(macroscopics[-1])
# Create Material and add to collection
material = openmc.Material(name=xsdata_name + '.' + xs_ids[i])
material.add_macroscopic(macroscopic)
materials.append(material)
materials_file = openmc.Materials(materials)
# Differentiate Geometry with new Material
if self.domain_type == 'material':
# Fill all appropriate Cells with new Material
for cell in all_cells:
if cell.fill.id == domain.id:
cell.fill = material
return (mgxs_file, materials_file)
elif self.domain_type == 'cell':
for cell in all_cells:
if cell.id == domain.id:
cell.fill = material
return mgxs_file, materials, geometry
def check_library_for_openmc_mgxs(self):
"""This routine will check the MGXS Types within a Library
@ -1092,8 +1113,9 @@ class Library(object):
needed to support tallies the user may wish to request.
- A nu-scatter matrix is required.
- Having a multiplicity matrix is preferred.
- Having both nu-scatter (of any order) and scatter
(at least isotropic) matrices is preferred
(at least isotropic) matrices is the second choice.
- If only nu-scatter, need total (not transport), to
be used in adjusting absorption
(i.e., reduced_abs = tot - nuscatt)
@ -1101,6 +1123,7 @@ class Library(object):
See also
--------
Library.create_mg_library()
Library.create_mg_mode()
"""
@ -1116,9 +1139,10 @@ class Library(object):
msg = '"nu-scatter matrix" MGXS type is required but not provided.'
warn(msg)
else:
# Ok, now see the status of scatter
if 'scatter matrix' not in self.mgxs_types:
# We dont have both nu-scatter and scatter, therefore
# Ok, now see the status of scatter and/or multiplicity
if ((('scatter matrix' not in self.mgxs_types) and
('multiplicity matrix' not in self.mgxs_types))):
# We dont have data needed for multiplicity matrix, therefore
# we need total, and not transport.
if 'total' not in self.mgxs_types:
error_flag = True

1195
openmc/mgxs/mgxs.py
View file

File diff suppressed because it is too large Load diff

88
openmc/mgxs_library.py
View file

@ -131,6 +131,8 @@ class XSdata(object):
num_polar : int
Number of equal width angular bins that the polar angular domain is
subdivided into. This only applies when ``representation`` is "angle".
use_chi : bool
Whether or not a chi vector or nu-fission matrix was used.
vector_shape : iterable of int
Dimensionality of vector multi-group cross sections (e.g., the total
cross section). The return result depends on the value of
@ -292,6 +294,10 @@ class XSdata(object):
def num_azimuthal(self):
return self._num_azimuthal
@property
def use_chi(self):
return self._use_chi
@property
def total(self):
return self._total
@ -461,6 +467,11 @@ class XSdata(object):
check_greater_than('num_azimuthal', num_azimuthal, 0)
self._num_azimuthal = num_azimuthal
@use_chi.setter
def use_chi(self, use_chi):
check_type('use_chi', use_chi, bool)
self._use_chi = use_chi
@total.setter
def total(self, total):
check_type('total', total, Iterable, expected_iter_type=Real)
@ -512,9 +523,9 @@ class XSdata(object):
@chi.setter
def chi(self, chi):
if self._use_chi is not None:
if not self._use_chi:
msg = 'Providing chi when nu_fission already provided as a' \
if self.use_chi is not None:
if not self.use_chi:
msg = 'Providing "chi" when "nu-fission" already provided as a' \
'matrix'
raise ValueError(msg)
@ -529,8 +540,8 @@ class XSdata(object):
self._chi = npchi
if self._use_chi is not None:
self._use_chi = True
if self.use_chi is not None:
self.use_chi = True
@scatter.setter
def scatter(self, scatter):
@ -574,8 +585,8 @@ class XSdata(object):
check_iterable_type('nu_fission', npnu_fission, Real,
max_depth=len(npnu_fission.shape))
if self._use_chi is not None:
if self._use_chi:
if self.use_chi is not None:
if self.use_chi:
check_value('nu_fission shape', npnu_fission.shape,
[self.vector_shape])
else:
@ -587,9 +598,9 @@ class XSdata(object):
# Find out if we have a nu-fission matrix or vector
# and set a flag to allow other methods to check this later.
if npnu_fission.shape == self.vector_shape:
self._use_chi = True
self.use_chi = True
else:
self._use_chi = False
self.use_chi = False
self._nu_fission = npnu_fission
if np.sum(self._nu_fission) > 0.0:
@ -728,10 +739,8 @@ class XSdata(object):
"""
# The NuFissionXS class does not have the capability to produce
# a fission matrix and therefore if this path is pursued, we know
# chi must be used.
check_type('nu_fission', nu_fission, openmc.mgxs.NuFissionXS)
check_type('nu_fission', nu_fission, (openmc.mgxs.NuFissionXS,
openmc.mgxs.NuFissionMatrixXS))
check_value('energy_groups', nu_fission.energy_groups,
[self.energy_groups])
check_value('domain_type', nu_fission.domain_type,
@ -744,7 +753,10 @@ class XSdata(object):
msg = 'Angular-Dependent MGXS have not yet been implemented'
raise ValueError(msg)
self._use_chi = True
if isinstance(nu_fission, openmc.mgxs.NuFissionMatrixXS):
self.use_chi = False
else:
self.use_chi = True
if np.sum(self._nu_fission) > 0.0:
self._fissionable = True
@ -809,8 +821,8 @@ class XSdata(object):
"""
if self._use_chi is not None:
if not self._use_chi:
if self.use_chi is not None:
if not self.use_chi:
msg = 'Providing chi when nu_fission already provided as a ' \
'matrix!'
raise ValueError(msg)
@ -827,8 +839,8 @@ class XSdata(object):
msg = 'Angular-Dependent MGXS have not yet been implemented'
raise ValueError(msg)
if self._use_chi is not None:
self._use_chi = True
if self.use_chi is not None:
self.use_chi = True
def set_scatter_mgxs(self, scatter, nuclide='total', xs_type='macro'):
"""This method allows for an openmc.mgxs.ScatterMatrixXS
@ -891,9 +903,10 @@ class XSdata(object):
msg = 'Angular-Dependent MGXS have not yet been implemented'
raise ValueError(msg)
def set_multiplicity_mgxs(self, nuscatter, scatter, nuclide='total',
def set_multiplicity_mgxs(self, nuscatter, scatter=None, nuclide='total',
xs_type='macro'):
"""This method allows for an openmc.mgxs.NuScatterMatrixXS and
"""This method allows for either the direct use of only an
openmc.mgxs.MultiplicityMatrixXS OR an openmc.mgxs.NuScatterMatrixXS and
openmc.mgxs.ScatterMatrixXS to be used to set the scattering
multiplicity for this XSdata object. Multiplicity,
in OpenMC parlance, is a factor used to account for the production
@ -903,9 +916,10 @@ class XSdata(object):
Parameters
----------
nuscatter: openmc.mgxs.NuScatterMatrixXS
MGXS Object containing the nu-scattering matrix cross section
for the domain of interest.
nuscatter: {openmc.mgxs.NuScatterMatrixXS,
openmc.mgxs.MultiplicityMatrixXS}
MGXS Object containing the matrix cross section for the domain
of interest.
scatter: openmc.mgxs.ScatterMatrixXS
MGXS Object containing the scattering matrix cross section
for the domain of interest.
@ -923,23 +937,33 @@ class XSdata(object):
"""
check_type('nuscatter', nuscatter, openmc.mgxs.NuScatterMatrixXS)
check_type('scatter', scatter, openmc.mgxs.ScatterMatrixXS)
check_type('nuscatter', nuscatter, (openmc.mgxs.NuScatterMatrixXS,
openmc.mgxs.MultiplicityMatrixXS))
check_value('energy_groups', nuscatter.energy_groups,
[self.energy_groups])
check_value('energy_groups', scatter.energy_groups,
[self.energy_groups])
check_value('domain_type', nuscatter.domain_type,
['universe', 'cell', 'material'])
check_value('domain_type', scatter.domain_type,
['universe', 'cell', 'material'])
if scatter is not None:
check_type('scatter', scatter, openmc.mgxs.ScatterMatrixXS)
if isinstance(nuscatter, openmc.mgxs.MultiplicityMatrixXS):
msg = 'Either an MultiplicityMatrixXS object must be passed ' \
'for "nuscatter" or the "scatter" argument must be ' \
'provided.'
raise ValueError(msg)
check_value('energy_groups', scatter.energy_groups,
[self.energy_groups])
check_value('domain_type', scatter.domain_type,
['universe', 'cell', 'material'])
if self.representation is 'isotropic':
nuscatt = nuscatter.get_xs(nuclides=nuclide,
xs_type=xs_type, moment=0)
scatt = scatter.get_xs(nuclides=nuclide,
xs_type=xs_type, moment=0)
self._multiplicity = np.divide(nuscatt, scatt)
if isinstance(nuscatter, openmc.mgxs.MultiplicityMatrixXS):
self._multiplicity = nuscatt
else:
scatt = scatter.get_xs(nuclides=nuclide,
xs_type=xs_type, moment=0)
self._multiplicity = np.divide(nuscatt, scatt)
elif self.representation is 'angle':
msg = 'Angular-Dependent MGXS have not yet been implemented'
raise ValueError(msg)

1
openmc/model/__init__.py Normal file
View file

@ -0,0 +1 @@
from .triso import *

155
openmc/model/triso.py Normal file
View file

@ -0,0 +1,155 @@
import copy
from collections import Iterable
from numbers import Real
import warnings
import numpy as np
import openmc
import openmc.checkvalue as cv
class TRISO(openmc.Cell):
"""Tristructural-isotopic (TRISO) micro fuel particle
Parameters
----------
outer_radius : float
Outer radius of TRISO particle
fill : openmc.Universe
Universe which contains all layers of the TRISO particle
center : Iterable of float
Cartesian coordinates of the center of the TRISO particle in cm
Attributes
----------
id : int
Unique identifier for the TRISO cell
name : str
Name of the TRISO cell
center : numpy.ndarray
Cartesian coordinates of the center of the TRISO particle in cm
fill : openmc.Universe
Universe that contains the TRISO layers
region : openmc.Region
Region of space within the TRISO particle
"""
def __init__(self, outer_radius, fill, center=(0., 0., 0.)):
self._surface = openmc.Sphere(R=outer_radius)
super(TRISO, self).__init__(fill=fill, region=-self._surface)
self.center = np.asarray(center)
@property
def center(self):
return self._center
@center.setter
def center(self, center):
cv.check_type('TRISO center', center, Iterable, Real)
self._surface.x0 = center[0]
self._surface.y0 = center[1]
self._surface.z0 = center[2]
self.translation = center
self._center = center
def classify(self, lattice):
"""Determine lattice element indices which might contain the TRISO particle.
Parameters
----------
lattice : openmc.RectLattice
Lattice to check
Returns
-------
list of tuple
(z,y,x) lattice element indices which might contain the TRISO
particle.
"""
ll, ur = self.region.bounding_box
if lattice.ndim == 2:
(i_min, j_min), p = lattice.find_element(ll)
(i_max, j_max), p = lattice.find_element(ur)
return list(np.broadcast(*np.ogrid[
j_min:j_max+1, i_min:i_max+1]))
else:
(i_min, j_min, k_min), p = lattice.find_element(ll)
(i_max, j_max, k_max), p = lattice.find_element(ur)
return list(np.broadcast(*np.ogrid[
k_min:k_max+1, j_min:j_max+1, i_min:i_max+1]))
def create_triso_lattice(trisos, lower_left, pitch, shape, background):
"""Create a lattice containing TRISO particles for optimized tracking.
Parameters
----------
trisos : list of openmc.model.TRISO
List of TRISO particles to put in lattice
lower_left : Iterable of float
Lower-left Cartesian coordinates of the lattice
pitch : Iterable of float
Pitch of the lattice elements in the x-, y-, and z-directions
shape : Iterable of float
Number of lattice elements in the x-, y-, and z-directions
background : openmc.Material
A background material that is used anywhere within the lattice but
outside a TRISO particle
Returns
-------
lattice : openmc.RectLattice
A lattice containing the TRISO particles
"""
lattice = openmc.RectLattice()
lattice.lower_left = lower_left
lattice.pitch = pitch
indices = list(np.broadcast(*np.ogrid[:shape[2], :shape[1], :shape[0]]))
triso_locations = {idx: [] for idx in indices}
for t in trisos:
for idx in t.classify(lattice):
if idx in sorted(triso_locations):
# Create copy of TRISO particle with materials preserved and
# different cell/surface IDs
t_copy = copy.deepcopy(t)
t_copy.id = None
t_copy.fill = t.fill
t_copy._surface.id = None
triso_locations[idx].append(t_copy)
else:
warnings.warn('TRISO particle is partially or completely '
'outside of the lattice.')
# Create universes
universes = np.empty(shape[::-1], dtype=openmc.Universe)
for idx, triso_list in sorted(triso_locations.items()):
if len(triso_list) > 0:
outside_trisos = openmc.Intersection(*[~t.region for t in triso_list])
background_cell = openmc.Cell(fill=background, region=outside_trisos)
else:
background_cell = openmc.Cell(fill=background)
u = openmc.Universe()
u.add_cell(background_cell)
for t in triso_list:
u.add_cell(t)
iz, iy, ix = idx
t.center = lattice.get_local_coordinates(t.center, (ix, iy, iz))
if len(shape) == 2:
universes[-1 - idx[0], idx[1]] = u
else:
universes[idx[0], -1 - idx[1], idx[2]] = u
lattice.universes = universes
# Set outer universe
background_cell = openmc.Cell(fill=background)
lattice.outer = openmc.Universe(cells=[background_cell])
return lattice

31
openmc/opencg_compatible.py
View file

@ -1,4 +1,5 @@
import copy
import operator
import numpy as np
@ -9,8 +10,6 @@ except ImportError:
raise ImportError(msg)
import openmc
from openmc.region import Intersection
from openmc.surface import Halfspace
import openmc.checkvalue as cv
@ -467,13 +466,13 @@ def get_opencg_cell(openmc_cell):
# half-spaces, i.e., no complex cells.
region = openmc_cell.region
if region is not None:
if isinstance(region, Halfspace):
if isinstance(region, openmc.Halfspace):
surface = region.surface
halfspace = -1 if region.side == '-' else 1
opencg_cell.add_surface(get_opencg_surface(surface), halfspace)
elif isinstance(region, Intersection):
elif isinstance(region, openmc.Intersection):
for node in region.nodes:
if not isinstance(node, Halfspace):
if not isinstance(node, openmc.Halfspace):
raise NotImplementedError("Complex cells not yet "
"supported in OpenCG.")
surface = node.surface
@ -697,12 +696,13 @@ def get_openmc_cell(opencg_cell):
translation = np.asarray(opencg_cell.translation, dtype=np.float64)
openmc_cell.translation = translation
surfaces = opencg_cell.surfaces
for surface_id in surfaces:
surface = surfaces[surface_id][0]
halfspace = surfaces[surface_id][1]
openmc_cell.add_surface(get_openmc_surface(surface), halfspace)
surfaces = []
operators = []
for surface, halfspace in opencg_cell.surfaces.values():
surfaces.append(get_openmc_surface(surface))
operators.append(operator.neg if halfspace == -1 else operator.pos)
openmc_cell.region = openmc.Intersection(
*[op(s) for op, s in zip(operators, surfaces)])
# Add the OpenMC Cell to the global collection of all OpenMC Cells
OPENMC_CELLS[cell_id] = openmc_cell
@ -861,8 +861,8 @@ def get_opencg_lattice(openmc_lattice):
universes = new_universes
# Initialize an empty array for the OpenCG nested Universes in this Lattice
universe_array = np.ndarray(tuple(np.array(dimension)[::-1]),
dtype=opencg.Universe)
universe_array = np.empty(tuple(np.array(dimension)[::-1]),
dtype=opencg.Universe)
# Create OpenCG Universes for each unique nested Universe in this Lattice
unique_universes = openmc_lattice.get_unique_universes()
@ -929,8 +929,8 @@ def get_openmc_lattice(opencg_lattice):
outer = opencg_lattice.outside
# Initialize an empty array for the OpenMC nested Universes in this Lattice
universe_array = np.ndarray(tuple(np.array(dimension)[::-1]),
dtype=openmc.Universe)
universe_array = np.empty(tuple(np.array(dimension)[::-1]),
dtype=openmc.Universe)
# Create OpenMC Universes for each unique nested Universe in this Lattice
unique_universes = opencg_lattice.get_unique_universes()
@ -953,7 +953,6 @@ def get_openmc_lattice(opencg_lattice):
np.array(dimension, dtype=np.float64))) / -2.0
openmc_lattice = openmc.RectLattice(lattice_id=lattice_id)
openmc_lattice.dimension = dimension
openmc_lattice.pitch = width
openmc_lattice.universes = universe_array
openmc_lattice.lower_left = lower_left

5
openmc/particle_restart.py
View file

@ -37,6 +37,11 @@ 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

62
openmc/region.py
View file

@ -28,6 +28,10 @@ class Region(object):
def __invert__(self):
return Complement(self)
@abstractmethod
def __contains__(self, point):
return False
@abstractmethod
def __str__(self):
return ''
@ -219,7 +223,7 @@ class Intersection(Region):
Attributes
----------
nodes : tuple of openmc.Region
nodes : list of openmc.Region
Regions to take the intersection of
bounding_box : tuple of numpy.array
Lower-left and upper-right coordinates of an axis-aligned bounding box
@ -229,6 +233,26 @@ class Intersection(Region):
def __init__(self, *nodes):
self.nodes = list(nodes)
def __iter__(self):
for n in self.nodes:
yield n
def __contains__(self, point):
"""Check whether a point is contained in the region.
Parameters
----------
point : 3-tuple of float
Cartesian coordinates, :math:`(x',y',z')`, of the point
Returns
-------
bool
Whether the point is in the region
"""
return all(point in n for n in self.nodes)
def __str__(self):
return '(' + ' '.join(map(str, self.nodes)) + ')'
@ -281,6 +305,26 @@ class Union(Region):
def __init__(self, *nodes):
self.nodes = list(nodes)
def __iter__(self):
for n in self.nodes:
yield n
def __contains__(self, point):
"""Check whether a point is contained in the region.
Parameters
----------
point : 3-tuple of float
Cartesian coordinates, :math:`(x',y',z')`, of the point
Returns
-------
bool
Whether the point is in the region
"""
return any(point in n for n in self.nodes)
def __str__(self):
return '(' + ' | '.join(map(str, self.nodes)) + ')'
@ -336,6 +380,22 @@ class Complement(Region):
def __init__(self, node):
self.node = node
def __contains__(self, point):
"""Check whether a point is contained in the region.
Parameters
----------
point : 3-tuple of float
Cartesian coordinates, :math:`(x',y',z')`, of the point
Returns
-------
bool
Whether the point is in the region
"""
return point not in self.node
def __str__(self):
return '~' + str(self.node)

38
openmc/statepoint.py
View file

@ -106,6 +106,11 @@ 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
@ -497,13 +502,18 @@ class StatePoint(object):
self.tallies[tally_id].sparse = self.sparse
def get_tally(self, scores=[], filters=[], nuclides=[],
name=None, id=None, estimator=None):
name=None, id=None, estimator=None, exact_filters=False,
exact_nuclides=False, exact_scores=False):
"""Finds and returns a Tally object with certain properties.
This routine searches the list of Tallies and returns the first Tally
found which satisfies all of the input parameters.
NOTE: The input parameters do not need to match the complete Tally
specification and may only represent a subset of the Tally's properties.
NOTE: If any of the "exact" parameters are False (default), the input
parameters do not need to match the complete Tally specification and
may only represent a subset of the Tally's properties. If an "exact"
parameter is True then number of scores, filters, or nuclides in the
parameters must precisely match those of any matching Tally.
Parameters
----------
@ -519,6 +529,18 @@ class StatePoint(object):
The id specified for the Tally (default is None).
estimator: str, optional
The type of estimator ('tracklength', 'analog'; default is None).
exact_filters : bool
If True, the number of filters in the parameters must be identical
to those in the matching Tally. If False (default), the filters in
the parameters may be a subset of those in the matching Tally.
exact_nuclides : bool
If True, the number of nuclides in the parameters must be identical
to those in the matching Tally. If False (default), the nuclides in
the parameters may be a subset of those in the matching Tally.
exact_scores : bool
If True, the number of scores in the parameters must be identical
to those in the matching Tally. If False (default), the scores
in the parameters may be a subset of those in the matching Tally.
Returns
-------
@ -547,7 +569,15 @@ class StatePoint(object):
continue
# Determine if Tally has queried estimator
if estimator and not estimator == test_tally.estimator:
if estimator and estimator != test_tally.estimator:
continue
# The number of filters, nuclides and scores must exactly match
if exact_scores and len(scores) != test_tally.num_scores:
continue
if exact_nuclides and len(nuclides) != test_tally.num_nuclides:
continue
if exact_filters and len(filters) != test_tally.num_filters:
continue
# Determine if Tally has the queried score(s)

4
openmc/stats/univariate.py
View file

@ -66,14 +66,14 @@ class Discrete(Univariate):
@x.setter
def x(self, x):
if cv._isinstance(x, Real):
if isinstance(x, Real):
x = [x]
cv.check_type('discrete values', x, Iterable, Real)
self._x = x
@p.setter
def p(self, p):
if cv._isinstance(p, Real):
if isinstance(p, Real):
p = [p]
cv.check_type('discrete probabilities', p, Iterable, Real)
for pk in p:

21
openmc/summary.py
View file

@ -26,6 +26,10 @@ class Summary(object):
# 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()
@ -362,7 +366,6 @@ class Summary(object):
# Create the Lattice
lattice = openmc.RectLattice(lattice_id=lattice_id, name=name)
lattice.dimension = tuple(dimension)
lattice.lower_left = lower_left
lattice.pitch = pitch
@ -372,7 +375,7 @@ class Summary(object):
# Build array of Universe pointers for the Lattice
universes = \
np.ndarray(tuple(universe_ids.shape), dtype=openmc.Universe)
np.empty(tuple(universe_ids.shape), dtype=openmc.Universe)
for z in range(universe_ids.shape[0]):
for y in range(universe_ids.shape[1]):
@ -407,8 +410,6 @@ class Summary(object):
# Create the Lattice
lattice = openmc.HexLattice(lattice_id=lattice_id, name=name)
lattice.num_rings = n_rings
lattice.num_axial = n_axial
lattice.center = center
lattice.pitch = pitch
@ -421,12 +422,12 @@ class Summary(object):
# (x, alpha, z) to the Python API's format of a ragged nested
# list of (z, ring, theta).
universes = []
for z in range(lattice.num_axial):
for z in range(n_axial):
# Add a list for this axial level.
universes.append([])
x = lattice.num_rings - 1
a = 2*lattice.num_rings - 2
for r in range(lattice.num_rings - 1, 0, -1):
x = n_rings - 1
a = 2*n_rings - 2
for r in range(n_rings - 1, 0, -1):
# Add a list for this ring.
universes[-1].append([])
@ -500,13 +501,13 @@ class Summary(object):
# Retrieve the object corresponding to the fill type and ID
if fill_type == 'normal':
if isinstance(fill_id, Iterable):
fill = [self.get_material_by_id(mat) if mat > 0 else 'void'
fill = [self.get_material_by_id(mat) if mat > 0 else None
for mat in fill_id]
else:
if fill_id > 0:
fill = self.get_material_by_id(fill_id)
else:
fill = 'void'
fill = None
elif fill_type == 'universe':
fill = self.get_universe_by_id(fill_id)
else:

331
openmc/surface.py
View file

@ -38,7 +38,9 @@ class Surface(object):
boundary_type : {'transmission, 'vacuum', 'reflective', 'periodic'}, optional
Boundary condition that defines the behavior for particles hitting the
surface. Defaults to transmissive boundary condition where particles
freely pass through the surface.
freely pass through the surface. Note that periodic boundary conditions
can only be applied to x-, y-, and z-planes, and only axis-aligned
periodicity is supported.
name : str, optional
Name of the surface. If not specified, the name will be the empty
string.
@ -193,7 +195,7 @@ class Plane(Surface):
surface_id : int, optional
Unique identifier for the surface. If not specified, an identifier will
automatically be assigned.
boundary_type : {'transmission, 'vacuum', 'reflective', 'periodic'}, optional
boundary_type : {'transmission, 'vacuum', 'reflective'}, optional
Boundary condition that defines the behavior for particles hitting the
surface. Defaults to transmissive boundary condition where particles
freely pass through the surface.
@ -218,9 +220,12 @@ class Plane(Surface):
The 'C' parameter for the plane
d : float
The 'D' parameter for the plane
boundary_type : {'transmission, 'vacuum', 'reflective', 'periodic'}
boundary_type : {'transmission, 'vacuum', 'reflective'}
Boundary condition that defines the behavior for particles hitting the
surface.
periodic_surface : openmc.Surface
If a periodic boundary condition is used, the surface with which this
one is periodic with
coefficients : dict
Dictionary of surface coefficients
id : int
@ -238,6 +243,7 @@ class Plane(Surface):
self._type = 'plane'
self._coeff_keys = ['A', 'B', 'C', 'D']
self._periodic_surface = None
self.a = A
self.b = B
self.c = C
@ -259,6 +265,10 @@ class Plane(Surface):
def d(self):
return self.coefficients['D']
@property
def periodic_surface(self):
return self._periodic_surface
@a.setter
def a(self, A):
check_type('A coefficient', A, Real)
@ -279,6 +289,40 @@ class Plane(Surface):
check_type('D coefficient', D, Real)
self._coefficients['D'] = D
@periodic_surface.setter
def periodic_surface(self, periodic_surface):
check_type('periodic surface', periodic_surface, Plane)
self._periodic_surface = periodic_surface
periodic_surface._periodic_surface = self
def evaluate(self, point):
"""Evaluate the surface equation at a given point.
Parameters
----------
point : 3-tuple of float
The Cartesian coordinates, :math:`(x',y',z')`, at which the surface
equation should be evaluated.
Returns
-------
float
:math:`Ax' + By' + Cz' - d`
"""
x, y, z = point
return self.a*x + self.b*y + self.c*z - self.d
def create_xml_subelement(self):
element = super(Plane, self).create_xml_subelement()
# Add periodic surface pair information
if self.boundary_type == 'periodic':
if self.periodic_surface is not None:
element.set("periodic_surface_id", str(self.periodic_surface.id))
return element
class XPlane(Plane):
"""A plane perpendicular to the x axis of the form :math:`x - x_0 = 0`
@ -291,7 +335,8 @@ class XPlane(Plane):
boundary_type : {'transmission, 'vacuum', 'reflective', 'periodic'}, optional
Boundary condition that defines the behavior for particles hitting the
surface. Defaults to transmissive boundary condition where particles
freely pass through the surface.
freely pass through the surface. Only axis-aligned periodicity is
supported, i.e., x-planes can only be paired with x-planes.
x0 : float, optional
Location of the plane. Defaults to 0.
name : str, optional
@ -304,6 +349,9 @@ class XPlane(Plane):
boundary_type : {'transmission, 'vacuum', 'reflective', 'periodic'}
Boundary condition that defines the behavior for particles hitting the
surface.
periodic_surface : openmc.Surface
If a periodic boundary condition is used, the surface with which this
one is periodic with
coefficients : dict
Dictionary of surface coefficients
id : int
@ -363,6 +411,23 @@ class XPlane(Plane):
return (np.array([self.x0, -np.inf, -np.inf]),
np.array([np.inf, np.inf, np.inf]))
def evaluate(self, point):
"""Evaluate the surface equation at a given point.
Parameters
----------
point : 3-tuple of float
The Cartesian coordinates, :math:`(x',y',z')`, at which the surface
equation should be evaluated.
Returns
-------
float
:math:`x' - x_0`
"""
return point[0] - self.x0
class YPlane(Plane):
"""A plane perpendicular to the y axis of the form :math:`y - y_0 = 0`
@ -375,7 +440,8 @@ class YPlane(Plane):
boundary_type : {'transmission, 'vacuum', 'reflective', 'periodic'}, optional
Boundary condition that defines the behavior for particles hitting the
surface. Defaults to transmissive boundary condition where particles
freely pass through the surface.
freely pass through the surface. Only axis-aligned periodicity is
supported, i.e., x-planes can only be paired with x-planes.
y0 : float, optional
Location of the plane
name : str, optional
@ -388,6 +454,9 @@ class YPlane(Plane):
boundary_type : {'transmission, 'vacuum', 'reflective', 'periodic'}
Boundary condition that defines the behavior for particles hitting the
surface.
periodic_surface : openmc.Surface
If a periodic boundary condition is used, the surface with which this
one is periodic with
coefficients : dict
Dictionary of surface coefficients
id : int
@ -448,6 +517,23 @@ class YPlane(Plane):
return (np.array([-np.inf, self.y0, -np.inf]),
np.array([np.inf, np.inf, np.inf]))
def evaluate(self, point):
"""Evaluate the surface equation at a given point.
Parameters
----------
point : 3-tuple of float
The Cartesian coordinates, :math:`(x',y',z')`, at which the surface
equation should be evaluated.
Returns
-------
float
:math:`y' - y_0`
"""
return point[1] - self.y0
class ZPlane(Plane):
"""A plane perpendicular to the z axis of the form :math:`z - z_0 = 0`
@ -460,7 +546,8 @@ class ZPlane(Plane):
boundary_type : {'transmission, 'vacuum', 'reflective', 'periodic'}, optional
Boundary condition that defines the behavior for particles hitting the
surface. Defaults to transmissive boundary condition where particles
freely pass through the surface.
freely pass through the surface. Only axis-aligned periodicity is
supported, i.e., x-planes can only be paired with x-planes.
z0 : float, optional
Location of the plane. Defaults to 0.
name : str, optional
@ -473,6 +560,9 @@ class ZPlane(Plane):
boundary_type : {'transmission, 'vacuum', 'reflective', 'periodic'}
Boundary condition that defines the behavior for particles hitting the
surface.
periodic_surface : openmc.Surface
If a periodic boundary condition is used, the surface with which this
one is periodic with
coefficients : dict
Dictionary of surface coefficients
id : int
@ -533,6 +623,23 @@ class ZPlane(Plane):
return (np.array([-np.inf, -np.inf, self.z0]),
np.array([np.inf, np.inf, np.inf]))
def evaluate(self, point):
"""Evaluate the surface equation at a given point.
Parameters
----------
point : 3-tuple of float
The Cartesian coordinates, :math:`(x',y',z')`, at which the surface
equation should be evaluated.
Returns
-------
float
:math:`z' - z_0`
"""
return point[2] - self.z0
class Cylinder(Surface):
"""A cylinder whose length is parallel to the x-, y-, or z-axis.
@ -542,7 +649,7 @@ class Cylinder(Surface):
surface_id : int, optional
Unique identifier for the surface. If not specified, an identifier will
automatically be assigned.
boundary_type : {'transmission, 'vacuum', 'reflective', 'periodic'}, optional
boundary_type : {'transmission, 'vacuum', 'reflective'}, optional
Boundary condition that defines the behavior for particles hitting the
surface. Defaults to transmissive boundary condition where particles
freely pass through the surface.
@ -556,7 +663,7 @@ class Cylinder(Surface):
----------
r : float
Radius of the cylinder
boundary_type : {'transmission, 'vacuum', 'reflective', 'periodic'}
boundary_type : {'transmission, 'vacuum', 'reflective'}
Boundary condition that defines the behavior for particles hitting the
surface.
coefficients : dict
@ -598,7 +705,7 @@ class XCylinder(Cylinder):
surface_id : int, optional
Unique identifier for the surface. If not specified, an identifier will
automatically be assigned.
boundary_type : {'transmission, 'vacuum', 'reflective', 'periodic'}, optional
boundary_type : {'transmission, 'vacuum', 'reflective'}, optional
Boundary condition that defines the behavior for particles hitting the
surface. Defaults to transmissive boundary condition where particles
freely pass through the surface.
@ -618,7 +725,7 @@ class XCylinder(Cylinder):
y-coordinate of the center of the cylinder
z0 : float
z-coordinate of the center of the cylinder
boundary_type : {'transmission, 'vacuum', 'reflective', 'periodic'}
boundary_type : {'transmission, 'vacuum', 'reflective'}
Boundary condition that defines the behavior for particles hitting the
surface.
coefficients : dict
@ -691,6 +798,25 @@ class XCylinder(Cylinder):
return (np.array([-np.inf, -np.inf, -np.inf]),
np.array([np.inf, np.inf, np.inf]))
def evaluate(self, point):
"""Evaluate the surface equation at a given point.
Parameters
----------
point : 3-tuple of float
The Cartesian coordinates, :math:`(x',y',z')`, at which the surface
equation should be evaluated.
Returns
-------
float
:math:`(y' - y_0)^2 + (z' - z_0)^2 - R^2`
"""
y = point[1] - self.y0
z = point[2] - self.z0
return y**2 + z**2 - self.r**2
class YCylinder(Cylinder):
"""An infinite cylinder whose length is parallel to the y-axis of the form
@ -701,7 +827,7 @@ class YCylinder(Cylinder):
surface_id : int, optional
Unique identifier for the surface. If not specified, an identifier will
automatically be assigned.
boundary_type : {'transmission, 'vacuum', 'reflective', 'periodic'}, optional
boundary_type : {'transmission, 'vacuum', 'reflective'}, optional
Boundary condition that defines the behavior for particles hitting the
surface. Defaults to transmissive boundary condition where particles
freely pass through the surface.
@ -721,7 +847,7 @@ class YCylinder(Cylinder):
x-coordinate of the center of the cylinder
z0 : float
z-coordinate of the center of the cylinder
boundary_type : {'transmission, 'vacuum', 'reflective', 'periodic'}
boundary_type : {'transmission, 'vacuum', 'reflective'}
Boundary condition that defines the behavior for particles hitting the
surface.
coefficients : dict
@ -794,6 +920,25 @@ class YCylinder(Cylinder):
return (np.array([-np.inf, -np.inf, -np.inf]),
np.array([np.inf, np.inf, np.inf]))
def evaluate(self, point):
"""Evaluate the surface equation at a given point.
Parameters
----------
point : 3-tuple of float
The Cartesian coordinates, :math:`(x',y',z')`, at which the surface
equation should be evaluated.
Returns
-------
float
:math:`(x' - x_0)^2 + (z' - z_0)^2 - R^2`
"""
x = point[0] - self.x0
z = point[2] - self.z0
return x**2 + z**2 - self.r**2
class ZCylinder(Cylinder):
"""An infinite cylinder whose length is parallel to the z-axis of the form
@ -804,7 +949,7 @@ class ZCylinder(Cylinder):
surface_id : int, optional
Unique identifier for the surface. If not specified, an identifier will
automatically be assigned.
boundary_type : {'transmission, 'vacuum', 'reflective', 'periodic'}, optional
boundary_type : {'transmission, 'vacuum', 'reflective'}, optional
Boundary condition that defines the behavior for particles hitting the
surface. Defaults to transmissive boundary condition where particles
freely pass through the surface.
@ -824,7 +969,7 @@ class ZCylinder(Cylinder):
x-coordinate of the center of the cylinder
y0 : float
y-coordinate of the center of the cylinder
boundary_type : {'transmission, 'vacuum', 'reflective', 'periodic'}
boundary_type : {'transmission, 'vacuum', 'reflective'}
Boundary condition that defines the behavior for particles hitting the
surface.
coefficients : dict
@ -897,6 +1042,25 @@ class ZCylinder(Cylinder):
return (np.array([-np.inf, -np.inf, -np.inf]),
np.array([np.inf, np.inf, np.inf]))
def evaluate(self, point):
"""Evaluate the surface equation at a given point.
Parameters
----------
point : 3-tuple of float
The Cartesian coordinates, :math:`(x',y',z')`, at which the surface
equation should be evaluated.
Returns
-------
float
:math:`(x' - x_0)^2 + (y' - y_0)^2 - R^2`
"""
x = point[0] - self.x0
y = point[1] - self.y0
return x**2 + y**2 - self.r**2
class Sphere(Surface):
"""A sphere of the form :math:`(x - x_0)^2 + (y - y_0)^2 + (z - z_0)^2 = R^2`.
@ -906,7 +1070,7 @@ class Sphere(Surface):
surface_id : int, optional
Unique identifier for the surface. If not specified, an identifier will
automatically be assigned.
boundary_type : {'transmission, 'vacuum', 'reflective', 'periodic'}, optional
boundary_type : {'transmission, 'vacuum', 'reflective'}, optional
Boundary condition that defines the behavior for particles hitting the
surface. Defaults to transmissive boundary condition where particles
freely pass through the surface.
@ -931,7 +1095,7 @@ class Sphere(Surface):
z-coordinate of the center of the sphere
R : float
Radius of the sphere
boundary_type : {'transmission, 'vacuum', 'reflective', 'periodic'}
boundary_type : {'transmission, 'vacuum', 'reflective'}
Boundary condition that defines the behavior for particles hitting the
surface.
coefficients : dict
@ -1025,6 +1189,26 @@ class Sphere(Surface):
return (np.array([-np.inf, -np.inf, -np.inf]),
np.array([np.inf, np.inf, np.inf]))
def evaluate(self, point):
"""Evaluate the surface equation at a given point.
Parameters
----------
point : 3-tuple of float
The Cartesian coordinates, :math:`(x',y',z')`, at which the surface
equation should be evaluated.
Returns
-------
float
:math:`(x' - x_0)^2 + (y' - y_0)^2 + (z' - z_0)^2 - R^2`
"""
x = point[0] - self.x0
y = point[1] - self.y0
z = point[2] - self.z0
return x**2 + y**2 + z**2 - self.r**2
class Cone(Surface):
"""A conical surface parallel to the x-, y-, or z-axis.
@ -1034,7 +1218,7 @@ class Cone(Surface):
surface_id : int, optional
Unique identifier for the surface. If not specified, an identifier will
automatically be assigned.
boundary_type : {'transmission, 'vacuum', 'reflective', 'periodic'}, optional
boundary_type : {'transmission, 'vacuum', 'reflective'}, optional
Boundary condition that defines the behavior for particles hitting the
surface. Defaults to transmissive boundary condition where particles
freely pass through the surface.
@ -1059,7 +1243,7 @@ class Cone(Surface):
z-coordinate of the apex
R2 : float
Parameter related to the aperature
boundary_type : {'transmission, 'vacuum', 'reflective', 'periodic'}
boundary_type : {'transmission, 'vacuum', 'reflective'}
Boundary condition that defines the behavior for particles hitting the
surface.
coefficients : dict
@ -1131,7 +1315,7 @@ class XCone(Cone):
surface_id : int, optional
Unique identifier for the surface. If not specified, an identifier will
automatically be assigned.
boundary_type : {'transmission, 'vacuum', 'reflective', 'periodic'}, optional
boundary_type : {'transmission, 'vacuum', 'reflective'}, optional
Boundary condition that defines the behavior for particles hitting the
surface. Defaults to transmissive boundary condition where particles
freely pass through the surface.
@ -1156,7 +1340,7 @@ class XCone(Cone):
z-coordinate of the apex
R2 : float
Parameter related to the aperature
boundary_type : {'transmission, 'vacuum', 'reflective', 'periodic'}
boundary_type : {'transmission, 'vacuum', 'reflective'}
Boundary condition that defines the behavior for particles hitting the
surface.
coefficients : dict
@ -1177,6 +1361,26 @@ class XCone(Cone):
self._type = 'x-cone'
def evaluate(self, point):
"""Evaluate the surface equation at a given point.
Parameters
----------
point : 3-tuple of float
The Cartesian coordinates, :math:`(x',y',z')`, at which the surface
equation should be evaluated.
Returns
-------
float
:math:`(y' - y_0)^2 + (z' - z_0)^2 - R^2(x' - x_0)^2`
"""
x = point[0] - self.x0
y = point[1] - self.y0
z = point[2] - self.z0
return y**2 + z**2 - self.r2*x**2
class YCone(Cone):
"""A cone parallel to the y-axis of the form :math:`(x - x_0)^2 + (z - z_0)^2 =
@ -1187,7 +1391,7 @@ class YCone(Cone):
surface_id : int, optional
Unique identifier for the surface. If not specified, an identifier will
automatically be assigned.
boundary_type : {'transmission, 'vacuum', 'reflective', 'periodic'}, optional
boundary_type : {'transmission, 'vacuum', 'reflective'}, optional
Boundary condition that defines the behavior for particles hitting the
surface. Defaults to transmissive boundary condition where particles
freely pass through the surface.
@ -1212,7 +1416,7 @@ class YCone(Cone):
z-coordinate of the apex
R2 : float
Parameter related to the aperature
boundary_type : {'transmission, 'vacuum', 'reflective', 'periodic'}
boundary_type : {'transmission, 'vacuum', 'reflective'}
Boundary condition that defines the behavior for particles hitting the
surface.
coefficients : dict
@ -1233,6 +1437,26 @@ class YCone(Cone):
self._type = 'y-cone'
def evaluate(self, point):
"""Evaluate the surface equation at a given point.
Parameters
----------
point : 3-tuple of float
The Cartesian coordinates, :math:`(x',y',z')`, at which the surface
equation should be evaluated.
Returns
-------
float
:math:`(x' - x_0)^2 + (z' - z_0)^2 - R^2(y' - y_0)^2`
"""
x = point[0] - self.x0
y = point[1] - self.y0
z = point[2] - self.z0
return x**2 + z**2 - self.r2*y**2
class ZCone(Cone):
"""A cone parallel to the x-axis of the form :math:`(x - x_0)^2 + (y - y_0)^2 =
@ -1243,7 +1467,7 @@ class ZCone(Cone):
surface_id : int, optional
Unique identifier for the surface. If not specified, an identifier will
automatically be assigned.
boundary_type : {'transmission, 'vacuum', 'reflective', 'periodic'}, optional
boundary_type : {'transmission, 'vacuum', 'reflective'}, optional
Boundary condition that defines the behavior for particles hitting the
surface. Defaults to transmissive boundary condition where particles
freely pass through the surface.
@ -1268,7 +1492,7 @@ class ZCone(Cone):
z-coordinate of the apex
R2 : float
Parameter related to the aperature
boundary_type : {'transmission, 'vacuum', 'reflective', 'periodic'}
boundary_type : {'transmission, 'vacuum', 'reflective'}
Boundary condition that defines the behavior for particles hitting the
surface.
coefficients : dict
@ -1289,6 +1513,26 @@ class ZCone(Cone):
self._type = 'z-cone'
def evaluate(self, point):
"""Evaluate the surface equation at a given point.
Parameters
----------
point : 3-tuple of float
The Cartesian coordinates, :math:`(x',y',z')`, at which the surface
equation should be evaluated.
Returns
-------
float
:math:`(x' - x_0)^2 + (y' - y_0)^2 - R^2(z' - z_0)^2`
"""
x = point[0] - self.x0
y = point[1] - self.y0
z = point[2] - self.z0
return x**2 + y**2 - self.r2*z**2
class Quadric(Surface):
"""A surface of the form :math:`Ax^2 + By^2 + Cz^2 + Dxy + Eyz + Fxz + Gx + Hy +
@ -1434,6 +1678,27 @@ class Quadric(Surface):
check_type('k coefficient', k, Real)
self._coefficients['k'] = k
def evaluate(self, point):
"""Evaluate the surface equation at a given point.
Parameters
----------
point : 3-tuple of float
The Cartesian coordinates, :math:`(x',y',z')`, at which the surface
equation should be evaluated.
Returns
-------
float
:math:`Ax'^2 + By'^2 + Cz'^2 + Dx'y' + Ey'z' + Fx'z' + Gx' + Hy' +
Jz' + K = 0`
"""
x, y, z = point
return x*(self.a*x + self.d*y + self.g) + \
y*(self.b*y + self.e*z + self.h) + \
z*(self.c*z + self.f*x + self.j) + self.k
class Halfspace(Region):
"""A positive or negative half-space region.
@ -1479,6 +1744,24 @@ class Halfspace(Region):
def __invert__(self):
return -self.surface if self.side == '+' else +self.surface
def __contains__(self, point):
"""Check whether a point is contained in the half-space.
Parameters
----------
point : 3-tuple of float
Cartesian coordinates, :math:`(x',y',z')`, of the point
Returns
-------
bool
Whether the point is in the half-space
"""
val = self.surface.evaluate(point)
return val >= 0. if self.side == '+' else val < 0.
@property
def surface(self):
return self._surface

62
openmc/tallies.py
View file

@ -129,51 +129,6 @@ class Tally(object):
self._sp_filename = None
self._results_read = False
def __deepcopy__(self, memo):
existing = memo.get(id(self))
# If this is the first time we have tried to copy this object, create a copy
if existing is None:
clone = type(self).__new__(type(self))
clone.id = self.id
clone.name = self.name
clone.estimator = self.estimator
clone.num_realizations = self.num_realizations
clone._sum = copy.deepcopy(self._sum, memo)
clone._sum_sq = copy.deepcopy(self._sum_sq, memo)
clone._mean = copy.deepcopy(self._mean, memo)
clone._std_dev = copy.deepcopy(self._std_dev, memo)
clone._with_summary = self.with_summary
clone._with_batch_statistics = self.with_batch_statistics
clone._derived = self.derived
clone._sparse = self.sparse
clone._sp_filename = self._sp_filename
clone._results_read = self._results_read
clone._filters = []
for self_filter in self.filters:
clone.filters.append(copy.deepcopy(self_filter, memo))
clone._nuclides = []
for nuclide in self.nuclides:
clone.nuclides.append(copy.deepcopy(nuclide, memo))
clone._scores = []
for score in self.scores:
clone.scores.append(score)
clone._triggers = []
for trigger in self.triggers:
clone.triggers.append(trigger)
memo[id(self)] = clone
return clone
# If this object has been copied before, return the first copy made
else:
return existing
def __eq__(self, other):
if not isinstance(other, Tally):
return False
@ -314,6 +269,10 @@ class Tally(object):
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')
@ -827,9 +786,7 @@ class Tally(object):
# Search for each of this tally's scores in the other tally
for score in self.scores:
if score not in other.scores:
all_scores_match = False
else:
if score in other.scores:
no_scores_match = False
# Search for each of the other tally's scores in this tally
@ -2877,7 +2834,7 @@ class Tally(object):
return other * self**-1
def __pos__(self):
def __abs__(self):
"""The absolute value of this tally.
Returns
@ -3471,7 +3428,8 @@ class Tallies(cv.CheckedList):
"""
if not isinstance(tally, Tally):
msg = 'Unable to add a non-Tally "{0}" to the Tallies instance'.format(tally)
msg = 'Unable to add a non-Tally "{0}" to the ' \
'Tallies instance'.format(tally)
raise TypeError(msg)
if merge:
@ -3482,13 +3440,13 @@ class Tallies(cv.CheckedList):
# If a mergeable tally is found
if tally2.can_merge(tally):
# Replace tally 2 with the merged tally
# Replace tally2 with the merged tally
merged_tally = tally2.merge(tally)
self[i] = merged_tally
merged = True
break
# If not mergeable tally was found, simply add this tally
# If no mergeable tally was found, simply add this tally
if not merged:
super(Tallies, self).append(tally)

19
openmc/trigger.py
View file

@ -39,25 +39,6 @@ class Trigger(object):
self.threshold = threshold
self._scores = []
def __deepcopy__(self, memo):
existing = memo.get(id(self))
# If this is first time we have tried to copy this object, create a copy
if existing is None:
clone = type(self).__new__(type(self))
clone._trigger_type = self._trigger_type
clone._threshold = self._threshold
clone.scores = self.scores
memo[id(self)] = clone
return clone
# If this object has been copied before, return the first copy made
else:
return existing
def __eq__(self, other):
if str(self) == str(other):
return True

145
openmc/universe.py
View file

@ -1,6 +1,7 @@
from collections import OrderedDict, Iterable
from numbers import Integral
from xml.etree import ElementTree as ET
import random
import sys
import warnings
@ -12,7 +13,6 @@ import openmc.checkvalue as cv
if sys.version_info[0] >= 3:
basestring = str
# A dictionary for storing IDs of cell elements that have already been written,
# used to optimize the writing process
WRITTEN_IDS = {}
@ -124,6 +124,149 @@ class Universe(object):
else:
self._name = ''
def find(self, point):
"""Find cells/universes/lattices which contain a given point
Parameters
----------
point : 3-tuple of float
Cartesian coordinates of the point
Returns
-------
list
Sequence of universes, cells, and lattices which are traversed to
find the given point
"""
p = np.asarray(point)
for cell in self._cells.values():
if p in cell:
if cell.fill_type in ('material', 'distribmat', 'void'):
return [self, cell]
elif cell.fill_type == 'universe':
if cell.translation is not None:
p -= cell.translation
if cell.rotation is not None:
p[:] = cell.rotation_matrix.dot(p)
return [self, cell] + cell.fill.find(p)
else:
return [self, cell] + cell.fill.find(p)
return []
def plot(self, center=(0., 0., 0.), width=(1., 1.), pixels=(200, 200),
basis='xy', color_by='cell', colors=None, filename=None, seed=None):
"""Display a slice plot of the universe.
Parameters
----------
center : Iterable of float
Coordinates at the center of the plot
width : Iterable of float
Width of the plot in each basis direction
pixels : Iterable of int
Number of pixels to use in each basis direction
basis : {'xy', 'xz', 'yz'}
The basis directions for the plot
color_by : {'cell', 'material'}
Indicate whether the plot should be colored by cell or by material
colors : dict
Assigns colors to specific materials or cells. Keys are instances of
:class:`Cell` or :class:`Material` and values are RGB 3-tuples or RGBA
4-tuples. Red, green, blue, and alpha should all be floats in the
range [0.0, 1.0], for example:
.. code-block:: python
# Make water blue
water = openmc.Cell(fill=h2o)
universe.plot(..., colors={water: (0., 0., 1.))
filename : str or None
Filename to save plot to. If no filename is given, the plot will be
displayed using the currently enabled matplotlib backend.
seed : hashable object or None
Hashable object which is used to seed the random number generator
used to select colors. If None, the generator is seeded from the
current time.
"""
import matplotlib.pyplot as plt
# Seed the random number generator
if seed is not None:
random.seed(seed)
if colors is None:
# Create default dictionary if none supplied
colors = {}
else:
# Convert to RGBA if necessary
for obj, rgb in colors.items():
if len(rgb) == 3:
colors[obj] = rgb + (1.0,)
if basis == 'xy':
x_min = center[0] - 0.5*width[0]
x_max = center[0] + 0.5*width[0]
y_min = center[1] - 0.5*width[1]
y_max = center[1] + 0.5*width[1]
elif basis == 'yz':
# The x-axis will correspond to physical y and the y-axis will correspond to physical z
x_min = center[1] - 0.5*width[0]
x_max = center[1] + 0.5*width[0]
y_min = center[2] - 0.5*width[1]
y_max = center[2] + 0.5*width[1]
elif basis == 'xz':
# The y-axis will correspond to physical z
x_min = center[0] - 0.5*width[0]
x_max = center[0] + 0.5*width[0]
y_min = center[2] - 0.5*width[1]
y_max = center[2] + 0.5*width[1]
# Determine locations to determine cells at
x_coords = np.linspace(x_min, x_max, pixels[0], endpoint=False) + \
0.5*(x_max - x_min)/pixels[0]
y_coords = np.linspace(y_max, y_min, pixels[1], endpoint=False) - \
0.5*(y_max - y_min)/pixels[1]
# Search for locations and assign colors
img = np.zeros(pixels + (4,)) # Use RGBA form
for i, x in enumerate(x_coords):
for j, y in enumerate(y_coords):
if basis == 'xy':
path = self.find((x, y, center[2]))
elif basis == 'yz':
path = self.find((center[0], x, y))
elif basis == 'xz':
path = self.find((x, center[1], y))
if len(path) > 0:
try:
if color_by == 'cell':
obj = path[-1]
elif color_by == 'material':
if path[-1].fill_type == 'material':
obj = path[-1].fill
else:
continue
except AttributeError:
continue
if obj not in colors:
colors[obj] = (random.random(), random.random(),
random.random(), 1.0)
img[j,i,:] = colors[obj]
# Display image
plt.imshow(img, extent=(x_min, x_max, y_min, y_max))
# Show or save the plot
if filename is None:
plt.show()
else:
plt.savefig(filename)
def add_cell(self, cell):
"""Add a cell to the universe.

3
setup.py
View file

@ -11,7 +11,8 @@ except ImportError:
kwargs = {'name': 'openmc',
'version': '0.7.1',
'packages': ['openmc', 'openmc.data', 'openmc.mgxs', 'openmc.stats'],
'packages': ['openmc', 'openmc.data', 'openmc.mgxs', 'openmc.model',
'openmc.stats'],
'scripts': glob.glob('scripts/openmc-*'),
# Metadata

28
src/constants.F90
View file

@ -281,7 +281,7 @@ module constants
EVENT_ABSORB = 2
! Tally score type
integer, parameter :: N_SCORE_TYPES = 22
integer, parameter :: N_SCORE_TYPES = 20
integer, parameter :: &
SCORE_FLUX = -1, & ! flux
SCORE_TOTAL = -2, & ! total reaction rate
@ -291,20 +291,18 @@ module constants
SCORE_SCATTER_PN = -6, & ! system for scoring 0th through nth moment
SCORE_NU_SCATTER_N = -7, & ! arbitrary nu-scattering moment
SCORE_NU_SCATTER_PN = -8, & ! system for scoring 0th through nth nu-scatter moment
SCORE_TRANSPORT = -9, & ! transport reaction rate
SCORE_N_1N = -10, & ! (n,1n) rate
SCORE_ABSORPTION = -11, & ! absorption rate
SCORE_FISSION = -12, & ! fission rate
SCORE_NU_FISSION = -13, & ! neutron production rate
SCORE_KAPPA_FISSION = -14, & ! fission energy production rate
SCORE_CURRENT = -15, & ! partial current
SCORE_FLUX_YN = -16, & ! angular moment of flux
SCORE_TOTAL_YN = -17, & ! angular moment of total reaction rate
SCORE_SCATTER_YN = -18, & ! angular flux-weighted scattering moment (0:N)
SCORE_NU_SCATTER_YN = -19, & ! angular flux-weighted nu-scattering moment (0:N)
SCORE_EVENTS = -20, & ! number of events
SCORE_DELAYED_NU_FISSION = -21, & ! delayed neutron production rate
SCORE_INVERSE_VELOCITY = -22 ! flux-weighted inverse velocity
SCORE_ABSORPTION = -9, & ! absorption rate
SCORE_FISSION = -10, & ! fission rate
SCORE_NU_FISSION = -11, & ! neutron production rate
SCORE_KAPPA_FISSION = -12, & ! fission energy production rate
SCORE_CURRENT = -13, & ! partial current
SCORE_FLUX_YN = -14, & ! angular moment of flux
SCORE_TOTAL_YN = -15, & ! angular moment of total reaction rate
SCORE_SCATTER_YN = -16, & ! angular flux-weighted scattering moment (0:N)
SCORE_NU_SCATTER_YN = -17, & ! angular flux-weighted nu-scattering moment (0:N)
SCORE_EVENTS = -18, & ! number of events
SCORE_DELAYED_NU_FISSION = -19, & ! delayed neutron production rate
SCORE_INVERSE_VELOCITY = -20 ! flux-weighted inverse velocity
! Maximum scattering order supported
integer, parameter :: MAX_ANG_ORDER = 10

4
src/endf.F90
View file

@ -34,10 +34,6 @@ contains
string = "nu-scatter-n"
case (SCORE_NU_SCATTER_PN)
string = "nu-scatter-pn"
case (SCORE_TRANSPORT)
string = "transport"
case (SCORE_N_1N)
string = "n1n"
case (SCORE_ABSORPTION)
string = "absorption"
case (SCORE_FISSION)

70
src/geometry.F90
View file

@ -408,6 +408,7 @@ contains
real(8) :: v ! y-component of direction
real(8) :: w ! z-component of direction
real(8) :: norm ! "norm" of surface normal
real(8) :: xyz(3) ! Saved global coordinate
integer :: i_surface ! index in surfaces
logical :: found ! particle found in universe?
class(Surface), pointer :: surf
@ -462,12 +463,13 @@ contains
! Score surface currents since reflection causes the direction of the
! particle to change -- artificially move the particle slightly back in
! case the surface crossing in coincident with a mesh boundary
! case the surface crossing is coincident with a mesh boundary
if (active_current_tallies % size() > 0) then
xyz = p % coord(1) % xyz
p % coord(1) % xyz = p % coord(1) % xyz - TINY_BIT * p % coord(1) % uvw
call score_surface_current(p)
p % coord(1) % xyz = p % coord(1) % xyz + TINY_BIT * p % coord(1) % uvw
p % coord(1) % xyz = xyz
end if
! Reflect particle off surface
@ -505,6 +507,70 @@ contains
&// trim(to_str(surf%id)))
end if
return
elseif (surf % bc == BC_PERIODIC .and. run_mode /= MODE_PLOTTING) then
! =======================================================================
! PERIODIC BOUNDARY
! Do not handle periodic boundary conditions on lower universes
if (p % n_coord /= 1) then
call handle_lost_particle(p, "Cannot transfer particle " &
// trim(to_str(p % id)) // " across surface in a lower universe.&
& Boundary conditions must be applied to universe 0.")
return
end if
! Score surface currents since reflection causes the direction of the
! particle to change -- artificially move the particle slightly back in
! case the surface crossing is coincident with a mesh boundary
if (active_current_tallies % size() > 0) then
xyz = p % coord(1) % xyz
p % coord(1) % xyz = p % coord(1) % xyz - TINY_BIT * p % coord(1) % uvw
call score_surface_current(p)
p % coord(1) % xyz = xyz
end if
select type (surf)
type is (SurfaceXPlane)
select type (opposite => surfaces(surf % i_periodic) % obj)
type is (SurfaceXPlane)
p % coord(1) % xyz(1) = opposite % x0
end select
type is (SurfaceYPlane)
select type (opposite => surfaces(surf % i_periodic) % obj)
type is (SurfaceYPlane)
p % coord(1) % xyz(2) = opposite % y0
end select
type is (SurfaceZPlane)
select type (opposite => surfaces(surf % i_periodic) % obj)
type is (SurfaceZPlane)
p % coord(1) % xyz(3) = opposite % z0
end select
end select
! Reassign particle's surface
p % surface = sign(surf % i_periodic, p % surface)
! Figure out what cell particle is in now
p % n_coord = 1
call find_cell(p, found)
if (.not. found) then
call handle_lost_particle(p, "Couldn't find particle after hitting &
&periodic boundary on surface " // trim(to_str(surf%id)) // ".")
return
end if
! Set previous coordinate going slightly past surface crossing
p % last_xyz = p % coord(1) % xyz + TINY_BIT * p % coord(1) % uvw
! Diagnostic message
if (verbosity >= 10 .or. trace) then
call write_message(" Hit periodic boundary on surface " &
// trim(to_str(surf%id)))
end if
return
end if
! ==========================================================================

47
src/geometry_header.F90
View file

@ -1,6 +1,6 @@
module geometry_header
use constants, only: HALF, TWO, THREE
use constants, only: HALF, TWO, THREE, INFINITY
implicit none
@ -204,20 +204,22 @@ contains
real(8), intent(in) :: global_xyz(3)
integer :: i_xyz(3)
real(8) :: xyz(3) ! global_xyz alias
real(8) :: xyz(3) ! global xyz relative to the center
real(8) :: alpha ! Skewed coord axis
real(8) :: xyz_t(3) ! Local xyz
real(8) :: dists(4) ! Squared distances from cell centers
real(8) :: d, d_min ! Squared distance from cell centers
integer :: i, j, k ! Iterators
integer :: loc(1) ! Minimum distance index
integer :: k_min ! Minimum distance index
xyz = global_xyz
xyz(1) = global_xyz(1) - this % center(1)
xyz(2) = global_xyz(2) - this % center(2)
! Index z direction.
if (this % is_3d) then
i_xyz(3) = ceiling((xyz(3) - this % center(3))/this % pitch(2) + HALF)&
+ this % n_axial/2
xyz(3) = global_xyz(3) - this % center(3)
i_xyz(3) = ceiling(xyz(3)/this % pitch(2) + HALF*this % n_axial)
else
xyz(3) = global_xyz(3)
i_xyz(3) = 1
end if
@ -236,28 +238,33 @@ contains
! the four possible cells. Regular hexagonal tiles form a centroidal
! Voronoi tessellation so the global xyz should be in the hexagonal cell
! that it is closest to the center of. This method is used over a
! method that uses the remainders of the floor divisions above becasue it
! method that uses the remainders of the floor divisions above because it
! provides better finite precision performance. Squared distances are
! used becasue they are more computationally efficient than normal
! distances.
k = 1
do i=0,1
do j=0,1
xyz_t = this % get_local_xyz(xyz, i_xyz + (/j, i, 0/))
dists(k) = xyz_t(1)**2 + xyz_t(2)**2
d_min = INFINITY
do i = 0, 1
do j = 0, 1
xyz_t = this % get_local_xyz(global_xyz, i_xyz + [j, i, 0])
d = xyz_t(1)**2 + xyz_t(2)**2
if (d < d_min) then
d_min = d
k_min = k
end if
k = k + 1
end do
end do
! Select the minimum squared distance which corresponds to the cell the
! coordinates are in.
loc = minloc(dists)
if (loc(1) == 2) then
i_xyz = i_xyz + (/1, 0, 0/)
else if (loc(1) == 3) then
i_xyz = i_xyz + (/0, 1, 0/)
else if (loc(1) == 4) then
i_xyz = i_xyz + (/1, 1, 0/)
if (k_min == 2) then
i_xyz(1) = i_xyz(1) + 1
else if (k_min == 3) then
i_xyz(2) = i_xyz(2) + 1
else if (k_min == 4) then
i_xyz(1) = i_xyz(1) + 1
i_xyz(2) = i_xyz(2) + 1
end if
end function get_inds_hex
@ -306,7 +313,7 @@ contains
(i_xyz(1) - this % n_rings) * this % pitch(1) / TWO)
if (this % is_3d) then
local_xyz(3) = xyz(3) - this % center(3) &
+ (this % n_axial/2 - i_xyz(3) + 1) * this % pitch(2)
+ (HALF*this % n_axial - i_xyz(3) + HALF) * this % pitch(2)
else
local_xyz(3) = xyz(3)
end if

116
src/input_xml.F90
View file

@ -1142,6 +1142,8 @@ contains
integer :: universe_num
integer :: n_cells_in_univ
integer :: coeffs_reqd
integer :: i_xmin, i_xmax, i_ymin, i_ymax, i_zmin, i_zmax
real(8) :: xmin, xmax, ymin, ymax, zmin, zmax
integer, allocatable :: temp_int_array(:)
real(8) :: phi, theta, psi
real(8), allocatable :: coeffs(:)
@ -1458,6 +1460,13 @@ contains
call fatal_error("No surfaces found in geometry.xml!")
end if
xmin = INFINITY
xmax = -INFINITY
ymin = INFINITY
ymax = -INFINITY
zmin = INFINITY
zmax = -INFINITY
! Allocate cells array
allocate(surfaces(n_surfaces))
@ -1549,10 +1558,28 @@ contains
select type(s)
type is (SurfaceXPlane)
s%x0 = coeffs(1)
! Determine outer surfaces
xmin = min(xmin, s % x0)
xmax = max(xmax, s % x0)
if (xmin == s % x0) i_xmin = i
if (xmax == s % x0) i_xmax = i
type is (SurfaceYPlane)
s%y0 = coeffs(1)
! Determine outer surfaces
ymin = min(ymin, s % y0)
ymax = max(ymax, s % y0)
if (ymin == s % y0) i_ymin = i
if (ymax == s % y0) i_ymax = i
type is (SurfaceZPlane)
s%z0 = coeffs(1)
! Determine outer surfaces
zmin = min(zmin, s % z0)
zmax = max(zmax, s % z0)
if (zmin == s % z0) i_zmin = i
if (zmax == s % z0) i_zmax = i
type is (SurfacePlane)
s%A = coeffs(1)
s%B = coeffs(2)
@ -1619,11 +1646,19 @@ contains
case ('reflective', 'reflect', 'reflecting')
s%bc = BC_REFLECT
boundary_exists = .true.
case ('periodic')
s%bc = BC_PERIODIC
boundary_exists = .true.
! Check for specification of periodic surface
if (check_for_node(node_surf, "periodic_surface_id")) then
call get_node_value(node_surf, "periodic_surface_id", &
s % i_periodic)
end if
case default
call fatal_error("Unknown boundary condition '" // trim(word) // &
&"' specified on surface " // trim(to_str(s%id)))
end select
! Add surface to dictionary
call surface_dict % add_key(s%id, i)
end do
@ -1634,6 +1669,67 @@ contains
call fatal_error("No boundary conditions were applied to any surfaces!")
end if
! Determine opposite side for periodic boundaries
do i = 1, size(surfaces)
if (surfaces(i) % obj % bc == BC_PERIODIC) then
select type (surf => surfaces(i) % obj)
type is (SurfaceXPlane)
if (surf % i_periodic == NONE) then
if (i == i_xmin) then
surf % i_periodic = i_xmax
elseif (i == i_xmax) then
surf % i_periodic = i_xmin
else
call fatal_error("Periodic boundary condition applied to &
&interior surface.")
end if
else
surf % i_periodic = surface_dict % get_key(surf % i_periodic)
end if
type is (SurfaceYPlane)
if (surf % i_periodic == NONE) then
if (i == i_ymin) then
surf % i_periodic = i_ymax
elseif (i == i_ymax) then
surf % i_periodic = i_ymin
else
call fatal_error("Periodic boundary condition applied to &
&interior surface.")
end if
else
surf % i_periodic = surface_dict % get_key(surf % i_periodic)
end if
type is (SurfaceZPlane)
if (surf % i_periodic == NONE) then
if (i == i_zmin) then
surf % i_periodic = i_zmax
elseif (i == i_zmax) then
surf % i_periodic = i_zmin
else
call fatal_error("Periodic boundary condition applied to &
&interior surface.")
end if
else
surf % i_periodic = surface_dict % get_key(surf % i_periodic)
end if
class default
call fatal_error("Periodic boundary condition applied to &
&non-planar surface.")
end select
! Make sure opposite surface is also periodic
associate (surf => surfaces(i) % obj)
if (surfaces(surf % i_periodic) % obj % bc /= BC_PERIODIC) then
call fatal_error("Could not find matching surface for periodic &
&boundary on surface " // trim(to_str(surf % id)) // ".")
end if
end associate
end if
end do
! ==========================================================================
! READ LATTICES FROM GEOMETRY.XML
@ -3464,22 +3560,12 @@ contains
j = j + n_bins - 1
case('transport')
t % score_bins(j) = SCORE_TRANSPORT
! Set tally estimator to analog
t % estimator = ESTIMATOR_ANALOG
case ('diffusion')
call fatal_error("Diffusion score no longer supported for tallies, &
call fatal_error("Transport score no longer supported for tallies, &
&please remove")
case ('n1n')
if (run_CE) then
t % score_bins(j) = SCORE_N_1N
! Set tally estimator to analog
t % estimator = ESTIMATOR_ANALOG
else
call fatal_error("Cannot tally n1n rate in multi-group mode!")
end if
case ('n1n')
call fatal_error("n1n score no longer supported for tallies, &
&please remove")
case ('n2n', '(n,2n)')
t % score_bins(j) = N_2N

23
src/mgxs_header.F90
View file

@ -1485,8 +1485,13 @@ module mgxs_header
nuc % scatter % energy(gin) % data(gout)
mult_num(gout, gin) = mult_num(gout, gin) + atom_density * &
nuscatt
mult_denom(gout, gin) = mult_denom(gout,gin) + atom_density * &
nuscatt / nuc % scatter % mult(gin) % data(gout)
if (nuc % scatter % mult(gin) % data(gout) > ZERO) then
mult_denom(gout, gin) = mult_denom(gout,gin) + atom_density * &
nuscatt / nuc % scatter % mult(gin) % data(gout)
else
! Avoid division by zero
mult_denom(gout, gin) = mult_denom(gout,gin) + atom_density
end if
end do
end do
@ -1722,10 +1727,16 @@ module mgxs_header
nuc % scatter(iazi, ipol) % obj % energy(gin) % data(gout)
mult_num(gout, gin, iazi, ipol) = mult_num(gout, gin, iazi, ipol) + &
atom_density * nuscatt
mult_denom(gout, gin, iazi, ipol) = &
mult_denom(gout, gin, iazi, ipol) + &
atom_density * nuscatt / &
nuc % scatter(iazi, ipol) % obj % mult(gin) % data(gout)
if (nuc % scatter(iazi, ipol) % obj % mult(gin) % data(gout) > ZERO) then
mult_denom(gout, gin, iazi, ipol) = &
mult_denom(gout, gin, iazi, ipol) + &
atom_density * nuscatt / &
nuc % scatter(iazi, ipol) % obj % mult(gin) % data(gout)
else
! Avoid division by zero
mult_denom(gout, gin, iazi, ipol) = &
mult_denom(gout,gin, iazi, ipol) + atom_density
end if
end do
end do
end do

2
src/output.F90
View file

@ -777,8 +777,6 @@ contains
score_names(abs(SCORE_TOTAL)) = "Total Reaction Rate"
score_names(abs(SCORE_SCATTER)) = "Scattering Rate"
score_names(abs(SCORE_NU_SCATTER)) = "Scattering Production Rate"
score_names(abs(SCORE_TRANSPORT)) = "Transport Rate"
score_names(abs(SCORE_N_1N)) = "(n,1n) Rate"
score_names(abs(SCORE_ABSORPTION)) = "Absorption Rate"
score_names(abs(SCORE_FISSION)) = "Fission Rate"
score_names(abs(SCORE_NU_FISSION)) = "Nu-Fission Rate"

5
src/relaxng/geometry.rnc
View file

@ -23,8 +23,9 @@ element geometry {
(element type { xsd:string { maxLength = "15" } } |
attribute type { xsd:string { maxLength = "15" } }) &
(element coeffs { list { xsd:double+ } } | attribute coeffs { list { xsd:double+ } }) &
(element boundary { ( "transmit" | "reflective" | "vacuum" ) } |
attribute boundary { ( "transmit" | "reflective" | "vacuum" ) })?
(element boundary { ( "transmit" | "reflective" | "vacuum" | "periodic" ) } |
attribute boundary { ( "transmit" | "reflective" | "vacuum" | "periodic" ) })? &
(element periodic_surface_id { xsd:int } | attribute periodic_surface_id { xsd:int })?
}*
& element lattice {

12
src/relaxng/geometry.rng
View file

@ -191,6 +191,7 @@
<value>transmit</value>
<value>reflective</value>
<value>vacuum</value>
<value>periodic</value>
</choice>
</element>
<attribute name="boundary">
@ -198,10 +199,21 @@
<value>transmit</value>
<value>reflective</value>
<value>vacuum</value>
<value>periodic</value>
</choice>
</attribute>
</choice>
</optional>
<optional>
<choice>
<element name="periodic_surface_id">
<data type="int"/>
</element>
<attribute name="periodic_surface_id">
<data type="int"/>
</attribute>
</choice>
</optional>
</interleave>
</element>
</zeroOrMore>

5
src/surface_header.F90
View file

@ -1,6 +1,6 @@
module surface_header
use constants, only: ONE, TWO, ZERO, HALF, INFINITY, FP_COINCIDENT
use constants, only: NONE, ONE, TWO, ZERO, HALF, INFINITY, FP_COINCIDENT
implicit none
@ -15,7 +15,8 @@ module surface_header
neighbor_pos(:), & ! List of cells on positive side
neighbor_neg(:) ! List of cells on negative side
integer :: bc ! Boundary condition
character(len=104) :: name = "" ! User-defined name
integer :: i_periodic = NONE ! Index of corresponding periodic surface
character(len=104) :: name = "" ! User-defined name
contains
procedure :: sense
procedure :: reflect

38
src/tally.F90
View file

@ -346,30 +346,6 @@ contains
end if
case (SCORE_TRANSPORT)
! Only analog estimators are available.
! Skip any event where the particle didn't scatter
if (p % event /= EVENT_SCATTER) cycle SCORE_LOOP
! get material macros
macro_total = material_xs % total
macro_scatt = material_xs % total - material_xs % absorption
! Score total rate - p1 scatter rate Note estimator needs to be
! adjusted since tallying is only occuring when a scatter has
! happened. Effectively this means multiplying the estimator by
! total/scatter macro
score = (macro_total - p % mu * macro_scatt) * (ONE / macro_scatt)
case (SCORE_N_1N)
! Only analog estimators are available.
! Skip any event where the particle didn't scatter
if (p % event /= EVENT_SCATTER) cycle SCORE_LOOP
! Skip any events where weight of particle changed
if (p % wgt /= p % last_wgt) cycle SCORE_LOOP
! All events that reach this point are (n,1n) reactions
score = p % last_wgt
case (SCORE_ABSORPTION)
if (t % estimator == ESTIMATOR_ANALOG) then
if (survival_biasing) then
@ -1021,20 +997,6 @@ contains
end if
case (SCORE_TRANSPORT)
! Only analog estimators are available.
! Skip any event where the particle didn't scatter
if (p % event /= EVENT_SCATTER) cycle SCORE_LOOP
! Score total rate - p1 scatter rate Note estimator needs to be
! adjusted since tallying is only occuring when a scatter has
! happened. Effectively this means multiplying the estimator by
! total/scatter macro
score = (material_xs % total - p % mu * material_xs % elastic)
if (material_xs % elastic /= ZERO) then
score = score / material_xs % elastic
end if
case (SCORE_ABSORPTION)
if (t % estimator == ESTIMATOR_ANALOG) then
if (survival_biasing) then

145
tests/input_set.py
View file

@ -15,8 +15,10 @@ class InputSet(object):
self.settings.export_to_xml()
self.materials.export_to_xml()
self.geometry.export_to_xml()
if self.tallies is not None: self.tallies.export_to_xml()
if self.plots is not None: self.plots.export_to_xml()
if self.tallies is not None:
self.tallies.export_to_xml()
if self.plots is not None:
self.plots.export_to_xml()
def build_default_materials_and_geometry(self):
# Define materials.
@ -82,7 +84,7 @@ class InputSet(object):
hot_water.add_s_alpha_beta('HH2O', '71t')
rpv_steel = openmc.Material(name='Reactor pressure vessel steel',
material_id=5)
material_id=5)
rpv_steel.set_density('g/cm3', 7.9)
rpv_steel.add_nuclide("Fe-54", 0.05437098, 'wo')
rpv_steel.add_nuclide("Fe-56", 0.88500663, 'wo')
@ -113,7 +115,7 @@ class InputSet(object):
rpv_steel.add_nuclide("Cu-65", 0.0006304, 'wo')
lower_rad_ref = openmc.Material(name='Lower radial reflector',
material_id=6)
material_id=6)
lower_rad_ref.set_density('g/cm3', 4.32)
lower_rad_ref.add_nuclide("H-1", 0.0095661, 'wo')
lower_rad_ref.add_nuclide("O-16", 0.0759107, 'wo')
@ -189,7 +191,8 @@ class InputSet(object):
bot_plate.add_nuclide("Cr-54", 0.004612692337, 'wo')
bot_plate.add_s_alpha_beta('HH2O', '71t')
bot_nozzle = openmc.Material(name='Bottom nozzle region', material_id=9)
bot_nozzle = openmc.Material(name='Bottom nozzle region',
material_id=9)
bot_nozzle.set_density('g/cm3', 2.53)
bot_nozzle.add_nuclide("H-1", 0.0245014, 'wo')
bot_nozzle.add_nuclide("O-16", 0.1944274, 'wo')
@ -252,7 +255,8 @@ class InputSet(object):
top_fa.add_nuclide("Zr-96", 0.02511169542, 'wo')
top_fa.add_s_alpha_beta('HH2O', '71t')
bot_fa = openmc.Material(name='Bottom of fuel assemblies', material_id=12)
bot_fa = openmc.Material(name='Bottom of fuel assemblies',
material_id=12)
bot_fa.set_density('g/cm3', 1.762)
bot_fa.add_nuclide("H-1", 0.0292856, 'wo')
bot_fa.add_nuclide("O-16", 0.2323919, 'wo')
@ -350,7 +354,6 @@ class InputSet(object):
# Define fuel lattices.
l100 = openmc.RectLattice(name='Fuel assembly (lower half)',
lattice_id=100)
l100.dimension = (17, 17)
l100.lower_left = (-10.71, -10.71)
l100.pitch = (1.26, 1.26)
l100.universes = [
@ -384,7 +387,6 @@ class InputSet(object):
l101 = openmc.RectLattice(name='Fuel assembly (upper half)',
lattice_id=101)
l101.dimension = (17, 17)
l101.lower_left = (-10.71, -10.71)
l101.pitch = (1.26, 1.26)
l101.universes = [
@ -444,7 +446,6 @@ class InputSet(object):
# Define core lattices
l200 = openmc.RectLattice(name='Core lattice (lower half)',
lattice_id=200)
l200.dimension = (21, 21)
l200.lower_left = (-224.91, -224.91)
l200.pitch = (21.42, 21.42)
l200.universes = [
@ -472,7 +473,6 @@ class InputSet(object):
l201 = openmc.RectLattice(name='Core lattice (lower half)',
lattice_id=201)
l201.dimension = (21, 21)
l201.lower_left = (-224.91, -224.91)
l201.pitch = (21.42, 21.42)
l201.universes = [
@ -570,6 +570,109 @@ class InputSet(object):
self.plots.add_plot(plot)
class PinCellInputSet(object):
def __init__(self):
self.settings = openmc.Settings()
self.materials = openmc.Materials()
self.geometry = openmc.Geometry()
self.tallies = None
self.plots = None
def export(self):
self.settings.export_to_xml()
self.materials.export_to_xml()
self.geometry.export_to_xml()
if self.tallies is not None:
self.tallies.export_to_xml()
if self.plots is not None:
self.plots.export_to_xml()
def build_default_materials_and_geometry(self):
# Define materials.
fuel = openmc.Material(name='Fuel')
fuel.set_density('g/cm3', 10.29769)
fuel.add_nuclide("U-234", 4.4843e-6)
fuel.add_nuclide("U-235", 5.5815e-4)
fuel.add_nuclide("U-238", 2.2408e-2)
fuel.add_nuclide("O-16", 4.5829e-2)
clad = openmc.Material(name='Cladding')
clad.set_density('g/cm3', 6.55)
clad.add_nuclide("Zr-90", 2.1827e-2)
clad.add_nuclide("Zr-91", 4.7600e-3)
clad.add_nuclide("Zr-92", 7.2758e-3)
clad.add_nuclide("Zr-94", 7.3734e-3)
clad.add_nuclide("Zr-96", 1.1879e-3)
hot_water = openmc.Material(name='Hot borated water')
hot_water.set_density('g/cm3', 0.740582)
hot_water.add_nuclide("H-1", 4.9457e-2)
hot_water.add_nuclide("O-16", 2.4672e-2)
hot_water.add_nuclide("B-10", 8.0042e-6)
hot_water.add_nuclide("B-11", 3.2218e-5)
hot_water.add_s_alpha_beta('HH2O', '71t')
# Define the materials file.
self.materials.default_xs = '71c'
self.materials += (fuel, clad, hot_water)
# Instantiate ZCylinder surfaces
fuel_or = openmc.ZCylinder(x0=0, y0=0, R=0.39218, name='Fuel OR')
clad_or = openmc.ZCylinder(x0=0, y0=0, R=0.45720, name='Clad OR')
left = openmc.XPlane(x0=-0.63, name='left')
right = openmc.XPlane(x0=0.63, name='right')
bottom = openmc.YPlane(y0=-0.63, name='bottom')
top = openmc.YPlane(y0=0.63, name='top')
left.boundary_type = 'reflective'
right.boundary_type = 'reflective'
top.boundary_type = 'reflective'
bottom.boundary_type = 'reflective'
# Instantiate Cells
fuel_pin = openmc.Cell(name='cell 1')
cladding = openmc.Cell(name='cell 3')
water = openmc.Cell(name='cell 2')
# Use surface half-spaces to define regions
fuel_pin.region = -fuel_or
cladding.region = +fuel_or & -clad_or
water.region = +clad_or & +left & -right & +bottom & -top
# Register Materials with Cells
fuel_pin.fill = fuel
cladding.fill = clad
water.fill = hot_water
# Instantiate Universe
root = openmc.Universe(universe_id=0, name='root universe')
# Register Cells with Universe
root.add_cells([fuel_pin, cladding, water])
# Instantiate a Geometry, register the root Universe, and export to XML
self.geometry.root_universe = root
def build_default_settings(self):
self.settings.batches = 10
self.settings.inactive = 5
self.settings.particles = 100
self.settings.source = Source(space=Box([-0.63, -0.63, -1],
[0.63, 0.63, 1],
only_fissionable=True))
def build_defualt_plots(self):
plot = openmc.Plot()
plot.filename = 'mat'
plot.origin = (0.0, 0.0, 0)
plot.width = (1.26, 1.26)
plot.pixels = (300, 300)
plot.color = 'mat'
self.plots.add_plot(plot)
class MGInputSet(InputSet):
def build_default_materials_and_geometry(self):
# Define materials needed for 1D/1G slab problem
@ -595,21 +698,21 @@ class MGInputSet(InputSet):
# Define surfaces.
# Assembly/Problem Boundary
left = openmc.XPlane(x0=0.0, surface_id=200,
boundary_type='reflective')
right = openmc.XPlane(x0=10.0, surface_id=201,
boundary_type='reflective')
left = openmc.XPlane(x0=0.0, surface_id=200,
boundary_type='reflective')
right = openmc.XPlane(x0=10.0, surface_id=201,
boundary_type='reflective')
bottom = openmc.YPlane(y0=0.0, surface_id=300,
boundary_type='reflective')
top = openmc.YPlane(y0=10.0, surface_id=301,
boundary_type='reflective')
top = openmc.YPlane(y0=10.0, surface_id=301,
boundary_type='reflective')
down = openmc.ZPlane(z0=0.0, surface_id=0,
boundary_type='reflective')
down = openmc.ZPlane(z0=0.0, surface_id=0,
boundary_type='reflective')
fuel_clad_intfc = openmc.ZPlane(z0=2.0, surface_id=1)
clad_lwtr_intfc = openmc.ZPlane(z0=2.4, surface_id=2)
up = openmc.ZPlane(z0=5.0, surface_id=3,
boundary_type='reflective')
up = openmc.ZPlane(z0=5.0, surface_id=3,
boundary_type='reflective')
# Define cells
c1 = openmc.Cell(cell_id=1)
@ -625,7 +728,7 @@ class MGInputSet(InputSet):
# Define root universe.
root = openmc.Universe(universe_id=0, name='root universe')
root.add_cells((c1,c2,c3))
root.add_cells((c1, c2, c3))
# Assign root universe to geometry
self.geometry.root_universe = root

1
tests/test_asymmetric_lattice/test_asymmetric_lattice.py
View file

@ -24,7 +24,6 @@ class AsymmetricLatticeTestHarness(PyAPITestHarness):
# Construct a 3x3 lattice of fuel assemblies
core_lat = openmc.RectLattice(name='3x3 Core Lattice', lattice_id=202)
core_lat.dimension = (3, 3)
core_lat.lower_left = (-32.13, -32.13)
core_lat.pitch = (21.42, 21.42)
core_lat.universes = [[fuel, water, water],

2
tests/test_distribmat/results_true.dat
View file

@ -3,7 +3,7 @@ k-combined:
Cell
ID = 11
Name =
Material = [2, 3, void, 2]
Fill = [2, 3, None, 2]
Region = -10000
Rotation = None
Translation = None

3
tests/test_distribmat/test_distribmat.py
View file

@ -45,7 +45,7 @@ class DistribmatTestHarness(PyAPITestHarness):
r0 = openmc.ZCylinder(R=0.3)
c11 = openmc.Cell(cell_id=11)
c11.region = -r0
c11.fill = [dense_fuel, light_fuel, 'void', dense_fuel]
c11.fill = [dense_fuel, light_fuel, None, dense_fuel]
c12 = openmc.Cell(cell_id=12)
c12.region = +r0
c12.fill = moderator
@ -53,7 +53,6 @@ class DistribmatTestHarness(PyAPITestHarness):
fuel_univ.add_cells((c11, c12))
lat = openmc.RectLattice(lattice_id=101)
lat.dimension = [2, 2]
lat.lower_left = [-2.0, -2.0]
lat.pitch = [2.0, 2.0]
lat.universes = [[fuel_univ]*2]*2

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

@ -0,0 +1 @@
34d5891f6f17c2d4b686b814ba61ba0045bc4289e278b1c3c47dbba59b83837fcfe15f2b8d58e7a2b07627b73d51e40348d70e9ed36dbb7cc94468d61c068c4c

2
tests/test_mgxs_library_ce_to_mg/results_true.dat Normal file
View file

@ -0,0 +1,2 @@
k-combined:
1.094839E+00 1.203524E-02

94
tests/test_mgxs_library_ce_to_mg/test_mgxs_library_ce_to_mg.py Normal file
View file

@ -0,0 +1,94 @@
#!/usr/bin/env python
import os
import sys
import glob
import hashlib
sys.path.insert(0, os.pardir)
from testing_harness import PyAPITestHarness
from input_set import PinCellInputSet
import openmc
import openmc.mgxs
class MGXSTestHarness(PyAPITestHarness):
def _build_inputs(self):
# Set the input set to use the pincell model
self._input_set = PinCellInputSet()
# Generate inputs using parent class routine
super(MGXSTestHarness, self)._build_inputs()
# Initialize a two-group structure
energy_groups = openmc.mgxs.EnergyGroups(group_edges=[0, 0.625e-6,
20.])
# Initialize MGXS Library for a few cross section types
self.mgxs_lib = openmc.mgxs.Library(self._input_set.geometry)
self.mgxs_lib.by_nuclide = False
self.mgxs_lib.mgxs_types = ['total', 'absorption', 'nu-fission matrix',
'nu-scatter matrix', 'multiplicity matrix']
self.mgxs_lib.energy_groups = energy_groups
self.mgxs_lib.correction = None
self.mgxs_lib.legendre_order = 3
self.mgxs_lib.domain_type = 'material'
self.mgxs_lib.build_library()
# Initialize a tallies file
self._input_set.tallies = openmc.Tallies()
self.mgxs_lib.add_to_tallies_file(self._input_set.tallies, merge=False)
self._input_set.tallies.export_to_xml()
def _run_openmc(self):
# Initial run
if self._opts.mpi_exec is not None:
returncode = openmc.run(mpi_procs=self._opts.mpi_np,
openmc_exec=self._opts.exe,
mpi_exec=self._opts.mpi_exec)
else:
returncode = openmc.run(openmc_exec=self._opts.exe)
assert returncode == 0, 'CE OpenMC calculation did not exit' \
'successfully.'
# Build MG Inputs
# Get data needed to execute Library calculations.
statepoint = glob.glob(os.path.join(os.getcwd(), self._sp_name))[0]
sp = openmc.StatePoint(statepoint)
self.mgxs_lib.load_from_statepoint(sp)
self._input_set.mgxs_file, self._input_set.materials, \
self._input_set.geometry = self.mgxs_lib.create_mg_mode()
# Modify settings so we can run in MG mode
self._input_set.settings.cross_sections = './mgxs.xml'
self._input_set.settings.energy_mode = 'multi-group'
# Write modified input files
self._input_set.settings.export_to_xml()
self._input_set.geometry.export_to_xml()
self._input_set.materials.export_to_xml()
self._input_set.mgxs_file.export_to_xml()
# Dont need tallies.xml, so remove the file
if os.path.exists('./tallies.xml'):
os.remove('./tallies.xml')
# Re-run MG mode.
if self._opts.mpi_exec is not None:
returncode = openmc.run(mpi_procs=self._opts.mpi_np,
openmc_exec=self._opts.exe,
mpi_exec=self._opts.mpi_exec)
else:
returncode = openmc.run(openmc_exec=self._opts.exe)
def _cleanup(self):
super(MGXSTestHarness, self)._cleanup()
f = os.path.join(os.getcwd(), 'mgxs.xml')
if os.path.exists(f):
os.remove(f)
if __name__ == '__main__':
harness = MGXSTestHarness('statepoint.10.*', False)
harness.main()

2
tests/test_mgxs_library_condense/inputs_true.dat
View file

@ -1 +1 @@
104e7fb527770ac5d3fc636da7716e8fb05d55761253d30516c899f466e6b38ffd881611a3d0cdf65c6af058c32f6f6758c68782be7a170d21024bdae751862f
317a63a9dd3bfd84e969667b00f46018e56c04c356461a75103f63569e6b70c84d0da7f5e611faaf1b2631330b05ab4346223d3d843018ce0ce8876671a450c0

191
tests/test_mgxs_library_condense/results_true.dat
View file

@ -1,85 +1,108 @@
material group in nuclide mean std. dev.
0 1 1 total 0.412084 0.02359 material group in nuclide mean std. dev.
0 1 1 total 0.076425 0.003691 material group in group out nuclide moment mean std. dev.
0 1 1 1 total P0 0.384780 0.022253
1 1 1 1 total P1 0.039277 0.004308
2 1 1 1 total P2 0.017574 0.002402
3 1 1 1 total P3 0.012203 0.002164 material group out nuclide mean std. dev.
0 1 1 total 1.0 0.055333 material group in nuclide mean std. dev.
0 2 1 total 0.241262 0.00841 material group in nuclide mean std. dev.
0 2 1 total 0.0 0.0 material group in group out nuclide moment mean std. dev.
0 2 1 1 total P0 0.272369 0.006872
1 2 1 1 total P1 0.031107 0.005483
2 2 1 1 total P2 0.025999 0.006151
3 2 1 1 total P3 0.003219 0.003312 material group out nuclide mean std. dev.
0 2 1 total 0.0 0.0 material group in nuclide mean std. dev.
0 3 1 total 0.400028 0.034667 material group in nuclide mean std. dev.
0 3 1 total 0.0 0.0 material group in group out nuclide moment mean std. dev.
0 3 1 1 total P0 0.794999 0.036548
1 3 1 1 total P1 0.401537 0.016175
2 3 1 1 total P2 0.143623 0.008719
3 3 1 1 total P3 0.001991 0.004433 material group out nuclide mean std. dev.
0 3 1 total 0.0 0.0 material group in nuclide mean std. dev.
0 4 1 total 0.377402 0.072937 material group in nuclide mean std. dev.
0 4 1 total 0.0 0.0 material group in group out nuclide moment mean std. dev.
0 4 1 1 total P0 0.727311 0.080096
1 4 1 1 total P1 0.355839 0.037901
2 4 1 1 total P2 0.124483 0.015823
3 4 1 1 total P3 0.012168 0.006224 material group out nuclide mean std. dev.
0 4 1 total 0.0 0.0 material group in nuclide mean std. dev.
0 5 1 total 0.0 0.0 material group in nuclide mean std. dev.
0 5 1 total 0.0 0.0 material group in group out nuclide moment mean std. dev.
0 5 1 1 total P0 0.0 0.0
1 5 1 1 total P1 0.0 0.0
2 5 1 1 total P2 0.0 0.0
3 5 1 1 total P3 0.0 0.0 material group out nuclide mean std. dev.
0 5 1 total 0.0 0.0 material group in nuclide mean std. dev.
0 6 1 total 0.0 0.0 material group in nuclide mean std. dev.
0 6 1 total 0.0 0.0 material group in group out nuclide moment mean std. dev.
0 6 1 1 total P0 0.0 0.0
1 6 1 1 total P1 0.0 0.0
2 6 1 1 total P2 0.0 0.0
3 6 1 1 total P3 0.0 0.0 material group out nuclide mean std. dev.
0 6 1 total 0.0 0.0 material group in nuclide mean std. dev.
0 7 1 total 0.0 0.0 material group in nuclide mean std. dev.
0 7 1 total 0.0 0.0 material group in group out nuclide moment mean std. dev.
0 7 1 1 total P0 0.0 0.0
1 7 1 1 total P1 0.0 0.0
2 7 1 1 total P2 0.0 0.0
3 7 1 1 total P3 0.0 0.0 material group out nuclide mean std. dev.
0 7 1 total 0.0 0.0 material group in nuclide mean std. dev.
0 8 1 total 0.0 0.0 material group in nuclide mean std. dev.
0 8 1 total 0.0 0.0 material group in group out nuclide moment mean std. dev.
0 8 1 1 total P0 0.0 0.0
1 8 1 1 total P1 0.0 0.0
2 8 1 1 total P2 0.0 0.0
3 8 1 1 total P3 0.0 0.0 material group out nuclide mean std. dev.
0 8 1 total 0.0 0.0 material group in nuclide mean std. dev.
0 9 1 total 0.600536 0.748875 material group in nuclide mean std. dev.
0 9 1 total 0.0 0.0 material group in group out nuclide moment mean std. dev.
0 9 1 1 total P0 0.720380 0.771015
1 9 1 1 total P1 0.119844 0.184691
2 9 1 1 total P2 0.038522 0.064485
3 9 1 1 total P3 0.056023 0.050595 material group out nuclide mean std. dev.
0 9 1 total 0.0 0.0 material group in nuclide mean std. dev.
0 10 1 total 0.235515 0.613974 material group in nuclide mean std. dev.
0 10 1 total 0.0 0.0 material group in group out nuclide moment mean std. dev.
0 10 1 1 total P0 0.501009 0.708534
1 10 1 1 total P1 0.265494 0.375465
2 10 1 1 total P2 0.141979 0.200788
3 10 1 1 total P3 0.074258 0.105017 material group out nuclide mean std. dev.
0 10 1 total 0.0 0.0 material group in nuclide mean std. dev.
0 11 1 total 0.510145 0.741941 material group in nuclide mean std. dev.
0 11 1 total 0.0 0.0 material group in group out nuclide moment mean std. dev.
0 11 1 1 total P0 0.804661 0.817658
1 11 1 1 total P1 0.312803 0.315315
2 11 1 1 total P2 0.168113 0.172935
3 11 1 1 total P3 0.003808 0.037911 material group out nuclide mean std. dev.
0 11 1 total 0.0 0.0 material group in nuclide mean std. dev.
0 12 1 total 0.73836 0.825631 material group in nuclide mean std. dev.
0 12 1 total 0.0 0.0 material group in group out nuclide moment mean std. dev.
0 12 1 1 total P0 0.943429 0.856119
1 12 1 1 total P1 0.220164 0.163180
2 12 1 1 total P2 0.052884 0.042440
3 12 1 1 total P3 0.039939 0.032867 material group out nuclide mean std. dev.
0 12 1 total 0.0 0.0
0 10000 1 total 0.453624 0.021053
material group in nuclide mean std. dev.
0 10000 1 total 0.400852 0.022858
material group in nuclide mean std. dev.
0 10000 1 total 0.400852 0.022858
material group in nuclide mean std. dev.
0 10000 1 total 0.064903 0.004313
material group in nuclide mean std. dev.
0 10000 1 total 0.028048 0.00458
material group in nuclide mean std. dev.
0 10000 1 total 0.036855 0.002622
material group in nuclide mean std. dev.
0 10000 1 total 0.090649 0.00641
material group in nuclide mean std. dev.
0 10000 1 total 7.137955 0.507364
material group in nuclide mean std. dev.
0 10000 1 total 0.388721 0.01783
material group in nuclide mean std. dev.
0 10000 1 total 0.389304 0.023076
material group in group out nuclide moment mean std. dev.
0 10000 1 1 total P0 0.389304 0.023146
1 10000 1 1 total P1 0.046224 0.005907
2 10000 1 1 total P2 0.017984 0.002883
3 10000 1 1 total P3 0.006628 0.002457
material group in group out nuclide moment mean std. dev.
0 10000 1 1 total P0 0.389304 0.023146
1 10000 1 1 total P1 0.046224 0.005907
2 10000 1 1 total P2 0.017984 0.002883
3 10000 1 1 total P3 0.006628 0.002457
material group in group out nuclide mean std. dev.
0 10000 1 1 total 1.0 0.066111
material group in group out nuclide mean std. dev.
0 10000 1 1 total 0.085835 0.005592
material group out nuclide mean std. dev.
0 10000 1 total 1.0 0.046071
material group in nuclide mean std. dev.
0 10001 1 total 0.311594 0.013793
material group in nuclide mean std. dev.
0 10001 1 total 0.279255 0.02919
material group in nuclide mean std. dev.
0 10001 1 total 0.279255 0.02919
material group in nuclide mean std. dev.
0 10001 1 total 0.00221 0.000286
material group in nuclide mean std. dev.
0 10001 1 total 0.00221 0.000286
material group in nuclide mean std. dev.
0 10001 1 total 0.0 0.0
material group in nuclide mean std. dev.
0 10001 1 total 0.0 0.0
material group in nuclide mean std. dev.
0 10001 1 total 0.0 0.0
material group in nuclide mean std. dev.
0 10001 1 total 0.309384 0.013551
material group in nuclide mean std. dev.
0 10001 1 total 0.307987 0.029308
material group in group out nuclide moment mean std. dev.
0 10001 1 1 total P0 0.307987 0.029308
1 10001 1 1 total P1 0.030617 0.007464
2 10001 1 1 total P2 0.018911 0.004323
3 10001 1 1 total P3 0.006235 0.003338
material group in group out nuclide moment mean std. dev.
0 10001 1 1 total P0 0.307987 0.029308
1 10001 1 1 total P1 0.030617 0.007464
2 10001 1 1 total P2 0.018911 0.004323
3 10001 1 1 total P3 0.006235 0.003338
material group in group out nuclide mean std. dev.
0 10001 1 1 total 1.0 0.095039
material group in group out nuclide mean std. dev.
0 10001 1 1 total 0.0 0.0
material group out nuclide mean std. dev.
0 10001 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.
0 10002 1 total 0.499184 0.040914
material group in nuclide mean std. dev.
0 10002 1 total 0.499184 0.040914
material group in nuclide mean std. dev.
0 10002 1 total 0.00606 0.000555
material group in nuclide mean std. dev.
0 10002 1 total 0.00606 0.000555
material group in nuclide mean std. dev.
0 10002 1 total 0.0 0.0
material group in nuclide mean std. dev.
0 10002 1 total 0.0 0.0
material group in nuclide mean std. dev.
0 10002 1 total 0.0 0.0
material group in nuclide mean std. dev.
0 10002 1 total 0.898938 0.043493
material group in nuclide mean std. dev.
0 10002 1 total 0.903415 0.043959
material group in group out nuclide moment mean std. dev.
0 10002 1 1 total P0 0.903415 0.043586
1 10002 1 1 total P1 0.410417 0.015877
2 10002 1 1 total P2 0.143301 0.007187
3 10002 1 1 total P3 0.008739 0.003571
material group in group out nuclide moment mean std. dev.
0 10002 1 1 total P0 0.903415 0.043586
1 10002 1 1 total P1 0.410417 0.015877
2 10002 1 1 total P2 0.143301 0.007187
3 10002 1 1 total P3 0.008739 0.003571
material group in group out nuclide mean std. dev.
0 10002 1 1 total 1.0 0.056867
material group in group out nuclide mean std. dev.
0 10002 1 1 total 0.0 0.0
material group out nuclide mean std. dev.
0 10002 1 total 0.0 0.0

15
tests/test_mgxs_library_condense/test_mgxs_library_condense.py
View file

@ -6,27 +6,28 @@ import glob
import hashlib
sys.path.insert(0, os.pardir)
from testing_harness import PyAPITestHarness
from input_set import PinCellInputSet
import openmc
import openmc.mgxs
class MGXSTestHarness(PyAPITestHarness):
def _build_inputs(self):
# The openmc.mgxs module needs a summary.h5 file
self._input_set.settings.output = {'summary': True}
# Set the input set to use the pincell model
self._input_set = PinCellInputSet()
# Generate inputs using parent class routine
super(MGXSTestHarness, self)._build_inputs()
# Initialize a two-group structure
energy_groups = openmc.mgxs.EnergyGroups(group_edges=[0, 0.625e-6, 20.])
energy_groups = openmc.mgxs.EnergyGroups(group_edges=[0, 0.625e-6,
20.])
# Initialize MGXS Library for a few cross section types
self.mgxs_lib = openmc.mgxs.Library(self._input_set.geometry)
self.mgxs_lib.by_nuclide = False
self.mgxs_lib.mgxs_types = ['transport', 'nu-fission',
'nu-scatter matrix', 'chi']
# Test all MGXS types
self.mgxs_lib.mgxs_types = openmc.mgxs.MGXS_TYPES
self.mgxs_lib.energy_groups = energy_groups
self.mgxs_lib.legendre_order = 3
self.mgxs_lib.domain_type = 'material'
@ -57,7 +58,7 @@ class MGXSTestHarness(PyAPITestHarness):
for mgxs_type in condense_lib.mgxs_types:
mgxs = condense_lib.get_mgxs(domain, mgxs_type)
df = mgxs.get_pandas_dataframe()
outstr += df.to_string()
outstr += df.to_string() + '\n'
# Hash the results if necessary
if hash_output:

2
tests/test_mgxs_library_distribcell/inputs_true.dat
View file

@ -1 +1 @@
018bbbc2099f7b94180b391e46e42fc9a82498c60b3f8f7f4c91480ea373427932d287fe571d53b2397f329e71485e7155d7644f0f995bbcb458ba3e872ab043
88849ac150f9c389e67de96356dfceb0bde08643f68ca25699e67d263995b95893d7340a2b08b2f0f5075fc5020f73553c5287ec6c56ace2f35ce0214961e123

36
tests/test_mgxs_library_distribcell/results_true.dat
View file

@ -1,8 +1,36 @@
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.718919 0.520644 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.0 0.0 avg(distribcell) group in group out nuclide moment mean std. dev.
0 (0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13,... 1 total 1.145934 0.553822
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.718919 0.520644
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.718919 0.520644
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.019762 0.010629
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.019762 0.010629
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.0 0.0
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.0 0.0
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.0 0.0
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 1.126172 0.54344
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 1.142547 0.570131
avg(distribcell) group in group out nuclide moment mean std. dev.
0 (0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13,... 1 1 total P0 1.142547 0.570131
1 (0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13,... 1 1 total P1 0.447381 0.216322
2 (0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13,... 1 1 total P2 0.141202 0.066504
3 (0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13,... 1 1 total P3 0.039228 0.024621 avg(distribcell) group out nuclide mean std. dev.
0 (0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13,... 1 total 0.0 0.0
3 (0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13,... 1 1 total P3 0.039228 0.024621
avg(distribcell) group in group out nuclide moment mean std. dev.
0 (0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13,... 1 1 total P0 1.142547 0.570131
1 (0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13,... 1 1 total P1 0.447381 0.216322
2 (0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13,... 1 1 total P2 0.141202 0.066504
3 (0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13,... 1 1 total P3 0.039228 0.024621
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 1.0 0.529717
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.0 0.0
avg(distribcell) group out nuclide mean std. dev.
0 (0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13,... 1 total 0.0 0.0

10
tests/test_mgxs_library_distribcell/test_mgxs_library_distribcell.py
View file

@ -12,10 +12,6 @@ import openmc.mgxs
class MGXSTestHarness(PyAPITestHarness):
def _build_inputs(self):
# The openmc.mgxs module needs a summary.h5 file
self._input_set.settings.output = {'summary': True}
# Generate inputs using parent class routine
super(MGXSTestHarness, self)._build_inputs()
@ -26,8 +22,8 @@ class MGXSTestHarness(PyAPITestHarness):
# for one material-filled cell in the geometry
self.mgxs_lib = openmc.mgxs.Library(self._input_set.geometry)
self.mgxs_lib.by_nuclide = False
self.mgxs_lib.mgxs_types = ['transport', 'nu-fission',
'nu-scatter matrix', 'chi']
# Test all MGXS types
self.mgxs_lib.mgxs_types = openmc.mgxs.MGXS_TYPES
self.mgxs_lib.energy_groups = energy_groups
self.mgxs_lib.legendre_order = 3
self.mgxs_lib.domain_type = 'distribcell'
@ -59,7 +55,7 @@ class MGXSTestHarness(PyAPITestHarness):
for mgxs_type in avg_lib.mgxs_types:
mgxs = avg_lib.get_mgxs(domain, mgxs_type)
df = mgxs.get_pandas_dataframe()
outstr += df.to_string()
outstr += df.to_string() + '\n'
# Hash the results if necessary
if hash_output:

2
tests/test_mgxs_library_hdf5/inputs_true.dat
View file

@ -1 +1 @@
104e7fb527770ac5d3fc636da7716e8fb05d55761253d30516c899f466e6b38ffd881611a3d0cdf65c6af058c32f6f6758c68782be7a170d21024bdae751862f
317a63a9dd3bfd84e969667b00f46018e56c04c356461a75103f63569e6b70c84d0da7f5e611faaf1b2631330b05ab4346223d3d843018ce0ce8876671a450c0

381
tests/test_mgxs_library_hdf5/results_true.dat
View file

@ -1,240 +1,195 @@
domain=1 type=transport
[ 0.37274472 0.86160691]
[ 0.02426918 0.03234902]
domain=1 type=nu-fission
[ 0.02178897 0.71407658]
[ 0.00118187 0.04055185]
domain=1 type=nu-scatter matrix
[[[ 3.81546297e-01 4.43012537e-02 2.06462886e-02 1.36952959e-02]
[ 1.55945353e-03 -5.97269486e-04 -2.38789528e-04 1.75508083e-04]]
domain=10000 type=total
[ 0.41482549 0.66016992]
[ 0.02279291 0.04751893]
domain=10000 type=transport
[ 0.35685964 0.64764766]
[ 0.0254936 0.02370374]
domain=10000 type=nu-transport
[ 0.35685964 0.64764766]
[ 0.0254936 0.02370374]
domain=10000 type=absorption
[ 0.02740784 0.26451074]
[ 0.0026925 0.02336708]
domain=10000 type=capture
[ 0.01984455 0.07171935]
[ 0.0026433 0.02520786]
domain=10000 type=fission
[ 0.00756329 0.19279139]
[ 0.00050848 0.01710592]
domain=10000 type=nu-fission
[ 0.01943174 0.46977478]
[ 0.00132298 0.041682 ]
domain=10000 type=kappa-fission
[ 1.47456982 37.28689641]
[ 0.09923532 3.30837772]
domain=10000 type=scatter
[ 0.38741765 0.39565918]
[ 0.02062573 0.02512506]
domain=10000 type=nu-scatter
[ 0.38518839 0.4123894 ]
[ 0.02694562 0.01542528]
domain=10000 type=scatter matrix
[[[ 3.84199458e-01 5.18702843e-02 2.00688453e-02 9.47771571e-03]
[ 9.88930393e-04 -2.07234596e-04 -1.03366181e-04 2.34290623e-04]]
[[ 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00]
[ 4.03915981e-01 -1.13103276e-02 -1.48065932e-02 -6.85505346e-03]]]
[[[ 0.02403322 0.00472203 0.00253903 0.00222437]
[ 0.00051015 0.00022485 0.00022157 0.00020939]]
[[ 9.24639909e-04 -7.67704968e-04 4.93788872e-04 -1.71497229e-04]
[ 4.11464759e-01 1.64817280e-02 6.37149049e-03 -1.04991221e-02]]]
[[[ 0.02700101 0.00698255 0.0028465 0.00223352]
[ 0.00048242 0.00014901 0.00018432 0.00012817]]
[[ 0. 0. 0. 0. ]
[ 0.01896646 0.00783919 0.00862908 0.00904704]]]
domain=1 type=chi
[[ 0.00092488 0.00076791 0.00049392 0.00017154]
[ 0.01524494 0.00450173 0.01055075 0.01043819]]]
domain=10000 type=nu-scatter matrix
[[[ 3.84199458e-01 5.18702843e-02 2.00688453e-02 9.47771571e-03]
[ 9.88930393e-04 -2.07234596e-04 -1.03366181e-04 2.34290623e-04]]
[[ 9.24639909e-04 -7.67704968e-04 4.93788872e-04 -1.71497229e-04]
[ 4.11464759e-01 1.64817280e-02 6.37149049e-03 -1.04991221e-02]]]
[[[ 0.02700101 0.00698255 0.0028465 0.00223352]
[ 0.00048242 0.00014901 0.00018432 0.00012817]]
[[ 0.00092488 0.00076791 0.00049392 0.00017154]
[ 0.01524494 0.00450173 0.01055075 0.01043819]]]
domain=10000 type=multiplicity matrix
[[ 1. 1.]
[ 1. 1.]]
[[ 0.07851646 0.68718427]
[ 1.41421356 0.04113035]]
domain=10000 type=nu-fission matrix
[[ 0.02014243 0. ]
[ 0.45436647 0. ]]
[[ 0.00314909 0. ]
[ 0.02742551 0. ]]
domain=10000 type=chi
[ 1. 0.]
[ 0.05533329 0. ]
domain=2 type=transport
[ 0.23725441 0.28593027]
[ 0.00818357 0.04879593]
domain=2 type=nu-fission
[ 0.04607052 0. ]
domain=10001 type=total
[ 0.31373767 0.3008214 ]
[ 0.0155819 0.02805245]
domain=10001 type=transport
[ 0.27322787 0.31237484]
[ 0.03311537 0.04960583]
domain=10001 type=nu-transport
[ 0.27322787 0.31237484]
[ 0.03311537 0.04960583]
domain=10001 type=absorption
[ 0.00157499 0.00540038]
[ 0.00032255 0.00061814]
domain=10001 type=capture
[ 0.00157499 0.00540038]
[ 0.00032255 0.00061814]
domain=10001 type=fission
[ 0. 0.]
[ 0. 0.]
domain=2 type=nu-scatter matrix
[[[ 0.27311543 0.03586102 0.02970389 0.00224892]
domain=10001 type=nu-fission
[ 0. 0.]
[ 0. 0.]
domain=10001 type=kappa-fission
[ 0. 0.]
[ 0. 0.]
domain=10001 type=scatter
[ 0.31216268 0.29542102]
[ 0.01532192 0.02744549]
domain=10001 type=nu-scatter
[ 0.31012074 0.29626427]
[ 0.03378811 0.04379223]
domain=10001 type=scatter matrix
[[[ 0.31012074 0.03822959 0.02074494 0.0079643 ]
[ 0. 0. 0. 0. ]]
[[ 0. 0. 0. 0. ]
[ 0.26405068 -0.02187959 -0.01529469 0.01403395]]]
[[[ 0.00625287 0.00587756 0.00664018 0.00337568]
[ 0.29626427 -0.01121364 0.00883657 -0.00327007]]]
[[[ 0.03378811 0.008484 0.00469561 0.00373162]
[ 0. 0. 0. 0. ]]
[[ 0. 0. 0. 0. ]
[ 0.04539742 0.01221814 0.01027609 0.01431818]]]
domain=2 type=chi
[ 0. 0.]
[ 0. 0.]
domain=3 type=transport
[ 0.28690578 1.41815062]
[ 0.02740142 0.26530756]
domain=3 type=nu-fission
[ 0. 0.]
[ 0. 0.]
domain=3 type=nu-scatter matrix
[[[ 0.64334557 0.38340871 0.15218526 0.00303724]
[ 0.02618721 0.00736219 -0.00273849 -0.00271989]]
[[ 0. 0. 0. 0. ]
[ 1.92421362 0.4984312 0.09120485 0.01705441]]]
[[[ 0.02837604 0.01644677 0.00957372 0.00464802]
[ 0.00166461 0.00093414 0.00075617 0.00055807]]
[[ 0. 0. 0. 0. ]
[ 0.28406198 0.06342067 0.01372628 0.01391602]]]
domain=3 type=chi
[ 0. 0.]
[ 0. 0.]
domain=4 type=transport
[ 0.24244686 1.25395921]
[ 0.06103082 0.38836257]
domain=4 type=nu-fission
[ 0. 0.]
[ 0. 0.]
domain=4 type=nu-scatter matrix
[[[ 0.54394096 0.32601136 0.13113269 0.01210477]
[ 0.023662 0.00752551 -0.00272975 -0.0031405 ]]
[[ 0. 0. 0. 0. ]
[ 1.76464845 0.50069481 0.09902596 0.03297543]]]
[[[ 0.06542705 0.03860196 0.0174751 0.00607268]
[ 0.00308328 0.00130111 0.00084112 0.00057761]]
[[ 0. 0. 0. 0. ]
[ 0.41620952 0.12217802 0.03871874 0.02510259]]]
domain=4 type=chi
[ 0. 0.]
[ 0. 0.]
domain=5 type=transport
[ 0. 0.]
[ 0. 0.]
domain=5 type=nu-fission
[ 0. 0.]
[ 0. 0.]
domain=5 type=nu-scatter matrix
[[[ 0. 0. 0. 0.]
[ 0. 0. 0. 0.]]
[[ 0. 0. 0. 0.]
[ 0. 0. 0. 0.]]]
[[[ 0. 0. 0. 0.]
[ 0. 0. 0. 0.]]
[[ 0. 0. 0. 0.]
[ 0. 0. 0. 0.]]]
domain=5 type=chi
[ 0. 0.]
[ 0. 0.]
domain=6 type=transport
[ 0. 0.]
[ 0. 0.]
domain=6 type=nu-fission
[ 0. 0.]
[ 0. 0.]
domain=6 type=nu-scatter matrix
[[[ 0. 0. 0. 0.]
[ 0. 0. 0. 0.]]
[[ 0. 0. 0. 0.]
[ 0. 0. 0. 0.]]]
[[[ 0. 0. 0. 0.]
[ 0. 0. 0. 0.]]
[[ 0. 0. 0. 0.]
[ 0. 0. 0. 0.]]]
domain=6 type=chi
[ 0. 0.]
[ 0. 0.]
domain=7 type=transport
[ 0. 0.]
[ 0. 0.]
domain=7 type=nu-fission
[ 0. 0.]
[ 0. 0.]
domain=7 type=nu-scatter matrix
[[[ 0. 0. 0. 0.]
[ 0. 0. 0. 0.]]
[[ 0. 0. 0. 0.]
[ 0. 0. 0. 0.]]]
[[[ 0. 0. 0. 0.]
[ 0. 0. 0. 0.]]
[[ 0. 0. 0. 0.]
[ 0. 0. 0. 0.]]]
domain=7 type=chi
[ 0. 0.]
[ 0. 0.]
domain=8 type=transport
[ 0. 0.]
[ 0. 0.]
domain=8 type=nu-fission
[ 0. 0.]
[ 0. 0.]
domain=8 type=nu-scatter matrix
[[[ 0. 0. 0. 0.]
[ 0. 0. 0. 0.]]
[[ 0. 0. 0. 0.]
[ 0. 0. 0. 0.]]]
[[[ 0. 0. 0. 0.]
[ 0. 0. 0. 0.]]
[[ 0. 0. 0. 0.]
[ 0. 0. 0. 0.]]]
domain=8 type=chi
[ 0. 0.]
[ 0. 0.]
domain=9 type=transport
[ 0.60053598 0. ]
[ 0.74887543 0. ]
domain=9 type=nu-fission
[ 0. 0.]
[ 0. 0.]
domain=9 type=nu-scatter matrix
[[[ 0.72037987 0.11984389 0.03852204 0.05602285]
[ 0.04379223 0.01618037 0.01150396 0.00732885]]]
domain=10001 type=nu-scatter matrix
[[[ 0.31012074 0.03822959 0.02074494 0.0079643 ]
[ 0. 0. 0. 0. ]]
[[ 0. 0. 0. 0. ]
[ 0. 0. 0. 0. ]]]
[[[ 0.77101455 0.18469083 0.06448453 0.05059534]
[ 0.29626427 -0.01121364 0.00883657 -0.00327007]]]
[[[ 0.03378811 0.008484 0.00469561 0.00373162]
[ 0. 0. 0. 0. ]]
[[ 0. 0. 0. 0. ]
[ 0. 0. 0. 0. ]]]
domain=9 type=chi
[ 0.04379223 0.01618037 0.01150396 0.00732885]]]
domain=10001 type=multiplicity matrix
[[ 1. 0.]
[ 0. 1.]]
[[ 0.1087787 0. ]
[ 0. 0.14242717]]
domain=10001 type=nu-fission matrix
[[ 0. 0.]
[ 0. 0.]]
[[ 0. 0.]
[ 0. 0.]]
domain=10001 type=chi
[ 0. 0.]
[ 0. 0.]
domain=10 type=transport
[ 0.23551495 0. ]
[ 0.61397415 0. ]
domain=10 type=nu-fission
domain=10002 type=total
[ 0.66457226 2.05238401]
[ 0.03121475 0.22434291]
domain=10002 type=transport
[ 0.29056526 1.51643801]
[ 0.02385185 0.23519727]
domain=10002 type=nu-transport
[ 0.29056526 1.51643801]
[ 0.02385185 0.23519727]
domain=10002 type=absorption
[ 0.0006904 0.03168726]
[ 4.41475687e-05 3.74655858e-03]
domain=10002 type=capture
[ 0.0006904 0.03168726]
[ 4.41475687e-05 3.74655858e-03]
domain=10002 type=fission
[ 0. 0.]
[ 0. 0.]
domain=10 type=nu-scatter matrix
[[[ 0.50100891 0.26549396 0.14197875 0.07425836]
[ 0. 0. 0. 0. ]]
domain=10002 type=nu-fission
[ 0. 0.]
[ 0. 0.]
domain=10002 type=kappa-fission
[ 0. 0.]
[ 0. 0.]
domain=10002 type=scatter
[ 0.66388186 2.02069676]
[ 0.03117268 0.22060445]
domain=10002 type=nu-scatter
[ 0.6712692 2.03538833]
[ 0.02618637 0.25806033]
domain=10002 type=scatter matrix
[[[ 6.39901485e-01 3.81167449e-01 1.52391898e-01 9.14802229e-03]
[ 3.13677198e-02 8.75772321e-03 -2.56790106e-03 -3.78480288e-03]]
[[ 0. 0. 0. 0. ]
[ 0. 0. 0. 0. ]]]
[[[ 0.70853359 0.37546516 0.20078827 0.10501718]
[ 0. 0. 0. 0. ]]
[[ 4.43343134e-04 3.99960414e-04 3.19562707e-04 2.13846969e-04]
[ 2.03494499e+00 5.09940513e-01 1.11174609e-01 2.49884357e-02]]]
[[[ 2.47091228e-02 1.62432649e-02 8.15627770e-03 3.88856214e-03]
[ 1.72811290e-03 9.25670501e-04 1.01398475e-03 8.17075571e-04]]
[[ 0. 0. 0. 0. ]
[ 0. 0. 0. 0. ]]]
domain=10 type=chi
[ 0. 0.]
[ 0. 0.]
domain=11 type=transport
[ 0.18632392 0.94598628]
[ 0.63212919 1.59113341]
domain=11 type=nu-fission
[ 0. 0.]
[ 0. 0.]
domain=11 type=nu-scatter matrix
[[[ 0.47812753 0.32367878 0.14337507 0.05400336]
[ 0.03187517 0.00858456 -0.01246962 -0.01132019]]
[[ 4.44850393e-04 4.01320183e-04 3.20649143e-04 2.14573997e-04]
[ 2.57799889e-01 5.12359063e-02 1.30198170e-02 8.31235256e-03]]]
domain=10002 type=nu-scatter matrix
[[[ 6.39901485e-01 3.81167449e-01 1.52391898e-01 9.14802229e-03]
[ 3.13677198e-02 8.75772321e-03 -2.56790106e-03 -3.78480288e-03]]
[[ 0. 0. 0. 0. ]
[ 1.20124973 0.28661101 0.21819147 -0.04851424]]]
[[[ 0.67617444 0.45775092 0.20276296 0.07637229]
[ 0.0450783 0.0121404 0.01763471 0.01600917]]
[[ 4.43343134e-04 3.99960414e-04 3.19562707e-04 2.13846969e-04]
[ 2.03494499e+00 5.09940513e-01 1.11174609e-01 2.49884357e-02]]]
[[[ 2.47091228e-02 1.62432649e-02 8.15627770e-03 3.88856214e-03]
[ 1.72811290e-03 9.25670501e-04 1.01398475e-03 8.17075571e-04]]
[[ 0. 0. 0. 0. ]
[ 1.69882367 0.40532917 0.30856933 0.0686095 ]]]
domain=11 type=chi
[ 0. 0.]
[ 0. 0.]
domain=12 type=transport
[ 0.21329208 1.3909745 ]
[ 0.27144387 2.13734565]
domain=12 type=nu-fission
[ 0. 0.]
[ 0. 0.]
domain=12 type=nu-scatter matrix
[[[ 0.40859392 0.22254143 0.0909719 0.03100368]
[ 0.02723959 -0.01008785 -0.00694631 0.00969231]]
[[ 0. 0. 0. 0. ]
[ 1.57432766 0.22974802 0.01417839 0.03899727]]]
[[[ 0.27812309 0.14577636 0.06962553 0.03598053]
[ 0.02955488 0.01094529 0.00753673 0.01051613]]
[[ 0. 0. 0. 0. ]
[ 2.22643553 0.32491277 0.02005128 0.05515046]]]
domain=12 type=chi
[[ 4.44850393e-04 4.01320183e-04 3.20649143e-04 2.14573997e-04]
[ 2.57799889e-01 5.12359063e-02 1.30198170e-02 8.31235256e-03]]]
domain=10002 type=multiplicity matrix
[[ 1. 1.]
[ 1. 1.]]
[[ 0.03860919 0.06766735]
[ 1.41421356 0.13592921]]
domain=10002 type=nu-fission matrix
[[ 0. 0.]
[ 0. 0.]]
[[ 0. 0.]
[ 0. 0.]]
domain=10002 type=chi
[ 0. 0.]
[ 0. 0.]

14
tests/test_mgxs_library_hdf5/test_mgxs_library_hdf5.py
View file

@ -7,27 +7,28 @@ import hashlib
import h5py
sys.path.insert(0, os.pardir)
from testing_harness import PyAPITestHarness
from input_set import PinCellInputSet
import openmc
import openmc.mgxs
class MGXSTestHarness(PyAPITestHarness):
def _build_inputs(self):
# The openmc.mgxs module needs a summary.h5 file
self._input_set.settings.output = {'summary': True}
# Set the input set to use the pincell model
self._input_set = PinCellInputSet()
# Generate inputs using parent class routine
super(MGXSTestHarness, self)._build_inputs()
# Initialize a two-group structure
energy_groups = openmc.mgxs.EnergyGroups(group_edges=[0, 0.625e-6, 20.])
energy_groups = openmc.mgxs.EnergyGroups(group_edges=[0, 0.625e-6,
20.])
# Initialize MGXS Library for a few cross section types
self.mgxs_lib = openmc.mgxs.Library(self._input_set.geometry)
self.mgxs_lib.by_nuclide = False
self.mgxs_lib.mgxs_types = ['transport', 'nu-fission',
'nu-scatter matrix', 'chi']
# Test all MGXS types
self.mgxs_lib.mgxs_types = openmc.mgxs.MGXS_TYPES
self.mgxs_lib.energy_groups = energy_groups
self.mgxs_lib.legendre_order = 3
self.mgxs_lib.domain_type = 'material'
@ -75,7 +76,6 @@ class MGXSTestHarness(PyAPITestHarness):
return outstr
def _cleanup(self):
super(MGXSTestHarness, self)._cleanup()
f = os.path.join(os.getcwd(), 'tallies.xml')

2
tests/test_mgxs_library_no_nuclides/inputs_true.dat
View file

@ -1 +1 @@
104e7fb527770ac5d3fc636da7716e8fb05d55761253d30516c899f466e6b38ffd881611a3d0cdf65c6af058c32f6f6758c68782be7a170d21024bdae751862f
317a63a9dd3bfd84e969667b00f46018e56c04c356461a75103f63569e6b70c84d0da7f5e611faaf1b2631330b05ab4346223d3d843018ce0ce8876671a450c0

494
tests/test_mgxs_library_no_nuclides/results_true.dat
View file

@ -1,265 +1,231 @@
material group in nuclide mean std. dev.
1 1 1 total 0.372745 0.024269
0 1 2 total 0.861607 0.032349 material group in nuclide mean std. dev.
1 1 1 total 0.021789 0.001182
0 1 2 total 0.714077 0.040552 material group in group out nuclide moment mean std. dev.
12 1 1 1 total P0 0.381546 0.024033
13 1 1 1 total P1 0.044301 0.004722
14 1 1 1 total P2 0.020646 0.002539
15 1 1 1 total P3 0.013695 0.002224
8 1 1 2 total P0 0.001559 0.000510
9 1 1 2 total P1 -0.000597 0.000225
10 1 1 2 total P2 -0.000239 0.000222
11 1 1 2 total P3 0.000176 0.000209
4 1 2 1 total P0 0.000000 0.000000
5 1 2 1 total P1 0.000000 0.000000
6 1 2 1 total P2 0.000000 0.000000
7 1 2 1 total P3 0.000000 0.000000
0 1 2 2 total P0 0.403916 0.018966
1 1 2 2 total P1 -0.011310 0.007839
2 1 2 2 total P2 -0.014807 0.008629
3 1 2 2 total P3 -0.006855 0.009047 material group out nuclide mean std. dev.
1 1 1 total 1.0 0.055333
0 1 2 total 0.0 0.000000 material group in nuclide mean std. dev.
1 2 1 total 0.237254 0.008184
0 2 2 total 0.285930 0.048796 material group in nuclide mean std. dev.
1 2 1 total 0.0 0.0
0 2 2 total 0.0 0.0 material group in group out nuclide moment mean std. dev.
12 2 1 1 total P0 0.273115 0.006253
13 2 1 1 total P1 0.035861 0.005878
14 2 1 1 total P2 0.029704 0.006640
15 2 1 1 total P3 0.002249 0.003376
8 2 1 2 total P0 0.000000 0.000000
9 2 1 2 total P1 0.000000 0.000000
10 2 1 2 total P2 0.000000 0.000000
11 2 1 2 total P3 0.000000 0.000000
4 2 2 1 total P0 0.000000 0.000000
5 2 2 1 total P1 0.000000 0.000000
6 2 2 1 total P2 0.000000 0.000000
7 2 2 1 total P3 0.000000 0.000000
0 2 2 2 total P0 0.264051 0.045397
1 2 2 2 total P1 -0.021880 0.012218
2 2 2 2 total P2 -0.015295 0.010276
3 2 2 2 total P3 0.014034 0.014318 material group out nuclide mean std. dev.
1 2 1 total 0.0 0.0
0 2 2 total 0.0 0.0 material group in nuclide mean std. dev.
1 3 1 total 0.286906 0.027401
0 3 2 total 1.418151 0.265308 material group in nuclide mean std. dev.
1 3 1 total 0.0 0.0
0 3 2 total 0.0 0.0 material group in group out nuclide moment mean std. dev.
12 3 1 1 total P0 0.643346 0.028376
13 3 1 1 total P1 0.383409 0.016447
14 3 1 1 total P2 0.152185 0.009574
15 3 1 1 total P3 0.003037 0.004648
8 3 1 2 total P0 0.026187 0.001665
9 3 1 2 total P1 0.007362 0.000934
10 3 1 2 total P2 -0.002738 0.000756
11 3 1 2 total P3 -0.002720 0.000558
4 3 2 1 total P0 0.000000 0.000000
5 3 2 1 total P1 0.000000 0.000000
6 3 2 1 total P2 0.000000 0.000000
7 3 2 1 total P3 0.000000 0.000000
0 3 2 2 total P0 1.924214 0.284062
1 3 2 2 total P1 0.498431 0.063421
2 3 2 2 total P2 0.091205 0.013726
3 3 2 2 total P3 0.017054 0.013916 material group out nuclide mean std. dev.
1 3 1 total 0.0 0.0
0 3 2 total 0.0 0.0 material group in nuclide mean std. dev.
1 4 1 total 0.242447 0.061031
0 4 2 total 1.253959 0.388363 material group in nuclide mean std. dev.
1 4 1 total 0.0 0.0
0 4 2 total 0.0 0.0 material group in group out nuclide moment mean std. dev.
12 4 1 1 total P0 0.543941 0.065427
13 4 1 1 total P1 0.326011 0.038602
14 4 1 1 total P2 0.131133 0.017475
15 4 1 1 total P3 0.012105 0.006073
8 4 1 2 total P0 0.023662 0.003083
9 4 1 2 total P1 0.007526 0.001301
10 4 1 2 total P2 -0.002730 0.000841
11 4 1 2 total P3 -0.003140 0.000578
4 4 2 1 total P0 0.000000 0.000000
5 4 2 1 total P1 0.000000 0.000000
6 4 2 1 total P2 0.000000 0.000000
7 4 2 1 total P3 0.000000 0.000000
0 4 2 2 total P0 1.764648 0.416210
1 4 2 2 total P1 0.500695 0.122178
2 4 2 2 total P2 0.099026 0.038719
3 4 2 2 total P3 0.032975 0.025103 material group out nuclide mean std. dev.
1 4 1 total 0.0 0.0
0 4 2 total 0.0 0.0 material group in nuclide mean std. dev.
1 5 1 total 0.0 0.0
0 5 2 total 0.0 0.0 material group in nuclide mean std. dev.
1 5 1 total 0.0 0.0
0 5 2 total 0.0 0.0 material group in group out nuclide moment mean std. dev.
12 5 1 1 total P0 0.0 0.0
13 5 1 1 total P1 0.0 0.0
14 5 1 1 total P2 0.0 0.0
15 5 1 1 total P3 0.0 0.0
8 5 1 2 total P0 0.0 0.0
9 5 1 2 total P1 0.0 0.0
10 5 1 2 total P2 0.0 0.0
11 5 1 2 total P3 0.0 0.0
4 5 2 1 total P0 0.0 0.0
5 5 2 1 total P1 0.0 0.0
6 5 2 1 total P2 0.0 0.0
7 5 2 1 total P3 0.0 0.0
0 5 2 2 total P0 0.0 0.0
1 5 2 2 total P1 0.0 0.0
2 5 2 2 total P2 0.0 0.0
3 5 2 2 total P3 0.0 0.0 material group out nuclide mean std. dev.
1 5 1 total 0.0 0.0
0 5 2 total 0.0 0.0 material group in nuclide mean std. dev.
1 6 1 total 0.0 0.0
0 6 2 total 0.0 0.0 material group in nuclide mean std. dev.
1 6 1 total 0.0 0.0
0 6 2 total 0.0 0.0 material group in group out nuclide moment mean std. dev.
12 6 1 1 total P0 0.0 0.0
13 6 1 1 total P1 0.0 0.0
14 6 1 1 total P2 0.0 0.0
15 6 1 1 total P3 0.0 0.0
8 6 1 2 total P0 0.0 0.0
9 6 1 2 total P1 0.0 0.0
10 6 1 2 total P2 0.0 0.0
11 6 1 2 total P3 0.0 0.0
4 6 2 1 total P0 0.0 0.0
5 6 2 1 total P1 0.0 0.0
6 6 2 1 total P2 0.0 0.0
7 6 2 1 total P3 0.0 0.0
0 6 2 2 total P0 0.0 0.0
1 6 2 2 total P1 0.0 0.0
2 6 2 2 total P2 0.0 0.0
3 6 2 2 total P3 0.0 0.0 material group out nuclide mean std. dev.
1 6 1 total 0.0 0.0
0 6 2 total 0.0 0.0 material group in nuclide mean std. dev.
1 7 1 total 0.0 0.0
0 7 2 total 0.0 0.0 material group in nuclide mean std. dev.
1 7 1 total 0.0 0.0
0 7 2 total 0.0 0.0 material group in group out nuclide moment mean std. dev.
12 7 1 1 total P0 0.0 0.0
13 7 1 1 total P1 0.0 0.0
14 7 1 1 total P2 0.0 0.0
15 7 1 1 total P3 0.0 0.0
8 7 1 2 total P0 0.0 0.0
9 7 1 2 total P1 0.0 0.0
10 7 1 2 total P2 0.0 0.0
11 7 1 2 total P3 0.0 0.0
4 7 2 1 total P0 0.0 0.0
5 7 2 1 total P1 0.0 0.0
6 7 2 1 total P2 0.0 0.0
7 7 2 1 total P3 0.0 0.0
0 7 2 2 total P0 0.0 0.0
1 7 2 2 total P1 0.0 0.0
2 7 2 2 total P2 0.0 0.0
3 7 2 2 total P3 0.0 0.0 material group out nuclide mean std. dev.
1 7 1 total 0.0 0.0
0 7 2 total 0.0 0.0 material group in nuclide mean std. dev.
1 8 1 total 0.0 0.0
0 8 2 total 0.0 0.0 material group in nuclide mean std. dev.
1 8 1 total 0.0 0.0
0 8 2 total 0.0 0.0 material group in group out nuclide moment mean std. dev.
12 8 1 1 total P0 0.0 0.0
13 8 1 1 total P1 0.0 0.0
14 8 1 1 total P2 0.0 0.0
15 8 1 1 total P3 0.0 0.0
8 8 1 2 total P0 0.0 0.0
9 8 1 2 total P1 0.0 0.0
10 8 1 2 total P2 0.0 0.0
11 8 1 2 total P3 0.0 0.0
4 8 2 1 total P0 0.0 0.0
5 8 2 1 total P1 0.0 0.0
6 8 2 1 total P2 0.0 0.0
7 8 2 1 total P3 0.0 0.0
0 8 2 2 total P0 0.0 0.0
1 8 2 2 total P1 0.0 0.0
2 8 2 2 total P2 0.0 0.0
3 8 2 2 total P3 0.0 0.0 material group out nuclide mean std. dev.
1 8 1 total 0.0 0.0
0 8 2 total 0.0 0.0 material group in nuclide mean std. dev.
1 9 1 total 0.600536 0.748875
0 9 2 total 0.000000 0.000000 material group in nuclide mean std. dev.
1 9 1 total 0.0 0.0
0 9 2 total 0.0 0.0 material group in group out nuclide moment mean std. dev.
12 9 1 1 total P0 0.720380 0.771015
13 9 1 1 total P1 0.119844 0.184691
14 9 1 1 total P2 0.038522 0.064485
15 9 1 1 total P3 0.056023 0.050595
8 9 1 2 total P0 0.000000 0.000000
9 9 1 2 total P1 0.000000 0.000000
10 9 1 2 total P2 0.000000 0.000000
11 9 1 2 total P3 0.000000 0.000000
4 9 2 1 total P0 0.000000 0.000000
5 9 2 1 total P1 0.000000 0.000000
6 9 2 1 total P2 0.000000 0.000000
7 9 2 1 total P3 0.000000 0.000000
0 9 2 2 total P0 0.000000 0.000000
1 9 2 2 total P1 0.000000 0.000000
2 9 2 2 total P2 0.000000 0.000000
3 9 2 2 total P3 0.000000 0.000000 material group out nuclide mean std. dev.
1 9 1 total 0.0 0.0
0 9 2 total 0.0 0.0 material group in nuclide mean std. dev.
1 10 1 total 0.235515 0.613974
0 10 2 total 0.000000 0.000000 material group in nuclide mean std. dev.
1 10 1 total 0.0 0.0
0 10 2 total 0.0 0.0 material group in group out nuclide moment mean std. dev.
12 10 1 1 total P0 0.501009 0.708534
13 10 1 1 total P1 0.265494 0.375465
14 10 1 1 total P2 0.141979 0.200788
15 10 1 1 total P3 0.074258 0.105017
8 10 1 2 total P0 0.000000 0.000000
9 10 1 2 total P1 0.000000 0.000000
10 10 1 2 total P2 0.000000 0.000000
11 10 1 2 total P3 0.000000 0.000000
4 10 2 1 total P0 0.000000 0.000000
5 10 2 1 total P1 0.000000 0.000000
6 10 2 1 total P2 0.000000 0.000000
7 10 2 1 total P3 0.000000 0.000000
0 10 2 2 total P0 0.000000 0.000000
1 10 2 2 total P1 0.000000 0.000000
2 10 2 2 total P2 0.000000 0.000000
3 10 2 2 total P3 0.000000 0.000000 material group out nuclide mean std. dev.
1 10 1 total 0.0 0.0
0 10 2 total 0.0 0.0 material group in nuclide mean std. dev.
1 11 1 total 0.186324 0.632129
0 11 2 total 0.945986 1.591133 material group in nuclide mean std. dev.
1 11 1 total 0.0 0.0
0 11 2 total 0.0 0.0 material group in group out nuclide moment mean std. dev.
12 11 1 1 total P0 0.478128 0.676174
13 11 1 1 total P1 0.323679 0.457751
14 11 1 1 total P2 0.143375 0.202763
15 11 1 1 total P3 0.054003 0.076372
8 11 1 2 total P0 0.031875 0.045078
9 11 1 2 total P1 0.008585 0.012140
10 11 1 2 total P2 -0.012470 0.017635
11 11 1 2 total P3 -0.011320 0.016009
4 11 2 1 total P0 0.000000 0.000000
5 11 2 1 total P1 0.000000 0.000000
6 11 2 1 total P2 0.000000 0.000000
7 11 2 1 total P3 0.000000 0.000000
0 11 2 2 total P0 1.201250 1.698824
1 11 2 2 total P1 0.286611 0.405329
2 11 2 2 total P2 0.218191 0.308569
3 11 2 2 total P3 -0.048514 0.068609 material group out nuclide mean std. dev.
1 11 1 total 0.0 0.0
0 11 2 total 0.0 0.0 material group in nuclide mean std. dev.
1 12 1 total 0.213292 0.271444
0 12 2 total 1.390975 2.137346 material group in nuclide mean std. dev.
1 12 1 total 0.0 0.0
0 12 2 total 0.0 0.0 material group in group out nuclide moment mean std. dev.
12 12 1 1 total P0 0.408594 0.278123
13 12 1 1 total P1 0.222541 0.145776
14 12 1 1 total P2 0.090972 0.069626
15 12 1 1 total P3 0.031004 0.035981
8 12 1 2 total P0 0.027240 0.029555
9 12 1 2 total P1 -0.010088 0.010945
10 12 1 2 total P2 -0.006946 0.007537
11 12 1 2 total P3 0.009692 0.010516
4 12 2 1 total P0 0.000000 0.000000
5 12 2 1 total P1 0.000000 0.000000
6 12 2 1 total P2 0.000000 0.000000
7 12 2 1 total P3 0.000000 0.000000
0 12 2 2 total P0 1.574328 2.226436
1 12 2 2 total P1 0.229748 0.324913
2 12 2 2 total P2 0.014178 0.020051
3 12 2 2 total P3 0.038997 0.055150 material group out nuclide mean std. dev.
1 12 1 total 0.0 0.0
0 12 2 total 0.0 0.0
1 10000 1 total 0.414825 0.022793
0 10000 2 total 0.660170 0.047519
material group in nuclide mean std. dev.
1 10000 1 total 0.356860 0.025494
0 10000 2 total 0.647648 0.023704
material group in nuclide mean std. dev.
1 10000 1 total 0.356860 0.025494
0 10000 2 total 0.647648 0.023704
material group in nuclide mean std. dev.
1 10000 1 total 0.027408 0.002692
0 10000 2 total 0.264511 0.023367
material group in nuclide mean std. dev.
1 10000 1 total 0.019845 0.002643
0 10000 2 total 0.071719 0.025208
material group in nuclide mean std. dev.
1 10000 1 total 0.007563 0.000508
0 10000 2 total 0.192791 0.017106
material group in nuclide mean std. dev.
1 10000 1 total 0.019432 0.001323
0 10000 2 total 0.469775 0.041682
material group in nuclide mean std. dev.
1 10000 1 total 1.474570 0.099235
0 10000 2 total 37.286896 3.308378
material group in nuclide mean std. dev.
1 10000 1 total 0.387418 0.020626
0 10000 2 total 0.395659 0.025125
material group in nuclide mean std. dev.
1 10000 1 total 0.385188 0.026946
0 10000 2 total 0.412389 0.015425
material group in group out nuclide moment mean std. dev.
12 10000 1 1 total P0 0.384199 0.027001
13 10000 1 1 total P1 0.051870 0.006983
14 10000 1 1 total P2 0.020069 0.002846
15 10000 1 1 total P3 0.009478 0.002234
8 10000 1 2 total P0 0.000989 0.000482
9 10000 1 2 total P1 -0.000207 0.000149
10 10000 1 2 total P2 -0.000103 0.000184
11 10000 1 2 total P3 0.000234 0.000128
4 10000 2 1 total P0 0.000925 0.000925
5 10000 2 1 total P1 -0.000768 0.000768
6 10000 2 1 total P2 0.000494 0.000494
7 10000 2 1 total P3 -0.000171 0.000172
0 10000 2 2 total P0 0.411465 0.015245
1 10000 2 2 total P1 0.016482 0.004502
2 10000 2 2 total P2 0.006371 0.010551
3 10000 2 2 total P3 -0.010499 0.010438
material group in group out nuclide moment mean std. dev.
12 10000 1 1 total P0 0.384199 0.027001
13 10000 1 1 total P1 0.051870 0.006983
14 10000 1 1 total P2 0.020069 0.002846
15 10000 1 1 total P3 0.009478 0.002234
8 10000 1 2 total P0 0.000989 0.000482
9 10000 1 2 total P1 -0.000207 0.000149
10 10000 1 2 total P2 -0.000103 0.000184
11 10000 1 2 total P3 0.000234 0.000128
4 10000 2 1 total P0 0.000925 0.000925
5 10000 2 1 total P1 -0.000768 0.000768
6 10000 2 1 total P2 0.000494 0.000494
7 10000 2 1 total P3 -0.000171 0.000172
0 10000 2 2 total P0 0.411465 0.015245
1 10000 2 2 total P1 0.016482 0.004502
2 10000 2 2 total P2 0.006371 0.010551
3 10000 2 2 total P3 -0.010499 0.010438
material group in group out nuclide mean std. dev.
3 10000 1 1 total 1.0 0.078516
2 10000 1 2 total 1.0 0.687184
1 10000 2 1 total 1.0 1.414214
0 10000 2 2 total 1.0 0.041130
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.454366 0.027426
0 10000 2 2 total 0.000000 0.000000
material group out nuclide mean std. dev.
1 10000 1 total 1.0 0.046071
0 10000 2 total 0.0 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
material group in nuclide mean std. dev.
1 10001 1 total 0.273228 0.033115
0 10001 2 total 0.312375 0.049606
material group in nuclide mean std. dev.
1 10001 1 total 0.273228 0.033115
0 10001 2 total 0.312375 0.049606
material group in nuclide mean std. dev.
1 10001 1 total 0.001575 0.000323
0 10001 2 total 0.005400 0.000618
material group in nuclide mean std. dev.
1 10001 1 total 0.001575 0.000323
0 10001 2 total 0.005400 0.000618
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 nuclide mean std. dev.
1 10001 1 total 0.0 0.0
0 10001 2 total 0.0 0.0
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 nuclide mean std. dev.
1 10001 1 total 0.312163 0.015322
0 10001 2 total 0.295421 0.027445
material group in nuclide mean std. dev.
1 10001 1 total 0.310121 0.033788
0 10001 2 total 0.296264 0.043792
material group in group out nuclide moment mean std. dev.
12 10001 1 1 total P0 0.310121 0.033788
13 10001 1 1 total P1 0.038230 0.008484
14 10001 1 1 total P2 0.020745 0.004696
15 10001 1 1 total P3 0.007964 0.003732
8 10001 1 2 total P0 0.000000 0.000000
9 10001 1 2 total P1 0.000000 0.000000
10 10001 1 2 total P2 0.000000 0.000000
11 10001 1 2 total P3 0.000000 0.000000
4 10001 2 1 total P0 0.000000 0.000000
5 10001 2 1 total P1 0.000000 0.000000
6 10001 2 1 total P2 0.000000 0.000000
7 10001 2 1 total P3 0.000000 0.000000
0 10001 2 2 total P0 0.296264 0.043792
1 10001 2 2 total P1 -0.011214 0.016180
2 10001 2 2 total P2 0.008837 0.011504
3 10001 2 2 total P3 -0.003270 0.007329
material group in group out nuclide moment mean std. dev.
12 10001 1 1 total P0 0.310121 0.033788
13 10001 1 1 total P1 0.038230 0.008484
14 10001 1 1 total P2 0.020745 0.004696
15 10001 1 1 total P3 0.007964 0.003732
8 10001 1 2 total P0 0.000000 0.000000
9 10001 1 2 total P1 0.000000 0.000000
10 10001 1 2 total P2 0.000000 0.000000
11 10001 1 2 total P3 0.000000 0.000000
4 10001 2 1 total P0 0.000000 0.000000
5 10001 2 1 total P1 0.000000 0.000000
6 10001 2 1 total P2 0.000000 0.000000
7 10001 2 1 total P3 0.000000 0.000000
0 10001 2 2 total P0 0.296264 0.043792
1 10001 2 2 total P1 -0.011214 0.016180
2 10001 2 2 total P2 0.008837 0.011504
3 10001 2 2 total P3 -0.003270 0.007329
material group in group out nuclide mean std. dev.
3 10001 1 1 total 1.0 0.108779
2 10001 1 2 total 0.0 0.000000
1 10001 2 1 total 0.0 0.000000
0 10001 2 2 total 1.0 0.142427
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 group out nuclide mean std. dev.
1 10001 1 total 0.0 0.0
0 10001 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
material group in nuclide mean std. dev.
1 10002 1 total 0.290565 0.023852
0 10002 2 total 1.516438 0.235197
material group in nuclide mean std. dev.
1 10002 1 total 0.290565 0.023852
0 10002 2 total 1.516438 0.235197
material group in nuclide mean std. dev.
1 10002 1 total 0.000690 0.000044
0 10002 2 total 0.031687 0.003747
material group in nuclide mean std. dev.
1 10002 1 total 0.000690 0.000044
0 10002 2 total 0.031687 0.003747
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 nuclide mean std. dev.
1 10002 1 total 0.0 0.0
0 10002 2 total 0.0 0.0
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 nuclide mean std. dev.
1 10002 1 total 0.663882 0.031173
0 10002 2 total 2.020697 0.220604
material group in nuclide mean std. dev.
1 10002 1 total 0.671269 0.026186
0 10002 2 total 2.035388 0.258060
material group in group out nuclide moment mean std. dev.
12 10002 1 1 total P0 0.639901 0.024709
13 10002 1 1 total P1 0.381167 0.016243
14 10002 1 1 total P2 0.152392 0.008156
15 10002 1 1 total P3 0.009148 0.003889
8 10002 1 2 total P0 0.031368 0.001728
9 10002 1 2 total P1 0.008758 0.000926
10 10002 1 2 total P2 -0.002568 0.001014
11 10002 1 2 total P3 -0.003785 0.000817
4 10002 2 1 total P0 0.000443 0.000445
5 10002 2 1 total P1 0.000400 0.000401
6 10002 2 1 total P2 0.000320 0.000321
7 10002 2 1 total P3 0.000214 0.000215
0 10002 2 2 total P0 2.034945 0.257800
1 10002 2 2 total P1 0.509941 0.051236
2 10002 2 2 total P2 0.111175 0.013020
3 10002 2 2 total P3 0.024988 0.008312
material group in group out nuclide moment mean std. dev.
12 10002 1 1 total P0 0.639901 0.024709
13 10002 1 1 total P1 0.381167 0.016243
14 10002 1 1 total P2 0.152392 0.008156
15 10002 1 1 total P3 0.009148 0.003889
8 10002 1 2 total P0 0.031368 0.001728
9 10002 1 2 total P1 0.008758 0.000926
10 10002 1 2 total P2 -0.002568 0.001014
11 10002 1 2 total P3 -0.003785 0.000817
4 10002 2 1 total P0 0.000443 0.000445
5 10002 2 1 total P1 0.000400 0.000401
6 10002 2 1 total P2 0.000320 0.000321
7 10002 2 1 total P3 0.000214 0.000215
0 10002 2 2 total P0 2.034945 0.257800
1 10002 2 2 total P1 0.509941 0.051236
2 10002 2 2 total P2 0.111175 0.013020
3 10002 2 2 total P3 0.024988 0.008312
material group in group out nuclide mean std. dev.
3 10002 1 1 total 1.0 0.038609
2 10002 1 2 total 1.0 0.067667
1 10002 2 1 total 1.0 1.414214
0 10002 2 2 total 1.0 0.135929
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 group out nuclide mean std. dev.
1 10002 1 total 0.0 0.0
0 10002 2 total 0.0 0.0

15
tests/test_mgxs_library_no_nuclides/test_mgxs_library_no_nuclides.py
View file

@ -6,27 +6,28 @@ import glob
import hashlib
sys.path.insert(0, os.pardir)
from testing_harness import PyAPITestHarness
from input_set import PinCellInputSet
import openmc
import openmc.mgxs
class MGXSTestHarness(PyAPITestHarness):
def _build_inputs(self):
# The openmc.mgxs module needs a summary.h5 file
self._input_set.settings.output = {'summary': True}
# Set the input set to use the pincell model
self._input_set = PinCellInputSet()
# Generate inputs using parent class routine
super(MGXSTestHarness, self)._build_inputs()
# Initialize a two-group structure
energy_groups = openmc.mgxs.EnergyGroups(group_edges=[0, 0.625e-6, 20.])
energy_groups = openmc.mgxs.EnergyGroups(group_edges=[0, 0.625e-6,
20.])
# Initialize MGXS Library for a few cross section types
self.mgxs_lib = openmc.mgxs.Library(self._input_set.geometry)
self.mgxs_lib.by_nuclide = False
self.mgxs_lib.mgxs_types = ['transport', 'nu-fission',
'nu-scatter matrix', 'chi']
# Test all MGXS types
self.mgxs_lib.mgxs_types = openmc.mgxs.MGXS_TYPES
self.mgxs_lib.energy_groups = energy_groups
self.mgxs_lib.legendre_order = 3
self.mgxs_lib.domain_type = 'material'
@ -53,7 +54,7 @@ class MGXSTestHarness(PyAPITestHarness):
for mgxs_type in self.mgxs_lib.mgxs_types:
mgxs = self.mgxs_lib.get_mgxs(domain, mgxs_type)
df = mgxs.get_pandas_dataframe()
outstr += df.to_string()
outstr += df.to_string() + '\n'
# Hash the results if necessary
if hash_output:

2
tests/test_mgxs_library_nuclides/inputs_true.dat
View file

@ -1 +1 @@
791a2bd647b8bae03aafc39e29ff1ce1ffc44063b0d757ccba4e1eda6eb73b8a275020f4f5774b17dede49fbf15549787279c8b2fc45caba0097155b32e56fa8
eebb1469278f470b5859ed83e9b6526e7c4e3fed503bd22e414c6dc13b19b8e4cb6a44e3c14269e6e173f43056eda78268f455662ae119280bc18ea6a071dac7

2
tests/test_mgxs_library_nuclides/results_true.dat
View file

@ -1 +1 @@
1ee58383dc8ac46c5e0d72321cbc34b0dba531435d5e0e632cbbf9572eb7d669c8c8ad9f370345325afa0bdeb2f818b0f5204b7c4a7c4aaf58ded7acbd715ef8
a631b8a347f344d822e6300ed2576caa7c05a74daedeb4aaaabfb89570942cff1bbd47ad7f81306e668e12266404f7abdcf680fdfeb5a4835579892e32bf57e8

15
tests/test_mgxs_library_nuclides/test_mgxs_library_nuclides.py
View file

@ -6,27 +6,28 @@ import glob
import hashlib
sys.path.insert(0, os.pardir)
from testing_harness import PyAPITestHarness
from input_set import PinCellInputSet
import openmc
import openmc.mgxs
class MGXSTestHarness(PyAPITestHarness):
def _build_inputs(self):
# The openmc.mgxs module needs a summary.h5 file
self._input_set.settings.output = {'summary': True}
# Set the input set to use the pincell model
self._input_set = PinCellInputSet()
# Generate inputs using parent class routine
super(MGXSTestHarness, self)._build_inputs()
# Initialize a two-group structure
energy_groups = openmc.mgxs.EnergyGroups(group_edges=[0, 0.625e-6, 20.])
energy_groups = openmc.mgxs.EnergyGroups(group_edges=[0, 0.625e-6,
20.])
# Initialize MGXS Library for a few cross section types
self.mgxs_lib = openmc.mgxs.Library(self._input_set.geometry)
self.mgxs_lib.by_nuclide = True
self.mgxs_lib.mgxs_types = ['transport', 'nu-fission',
'nu-scatter matrix', 'chi']
# Test all MGXS types
self.mgxs_lib.mgxs_types = openmc.mgxs.MGXS_TYPES
self.mgxs_lib.energy_groups = energy_groups
self.mgxs_lib.legendre_order = 3
self.mgxs_lib.domain_type = 'material'
@ -53,7 +54,7 @@ class MGXSTestHarness(PyAPITestHarness):
for mgxs_type in self.mgxs_lib.mgxs_types:
mgxs = self.mgxs_lib.get_mgxs(domain, mgxs_type)
df = mgxs.get_pandas_dataframe()
outstr += df.to_string()
outstr += df.to_string() + '\n'
# Hash the results if necessary
if hash_output:

2
tests/test_multipole/results_true.dat
View file

@ -3,7 +3,7 @@ k-combined:
Cell
ID = 11
Name =
Material = 2
Fill = Material 2
Region = -10000
Rotation = None
Temperature = [ 500. 0. 700. 800.]

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

@ -0,0 +1 @@
af589996f2930337afe34ba9894098ff5efe3b29b6e927117220b718bf29b630ffdbc931754d465a8e8100125a8aa997dbe10aab322b43f69d59710573996a6d

2
tests/test_periodic/results_true.dat Normal file
View file

@ -0,0 +1,2 @@
k-combined:
1.040109E+00 6.527490E-02

14
tests/test_periodic/tallies.xml Normal file
View file

@ -0,0 +1,14 @@
<?xml version="1.0"?>
<tallies>
<mesh id="1">
<type>regular</type>
<lower_left>-200. -1e50</lower_left>
<upper_right>200. 1e50</upper_right>
<dimension>50 1</dimension>
</mesh>
<tally id="1">
<estimator>collision</estimator>
<filter type="mesh" bins="1" />
<scores>fission</scores>
</tally>
</tallies>

60
tests/test_periodic/test_periodic.py Normal file
View file

@ -0,0 +1,60 @@
#!/usr/bin/env python
import os
import sys
sys.path.insert(0, os.pardir)
from testing_harness import PyAPITestHarness
import openmc
class PeriodicTest(PyAPITestHarness):
def _build_inputs(self):
# Define materials
water = openmc.Material(1)
water.add_nuclide('H-1', 2.0)
water.add_nuclide('O-16', 1.0)
water.add_s_alpha_beta('HH2O', '71t')
water.set_density('g/cc', 1.0)
fuel = openmc.Material(2)
fuel.add_nuclide('U-235', 1.0)
fuel.set_density('g/cc', 4.5)
materials = openmc.Materials((water, fuel))
materials.default_xs = '71c'
materials.export_to_xml()
# Define geometry
x_min = openmc.XPlane(1, x0=-5., boundary_type='periodic')
x_max = openmc.XPlane(2, x0=5., boundary_type='periodic')
x_max.periodic_surface = x_min
y_min = openmc.YPlane(3, y0=-5., boundary_type='periodic')
y_max = openmc.YPlane(4, y0=5., boundary_type='periodic')
z_min = openmc.ZPlane(5, z0=-5., boundary_type='reflective')
z_max = openmc.ZPlane(6, z0=5., boundary_type='reflective')
z_cyl = openmc.ZCylinder(7, x0=-2.5, y0=2.5, R=2.0)
outside_cyl = openmc.Cell(1, fill=water, region=(
+x_min & -x_max & +y_min & -y_max & +z_min & -z_max & +z_cyl))
inside_cyl = openmc.Cell(2, fill=fuel, region=+z_min & -z_max & -z_cyl)
root_universe = openmc.Universe(0, cells=(outside_cyl, inside_cyl))
geometry = openmc.Geometry()
geometry.root_universe = root_universe
geometry.export_to_xml()
# Define settings
settings = openmc.Settings()
settings.particles = 1000
settings.batches = 4
settings.inactive = 0
settings.source = openmc.Source(space=openmc.stats.Box(
*outside_cyl.region.bounding_box))
settings.export_to_xml()
if __name__ == '__main__':
harness = PeriodicTest('statepoint.4.h5')
harness.main()

2
tests/test_tallies/inputs_true.dat
View file

@ -1 +1 @@
0597eff3fddbc45a09b5b324c9704e540b694b07c136f2040426fdcfe5ec544f036073e4afa34a5fb0fbd721a4c0a609b9b68bf17ce4ec78302023b46b71930c
ea09926d8f5c6c96529bf5529f4deb3be78eda2da80adbbf3440147c337587358c2b1823bc72df9463676135573eb481dcd361b735f18365216645ee81092f1e

2
tests/test_tallies/results_true.dat
View file

@ -1 +1 @@
9f14aaa1694489032b3ce193ad29ecf6ac8976c88c2dd6b26d4c30ae88348e249a9b702b1d39c22204350b8f3bd689800c1b6a6003f19c7bdaf64084a209a2cc
a0c7d6ca246ecd7dd5fed06373af142390971401c4e97744f29e55810ab9c231c97c4d8947cdf0b3d2df0ae829a9ddf768e5b2d889bbea34f2b6db0e567db884

4
tests/test_tallies/test_tallies.py
View file

@ -160,9 +160,6 @@ class TalliesTestHarness(PyAPITestHarness):
total_tallies[2].estimator = 'analog'
total_tallies[3].estimator = 'collision'
questionable_tally = Tally()
questionable_tally.scores = ['transport', 'n1n']
all_nuclide_tallies = [Tally(), Tally()]
for t in all_nuclide_tallies:
t.filters = [cell_filter]
@ -182,7 +179,6 @@ class TalliesTestHarness(PyAPITestHarness):
self._input_set.tallies += flux_tallies
self._input_set.tallies += (scatter_tally1, scatter_tally2)
self._input_set.tallies += total_tallies
self._input_set.tallies.append(questionable_tally)
self._input_set.tallies += all_nuclide_tallies
self._input_set.export()

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

@ -0,0 +1 @@
2dcfd1a17cba671874e60192a7355deb57e2e51467a474fd168c8b51e454a977edb34df07ae11625c0a43906112152c75113e442a9a8f240a4c9d1a11ee4771d

6
tests/test_triso/plots.xml Normal file
View file

@ -0,0 +1,6 @@
<?xml version="1.0"?>
<plots>
<plot id="1" type="slice" basis="xy" color="material"
origin="0.0 0.0 0.0" width="1.0 1.0" pixels="400 400">
</plot>
</plots>

2
tests/test_triso/results_true.dat Normal file
View file

@ -0,0 +1,2 @@
k-combined:
1.662675E+00 1.475968E-02

117
tests/test_triso/test_triso.py Normal file
View file

@ -0,0 +1,117 @@
#!/usr/bin/env python
import os
import sys
import glob
import random
from math import sqrt
import numpy as np
sys.path.insert(0, os.pardir)
from testing_harness import PyAPITestHarness
import openmc
import openmc.model
class TRISOTestHarness(PyAPITestHarness):
def _build_inputs(self):
# Define TRISO matrials
fuel = openmc.Material()
fuel.set_density('g/cm3', 10.5)
fuel.add_nuclide('U-235', 0.14154)
fuel.add_nuclide('U-238', 0.85846)
fuel.add_nuclide('C-Nat', 0.5)
fuel.add_nuclide('O-16', 1.5)
porous_carbon = openmc.Material()
porous_carbon.set_density('g/cm3', 1.0)
porous_carbon.add_nuclide('C-Nat', 1.0)
porous_carbon.add_s_alpha_beta('Graph', '71t')
ipyc = openmc.Material()
ipyc.set_density('g/cm3', 1.90)
ipyc.add_nuclide('C-Nat', 1.0)
ipyc.add_s_alpha_beta('Graph', '71t')
sic = openmc.Material()
sic.set_density('g/cm3', 3.20)
sic.add_element('Si', 1.0)
sic.add_nuclide('C-Nat', 1.0)
opyc = openmc.Material()
opyc.set_density('g/cm3', 1.87)
opyc.add_nuclide('C-Nat', 1.0)
opyc.add_s_alpha_beta('Graph', '71t')
graphite = openmc.Material()
graphite.set_density('g/cm3', 1.1995)
graphite.add_nuclide('C-Nat', 1.0)
graphite.add_s_alpha_beta('Graph', '71t')
# Create TRISO particles
spheres = [openmc.Sphere(R=r*1e-4)
for r in [212.5, 312.5, 347.5, 382.5]]
c1 = openmc.Cell(fill=fuel, region=-spheres[0])
c2 = openmc.Cell(fill=porous_carbon, region=+spheres[0] & -spheres[1])
c3 = openmc.Cell(fill=ipyc, region=+spheres[1] & -spheres[2])
c4 = openmc.Cell(fill=sic, region=+spheres[2] & -spheres[3])
c5 = openmc.Cell(fill=opyc, region=+spheres[3])
inner_univ = openmc.Universe(cells=[c1, c2, c3, c4, c5])
outer_radius = 422.5*1e-4
trisos = []
random.seed(1)
for i in range(100):
# Randomly sample location
lim = 0.5 - outer_radius*1.001
x = random.uniform(-lim, lim)
y = random.uniform(-lim, lim)
z = random.uniform(-lim, lim)
t = openmc.model.TRISO(outer_radius, inner_univ, (x, y, z))
# Make sure TRISO doesn't overlap with another
for tp in trisos:
xp, yp, zp = tp.center
distance = sqrt((x - xp)**2 + (y - yp)**2 + (z - zp)**2)
if distance <= 2*outer_radius:
break
else:
trisos.append(t)
# Define box to contain lattice
min_x = openmc.XPlane(x0=-0.5, boundary_type='reflective')
max_x = openmc.XPlane(x0=0.5, boundary_type='reflective')
min_y = openmc.YPlane(y0=-0.5, boundary_type='reflective')
max_y = openmc.YPlane(y0=0.5, boundary_type='reflective')
min_z = openmc.ZPlane(z0=-0.5, boundary_type='reflective')
max_z = openmc.ZPlane(z0=0.5, boundary_type='reflective')
box = openmc.Cell(region=+min_x & -max_x & +min_y & -max_y & +min_z & -max_z)
# Create lattice
ll, ur = box.region.bounding_box
shape = (3, 3, 3)
pitch = (ur - ll) / shape
lattice = openmc.model.create_triso_lattice(
trisos, ll, pitch, shape, graphite)
box.fill = lattice
root = openmc.Universe(0, cells=[box])
geom = openmc.Geometry(root)
geom.export_to_xml()
settings = openmc.Settings()
settings.batches = 5
settings.inactive = 0
settings.particles = 100
settings.source = openmc.Source(space=openmc.stats.Point())
settings.export_to_xml()
mats = openmc.Materials([fuel, porous_carbon, ipyc, sic, opyc, graphite])
mats.default_xs = '71c'
mats.export_to_xml()
if __name__ == '__main__':
harness = TRISOTestHarness('statepoint.5.h5')
harness.main()

33
tests/testing_harness.py
View file

@ -6,7 +6,6 @@ import hashlib
from optparse import OptionParser
import os
import shutil
from subprocess import Popen, STDOUT, PIPE, call
import sys
import numpy as np
@ -18,6 +17,7 @@ import openmc
class TestHarness(object):
"""General class for running OpenMC regression tests."""
def __init__(self, statepoint_name, tallies_present=False):
self._sp_name = statepoint_name
self._tallies = tallies_present
@ -74,13 +74,13 @@ class TestHarness(object):
def _test_output_created(self):
"""Make sure statepoint.* and tallies.out have been created."""
statepoint = glob.glob(os.path.join(os.getcwd(), self._sp_name))
assert len(statepoint) == 1, 'Either multiple or no statepoint files ' \
'exist.'
assert len(statepoint) == 1, 'Either multiple or no statepoint files' \
' exist.'
assert statepoint[0].endswith('h5'), \
'Statepoint file is not a HDF5 file.'
'Statepoint file is not a HDF5 file.'
if self._tallies:
assert os.path.exists(os.path.join(os.getcwd(), 'tallies.out')), \
'Tally output file does not exist.'
'Tally output file does not exist.'
def _get_results(self, hash_output=False):
"""Digest info in the statepoint and return as a string."""
@ -98,7 +98,7 @@ class TestHarness(object):
tally_num = 1
for tally_ind in sp.tallies:
tally = sp.tallies[tally_ind]
results = np.zeros((tally.sum.size*2, ))
results = np.zeros((tally.sum.size * 2, ))
results[0::2] = tally.sum.ravel()
results[1::2] = tally.sum_sq.ravel()
results = ['{0:12.6E}'.format(x) for x in results]
@ -133,9 +133,10 @@ class TestHarness(object):
def _cleanup(self):
"""Delete statepoints, tally, and test files."""
output = glob.glob(os.path.join(os.getcwd(), 'statepoint.*.*'))
output = glob.glob(os.path.join(os.getcwd(), 'statepoint.*.h5'))
output.append(os.path.join(os.getcwd(), 'tallies.out'))
output.append(os.path.join(os.getcwd(), 'results_test.dat'))
output.append(os.path.join(os.getcwd(), 'summary.h5'))
for f in output:
if os.path.exists(f):
os.remove(f)
@ -143,6 +144,7 @@ class TestHarness(object):
class HashedTestHarness(TestHarness):
"""Specialized TestHarness that hashes the results."""
def _get_results(self):
"""Digest info in the statepoint and return as a string."""
return super(HashedTestHarness, self)._get_results(True)
@ -150,6 +152,7 @@ class HashedTestHarness(TestHarness):
class CMFDTestHarness(TestHarness):
"""Specialized TestHarness for running OpenMC CMFD tests."""
def _get_results(self):
"""Digest info in the statepoint and return as a string."""
# Read the statepoint file.
@ -183,6 +186,7 @@ class CMFDTestHarness(TestHarness):
class ParticleRestartTestHarness(TestHarness):
"""Specialized TestHarness for running OpenMC particle restart tests."""
def _run_openmc(self):
# Set arguments
args = {'openmc_exec': self._opts.exe}
@ -203,9 +207,9 @@ class ParticleRestartTestHarness(TestHarness):
"""Make sure the restart file has been created."""
particle = glob.glob(os.path.join(os.getcwd(), self._sp_name))
assert len(particle) == 1, 'Either multiple or no particle restart ' \
'files exist.'
'files exist.'
assert particle[0].endswith('h5'), \
'Particle restart file is not a HDF5 file.'
'Particle restart file is not a HDF5 file.'
def _get_results(self):
"""Digest info in the statepoint and return as a string."""
@ -228,10 +232,10 @@ class ParticleRestartTestHarness(TestHarness):
outstr += 'particle energy:\n'
outstr += "{0:12.6E}\n".format(p.energy)
outstr += 'particle xyz:\n'
outstr += "{0:12.6E} {1:12.6E} {2:12.6E}\n".format(p.xyz[0],p.xyz[1],
outstr += "{0:12.6E} {1:12.6E} {2:12.6E}\n".format(p.xyz[0], p.xyz[1],
p.xyz[2])
outstr += 'particle uvw:\n'
outstr += "{0:12.6E} {1:12.6E} {2:12.6E}\n".format(p.uvw[0],p.uvw[1],
outstr += "{0:12.6E} {1:12.6E} {2:12.6E}\n".format(p.uvw[0], p.uvw[1],
p.uvw[2])
return outstr
@ -239,13 +243,15 @@ class ParticleRestartTestHarness(TestHarness):
class PyAPITestHarness(TestHarness):
def __init__(self, statepoint_name, tallies_present=False, mg=False):
super(PyAPITestHarness, self).__init__(statepoint_name, tallies_present)
super(PyAPITestHarness, self).__init__(statepoint_name,
tallies_present)
self.parser.add_option('--build-inputs', dest='build_only',
action='store_true', default=False)
if mg:
self._input_set = MGInputSet()
else:
self._input_set = InputSet()
def main(self):
"""Accept commandline arguments and either run or update tests."""
(self._opts, self._args) = self.parser.parse_args()
@ -320,7 +326,8 @@ class PyAPITestHarness(TestHarness):
compare = filecmp.cmp('inputs_test.dat', 'inputs_true.dat')
if not compare:
f = open('inputs_test.dat')
for line in f.readlines(): print(line)
for line in f.readlines():
print(line)
f.close()
os.rename('inputs_test.dat', 'inputs_error.dat')
assert compare, 'Input files are broken.'