Merge remote-tracking branch 'upstream/develop' into simplified_tests

Conflicts:
	openmc/executor.py
This commit is contained in:
Sterling Harper 2015-06-29 17:35:07 -06:00
commit 6cb48baff1
42 changed files with 3069 additions and 1950 deletions

View file

@ -17,13 +17,18 @@ import sys, os
# add these directories to sys.path here. If the directory is relative to the
# documentation root, use os.path.abspath to make it absolute, like shown here.
sys.path.insert(0, os.path.abspath('../sphinxext'))
sys.path.insert(0, os.path.abspath('../..'))
# -- General configuration -----------------------------------------------------
# Add any Sphinx extension module names here, as strings. They can be extensions
# coming with Sphinx (named 'sphinx.ext.*') or your custom ones.
extensions = ['sphinx.ext.pngmath', 'sphinxcontrib.tikz', 'sphinx_numfig']
extensions = ['sphinx.ext.autodoc',
'sphinx.ext.napoleon',
'sphinx.ext.pngmath',
'sphinxcontrib.tikz',
'sphinx_numfig']
# Add any paths that contain templates here, relative to this directory.
templates_path = ['_templates']

View file

@ -18,9 +18,7 @@ from the git repository as such:
.. code-block:: sh
git clone https://bitbucket.org/philexander/tikz.git
cd tikz
sudo python setup.py install
sudo pip install https://bitbucket.org/philexander/tikz/get/HEAD.tar.gz
The Numfig_ package can be installed directly with pip:

View file

@ -32,6 +32,7 @@ free to send a message to the User's Group `mailing list`_.
methods/index
usersguide/index
devguide/index
pythonapi/index
publications
license
developers

View file

@ -0,0 +1,8 @@
.. _pythonapi_ace:
==========
ACE Format
==========
.. automodule:: openmc.ace
:members:

View file

@ -0,0 +1,8 @@
.. _pythonapi_cmfd:
====
CMFD
====
.. automodule:: openmc.cmfd
:members:

View file

@ -0,0 +1,8 @@
.. _pythonapi_element:
=======
Element
=======
.. automodule:: openmc.element
:members:

View file

@ -0,0 +1,8 @@
.. _pythonapi_executor:
========
Executor
========
.. automodule:: openmc.executor
:members:

View file

@ -0,0 +1,8 @@
.. _pythonapi_filter:
======
Filter
======
.. automodule:: openmc.filter
:members:

View file

@ -0,0 +1,8 @@
.. _pythonapi_geometry:
========
Geometry
========
.. automodule:: openmc.geometry
:members:

View file

@ -0,0 +1,32 @@
.. _pythonapi:
==========
Python API
==========
--------
Contents
--------
.. toctree::
:maxdepth: 1
ace
cmfd
element
executor
filter
geometry
material
mesh
nuclide
opencg_compatible
particle_restart
plots
settings
statepoint
summary
surface
tallies
trigger
universe

View file

@ -0,0 +1,8 @@
.. _pythonapi_material:
=========
Materials
=========
.. automodule:: openmc.material
:members:

View file

@ -0,0 +1,8 @@
.. _pythonapi_mesh:
====
Mesh
====
.. automodule:: openmc.mesh
:members:

View file

@ -0,0 +1,8 @@
.. _pythonapi_nuclide:
=======
Nuclide
=======
.. automodule:: openmc.nuclide
:members:

View file

@ -0,0 +1,8 @@
.. _pythonapi_opencg_compatible:
====================
OpenCG Compatibility
====================
.. automodule:: openmc.opencg_compatible
:members:

View file

@ -0,0 +1,8 @@
.. _pythonapi_particle_restart:
================
Particle Restart
================
.. automodule:: openmc.particle_restart
:members:

View file

@ -0,0 +1,8 @@
.. _pythonapi_plots:
=====
Plots
=====
.. automodule:: openmc.plots
:members:

View file

@ -0,0 +1,8 @@
.. _pythonapi_settings:
========
Settings
========
.. automodule:: openmc.settings
:members:

View file

@ -0,0 +1,8 @@
.. _pythonapi_statepoint:
==========
Statepoint
==========
.. automodule:: openmc.statepoint
:members:

View file

@ -0,0 +1,8 @@
.. _pythonapi_summary:
=======
Summary
=======
.. automodule:: openmc.summary
:members:

View file

@ -0,0 +1,8 @@
.. _pythonapi_surface:
=======
Surface
=======
.. automodule:: openmc.surface
:members:

View file

@ -0,0 +1,8 @@
.. _pythonapi_tallies:
=======
Tallies
=======
.. automodule:: openmc.tallies
:members:

View file

@ -0,0 +1,8 @@
.. _pythonapi_trigger:
=======
Trigger
=======
.. automodule:: openmc.trigger
:members:

View file

@ -0,0 +1,8 @@
.. _pythonapi_universe:
========
Universe
========
.. automodule:: openmc.universe
:members:

View file

@ -441,24 +441,22 @@ attributes/sub-elements:
has the following attributes:
:type:
The type of spatial distribution. Valid options are "box" and "point". A
"box" spatial distribution has coordinates sampled uniformly in a
parallelepiped. A "point" spatial distribution has coordinates specified
by a triplet.
The type of spatial distribution. Valid options are "box", "fission", and
"point". A "box" spatial distribution has coordinates sampled uniformly in
a parallelepiped. A "fission" spatial distribution samples locations from
a "box" distribution but only locations in fissionable materials are
accepted. A "point" spatial distribution has coordinates specified by a
triplet.
*Default*: None
:parameters:
For a "box" spatial distribution, ``parameters`` should be given as six
real numbers, the first three of which specify the lower-left corner of a
parallelepiped and the last three of which specify the upper-right
corner. Source sites are sampled uniformly through that parallelepiped.
To filter a "box" spatial distribution by fissionable material, specify
"fission" tag instead of "box". The ``parameters`` should be given as six
real numbers, the first three of which specify the lower-left corner of a
parallelepiped and the last three of which specify the upper-right
corner. Source sites are sampled uniformly through that parallelepiped.
For a "box" or "fission" spatial distribution, ``parameters`` should be
given as six real numbers, the first three of which specify the lower-left
corner of a parallelepiped and the last three of which specify the
upper-right corner. Source sites are sampled uniformly through that
parallelepiped.
For a "point" spatial distribution, ``parameters`` should be given as
three real numbers which specify the (x,y,z) location of an isotropic
@ -481,7 +479,7 @@ attributes/sub-elements:
:parameters:
For an "isotropic" angular distribution, ``parameters`` should not be
specified
specified.
For a "monodirectional" angular distribution, ``parameters`` should be
given as three real numbers which specify the angular cosines with respect
@ -499,7 +497,7 @@ attributes/sub-elements:
"watt", and "maxwell". The "monoenergetic" option produces source sites at
a single energy. The "watt" option produces source sites whose energy is
sampled from a Watt fission spectrum. The "maxwell" option produce source
sites whose energy is sampled from a Maxwell fission spectrum
sites whose energy is sampled from a Maxwell fission spectrum.
*Default*: watt
@ -605,8 +603,6 @@ survival biasing, otherwise known as implicit capture or absorption.
*Default*: false
.. _trace:
``<threads>`` Element
---------------------
@ -615,6 +611,8 @@ a simulation. It has no attributes and accepts a positive integer value.
*Default*: None (Determined by environment variable :envvar:`OMP_NUM_THREADS`)
.. _trace:
``<trace>`` Element
-------------------
@ -631,9 +629,9 @@ integers: the batch number, generation number, and particle number.
The ``<track>`` element specifies particles for which OpenMC will output binary
files describing particle position at every step of its transport. This element
should be followed by triplets of integers. Each triplet describes one particle
. The integers in each triplet specify the batch number, generation number, and
particle number, respectively.
should be followed by triplets of integers. Each triplet describes one
particle. The integers in each triplet specify the batch number, generation
number, and particle number, respectively.
*Default*: None
@ -1215,8 +1213,8 @@ The ``<tally>`` element accepts the following sub-elements:
The ``filter`` element has the following attributes/sub-elements:
:type:
The type of the filter. Accepted options are "cell", "cellborn",
"material", "universe", "energy", "energyout", "mesh", and
The type of the filter. Accepted options are "cell", "cellborn",
"material", "universe", "energy", "energyout", "mesh", and
"distribcell".
:bins:
@ -1666,6 +1664,15 @@ The ``<begin>`` element controls what batch CMFD calculations should begin.
*Default*: 1
``<dhat_reset>`` Element
------------------------
The ``<dhat_reset>`` element controls whether :math:`\widehat{D}` nonlinear
CMFD parameters should be reset to zero before solving CMFD eigenproblem.
It can be turned on with "true" and off with "false".
*Default*: false
``<display>`` Element
---------------------
@ -1682,15 +1689,6 @@ The ``<display>`` element sets one additional CMFD output column. Options are:
*Default*: balance
``<dhat_reset>`` Element
------------------------
The ``<dhat_reset>`` element controls whether :math:`\widehat{D}` nonlinear
CMFD parameters should be reset to zero before solving CMFD eigenproblem.
It can be turned on with "true" and off with "false".
*Default*: false
``<downscatter>`` Element
-------------------------
@ -1734,11 +1732,11 @@ The CMFD mesh is a structured Cartesian mesh. This element has the following
attributes/sub-elements:
:lower_left:
The lower-left corner of the structured mesh. If only two coordinate are
The lower-left corner of the structured mesh. If only two coordinates are
given, it is assumed that the mesh is an x-y mesh.
:upper_right:
The upper-right corner of the structrued mesh. If only two coordinate are
The upper-right corner of the structrued mesh. If only two coordinates are
given, it is assumed that the mesh is an x-y mesh.
:dimension:
@ -1761,7 +1759,7 @@ attributes/sub-elements:
:map:
An optional acceleration map can be specified to overlay on the coarse
mesh spatial grid. If this option is used a ``1`` is used for a
mesh spatial grid. If this option is used, a ``1`` is used for a
non-accelerated region and a ``2`` is used for an accelerated region.
For a simple 4x4 coarse mesh with a 2x2 fuel lattice surrounded by
reflector, the map is:

View file

@ -1,3 +1,15 @@
"""This module can be used to specify parameters used for coarse mesh finite
difference (CMFD) acceleration in OpenMC. CMFD was first proposed by [Smith]_
and is widely used in accelerating neutron transport problems.
References
----------
.. [Smith] K. Smith, "Nodal method storage reduction by non-linear
iteration", *Trans. Am. Nucl. Soc.*, **44**, 265 (1983).
"""
from xml.etree import ElementTree as ET
import numpy as np
@ -7,9 +19,52 @@ from openmc.clean_xml import *
class CMFDMesh(object):
"""A structured Cartesian mesh used for Coarse Mesh Finite Difference (CMFD)
acceleration.
Attributes
----------
lower_left : tuple or list or ndarray
The lower-left corner of the structured mesh. If only two coordinates are
given, it is assumed that the mesh is an x-y mesh.
upper_right : tuple or list or ndarray
The upper-right corner of the structrued mesh. If only two coordinates
are given, it is assumed that the mesh is an x-y mesh.
dimension : tuple or list or ndarray
The number of mesh cells in each direction.
width : tuple or list or ndarray
The width of mesh cells in each direction.
energy : tuple or list or ndarray
Energy bins in MeV, listed in ascending order (e.g. [0.0, 0.625e-7,
20.0]) for CMFD tallies and acceleration. If no energy bins are listed,
OpenMC automatically assumes a one energy group calculation over the
entire energy range.
albedo : tuple or list or ndarray
Surface ratio of incoming to outgoing partial currents on global
boundary conditions. They are listed in the following order: -x +x -y +y
-z +z.
map : tuple or list or ndarray
An optional acceleration map can be specified to overlay on the coarse
mesh spatial grid. If this option is used, a ``1`` is used for a
non-accelerated region and a ``2`` is used for an accelerated region.
For a simple 4x4 coarse mesh with a 2x2 fuel lattice surrounded by
reflector, the map is:
::
[1, 1, 1, 1,
1, 2, 2, 1,
1, 2, 2, 1,
1, 1, 1, 1]
Therefore a 2x2 system of equations is solved rather than a 4x4. This is
extremely important to use in reflectors as neutrons will not contribute
to any tallies far away from fission source neutron regions. A ``2``
must be used to identify any fission source region.
"""
def __init__(self):
self._lower_left = None
self._upper_right = None
self._dimension = None
@ -18,45 +73,36 @@ class CMFDMesh(object):
self._albedo = None
self._map = None
@property
def lower_left(self):
return self._lower_left
@property
def upper_right(self):
return self._upper_right
@property
def dimension(self):
return self._dimension
@property
def width(self):
return self._width
@property
def energy(self):
return self._energy
@property
def albedo(self):
return self._albedo
@property
def map(self):
return self._mape
return self._map
@lower_left.setter
def lower_left(self, lower_left):
if not isinstance(lower_left, (tuple, list, np.ndarray)):
msg = 'Unable to set CMFD Mesh with lower_left {0} which is ' \
'not a Python list, tuple or NumPy array'.format(lower_left)
@ -64,11 +110,10 @@ class CMFDMesh(object):
elif len(lower_left) != 2 and len(lower_left) != 3:
msg = 'Unable to set CMFD Mesh with lower_left {0} since it ' \
'must include 2 or 3 dimensions'.format(lower_left)
'must include 2 or 3 dimensions'.format(lower_left)
raise ValueError(msg)
for coord in lower_left:
if not is_integer(coord) and not is_float(coord):
msg = 'Unable to set CMFD Mesh with lower_left {0} which is ' \
'not an integer or a floating point value'.format(coord)
@ -76,7 +121,6 @@ class CMFDMesh(object):
self._lower_left = lower_left
@upper_right.setter
def upper_right(self, upper_right):
@ -91,7 +135,6 @@ class CMFDMesh(object):
raise ValueError(msg)
for coord in upper_right:
if not is_integer(coord) and not is_float(coord):
msg = 'Unable to set CMFD Mesh with upper_right {0} which ' \
'is not an integer or floating point value'.format(coord)
@ -99,10 +142,8 @@ class CMFDMesh(object):
self._upper_right = upper_right
@dimension.setter
def dimension(self, dimension):
if not isinstance(dimension, (tuple, list, np.ndarray)):
msg = 'Unable to set CMFD Mesh with dimension {0} which is ' \
'not a Python list, tuple or NumPy array'.format(dimension)
@ -114,7 +155,6 @@ class CMFDMesh(object):
raise ValueError(msg)
for dim in dimension:
if not is_integer(dim):
msg = 'Unable to set CMFD Mesh with dimension {0} which ' \
'is a non-integer'.format(dim)
@ -122,12 +162,9 @@ class CMFDMesh(object):
self._dimension = dimension
@width.setter
def width(self, width):
if not width is None:
if width is not None:
if not isinstance(width, (tuple, list, np.ndarray)):
msg = 'Unable to set CMFD Mesh with width {0} which ' \
'is not a Python list, tuple or NumPy array'.format(width)
@ -139,7 +176,6 @@ class CMFDMesh(object):
raise ValueError(msg)
for dim in width:
if not is_integer(dim) and not is_float(dim):
msg = 'Unable to set CMFD Mesh with width {0} which is ' \
'not an integer or floating point value'.format(width)
@ -147,22 +183,18 @@ class CMFDMesh(object):
self._width = width
@energy.setter
def energy(self, energy):
if not isinstance(energy, (tuple, list, np.ndarray)):
msg = 'Unable to set CMFD Mesh energy to {0} which is not ' \
'a Python tuple/list or NumPy array'.format(energy)
raise ValueError(msg)
for e in energy:
if not is_integer(e) and not is_float(e):
msg = 'Unable to set CMFD Mesh energy to {0} which is not ' \
'an integer or floating point value'.format(e)
raise ValueError(msg)
elif e < 0:
msg = 'Unable to set CMFD Mesh energy to {0} which is ' \
'is a negative integer'.format(e)
@ -170,10 +202,8 @@ class CMFDMesh(object):
self._energy = energy
@albedo.setter
def albedo(self, albedo):
if not isinstance(albedo, (tuple, list, np.ndarray)):
msg = 'Unable to set CMFD Mesh albedo to {0} which is not ' \
'a Python tuple/list or NumPy array'.format(albedo)
@ -185,12 +215,10 @@ class CMFDMesh(object):
raise ValueError(msg)
for a in albedo:
if not is_integer(a) and not is_float(a):
msg = 'Unable to set CMFD Mesh albedo to {0} which is not ' \
'an integer or floating point value'.format(a)
raise ValueError(msg)
elif a < 0 or a > 1:
msg = 'Unable to set CMFD Mesh albedo to {0} which is ' \
'is not in [0,1]'.format(a)
@ -198,7 +226,6 @@ class CMFDMesh(object):
self._albedo = albedo
@map.setter
def map(self, map):
@ -208,7 +235,6 @@ class CMFDMesh(object):
raise ValueError(msg)
for m in map:
if m != 1 and m != 2:
msg = 'Unable to set CMFD Mesh map to {0} which is ' \
'is not 1 or 2'.format(m)
@ -216,191 +242,181 @@ class CMFDMesh(object):
self._map = map
def get_mesh_xml(self):
def _get_xml_element(self):
element = ET.Element("mesh")
if len(self._lower_left) == 2:
subelement = ET.SubElement(element, "lower_left")
subelement.text = '{0} {1}'.format(self._lower_left[0],
self._lower_left[1])
else:
subelement = ET.SubElement(element, "lower_left")
subelement.text = '{0} {1} {2}'.format(self._lower_left[0],
self._lower_left[1],
self._lower_left[2])
subelement = ET.SubElement(element, "lower_left")
subelement.text = ' '.join(map(str, self._lower_left))
if not self._upper_right is None:
if len(self._upper_right) == 2:
subelement = ET.SubElement(element, "upper_right")
subelement.text = '{0} {1}'.format(self._upper_right[0],
self._upper_right[1])
else:
subelement = ET.SubElement(element, "upper_right")
subelement.text = '{0} {1} {2}'.format(self._upper_right[0],
self._upper_right[1],
self._upper_right[2])
if self.upper_right is not None:
subelement = ET.SubElement(element, "upper_right")
subelement.text = ' '.join(map(str, self.upper_right))
if len(self._dimension) == 2:
subelement = ET.SubElement(element, "dimension")
subelement.text = '{0} {1}'.format(self._dimension[0],
self._dimension[1])
else:
subelement = ET.SubElement(element, "dimension")
subelement.text = '{0} {1} {2}'.format(self._dimension[0],
self._dimension[1],
self._dimension[2])
subelement = ET.SubElement(element, "dimension")
subelement.text = ' '.join(map(str, self.dimension))
if not self._width is None:
if len(self._width) == 2:
subelement = ET.SubElement(element, "width")
subelement.text = '{0} {1}'.format(self._width[0],
self._width[1])
else:
subelement = ET.SubElement(element, "width")
subelement.text = '{0} {1} {2}'.format(self._width[0],
self._width[1],
self._width[2])
if not self._energy is None:
if self.width is not None:
subelement = ET.SubElement(element, "width")
subelement.text = ' '.join(map(str, self.width))
if self.energy is not None:
subelement = ET.SubElement(element, "energy")
subelement.text = ' '.join(map(str, self.energy))
energy = ''
for e in self._energy:
energy += '{0} '.format(e)
subelement.set("energy", energy.rstrip(' '))
if not self._albedo is None:
if self.albedo is not None:
subelement = ET.SubElement(element, "albedo")
subelement.text = ' '.join(map(str, self.albedo))
albedo = ''
for a in self._albedo:
albedo += '{0} '.format(a)
subelement.set("albedo", albedo.rstrip(' '))
if not self._map is None:
if self.map is not None:
subelement = ET.SubElement(element, "map")
map = ''
for m in self._map:
map += '{0} '.format(m)
subelement.set("map", map.rstrip(' '))
subelement.text = ' '.join(map(str, self.map))
return element
class CMFDFile(object):
"""Parameters that control the use of coarse-mesh finite difference acceleration
in OpenMC. This corresponds directly to the cmfd.xml input file.
Attributes
----------
begin : int
Batch number at which CMFD calculations should begin
dhat_reset : bool
Indicate whether :math:`\widehat{D}` nonlinear CMFD parameters should be
reset to zero before solving CMFD eigenproblem.
display : {'balance', 'dominance', 'entropy', 'source'}
Set one additional CMFD output column. Options are:
* "balance" - prints the RMS [%] of the resdiual from the neutron balance
equation on CMFD tallies.
* "dominance" - prints the estimated dominance ratio from the CMFD
iterations.
* "entropy" - prints the *entropy* of the CMFD predicted fission source.
* "source" - prints the RMS [%] between the OpenMC fission source and
CMFD fission source.
downscatter : bool
Indicate whether an effective downscatter cross section should be used
when using 2-group CMFD.
feedback : bool
Indicate or not the CMFD diffusion result is used to adjust the weight
of fission source neutrons on the next OpenMC batch. Defaults to False.
gauss_seidel_tolerance : tuple or list or ndarray of float
Two parameters specifying the absolute inner tolerance and the relative
inner tolerance for Gauss-Seidel iterations when performing CMFD.
ktol : float
Tolerance on the eigenvalue when performing CMFD power iteration
cmfd_mesh : CMFDMesh
Structured mesh to be used for acceleration
norm : float
Normalization factor applied to the CMFD fission source distribution
power_monitor : bool
View convergence of power iteration during CMFD acceleration
run_adjoint : bool
Perform adjoint calculation on the last batch
shift : float
Optional Wielandt shift parameter for accelerating power iterations. By
default, it is very large so there is effectively no impact.
spectral : float
Optional spectral radius that can be used to accelerate the convergence
of Gauss-Seidel iterations during CMFD power iteration.
stol : float
Tolerance on the fission source when performing CMFD power iteration
tally_reset : list of int
List of batch numbers at which CMFD tallies should be reset
write_matrices : bool
Write sparse matrices that are used during CMFD acceleration (loss,
production) to file
"""
def __init__(self):
self._active_flush = None
self._begin = None
self._dhat_reset = None
self._display = None
self._downscatter = None
self._feedback = None
self._inactive = None
self._inactive_flush = None
self._gauss_seidel_tolerance = None
self._ktol = None
self._cmfd_mesh = None
self._norm = None
self._num_flushes = None
self._power_monitor = None
self._run_adjoint = None
self._shift = None
self._spectral = None
self._stol = None
self._tally_reset = None
self._write_matrices = None
self._cmfd_file = ET.Element("cmfd")
self._cmfd_mesh_element = None
@property
def active_flush(self):
return self._active_flush
@property
def begin(self):
return self._begin
@property
def dhat_reset(self):
return self._dhat_reset
@property
def display(self):
return self._display
@property
def downscatter(self):
return self._downscatter
@property
def feedback(self):
return self._feedback
@property
def gauss_seidel_tolerance(self):
return self._gauss_seidel_tolerance
@property
def inactive(self):
return self._inactive
@property
def inactive_flush(self):
return self._inactive_flush
def ktol(self):
return self._ktol
@property
def cmfd_mesh(self):
return self._cmfd_mesh
@property
def norm(self):
return self._norm
@property
def num_flushes(self):
return self._num_flushes
@property
def power_monitor(self):
return self._power_monitor
@property
def run_adjoint(self):
return self._run_adjoint
@property
def shift(self):
return self._shift
@property
def solver(self):
return self._solver
def spectral(self):
return self._spectral
@property
def stol(self):
return self._stol
@property
def tally_reset(self):
return self._tally_reset
@property
def write_matrices(self):
return self._write_matrices
@active_flush.setter
def active_flush(self, active_flush):
if not is_integer(active_flush):
msg = 'Unable to set CMFD active flush batch to a non-integer ' \
'value {0}'.format(active_flush)
raise ValueError(msg)
if active_flush < 0:
msg = 'Unable to set CMFD active flush batch to a negative ' \
'value {0}'.format(active_flush)
raise ValueError(msg)
self._active_flush = active_flush
@begin.setter
def begin(self, begin):
if not is_integer(begin):
msg = 'Unable to set CMFD begin batch to a non-integer ' \
'value {0}'.format(begin)
@ -413,10 +429,17 @@ class CMFDFile(object):
self._begin = begin
@dhat_reset.setter
def dhat_reset(self, dhat_reset):
if not isinstance(dhat_reset, bool):
msg = 'Unable to set Dhat reset to {0} which is ' \
'a non-boolean value'.format(dhat_reset)
raise ValueError(msg)
self._dhat_reset = dhat_reset
@display.setter
def display(self, display):
if not is_string(display):
msg = 'Unable to set CMFD display to a non-string ' \
'value'.format(display)
@ -429,10 +452,17 @@ class CMFDFile(object):
self._display = display
@downscatter.setter
def downscatter(self, downscatter):
if not isinstance(downscatter, bool):
msg = 'Unable to set downscatter to {0} which is ' \
'a non-boolean value'.format(downscatter)
raise ValueError(msg)
self._downscatter = downscatter
@feedback.setter
def feedback(self, feedback):
if not isinstance(feedback, bool):
msg = 'Unable to set CMFD feedback to {0} which is ' \
'a non-boolean value'.format(feedback)
@ -440,37 +470,38 @@ class CMFDFile(object):
self._feedback = feedback
@inactive.setter
def inactive(self, inactive):
if not isinstance(inactive, bool):
msg = 'Unable to set CMFD inactive batch to {0} which is ' \
' a non-boolean value'.format(inactive)
@gauss_seidel_tolerance.setter
def gauss_seidel_tolerance(self, gauss_seidel_tolerance):
if not isinstance(gauss_seidel_tolerance, (float, list, np.ndarray)):
msg = 'Unable to set Gauss-Seidel tolerance to {0} which is ' \
'not a Python tuple/list or NumPy array'.format(
gauss_seidel_tolerance)
raise ValueError(msg)
self._inactive = inactive
@inactive_flush.setter
def inactive_flush(self, inactive_flush):
if not is_integer(inactive_flush):
msg = 'Unable to set CMFD inactive flush batch to {0} which is ' \
'a non-integer value'.format(inactive_flush)
if len(gauss_seidel_tolerance) != 2:
msg = 'Unable to set Gauss-Seidel tolerance with {0} since ' \
'it must be of length 2'.format(width)
raise ValueError(msg)
if inactive_flush <= 0:
msg = 'Unable to set CMFD inactive flush batch to {0} which is ' \
'a negative value {0}'.format(inactive_flush)
for t in gauss_seidel_tolerance:
if not is_integer(t) and not is_float(t):
msg = 'Unable to set Gauss-Seidel tolerance with {0} which ' \
'is not an integer or floating point value'.format(t)
raise ValueError(msg)
self._gauss_seidel_tolerance = gauss_seidel_tolerance
@ktol.setter
def ktol(self, ktol):
if not is_integer(ktol) and not is_float(ktol):
msg = 'Unable to set the eigenvalue tolerance to {0} which is ' \
'not an integer or floating point value'.format(ktol)
raise ValueError(msg)
self._inactive_flush = inactive_flush
self._ktol = ktol
@cmfd_mesh.setter
def cmfd_mesh(self, mesh):
if not isinstance(mesh, CMFDMesh):
msg = 'Unable to set CMFD mesh to {0} which is not a ' \
'CMFDMesh object'.format(mesh)
@ -478,10 +509,8 @@ class CMFDFile(object):
self._mesh = mesh
@norm.setter
def norm(self, norm):
if not is_integer(norm) and not is_float(norm):
msg = 'Unable to set the CMFD norm to {0} which is not ' \
'an integer or floating point value'.format(norm)
@ -489,26 +518,8 @@ class CMFDFile(object):
self._norm = norm
@num_flushes.setter
def num_flushes(self, num_flushes):
if not is_integer(num_flushes):
msg = 'Unable to set the CMFD number of flushes to {0} ' \
'which is not an integer value'.format(num_flushes)
raise ValueError(msg)
if num_flushes < 0:
msg = 'Unable to set CMFD number of flushes to a negative ' \
'value {0}'.format(num_flushes)
raise ValueError(msg)
self._num_flushes = num_flushes
@power_monitor.setter
def power_monitor(self, power_monitor):
if not isinstance(power_monitor, bool):
msg = 'Unable to set CMFD power monitor to {0} which is a ' \
'non-boolean value'.format(power_monitor)
@ -516,10 +527,8 @@ class CMFDFile(object):
self._power_monitor = power_monitor
@run_adjoint.setter
def run_adjoint(self, run_adjoint):
if not isinstance(run_adjoint, bool):
msg = 'Unable to set CMFD run adjoint to {0} which is a ' \
'non-boolean value'.format(run_adjoint)
@ -527,10 +536,51 @@ class CMFDFile(object):
self._run_adjoint = run_adjoint
@shift.setter
def shift(self, shift):
if not is_integer(shift) and not is_float(shift):
msg = 'Unable to set the Wielandt shift to {0} which is ' \
'not an integer or floating point value'.format(shift)
raise ValueError(msg)
self._shift = shift
@spectral.setter
def spectral(self, spectral):
if not is_integer(spectral) and not is_float(spectral):
msg = 'Unable to set the spectral radius to {0} which is ' \
'not an integer or floating point value'.format(spectral)
raise ValueError(msg)
self._spectral = spectral
@stol.setter
def stol(self, stol):
if not is_integer(stol) and not is_float(stol):
msg = 'Unable to set the fission source tolerance to {0} which ' \
'is not an integer or floating point value'.format(stol)
raise ValueError(msg)
self._stol = stol
@tally_reset.setter
def tally_reset(self, tally_reset):
if not isinstance(tally_reset, (tuple, list, np.ndarray)):
msg = 'Unable to set tally reset batches to {0} which is ' \
'not a Python tuple/list or NumPy array'.format(
tally_reset)
raise ValueError(msg)
for t in tally_reset:
if not is_integer(t):
msg = 'Unable to set tally reset batch to {0} which ' \
'is not an integer'.format(t)
raise ValueError(msg)
self._tally_reset = tally_reset
@write_matrices.setter
def write_matrices(self, write_matrices):
if not isinstance(write_matrices, bool):
msg = 'Unable to set CMFD write matrices to {0} which is a ' \
'non-boolean value'.format(write_matrices)
@ -538,105 +588,108 @@ class CMFDFile(object):
self._write_matrices = write_matrices
def create_active_flush_subelement(self):
if not self._active_flush is None:
element = ET.SubElement(self._cmfd_file, "active_flush")
element.text = '{0}'.format(str(self._active_flush))
def create_begin_subelement(self):
if not self._begin is None:
def _create_begin_subelement(self):
if self._begin is not None:
element = ET.SubElement(self._cmfd_file, "begin")
element.text = '{0}'.format(str(self._begin))
element.text = str(self._begin)
def _create_dhat_reset_subelement(self):
if self._dhat_reset is not None:
element = ET.SubElement(self._cmfd_file, "dhat_reset")
element.text = str(self._dhat_reset).lower()
def create_display_subelement(self):
if not self._display is None:
def _create_display_subelement(self):
if self._display is not None:
element = ET.SubElement(self._cmfd_file, "display")
element.text = '{0}'.format(str(self._display))
element.text = str(self._display)
def _create_downscatter_subelement(self):
if self._downscatter is not None:
element = ET.SubElement(self._cmfd_file, "downscatter")
element.text = str(self._downscatter).lower()
def create_feedback_subelement(self):
if not self._feedback is None:
def _create_feedback_subelement(self):
if self._feedback is not None:
element = ET.SubElement(self._cmfd_file, "feeback")
element.text = '{0}'.format(str(self._feedback).lower())
element.text = str(self._feedback).lower()
def _create_gauss_seidel_tolerance_subelement(self):
if self._gauss_seidel_tolerance is not None:
element = ET.SubElement(self._cmfd_file, "gauss_seidel_tolerance")
element.text = ' '.join(map(str, self._gauss_seidel_tolerance))
def create_inactive_subelement(self):
def _create_ktol_subelement(self):
if self._ktol is not None:
element = ET.SubElement(self._ktol, "ktol")
element.text = str(self._ktol)
if not self._inactive is None:
element = ET.SubElement(self._cmfd_file, "inactive")
element.text = '{0}'.format(str(self._inactive).lower())
def create_inactive_flush_subelement(self):
if not self._inactive_flush is None:
element = ET.SubElement(self._cmfd_file, "inactive_flush")
element.text = '{0}'.format(str(self._inactive_flush))
def create_mesh_subelement(self):
if not self._mesh is None:
xml_element = self._mesh.get_mesh_xml()
def _create_mesh_subelement(self):
if self._mesh is not None:
xml_element = self._mesh._get_xml_element()
self._cmfd_file.append(xml_element)
def create_norm_subelement(self):
if not self._num_flushes is None:
def _create_norm_subelement(self):
if self._norm is not None:
element = ET.SubElement(self._cmfd_file, "norm")
element.text = '{0}'.format(str(self._norm))
element.text = str(self._norm)
def create_num_flushes_subelement(self):
if not self._num_flushes is None:
element = ET.SubElement(self._cmfd_file, "num_flushes")
element.text = '{0}'.format(str(self._num_flushes))
def create_power_monitor_subelement(self):
if not self._power_monitor is None:
def _create_power_monitor_subelement(self):
if self._power_monitor is not None:
element = ET.SubElement(self._cmfd_file, "power_monitor")
element.text = '{0}'.format(str(self._power_monitor).lower())
element.text = str(self._power_monitor).lower()
def create_run_adjoint_subelement(self):
if not self._run_adjoint is None:
def _create_run_adjoint_subelement(self):
if self._run_adjoint is not None:
element = ET.SubElement(self._cmfd_file, "run_adjoint")
element.text = '{0}'.format(str(self._run_adjoint).lower())
element.text = str(self._run_adjoint).lower()
def _create_shift_subelement(self):
if self._shift is not None:
element = ET.SubElement(self._shift, "shift")
element.text = str(self._shift)
def create_write_matrices_subelement(self):
def _create_spectral_subelement(self):
if self._spectral is not None:
element = ET.SubElement(self._spectral, "spectral")
element.text = str(self._spectral)
if not self._write_matrices is None:
def _create_stol_subelement(self):
if self._stol is not None:
element = ET.SubElement(self._stol, "stol")
element.text = str(self._stol)
def _create_tally_reset_subelement(self):
if self._tally_reset is not None:
element = ET.SubElement(self._tally_reset, "tally_reset")
element.text = ' '.join(map(str, self._tally_reset))
def _create_write_matrices_subelement(self):
if self._write_matrices is not None:
element = ET.SubElement(self._cmfd_file, "write_matrices")
element.text = '{0}'.format(str(self._write_matrices).lower())
element.text = str(self._write_matrices).lower()
def export_to_xml(self):
"""Create a cmfd.xml file using the class data that can be used for an OpenMC
simulation.
self.create_active_flush_subelement()
self.create_begin_subelement()
self.create_display_subelement()
self.create_feedback_subelement()
self.create_inactive_subelement()
self.create_inactive_flush_subelement()
self.create_mesh_subelement()
self.create_norm_subelement()
self.create_num_flushes_subelement()
self.create_power_monitor_subelement()
self.create_run_adjoint_subelement()
self.create_write_matrices_subelement()
"""
self._create_begin_subelement()
self._create_dhat_reset_subelement()
self._create_display_subelement()
self._create_downscatter_subelement()
self._create_feedback_subelement()
self._create_gauss_seidel_tolerance_subelement()
self._create_ktol_subelement()
self._create_mesh_subelement()
self._create_norm_subelement()
self._create_power_monitor_subelement()
self._create_run_adjoint_subelement()
self._create_shift_subelement()
self._create_spectral_subelement()
self._create_stol_subelement()
self._create_tally_reset_subelement()
self._create_write_matrices_subelement()
# Clean the indentation in the file to be user-readable
clean_xml_indentation(self._cmfd_file)

View file

@ -1,58 +1,64 @@
from openmc.checkvalue import *
class Element(object):
"""A natural element used in a material via <element>. Internally, OpenMC will
expand the natural element into isotopes based on the known natural
abundances.
Parameters
----------
name : str
Chemical symbol of the element, e.g. Pu
xs : str
Cross section identifier, e.g. 71c
Attributes
----------
name : str
Chemical symbol of the element, e.g. Pu
xs : str
Cross section identifier, e.g. 71c
"""
def __init__(self, name='', xs=None):
# Initialize class attributes
self._name = ''
self._xs = None
# Set the Material class attributes
# Set class attributes
self.name = name
if not xs is None:
if xs is not None:
self.xs = xs
def __eq__(self, element2):
# Check type
if not isinstance(element2, Element):
return False
# Check name
# Check name and xs
if self._name != element2._name:
return False
# Check xs
elif self._xs != element2._xs:
return False
else:
return True
def __hash__(self):
hashable = []
hashable.append(self._name)
hashable.append(self._xs)
return hash(tuple(hashable))
return hash((self._name, self._xs))
@property
def xs(self):
return self._xs
@property
def name(self):
return self._name
@xs.setter
def xs(self, xs):
if not is_string(xs):
msg = 'Unable to set cross-section identifier xs for Element ' \
'with a non-string value {0}'.format(xs)
@ -60,10 +66,8 @@ class Element(object):
self._xs = xs
@name.setter
def name(self, name):
if not is_string(name):
msg = 'Unable to set name for Element with a non-string ' \
'value {0}'.format(name)
@ -71,9 +75,7 @@ class Element(object):
self._name = name
def __repr__(self):
string = 'Element - {0}\n'.format(self._name)
string += '{0: <16}{1}{2}\n'.format('\tXS', '=\t', self._xs)
return string

View file

@ -1,3 +1,4 @@
from __future__ import print_function
import subprocess
import os
@ -5,23 +6,28 @@ from openmc.checkvalue import *
class Executor(object):
"""Control execution of OpenMC
Attributes
----------
working_directory : str
Path to working directory to run in
"""
def __init__(self):
self._working_directory = '.'
def _run_openmc(self, command, output):
# Launch a subprocess to run OpenMC
p = subprocess.Popen(command, shell=True,
p = subprocess.Popen(command, shell=True,
cwd=self._working_directory,
stdout=subprocess.PIPE)
# Capture and re-print OpenMC output in real-time
while (True and output):
line = p.stdout.readline()
print(line),
print(line, end='')
# If OpenMC is finished, break loop
if line == '' and p.poll() != None:
@ -30,15 +36,12 @@ class Executor(object):
# Return the returncode (integer, zero if no problems encountered)
return p.returncode
@property
def working_directory(self):
return self._working_directory
@working_directory.setter
def working_directory(self, working_directory):
if not is_string(working_directory):
msg = 'Unable to set Executor\'s working directory to {0} ' \
'since it is not a string'.format(working_directory)
@ -51,14 +54,34 @@ class Executor(object):
self._working_directory = working_directory
def plot_geometry(self, output=True):
return self._run_openmc('openmc -p', output)
"""Run OpenMC in plotting mode"""
return self._run_openmc('openmc -p', output)
def run_simulation(self, particles=None, threads=None,
geometry_debug=False, restart_file=None,
tracks=False, mpi_procs=1, output=True):
"""Run an OpenMC simulation.
Parameters
----------
particles : int
Number of particles to simulate per generation
threads : int
Number of OpenMP threads
geometry_debug : bool
Turn on geometry debugging during simulation
restart_file : str
Path to restart file to use
tracks : bool
Write tracks for all particles
mpi_procs : int
Number of MPI processes
output : bool
Capture OpenMC output from standard out
"""
post_args = ' '
pre_args = ''

View file

@ -6,10 +6,30 @@ from openmc.constants import *
class Filter(object):
"""A filter used to constrain a tally to a specific criterion, e.g. only tally
events when the particle is in a certain cell and energy range.
Parameters
----------
type : str
The type of the tally filter. Acceptable values are "universe",
"material", "cell", "cellborn", "surface", "mesh", "energy",
"energyout", and "distribcell".
bins : int or list of int or list of float or ndarray
The bins for the filter. This takes on different meaning for different
filters.
Attributes
----------
type : str
The type of the tally filter.
bins : int or list of int or list of float or ndarray
The bins for the filter
"""
# Initialize Filter class attributes
def __init__(self, type=None, bins=None):
self.type = type
self._num_bins = 0
self.bins = bins
@ -17,9 +37,7 @@ class Filter(object):
self._offset = -1
self._stride = None
def __eq__(self, filter2):
# Check type
if self._type != filter2._type:
return False
@ -35,21 +53,14 @@ class Filter(object):
else:
return True
def __hash__(self):
hashable = []
hashable.append(self._type)
hashable.append(self._bins)
return hash(tuple(hashable))
return hash((self._type, self._bins))
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)
@ -66,57 +77,45 @@ class Filter(object):
else:
return existing
@property
def type(self):
return self._type
@property
def bins(self):
return self._bins
@property
def num_bins(self):
return self._num_bins
@property
def mesh(self):
return self._mesh
@property
def offset(self):
return self._offset
@property
def stride(self):
return self._stride
@type.setter
def type(self, type):
if type is None:
self._type = type
elif not type in FILTER_TYPES.values():
elif type not in FILTER_TYPES.values():
msg = 'Unable to set Filter type to "{0}" since it is not one ' \
'of the supported types'.format(type)
raise ValueError(msg)
self._type = type
@bins.setter
def bins(self, bins):
if bins is None:
self.num_bins = 0
elif self._type is None:
msg = 'Unable to set bins for Filter to "{0}" since ' \
'the Filter type has not yet been set'.format(bins)
@ -132,30 +131,23 @@ class Filter(object):
if self._type in ['cell', 'cellborn', 'surface', 'material',
'universe', 'distribcell']:
for edge in bins:
if not is_integer(edge):
msg = 'Unable to add bin "{0}" to a {1} Filter since ' \
'it is a non-integer'.format(edge, self._type)
raise ValueError(msg)
elif edge < 0:
msg = 'Unable to add bin "{0}" to a {1} Filter since ' \
'it is a negative integer'.format(edge, self._type)
raise ValueError(msg)
elif self._type in ['energy', 'energyout']:
for edge in bins:
if not is_integer(edge) and not is_float(edge):
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)
raise ValueError(msg)
elif edge < 0.:
msg = 'Unable to add bin edge "{0}" to a {1} Filter ' \
'since it is a negative value'.format(edge, self._type)
@ -163,27 +155,22 @@ class Filter(object):
# Check that bin edges are monotonically increasing
for index in range(len(bins)):
if index > 0 and bins[index] < bins[index-1]:
msg = 'Unable to add bin edges "{0}" to a {1} Filter ' \
'since they are not monotonically ' \
'increasing'.format(bins, self._type)
raise ValueError(msg)
# mesh filters
elif self._type == 'mesh':
if not len(bins) == 1:
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 is_integer(bins[0]):
msg = 'Unable to add bin "{0}" to mesh Filter since it ' \
'is a non-integer'.format(bins[0])
raise ValueError(msg)
elif bins[0] < 0:
msg = 'Unable to add bin "{0}" to mesh Filter since it ' \
'is a negative integer'.format(bins[0])
@ -192,11 +179,9 @@ class Filter(object):
# If all error checks passed, add bin edges
self._bins = bins
# FIXME
@num_bins.setter
def num_bins(self, num_bins):
if not is_integer(num_bins) or num_bins < 0:
msg = 'Unable to set the number of bins "{0}" for a {1} Filter ' \
'since it is not a positive ' \
@ -205,10 +190,8 @@ class Filter(object):
self._num_bins = num_bins
@mesh.setter
def mesh(self, mesh):
if not isinstance(mesh, Mesh):
msg = 'Unable to set Mesh to "{0}" for Filter since it is not a ' \
'Mesh object'.format(mesh)
@ -218,10 +201,8 @@ class Filter(object):
self.type = 'mesh'
self.bins = self._mesh._id
@offset.setter
def offset(self, offset):
if not is_integer(offset):
msg = 'Unable to set offset "{0}" for a {1} Filter since it is a ' \
'non-integer value'.format(offset, self._type)
@ -229,10 +210,8 @@ class Filter(object):
self._offset = offset
@stride.setter
def stride(self, stride):
if not is_integer(stride):
msg = 'Unable to set stride "{0}" for a {1} Filter since it is a ' \
'non-integer value'.format(stride, self._type)
@ -245,8 +224,20 @@ class Filter(object):
self._stride = stride
def can_merge(self, filter):
"""Determine if filter can be merged with another.
Parameters
----------
filter : Filter
Filter to compare with
Returns
-------
bool
Whether the filter can be merged
"""
if not isinstance(filter, Filter):
return False
@ -270,8 +261,20 @@ class Filter(object):
else:
return True
def merge(self, filter):
"""Merge this filter with another.
Parameters
----------
filter : Filter
Filter to merge with
Returns
-------
merged_filter : Filter
Filter resulting from the merge
"""
if not self.can_merge(filter):
msg = 'Unable to merge {0} with {1} filters'.format(self._type, filter._type)
@ -287,33 +290,30 @@ class Filter(object):
return merged_filter
def get_bin_index(self, filter_bin):
"""Returns the index in the Filter for some bin.
Parameters
----------
filter_bin : int, tuple
The bin is the integer ID for 'material', 'surface', 'cell',
'cellborn', and 'universe' Filters. The bin is an integer for
the cell instance ID for 'distribcell' Filters. The bin is
a 2-tuple of floats for 'energy' and 'energyout' filters
corresponding to the energy boundaries of the bin of interest.
The bin is a (x,y,z) 3-tuple for 'mesh' filters corresponding to
the mesh cell of interest.
filter_bin : int or tuple
The bin is the integer ID for 'material', 'surface', 'cell',
'cellborn', and 'universe' Filters. The bin is an integer for the
cell instance ID for 'distribcell' Filters. The bin is a 2-tuple of
floats for 'energy' and 'energyout' filters corresponding to the
energy boundaries of the bin of interest. The bin is a (x,y,z)
3-tuple for 'mesh' filters corresponding to the mesh cell of
interest.
Returns
-------
filter_index : int
The index in the Tally data array for this filter bin.
"""
# FIXME: This does not work for distribcells!!!
try:
# Filter bins for a mesh are an (x,y,z) tuple
if self.type == 'mesh':
# Convert (x,y,z) to a single bin -- this is similar to
# subroutine mesh_indices_to_bin in openmc/src/mesh.F90.
if (len(self.mesh.dimension) == 3):
@ -350,9 +350,7 @@ class Filter(object):
return filter_index
def __repr__(self):
string = 'Filter\n'
string += '{0: <16}{1}{2}\n'.format('\tType', '=\t', self._type)
string += '{0: <16}{1}{2}\n'.format('\tBins', '=\t', self._bins)

View file

@ -12,19 +12,24 @@ def reset_auto_ids():
class Geometry(object):
"""Geometry representing a collection of surfaces, cells, and universes.
Attributes
----------
root_universe : openmc.universe.Universe
Root universe which contains all others
"""
def __init__(self):
# Initialize Geometry class attributes
self._root_universe = None
self._offsets = {}
@property
def root_universe(self):
return self._root_universe
@root_universe.setter
def root_universe(self, root_universe):
@ -41,24 +46,26 @@ class Geometry(object):
self._root_universe = root_universe
def get_offset(self, path, filter_offset):
"""
Returns the corresponding location in the results array for a given
path and filter number.
"""Returns the corresponding location in the results array for a given path and
filter number. This is primarily intended to post-processing result when
a distribcell filter is used.
Parameters
----------
path : list
A list of IDs that form the path to the target. It should begin
with 0 for the base universe, and should cover every universe,
cell, and lattice passed through. For the case of the lattice,
a tuple should be provided to indicate which coordinates in the
lattice should be entered. This should be in the
form: (lat_id, i_x, i_y, i_z)
A list of IDs that form the path to the target. It should begin with
0 for the base universe, and should cover every universe, cell, and
lattice passed through. For the case of the lattice, a tuple should
be provided to indicate which coordinates in the lattice should be
entered. This should be in the form: (lat_id, i_x, i_y, i_z)
filter_offset : int
An integer that specifies which offset map the filter is using
An integer that specifies which offset map the filter is using
Returns
-------
offset : int
Location in the results array for the path and filter
"""
@ -74,16 +81,39 @@ class Geometry(object):
# Return the final offset
return offset
def get_all_cells(self):
"""Return all cells defined
Returns
-------
list of openmc.universe.Cell
Cells in the geometry
"""
return self._root_universe.get_all_cells()
def get_all_universes(self):
"""Return all universes defined
Returns
-------
list of openmc.universe.Universe
Universes in the geometry
"""
return self._root_universe.get_all_universes()
def get_all_nuclides(self):
"""Return all nuclides assigned to a material in the geometry
Returns
-------
list of openmc.nuclide.Nuclide
Nuclides in the geometry
"""
nuclides = {}
materials = self.get_all_materials()
@ -93,8 +123,15 @@ class Geometry(object):
return nuclides
def get_all_materials(self):
"""Return all materials assigned to a cell
Returns
-------
list of openmc.material.Material
Materials in the geometry
"""
material_cells = self.get_all_material_cells()
materials = set()
@ -104,9 +141,7 @@ class Geometry(object):
return list(materials)
def get_all_material_cells(self):
all_cells = self.get_all_cells()
material_cells = set()
@ -116,16 +151,21 @@ class Geometry(object):
return list(material_cells)
def get_all_material_universes(self):
"""Return all universes composed of at least one non-fill cell
Returns
-------
list of openmc.universe.Universe
Universes with non-fill cells
"""
all_universes = self.get_all_universes()
material_universes = set()
for universe_id, universe in all_universes.items():
cells = universe._cells
for cell_id, cell in cells.items():
if cell._type == 'normal':
material_universes.add(universe)
@ -133,24 +173,28 @@ class Geometry(object):
return list(material_universes)
class GeometryFile(object):
"""Geometry file used for an OpenMC simulation. Corresponds directly to the
geometry.xml input file.
Attributes
----------
geometry : Geometry
The geometry to be used
"""
def __init__(self):
# Initialize GeometryFile class attributes
self._geometry = None
self._geometry_file = ET.Element("geometry")
@property
def geometry(self):
return self._geometry
@geometry.setter
def geometry(self, geometry):
if not isinstance(geometry, Geometry):
msg = 'Unable to set the Geometry to {0} for the GeometryFile ' \
'since it is not a Geometry object'.format(geometry)
@ -158,8 +202,10 @@ class GeometryFile(object):
self._geometry = geometry
def export_to_xml(self):
"""Create a geometry.xml file that can be used for a simulation.
"""
root_universe = self._geometry._root_universe
root_universe.create_xml_subelement(self._geometry_file)

View file

@ -14,6 +14,7 @@ MATERIAL_IDS = []
# A static variable for auto-generated Material IDs
AUTO_MATERIAL_ID = 10000
def reset_auto_material_id():
global AUTO_MATERIAL_ID, MATERIAL_IDS
AUTO_MATERIAL_ID = 10000
@ -23,20 +24,36 @@ def reset_auto_material_id():
# Units for density supported by OpenMC
DENSITY_UNITS = ['g/cm3', 'g/cc', 'kg/cm3', 'at/b-cm', 'at/cm3', 'sum']
# ENDF temperatures
ENDF_TEMPS = np.array([300, 600, 700, 900, 1200, 1500])
# ENDF ZAIDs
ENDF_ZAIDS = np.array(['70c', '71c', '72c', '73c', '74c'])
# Constant for density when not needed
NO_DENSITY = 99999.
class Material(object):
"""A material composed of a collection of nuclides/elements that can be assigned
to a region of space.
Parameters
----------
material_id : int, optional
Unique identifier for the material. If not specified, an identifier will
automatically be assigned.
name : str, optional
Name of the material. If not specified, the name will be the empty
string.
Attributes
----------
id : int
Unique identifier for the material
density : float
Density of the material (units defined separately)
density_units : str
Units used for `density`. Can be one of 'g/cm3', 'g/cc', 'kg/cm3',
'atom/b-cm', 'atom/cm3', or 'sum'.
"""
def __init__(self, material_id=None, name=''):
# Initialize class attributes
self.id = material_id
self.name = name
@ -62,44 +79,36 @@ class Material(object):
# If specified, this file will be used instead of composition values
self._distrib_otf_file = None
@property
def id(self):
return self._id
@property
def name(self):
return self._name
@property
def density(self):
return self._density
@property
def density_units(self):
return self._density_units
@property
def convert_to_distrib_comps(self):
return self._convert_to_distrib_comps
@property
def distrib_otf_file(self):
return self._distrib_otf_file
@id.setter
def id(self, material_id):
global AUTO_MATERIAL_ID, MATERIAL_IDS
# If the Material already has an ID, remove it from global list
if hasattr(self, '_id') and not self._id is None:
if hasattr(self, '_id') and self._id is not None:
MATERIAL_IDS.remove(self._id)
if material_id is None:
@ -127,27 +136,34 @@ class Material(object):
self._id = material_id
MATERIAL_IDS.append(material_id)
@name.setter
def name(self, name):
if not is_string(name):
msg = 'Unable to set name for Material ID={0} with a non-string ' \
'value {1}'.format(self._id, name)
raise ValueError(msg)
else:
self._name = name
def set_density(self, units, density=NO_DENSITY):
"""Set the density of the material
Parameters
----------
units : str
Physical units of density
density : float, optional
Value of the density. Must be specified unless units is given as
'sum'.
"""
if not is_float(density):
msg = 'Unable to set the density for Material ID={0} to a ' \
'non-floating point value {1}'.format(self._id, density)
raise ValueError(msg)
elif not units in DENSITY_UNITS:
elif units not in DENSITY_UNITS:
msg = 'Unable to set the density for Material ID={0} with ' \
'units {1}'.format(self._id, units)
raise ValueError(msg)
@ -161,33 +177,40 @@ class Material(object):
self._density = density
self._density_units = units
@distrib_otf_file.setter
def distrib_otf_file(self, filename):
# TODO: remove this when distributed materials are merged
warnings.warn('This feature is not yet implemented in a release ' \
warnings.warn('This feature is not yet implemented in a release '
'version of openmc')
if not is_string(filename) and not filename is None:
if not is_string(filename) and filename is not None:
msg = 'Unable to add OTF material file to Material ID={0} with a ' \
'non-string name {1}'.format(self._id, filename)
raise ValueError(msg)
self._distrib_otf_file = filename
@convert_to_distrib_comps.setter
def convert_to_distrib_comps(self):
# TODO: remove this when distributed materials are merged
warnings.warn('This feature is not yet implemented in a release ' \
warnings.warn('This feature is not yet implemented in a release '
'version of openmc')
self._convert_to_distrib_comps = True
def add_nuclide(self, nuclide, percent, percent_type='ao'):
"""Add a nuclide to the material
Parameters
----------
nuclide : str or openmc.nuclide.Nuclide
Nuclide to add
percent : float
Atom or weight percent
percent_type : str
'ao' for atom percent and 'wo' for weight percent
"""
if not isinstance(nuclide, (openmc.Nuclide, str)):
msg = 'Unable to add a Nuclide to Material ID={0} with a ' \
@ -199,7 +222,7 @@ class Material(object):
'non-floating point value {1}'.format(self._id, percent)
raise ValueError(msg)
elif not percent_type in ['ao', 'wo', 'at/g-cm']:
elif percent_type not in ['ao', 'wo', 'at/g-cm']:
msg = 'Unable to add a Nuclide to Material ID={0} with a ' \
'percent type {1}'.format(self._id, percent_type)
raise ValueError(msg)
@ -213,8 +236,15 @@ class Material(object):
self._nuclides[nuclide._name] = (nuclide, percent, percent_type)
def remove_nuclide(self, nuclide):
"""Remove a nuclide from the material
Parameters
----------
nuclide : openmc.nuclide.Nuclide
Nuclide to remove
"""
if not isinstance(nuclide, openmc.Nuclide):
msg = 'Unable to remove a Nuclide {0} in Material ID={1} ' \
@ -225,8 +255,19 @@ class Material(object):
if nuclide._name in self._nuclides:
del self._nuclides[nuclide._name]
def add_element(self, element, percent, percent_type='ao'):
"""Add a natural element to the material
Parameters
----------
element : openmc.element.Element
Element to add
percent : float
Atom or weight percent
percent_type : str
'ao' for atom percent and 'wo' for weight percent
"""
if not isinstance(element, openmc.Element):
msg = 'Unable to add an Element to Material ID={0} with a ' \
@ -238,7 +279,7 @@ class Material(object):
'non-floating point value {1}'.format(self._id, percent)
raise ValueError(msg)
if not percent_type in ['ao', 'wo']:
if percent_type not in ['ao', 'wo']:
msg = 'Unable to add an Element to Material ID={0} with a ' \
'percent type {1}'.format(self._id, percent_type)
raise ValueError(msg)
@ -248,15 +289,31 @@ class Material(object):
self._elements[element._name] = (element, percent, percent_type)
def remove_element(self, element):
"""Remove a natural element from the material
Parameters
----------
element : openmc.element.Element
Element to remove
"""
# If the Material contains the Element, delete it
if element._name in self._elements:
del self._elements[element._name]
def add_s_alpha_beta(self, name, xs):
r"""Add an :math:`S(\alpha,\beta)` table to the material
Parameters
----------
name : str
Name of the :math:`S(\alpha,\beta)` table
xs : str
Cross section identifier, e.g. '71t'
"""
if not is_string(name):
msg = 'Unable to add an S(a,b) table to Material ID={0} with a ' \
@ -270,8 +327,16 @@ class Material(object):
self._sab.append((name, xs))
def get_all_nuclides(self):
"""Returns all nuclides in the material
Returns
-------
nuclides : dict
Dictionary whose keys are nuclide names and values are 2-tuples of
(nuclide, density)
"""
nuclides = {}
@ -282,9 +347,7 @@ class Material(object):
return nuclides
def _repr__(self):
string = 'Material\n'
string += '{0: <16}{1}{2}\n'.format('\tID', '=\t', self._id)
string += '{0: <16}{1}{2}\n'.format('\tName', '=\t', self._name)
@ -316,9 +379,7 @@ class Material(object):
return string
def get_nuclide_xml(self, nuclide, distrib=False):
def _get_nuclide_xml(self, nuclide, distrib=False):
xml_element = ET.Element("nuclide")
xml_element.set("name", nuclide[0]._name)
@ -328,14 +389,12 @@ class Material(object):
else:
xml_element.set("wo", str(nuclide[1]))
if not nuclide[0]._xs is None:
if nuclide[0]._xs is not None:
xml_element.set("xs", nuclide[0]._xs)
return xml_element
def get_element_xml(self, element, distrib=False):
def _get_element_xml(self, element, distrib=False):
xml_element = ET.Element("element")
xml_element.set("name", str(element[0]._name))
@ -347,28 +406,31 @@ class Material(object):
return xml_element
def get_nuclides_xml(self, nuclides, distrib=False):
def _get_nuclides_xml(self, nuclides, distrib=False):
xml_elements = []
for nuclide in nuclides.values():
xml_elements.append(self.get_nuclide_xml(nuclide, distrib))
xml_elements.append(self._get_nuclide_xml(nuclide, distrib))
return xml_elements
def get_elements_xml(self, elements, distrib=False):
def _get_elements_xml(self, elements, distrib=False):
xml_elements = []
for element in elements.values():
xml_elements.append(self.get_element_xml(element, distrib))
xml_elements.append(self._get_element_xml(element, distrib))
return xml_elements
def get_material_xml(self):
"""Return XML representation of the material
Returns
-------
element : xml.etree.ElementTree.Element
XML element containing material data
"""
# Create Material XML element
element = ET.Element("material")
@ -384,19 +446,17 @@ class Material(object):
subelement.set("units", self._density_units)
if not self._convert_to_distrib_comps:
# Create nuclide XML subelements
subelements = self.get_nuclides_xml(self._nuclides)
subelements = self._get_nuclides_xml(self._nuclides)
for subelement in subelements:
element.append(subelement)
# Create element XML subelements
subelements = self.get_elements_xml(self._elements)
subelements = self._get_elements_xml(self._elements)
for subelement in subelements:
element.append(subelement)
else:
subelement = ET.SubElement(element, "compositions")
comps = []
@ -409,29 +469,24 @@ class Material(object):
raise ValueError(msg)
comps.append(per)
if self._distrib_otf_file is None:
# Create values and units subelements
subsubelement = ET.SubElement(subelement, "values")
subsubelement.text = ' '.join([str(c) for c in comps])
subsubelement = ET.SubElement(subelement, "units")
subsubelement.text = dist_per_type
else:
# Specify the materials file
subsubelement = ET.SubElement(subelement, "otf_file_path")
subsubelement.text = self._distrib_otf_file
# Create nuclide XML subelements
subelements = self.get_nuclides_xml(self._nuclides, distrib=True)
for subelement_nuc in subelements:
subelement.append(subelement_nuc)
# Create element XML subelements
subelements = self.get_elements_xml(self._elements, distrib=True)
subelements = self._get_elements_xml(self._elements, distrib=True)
for subelement_ele in subelements:
subelement.append(subelement_ele)
@ -445,31 +500,44 @@ class Material(object):
class MaterialsFile(object):
"""Materials file used for an OpenMC simulation. Corresponds directly to the
materials.xml input file.
Attributes
----------
default_xs : str
The default cross section identifier applied to a nuclide when none is
specified
"""
def __init__(self):
# Initialize MaterialsFile class attributes
self._materials = []
self._default_xs = None
self._materials_file = ET.Element("materials")
@property
def default_xs(self):
return self._default_xs
@default_xs.setter
def default_xs(self, xs):
if not is_string(xs):
msg = 'Unable to set default xs to a non-string value'.format(xs)
raise ValueError(msg)
self._default_xs = xs
def add_material(self, material):
"""Add a material to the file.
Parameters
----------
material : Material
Material to add
"""
if not isinstance(material, Material):
msg = 'Unable to add a non-Material {0} to the ' \
@ -478,8 +546,15 @@ class MaterialsFile(object):
self._materials.append(material)
def add_materials(self, materials):
"""Add multiple materials to the file.
Parameters
----------
materials : tuple or list of Material
Materials to add
"""
if not isinstance(materials, (tuple, list, MappingView)):
msg = 'Unable to create OpenMC materials.xml file from {0} which ' \
@ -489,8 +564,15 @@ class MaterialsFile(object):
for material in materials:
self.add_material(material)
def remove_material(self, material):
"""Remove a material from the file
def remove_materials(self, material):
Parameters
----------
material : Material
Material to remove
"""
if not isinstance(material, Material):
msg = 'Unable to remove a non-Material {0} from the ' \
@ -499,22 +581,22 @@ class MaterialsFile(object):
self._materials.remove(material)
def create_material_subelements(self):
def _create_material_subelements(self):
subelement = ET.SubElement(self._materials_file, "default_xs")
if not self._default_xs is None:
if self._default_xs is not None:
subelement.text = self._default_xs
for material in self._materials:
xml_element = material.get_material_xml()
self._materials_file.append(xml_element)
def export_to_xml(self):
"""Create a materials.xml file that can be used for a simulation.
self.create_material_subelements()
"""
self._create_material_subelements()
# Clean the indentation in the file to be user-readable
sort_xml_elements(self._materials_file)

View file

@ -14,9 +14,37 @@ def reset_auto_mesh_id():
class Mesh(object):
"""A structured Cartesian mesh in two or three dimensions
Parameters
----------
mesh_id : int
Unique identifier for the mesh
name : str
Name of the mesh
Attributes
----------
id : int
Unique identifier for the mesh
name : str
Name of the mesh
type : str
Type of the mesh
dimension : tuple or list or ndarray
The number of mesh cells in each direction.
lower_left : tuple or list or ndarray
The lower-left corner of the structured mesh. If only two coordinate are
given, it is assumed that the mesh is an x-y mesh.
upper_right : tuple or list or ndarray
The upper-right corner of the structrued mesh. If only two coordinate
are given, it is assumed that the mesh is an x-y mesh.
width : tuple or list or ndarray
The width of mesh cells in each direction.
"""
def __init__(self, mesh_id=None, name=''):
# Initialize Mesh class attributes
self.id = mesh_id
self.name = name
@ -26,9 +54,7 @@ class Mesh(object):
self._upper_right = None
self._width = None
def __eq__(self, mesh2):
# Check type
if self._type != mesh2._type:
return False
@ -49,14 +75,11 @@ 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
@ -74,50 +97,40 @@ class Mesh(object):
else:
return existing
@property
def id(self):
return self._id
@property
def name(self):
return self._name
@property
def type(self):
return self._type
@property
def dimension(self):
return self._dimension
@property
def lower_left(self):
return self._lower_left
@property
def upper_right(self):
return self._upper_right
@property
def width(self):
return self._width
@property
def num_mesh_cells(self):
return np.prod(self._dimension)
@id.setter
def id(self, mesh_id):
if mesh_id is None:
global AUTO_MESH_ID
self._id = AUTO_MESH_ID
@ -136,28 +149,22 @@ class Mesh(object):
else:
self._id = mesh_id
@name.setter
def name(self, name):
if not is_string(name):
msg = 'Unable to set name for Mesh ID={0} with a non-string ' \
'value {1}'.format(self._id, name)
raise ValueError(msg)
else:
self._name = name
@type.setter
def type(self, type):
if not is_string(type):
msg = 'Unable to set Mesh ID={0} for type {1} which is not ' \
'a string'.format(self._id, type)
raise ValueError(msg)
elif not type in ['rectangular', 'hexagonal']:
elif type not in ['rectangular', 'hexagonal']:
msg = 'Unable to set Mesh ID={0} for type {1} which since ' \
'only rectangular and hexagonal meshes are ' \
'supported '.format(self._id, type)
@ -165,23 +172,19 @@ class Mesh(object):
self._type = type
@dimension.setter
def dimension(self, dimension):
if not isinstance(dimension, (tuple, list, np.ndarray)):
msg = 'Unable to set Mesh ID={0} with dimension {1} which is ' \
'not a Python list, tuple or NumPy ' \
'array'.format(self._id, dimension)
raise ValueError(msg)
elif len(dimension) != 2 and len(dimension) != 3:
msg = 'Unable to set Mesh ID={0} with dimension {1} since it ' \
'must include 2 or 3 dimensions'.format(self._id, dimension)
raise ValueError(msg)
for dim in dimension:
if not is_integer(dim):
msg = 'Unable to set Mesh ID={0} with dimension {1} which ' \
'is a non-integer'.format(self._id, dim)
@ -189,10 +192,8 @@ class Mesh(object):
self._dimension = dimension
@lower_left.setter
def lower_left(self, lower_left):
if not isinstance(lower_left, (tuple, list, np.ndarray)):
msg = 'Unable to set Mesh ID={0} with lower_left {1} which is ' \
'not a Python list, tuple or NumPy ' \
@ -205,7 +206,6 @@ class Mesh(object):
raise ValueError(msg)
for coord in lower_left:
if not is_integer(coord) and not is_float(coord):
msg = 'Unable to set Mesh ID={0} with lower_left {1} which ' \
'is neither neither an integer nor a floating point ' \
@ -214,10 +214,8 @@ class Mesh(object):
self._lower_left = lower_left
@upper_right.setter
def upper_right(self, upper_right):
if not isinstance(upper_right, (tuple, list, np.ndarray)):
msg = 'Unable to set Mesh ID={0} with upper_right {1} which ' \
'is not a Python list, tuple or NumPy ' \
@ -230,7 +228,6 @@ class Mesh(object):
raise ValueError(msg)
for coord in upper_right:
if not is_integer(coord) and not is_float(coord):
msg = 'Unable to set Mesh ID={0} with upper_right {1} which ' \
'is neither an integer nor a floating point ' \
@ -239,12 +236,9 @@ class Mesh(object):
self._upper_right = upper_right
@width.setter
def width(self, width):
if not width is None:
if width is not None:
if not isinstance(width, (tuple, list, np.ndarray)):
msg = 'Unable to set Mesh ID={0} with width {1} which ' \
'is not a Python list, tuple or NumPy ' \
@ -257,7 +251,6 @@ class Mesh(object):
raise ValueError(msg)
for dim in width:
if not is_integer(dim) and not is_float(dim):
msg = 'Unable to set Mesh ID={0} with width {1} which is ' \
'neither an integer nor a floating point ' \
@ -266,9 +259,7 @@ class Mesh(object):
self._width = width
def __repr__(self):
string = 'Mesh\n'
string += '{0: <16}{1}{2}\n'.format('\tID', '=\t', self._id)
string += '{0: <16}{1}{2}\n'.format('\tName', '=\t', self._name)
@ -279,53 +270,32 @@ class Mesh(object):
string += '{0: <16}{1}{2}\n'.format('\tPixels', '=\t', self._width)
return string
def get_mesh_xml(self):
"""Return XML representation of the mesh
Returns
-------
element : xml.etree.ElementTree.Element
XML element containing mesh data
"""
element = ET.Element("mesh")
element.set("id", str(self._id))
element.set("type", self._type)
if len(self._dimension) == 2:
subelement = ET.SubElement(element, "dimension")
subelement.text = '{0} {1}'.format(self._dimension[0],
self._dimension[1])
else:
subelement = ET.SubElement(element, "dimension")
subelement.text = '{0} {1} {2}'.format(self._dimension[0],
self._dimension[1],
self._dimension[2])
subelement = ET.SubElement(element, "dimension")
subelement.text = ' '.join(map(str, self._dimension))
if len(self._lower_left) == 2:
subelement = ET.SubElement(element, "lower_left")
subelement.text = '{0} {1}'.format(self._lower_left[0],
self._lower_left[1])
else:
subelement = ET.SubElement(element, "lower_left")
subelement.text = '{0} {1} {2}'.format(self._lower_left[0],
self._lower_left[1],
self._lower_left[2])
subelement = ET.SubElement(element, "lower_left")
subelement.text = ' '.join(map(str, self._lower_left))
if not self._upper_right is None:
if len(self._upper_right) == 2:
subelement = ET.SubElement(element, "upper_right")
subelement.text = '{0} {1}'.format(self._upper_right[0],
self._upper_right[1])
else:
subelement = ET.SubElement(element, "upper_right")
subelement.text = '{0} {1} {2}'.format(self._upper_right[0],
self._upper_right[1],
self._upper_right[2])
if self._upper_right is not None:
subelement = ET.SubElement(element, "upper_right")
subelement.text = ' '.join(map(str, self._upper_right))
if not self._width is None:
if len(self._width) == 2:
subelement = ET.SubElement(element, "width")
subelement.text = '{0} {1}'.format(self._width[0],
self._width[1])
else:
subelement = ET.SubElement(element, "width")
subelement.text = '{0} {1} {2}'.format(self._width[0],
self._width[1],
self._width[2])
if self._width is not None:
subelement = ET.SubElement(element, "width")
subelement.text = ' '.join(map(str, self._width))
return element
return element

View file

@ -2,9 +2,28 @@ from openmc.checkvalue import *
class Nuclide(object):
"""A nuclide that can be used in a material.
Parameters
----------
name : str
Name of the nuclide, e.g. U-235
xs : str
Cross section identifier, e.g. 71c
Attributes
----------
name : str
Name of the nuclide, e.g. U-235
xs : str
Cross section identifier, e.g. 71c
zaid : int
1000*(atomic number) + mass number. As an example, the zaid of U-235
would be 92235.
"""
def __init__(self, name='', xs=None):
# Initialize class attributes
self._name = ''
self._xs = None
@ -13,12 +32,10 @@ class Nuclide(object):
# Set the Material class attributes
self.name = name
if not xs is None:
if xs is not None:
self.xs = xs
def __eq__(self, nuclide2):
# Check type
if not isinstance(nuclide2, Nuclide):
return False
@ -34,29 +51,23 @@ class Nuclide(object):
else:
return True
def __hash__(self):
return hash((self._name, self._xs))
@property
def name(self):
return self._name
@property
def xs(self):
return self._xs
@property
def zaid(self):
return self._zaid
@name.setter
def name(self, name):
if not is_string(name):
msg = 'Unable to set name for Nuclide with a non-string ' \
'value {0}'.format(name)
@ -64,10 +75,8 @@ class Nuclide(object):
self._name = name
@xs.setter
def xs(self, xs):
if not is_string(xs):
msg = 'Unable to set cross-section identifier xs for Nuclide ' \
'with a non-string value {0}'.format(xs)
@ -75,10 +84,8 @@ class Nuclide(object):
self._xs = xs
@zaid.setter
def zaid(self, zaid):
if not is_integer(zaid):
msg = 'Unable to set zaid for Nuclide ' \
'with a non-integer {0}'.format(zaid)
@ -86,11 +93,9 @@ class Nuclide(object):
self._zaid = zaid
def __repr__(self):
string = 'Nuclide - {0}\n'.format(self._name)
string += '{0: <16}{1}{2}\n'.format('\tXS', '=\t', self._xs)
if self._zaid is not None:
string += '{0: <16}{1}{2}\n'.format('\tZAID', '=\t', self._zaid)
return string
string += '{0: <16}{1}{2}\n'.format('\tZAID', '=\t', self._zaid)
return string

View file

@ -62,8 +62,20 @@ OPENMC_LATTICES = {}
OPENCG_LATTICES = {}
def get_opencg_material(openmc_material):
"""Return an OpenCG material corresponding to an OpenMC material.
Parameters
----------
openmc_material : openmc.material.Material
OpenMC material
Returns
-------
opencg_material : opencg.Material
Equivalent OpenCG material
"""
if not isinstance(openmc_material, openmc.Material):
msg = 'Unable to create an OpenCG Material from {0} ' \
@ -91,6 +103,19 @@ def get_opencg_material(openmc_material):
def get_openmc_material(opencg_material):
"""Return an OpenMC material corresponding to an OpenCG material.
Parameters
----------
opencg_material : opencg.Material
OpenCG material
Returns
-------
openmc_material : openmc.material.Material
Equivalent OpenMC material
"""
if not isinstance(opencg_material, opencg.Material):
msg = 'Unable to create an OpenMC Material from {0} ' \
@ -118,6 +143,25 @@ def get_openmc_material(opencg_material):
def is_opencg_surface_compatible(opencg_surface):
"""Determine whether OpenCG surface is compatible with OpenMC geometry.
A surface is considered compatible if there is a one-to-one correspondence
between OpenMC and OpenCG surface types. Note that some OpenCG surfaces,
e.g. SquarePrism, do not have a one-to-one correspondence with OpenMC
surfaces but can still be converted into an equivalent collection of OpenMC
surfaces.
Parameters
----------
opencg_surface : opencg.Surface
OpenCG surface
Returns
-------
bool
Whether OpenCG surface is compatible with OpenMC
"""
if not isinstance(opencg_surface, opencg.Surface):
msg = 'Unable to check if OpenCG Surface is compatible' \
@ -132,6 +176,19 @@ def is_opencg_surface_compatible(opencg_surface):
def get_opencg_surface(openmc_surface):
"""Return an OpenCG surface corresponding to an OpenMC surface.
Parameters
----------
openmc_surface : openmc.surface.Surface
OpenMC surface
Returns
-------
opencg_surface : opencg.Surface
Equivalent OpenCG surface
"""
if not isinstance(openmc_surface, openmc.Surface):
msg = 'Unable to create an OpenCG Surface from {0} ' \
@ -205,6 +262,19 @@ def get_opencg_surface(openmc_surface):
def get_openmc_surface(opencg_surface):
"""Return an OpenMC surface corresponding to an OpenCG surface.
Parameters
----------
opencg_surface : opencg.Surface
OpenCG surface
Returns
-------
openmc_surface : openmc.surface.Surface
Equivalent OpenMC surface
"""
if not isinstance(opencg_surface, opencg.Surface):
msg = 'Unable to create an OpenMC Surface from {0} which ' \
@ -269,7 +339,6 @@ def get_openmc_surface(opencg_surface):
'Surface type in OpenMC'.format(opencg_surface._type)
raise ValueError(msg)
# Add the OpenMC Surface to the global collection of all OpenMC Surfaces
OPENMC_SURFACES[surface_id] = openmc_surface
@ -280,6 +349,23 @@ def get_openmc_surface(opencg_surface):
def get_compatible_opencg_surfaces(opencg_surface):
"""Generate OpenCG surfaces that are compatible with OpenMC equivalent to an
OpenCG surface that is not compatible. For example, this method may be used
to convert a ZSquarePrism OpenCG surface into a collection of equivalent
XPlane and YPlane OpenCG surfaces.
Parameters
----------
opencg_surface : opencg.Surface
OpenCG surface that is incompatible with OpenMC
Returns
-------
surfaces : list of opencg.Surface
Collection of surfaces equivalent to the original one but compatible
with OpenMC
"""
if not isinstance(opencg_surface, opencg.Surface):
msg = 'Unable to create an OpenMC Surface from {0} which ' \
@ -349,6 +435,19 @@ def get_compatible_opencg_surfaces(opencg_surface):
def get_opencg_cell(openmc_cell):
"""Return an OpenCG cell corresponding to an OpenMC cell.
Parameters
----------
openmc_cell : openmc.universe.Cell
OpenMC cell
Returns
-------
opencg_cell : opencg.Cell
Equivalent OpenCG cell
"""
if not isinstance(openmc_cell, openmc.Cell):
msg = 'Unable to create an OpenCG Cell from {0} which ' \
@ -398,7 +497,26 @@ def get_opencg_cell(openmc_cell):
def get_compatible_opencg_cells(opencg_cell, opencg_surface, halfspace):
"""Generate OpenCG cells that are compatible with OpenMC equivalent to an OpenCG
cell that is not compatible.
Parameters
----------
opencg_cell : opencg.Cell
OpenCG cell
opencg_surface : opencg.Surface
OpenCG surface that causes the incompatibility, e.g. an instance of
XSquarePrism
halfspace : {-1, 1}
Which halfspace defined by the surface is contained in the cell
Returns
-------
compatible_cells : list of opencg.Cell
Collection of cells equivalent to the original one but compatible with
OpenMC
"""
if not isinstance(opencg_cell, opencg.Cell):
msg = 'Unable to create compatible OpenMC Cell from {0} which ' \
'is not an OpenCG Cell'.format(opencg_cell)
@ -409,7 +527,7 @@ def get_compatible_opencg_cells(opencg_cell, opencg_surface, halfspace):
'not an OpenCG Surface'.format(opencg_surface)
raise ValueError(msg)
elif not halfspace in [-1, +1]:
elif halfspace not in [-1, +1]:
msg = 'Unable to create compatible Cell since {0}' \
'is not a +/-1 halfspace'.format(halfspace)
raise ValueError(msg)
@ -418,8 +536,8 @@ def get_compatible_opencg_cells(opencg_cell, opencg_surface, halfspace):
compatible_cells = []
# SquarePrism Surfaces
if opencg_surface._type in ['x-squareprism',
'y-squareprism', 'z-squareprism']:
if opencg_surface._type in ['x-squareprism', 'y-squareprism',
'z-squareprism']:
# Get the compatible Surfaces (XPlanes and YPlanes)
compatible_surfaces = get_compatible_opencg_surfaces(opencg_surface)
@ -436,13 +554,11 @@ def get_compatible_opencg_cells(opencg_cell, opencg_surface, halfspace):
# If Cell is outside SquarePrism, add "outside" of Surface halfspaces
else:
# Create 8 Cell clones to represent each of the disjoint planar
# Surface halfspace intersections
num_clones = 8
for clone_id in range(num_clones):
# Create a cloned OpenCG Cell with Surfaces compatible with OpenMC
clone = opencg_cell.clone()
compatible_cells.append(clone)
@ -500,6 +616,14 @@ def get_compatible_opencg_cells(opencg_cell, opencg_surface, halfspace):
def make_opencg_cells_compatible(opencg_universe):
"""Make all cells in an OpenCG universe compatible with OpenMC.
Parameters
----------
opencg_universe : opencg.Universe
Universe to check
"""
if not isinstance(opencg_universe, opencg.Universe):
msg = 'Unable to make compatible OpenCG Cells for {0} which ' \
@ -545,8 +669,20 @@ def make_opencg_cells_compatible(opencg_universe):
return
def get_openmc_cell(opencg_cell):
"""Return an OpenMC cell corresponding to an OpenCG cell.
Parameters
----------
opencg_cell : opencg.Cell
OpenCG cell
Returns
-------
openmc_cell : openmc.universe.Cell
Equivalent OpenMC cell
"""
if not isinstance(opencg_cell, opencg.Cell):
msg = 'Unable to create an OpenMC Cell from {0} which ' \
@ -597,8 +733,20 @@ def get_openmc_cell(opencg_cell):
return openmc_cell
def get_opencg_universe(openmc_universe):
"""Return an OpenCG universe corresponding to an OpenMC universe.
Parameters
----------
openmc_universe : openmc.universe.Universe
OpenMC universe
Returns
-------
opencg_universe : opencg.Universe
Equivalent OpenCG universe
"""
if not isinstance(openmc_universe, openmc.Universe):
msg = 'Unable to create an OpenCG Universe from {0} which ' \
@ -633,6 +781,19 @@ def get_opencg_universe(openmc_universe):
def get_openmc_universe(opencg_universe):
"""Return an OpenMC universe corresponding to an OpenCG universe.
Parameters
----------
opencg_universe : opencg.Universe
OpenCG universe
Returns
-------
openmc_universe : openmc.universe.Universe
Equivalent OpenMC universe
"""
if not isinstance(opencg_universe, opencg.Universe):
msg = 'Unable to create an OpenMC Universe from {0} which ' \
@ -670,6 +831,19 @@ def get_openmc_universe(opencg_universe):
def get_opencg_lattice(openmc_lattice):
"""Return an OpenCG lattice corresponding to an OpenMC lattice.
Parameters
----------
openmc_lattice : openmc.universe.Lattice
OpenMC lattice
Returns
-------
opencg_lattice : opencg.Lattice
Equivalent OpenCG lattice
"""
if not isinstance(openmc_lattice, openmc.Lattice):
msg = 'Unable to create an OpenCG Lattice from {0} which ' \
@ -701,7 +875,7 @@ def get_opencg_lattice(openmc_lattice):
lower_left = new_lower_left
# Initialize an empty array for the OpenCG nested Universes in this Lattice
universe_array = np.ndarray(tuple(np.array(dimension)[::-1]), \
universe_array = np.ndarray(tuple(np.array(dimension)[::-1]),
dtype=opencg.Universe)
# Create OpenCG Universes for each unique nested Universe in this Lattice
@ -723,8 +897,8 @@ def get_opencg_lattice(openmc_lattice):
opencg_lattice.setUniverses(universe_array)
offset = np.array(lower_left, dtype=np.float64) - \
((np.array(pitch, dtype=np.float64) * \
np.array(dimension, dtype=np.float64))) / -2.0
((np.array(pitch, dtype=np.float64) *
np.array(dimension, dtype=np.float64))) / -2.0
opencg_lattice.setOffset(offset)
# Add the OpenMC Lattice to the global collection of all OpenMC Lattices
@ -737,6 +911,19 @@ def get_opencg_lattice(openmc_lattice):
def get_openmc_lattice(opencg_lattice):
"""Return an OpenMC lattice corresponding to an OpenCG lattice.
Parameters
----------
opencg_lattice : opencg.Lattice
OpenCG lattice
Returns
-------
openmc_lattice : openmc.universe.Lattice
Equivalent OpenMC lattice
"""
if not isinstance(opencg_lattice, opencg.Lattice):
msg = 'Unable to create an OpenMC Lattice from {0} which ' \
@ -756,7 +943,7 @@ def get_openmc_lattice(opencg_lattice):
universes = opencg_lattice._universes
# Initialize an empty array for the OpenMC nested Universes in this Lattice
universe_array = np.ndarray(tuple(np.array(dimension)), \
universe_array = np.ndarray(tuple(np.array(dimension)),
dtype=openmc.Universe)
# Create OpenMC Universes for each unique nested Universe in this Lattice
@ -773,11 +960,11 @@ def get_openmc_lattice(opencg_lattice):
universe_array[x][y][z] = unique_universes[universe_id]
# Reverse y-dimension in array to match ordering in OpenCG
universe_array = universe_array[:,::-1,:]
universe_array = universe_array[:, ::-1, :]
lower_left = np.array(offset, dtype=np.float64) + \
((np.array(width, dtype=np.float64) * \
np.array(dimension, dtype=np.float64))) / -2.0
((np.array(width, dtype=np.float64) *
np.array(dimension, dtype=np.float64))) / -2.0
openmc_lattice = openmc.RectLattice(lattice_id=lattice_id)
openmc_lattice.dimension = dimension
@ -795,6 +982,19 @@ def get_openmc_lattice(opencg_lattice):
def get_opencg_geometry(openmc_geometry):
"""Return an OpenCG geometry corresponding to an OpenMC geometry.
Parameters
----------
openmc_geometry : openmc.universe.Geometry
OpenMC geometry
Returns
-------
opencg_geometry : opencg.Geometry
Equivalent OpenCG geometry
"""
if not isinstance(openmc_geometry, openmc.Geometry):
msg = 'Unable to get OpenCG geometry from {0} which is ' \
@ -822,6 +1022,19 @@ def get_opencg_geometry(openmc_geometry):
def get_openmc_geometry(opencg_geometry):
"""Return an OpenMC geometry corresponding to an OpenCG geometry.
Parameters
----------
opencg_geometry : opencg.Geometry
OpenCG geometry
Returns
-------
openmc_geometry : openmc.universe.Geometry
Equivalent OpenMC geometry
"""
if not isinstance(opencg_geometry, opencg.Geometry):
msg = 'Unable to get OpenMC geometry from {0} which is ' \
@ -849,8 +1062,8 @@ def get_openmc_geometry(opencg_geometry):
# Make the entire geometry "compatible" before assigning auto IDs
universes = opencg_geometry.getAllUniverses()
for universe_id, universe in universes.items():
if not isinstance(universe, opencg.Lattice):
make_opencg_cells_compatible(universe)
if not isinstance(universe, opencg.Lattice):
make_opencg_cells_compatible(universe)
opencg_geometry.assignAutoIds()

View file

@ -1,9 +1,44 @@
#!/usr/bin/env python
import struct
class Particle(object):
"""Information used to restart a specific particle that caused a simulation to
fail.
Parameters
----------
filename : str
Path to the particle restart file
Attributes
----------
filetype : int
Integer indicating the file type
revision : int
Revision of the particle restart format
current_batch : int
The batch containing the particle
gen_per_batch : int
Number of generations per batch
current_gen : int
The generation containing the particle
n_particles : int
Number of particles per generation
run_mode : int
Type of simulation (criticality or fixed source)
id : long
Identifier of the particle
weight : float
Weight of the particle
energy : float
Energy of the particle in MeV
xyz : list of float
Position of the particle
uvw : list of float
Directional cosines of the particle
"""
def __init__(self, filename):
if filename.endswith('.h5'):
import h5py

View file

@ -9,6 +9,7 @@ from openmc.clean_xml import *
# A static variable for auto-generated Plot IDs
AUTO_PLOT_ID = 10000
def reset_auto_plot_id():
global AUTO_PLOT_ID
AUTO_PLOT_ID = 10000
@ -18,9 +19,50 @@ BASES = ['xy', 'xz', 'yz']
class Plot(object):
"""Definition of a finite region of space to be plotted, either as a slice plot
in two dimensions or as a voxel plot in three dimensions.
Parameters
----------
plot_id : int
Unique identifier for the plot
name : str
Name of the plot
Attributes
----------
id : int
Unique identifier
name : str
Name of the plot
width : tuple or list or ndarray
Width of the plot in each basis direction
pixels : tuple or list or ndarray
Number of pixels to use in each basis direction
origin : tuple or list of ndarray
Origin (center) of the plot
filename :
Path to write the plot to
color : {'cell', 'mat'}
Indicate whether the plot should be colored by cell or by material
type : {'slice', 'voxel'}
The type of the plot
basis : {'xy', 'xz', 'yz'}
The basis directions for the plot
background : tuple or list of ndarray
Color of the background defined by RGB
mask_components : tuple or list or ndarray
Unique id numbers of the cells or materials to plot
mask_background : tuple or list or ndarray
Color to apply to all cells/materials not listed in mask_components
defined by RGB
col_spec : dict
Dictionary indicating that certain cells/materials (keys) should be
colored with a specific RGB (values)
"""
def __init__(self, plot_id=None, name=''):
# Initialize Plot class attributes
self.id = plot_id
self.name = name
@ -36,75 +78,60 @@ class Plot(object):
self._mask_background = None
self._col_spec = None
@property
def id(self):
return self._id
@property
def name(self):
return self._name
@property
def width(self):
return self._width
@property
def pixels(self):
return self._pixels
@property
def origin(self):
return self._origin
@property
def filename(self):
return self._filename
@property
def color(self):
return self._color
@property
def type(self):
return self._type
@property
def basis(self):
return self._basis
@property
def background(self):
return self._background
@property
def mask_componenets(self):
return self._mask_components
@property
def mask_background(self):
return self._mask_background
@property
def col_spec(self):
return self._col_spec
@id.setter
def id(self, plot_id):
if plot_id is None:
global AUTO_PLOT_ID
self._id = AUTO_PLOT_ID
@ -123,7 +150,6 @@ class Plot(object):
else:
self._id = plot_id
@name.setter
def name(self, name):
@ -135,10 +161,8 @@ class Plot(object):
else:
self._name = name
@width.setter
def width(self, width):
if not isinstance(width, (tuple, list, np.ndarray)):
msg = 'Unable to create Plot ID={0} with width {1} which is not ' \
'a Python tuple/list or NumPy array'.format(self._id, width)
@ -158,10 +182,8 @@ class Plot(object):
self._width = width
@origin.setter
def origin(self, origin):
if not isinstance(origin, (tuple, list, np.ndarray)):
msg = 'Unable to create Plot ID={0} with origin {1} which is not ' \
'a Python tuple/list or NumPy array'.format(self._id, origin)
@ -172,7 +194,6 @@ class Plot(object):
'a 3D coordinate must be input'.format(self._id, origin)
raise ValueError(msg)
for dim in origin:
if not is_integer(dim) and not is_float(dim):
msg = 'Unable to create Plot ID={0} with origin {1} since ' \
@ -182,10 +203,8 @@ class Plot(object):
self._origin = origin
@pixels.setter
def pixels(self, pixels):
if not isinstance(pixels, (tuple, list, np.ndarray)):
msg = 'Unable to create Plot ID={0} with pixels {1} which is not ' \
'a Python tuple/list or NumPy array'.format(self._id, pixels)
@ -197,7 +216,6 @@ class Plot(object):
raise ValueError(msg)
for dim in pixels:
if not is_integer(dim):
msg = 'Unable to create Plot ID={0} with pixel value {1} ' \
'which is not an integer'.format(self._id, dim)
@ -210,10 +228,8 @@ class Plot(object):
self._pixels = pixels
@filename.setter
def filename(self, filename):
if not is_string(filename):
msg = 'Unable to create Plot ID={0} with filename {1} which is ' \
'not a string'.format(self._id, filename)
@ -221,58 +237,50 @@ class Plot(object):
self._filename = filename
@color.setter
def color(self, color):
if not is_string(color):
msg = 'Unable to create Plot ID={0} with color {1} which is not ' \
'a string'.format(self._id, color)
raise ValueError(msg)
elif not color in ['cell', 'mat']:
elif color not in ['cell', 'mat']:
msg = 'Unable to create Plot ID={0} with color {1} which is not ' \
'a cell or mat'.format(self._id, color)
raise ValueError(msg)
self._color = color
@type.setter
def type(self, type):
if not is_string(type):
msg = 'Unable to create Plot ID={0} with type {1} which is not ' \
'a string'.format(self._id, type)
raise ValueError(msg)
elif not type in ['slice', 'voxel']:
elif type not in ['slice', 'voxel']:
msg = 'Unable to create Plot ID={0} with type {1} which is not ' \
'slice or voxel'.format(self._id, type)
raise ValueError(msg)
self._type = type
@basis.setter
def basis(self, basis):
if not is_string(basis):
msg = 'Unable to create Plot ID={0} with basis {1} which is not ' \
'a string'.format(self._id, basis)
raise ValueError(msg)
elif not basis in ['xy', 'xz', 'yz']:
elif basis not in ['xy', 'xz', 'yz']:
msg = 'Unable to create Plot ID={0} with basis {1} which is not ' \
'xy, xz, or yz'.format(self._id, basis)
raise ValueError(msg)
self._basis = basis
@background.setter
def background(self, background):
if not isinstance(background, (tuple, list, np.ndarray)):
msg = 'Unable to create Plot ID={0} with background {1} ' \
'which is not a Python tuple/list or NumPy ' \
@ -286,7 +294,6 @@ class Plot(object):
raise ValueError(msg)
for rgb in background:
if not is_integer(rgb):
msg = 'Unable to create Plot ID={0} with background RGB ' \
'value {1} which is not an integer'.format(self._id, rgb)
@ -299,18 +306,15 @@ class Plot(object):
self._background = background
@col_spec.setter
def col_spec(self, col_spec):
if not isinstance(col_spec, dict):
msg= 'Unable to create Plot ID={0} with col_spec parameter {1} ' \
'which is not a Python dictionary of IDs to ' \
'pixels'.format(self._id, col_spec)
msg = 'Unable to create Plot ID={0} with col_spec parameter {1} ' \
'which is not a Python dictionary of IDs to ' \
'pixels'.format(self._id, col_spec)
raise ValueError(msg)
for key in col_spec:
if not is_integer(key):
msg = 'Unable to create Plot ID={0} with col_spec ID {1} ' \
'which is not an integer'.format(self._id, key)
@ -335,10 +339,8 @@ class Plot(object):
self._col_spec = col_spec
@mask_componenets.setter
def mask_components(self, mask_components):
if not isinstance(mask_components, (list, tuple, np.ndarray)):
msg = 'Unable to create Plot ID={0} with mask components {1} ' \
'which is not a Python tuple/list or NumPy ' \
@ -358,10 +360,8 @@ class Plot(object):
self._mask_components = mask_components
@mask_background.setter
def mask_background(self, mask_background):
if not isinstance(mask_background, (list, tuple, np.ndarray)):
msg = 'Unable to create Plot ID={0} with mask background {1} ' \
'which is not a Python tuple/list or NumPy ' \
@ -389,9 +389,7 @@ class Plot(object):
self._mask_background = mask_background
def __repr__(self):
string = 'Plot\n'
string += '{0: <16}{1}{2}\n'.format('\tID', '=\t', self._id)
string += '{0: <16}{1}{2}\n'.format('\tName', '=\t', self._name)
@ -409,8 +407,15 @@ class Plot(object):
string += '{0: <16}{1}{2}\n'.format('\tCol Spec', '=\t', self._col_spec)
return string
def get_plot_xml(self):
"""Return XML representation of the plot
Returns
-------
element : xml.etree.ElementTree.Element
XML element containing plot data
"""
element = ET.Element("plot")
element.set("id", str(self._id))
@ -422,64 +427,55 @@ class Plot(object):
element.set("basis", self._basis)
subelement = ET.SubElement(element, "origin")
text = ''
for coord in self._origin:
text += str(coord) + ' '
subelement.text = text.rstrip(' ')
subelement.text = ' '.join(map(str, self._origin))
subelement = ET.SubElement(element, "width")
text = ''
for dim in self._width:
text += str(dim) + ' '
subelement.text = text.rstrip(' ')
subelement.text = ' '.join(map(str, self._width))
subelement = ET.SubElement(element, "pixels")
text = ''
for dim in self._pixels:
text += str(dim) + ' '
subelement.text = text.rstrip(' ')
subelement.text = ' '.join(map(str, self._pixels))
if not self._mask_background is None:
if self._mask_background is not None:
subelement = ET.SubElement(element, "background")
text = ''
for rgb in self._background:
text += str(rgb) + ' '
subelement.text = text.rstrip(' ')
if not self._col_spec is None:
subelement.text = ' '.join(map(str, self._background))
if self._col_spec is not None:
for key in self._col_spec:
subelement = ET.SubElement(element, "col_spec")
subelement.set("id", '{0}'.format(key))
value = self._col_spec[key]
subelement.set("rgb",'{0} {1} {2}'.format(value[0],
value[1], value[2]))
subelement.set("id", str(key))
subelement.set("rgb", ' '.join(map(
str, self._col_spec[key])))
if not self._mask_components is None:
if self._mask_components is not None:
subelement = ET.SubElement(element, "mask")
text = ''
for id in self._mask_components:
text += str(id) + ' '
subelement.set("components", text.rstrip(' '))
rgb = self._mask_background
subelement.set("background", '{0} {1} {2}'.format(rgb[0],
rgb[1], rgb[2]))
subelement.set("components", ' '.join(map(
str, self._mask_components)))
subelement.set("background", ' '.join(map(
str, self._mask_background)))
return element
class PlotsFile(object):
"""Plots file used for an OpenMC simulation. Corresponds directly to the
plots.xml input file.
"""
def __init__(self):
# Initialize PlotsFile class attributes
self._plots = []
self._plots_file = ET.Element("plots")
def add_plot(self, plot):
"""Add a plot to the file.
Parameters
----------
plot : Plot
Plot to add
"""
if not isinstance(plot, Plot):
msg = 'Unable to add a non-Plot {0} to the PlotsFile'.format(plot)
@ -487,15 +483,20 @@ class PlotsFile(object):
self._plots.append(plot)
def remove_plot(self, plot):
"""Remove a plot from the file.
Parameters
----------
plot : Plot
Plot to remove
"""
self._plots.remove(plot)
def create_plot_subelements(self):
def _create_plot_subelements(self):
for plot in self._plots:
xml_element = plot.get_plot_xml()
if len(plot._name) > 0:
@ -503,10 +504,12 @@ class PlotsFile(object):
self._plots_file.append(xml_element)
def export_to_xml(self):
"""Create a plots.xml file that can be used by OpenMC.
self.create_plot_subelements()
"""
self._create_plot_subelements()
# Clean the indentation in the file to be user-readable
clean_xml_indentation(self._plots_file)

File diff suppressed because it is too large Load diff

View file

@ -13,17 +13,28 @@ if sys.version > '3':
class SourceSite(object):
"""A single source site produced from fission.
Attributes
----------
weight : float
Weight of the particle arising from the site
xyz : list of float
Cartesian coordinates of the site
uvw : list of float
Directional cosines for particles emerging from the site
E : float
Energy of the emerging particle in MeV
"""
def __init__(self):
self._weight = None
self._xyz = None
self._uvw = None
self._E = None
def __repr__(self):
string = 'SourceSite\n'
string += '{0: <16}{1}{2}\n'.format('\tweight', '=\t', self._weight)
string += '{0: <16}{1}{2}\n'.format('\tE', '=\t', self._E)
@ -31,32 +42,56 @@ class SourceSite(object):
string += '{0: <16}{1}{2}\n'.format('\t(u,v,w)', '=\t', self._uvw)
return string
@property
def weight(self):
return self._weight
@property
def xyz(self):
return self._xyz
@property
def uvw(self):
return self._uvw
@property
def E(self):
return self._E
class StatePoint(object):
"""State information on a simulation at a certain point in time (at the end of a
given batch). Statepoints can be used to analyze tally results as well as
restart a simulation.
Attributes
----------
k_combined : list
Combined estimator for k-effective and its uncertainty
n_particles : int
Number of particles per generation
n_batches : int
Number of batches
current_batch :
Number of batches simulated
results : bool
Indicate whether tally results have been read
source : ndarray of SourceSite
Array of source sites
with_summary : bool
Indicate whether statepoint data has been linked against a summary file
tallies : dict
Dictionary whose keys are tally IDs and whose values are Tally objects
tallies_present : bool
Indicate whether user-defined tallies are present
global_tallies : ndarray
Global tallies and their uncertainties
n_realizations : int
Number of tally realizations
"""
def __init__(self, filename):
if filename.endswith('.h5'):
import h5py
self._f = h5py.File(filename, 'r')
@ -79,7 +114,6 @@ class StatePoint(object):
# Read tally metadata
self._read_tallies()
def close(self):
self._f.close()
@ -128,7 +162,6 @@ class StatePoint(object):
return self._n_realizations
def _read_metadata(self):
# Read filetype
self._filetype = self._get_int(path='filetype')[0]
@ -169,9 +202,7 @@ class StatePoint(object):
if self._run_mode == 2:
self._read_criticality()
def _read_criticality(self):
# Read criticality information
if self._run_mode == 2:
@ -191,9 +222,7 @@ class StatePoint(object):
# Read CMFD information (if used)
self._read_cmfd()
def _read_cmfd(self):
base = 'cmfd'
# Read CMFD information
@ -217,9 +246,7 @@ class StatePoint(object):
self._cmfd_srccmp = self._get_double(self._current_batch,
path='{0}/cmfd_srccmp'.format(base))
def _read_meshes(self):
# Initialize dictionaries for the Meshes
# Keys - Mesh IDs
# Values - Mesh objects
@ -282,9 +309,7 @@ class StatePoint(object):
# Add mesh to the global dictionary of all Meshes
self._meshes[mesh_id] = mesh
def _read_tallies(self):
# Initialize dictionaries for the Tallies
# Keys - Tally IDs
# Values - Tally objects
@ -423,7 +448,6 @@ class StatePoint(object):
# Add the scores to the Tally
for j, score in enumerate(scores):
# If this is a scattering moment, insert the scattering order
if '-n' in score:
score = score.replace('-n', '-' + str(moments[j]))
@ -437,10 +461,13 @@ class StatePoint(object):
# Add Tally to the global dictionary of all Tallies
self.tallies[tally_key] = tally
def read_results(self):
"""Read tally results and store them in the ``tallies`` attribute. No results
are read when the statepoint is instantiated.
# Number of realizations for global Tallies
"""
# Number of realizations for global Tallies
self._n_realizations = self._get_int(path='n_realizations')[0]
# Read global Tallies
@ -483,7 +510,8 @@ class StatePoint(object):
sum_sq = results[1::2]
# Define a routine to convert 0 to 1
nonzero = lambda val: 1 if not val else val
def nonzero(val):
return 1 if not val else val
# Reshape the results arrays
new_shape = (nonzero(tally.num_filter_bins),
@ -494,13 +522,17 @@ class StatePoint(object):
sum_sq = np.reshape(sum_sq, new_shape)
# Set the data for this Tally
tally.set_results(sum=sum, sum_sq=sum_sq)
tally.sum = sum
tally.sum_sq = sum_sq
# Indicate that Tally results have been read
self._results = True
def read_source(self):
"""Read and store source sites from the statepoint file. By default, source
sites are not loaded upon initialization.
"""
# Check whether Tally results have been read
if not self._results:
@ -520,7 +552,6 @@ class StatePoint(object):
# Initialize SourceSite object for each particle
for i in range(self._n_particles):
# Initialize new source site
site = SourceSite()
@ -536,9 +567,18 @@ class StatePoint(object):
# Store the source site in the NumPy array
self._source[i] = site
def compute_ci(self, confidence=0.95):
"""Computes confidence intervals for each Tally bin."""
"""Computes confidence intervals for each Tally bin.
This method is equivalent to calling compute_stdev(...) when the
confidence is known as opposed to its corresponding t value.
Parameters
----------
confidence : float, optional
Confidence level. Defaults to 0.95.
"""
# Determine significance level and percentile for two-sided CI
alpha = 1 - confidence
@ -548,11 +588,16 @@ class StatePoint(object):
t_value = scipy.stats.t.ppf(percentile, self._n_realizations - 1)
self.compute_stdev(t_value)
def compute_stdev(self, t_value=1.0):
"""
Computes the sample mean and the standard deviation of the mean
"""Computes the sample mean and the standard deviation of the mean
for each Tally bin.
Parameters
----------
t_value : float, optional
Student's t-value applied to the uncertainty. Defaults to 1.0,
meaning the reported value is the sample standard deviation.
"""
# Determine number of realizations
@ -572,49 +617,45 @@ class StatePoint(object):
if s != 0.0:
self._global_tallies[i, 1] = t_value * np.sqrt((s2 / n - s**2) / (n-1))
# Calculate sample mean and standard deviation for user-defined Tallies
for tally_id, tally in self.tallies.items():
tally.compute_std_dev(t_value)
def get_tally(self, scores=[], filters=[], nuclides=[],
name=None, id=None, estimator=None):
"""Finds and returns a Tally object with certain properties.
This routine searches the list of Tallies and returns the first Tally
found it finds which satisfieds all of the input parameters.
found it finds 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 Tallies properties.
specification and may only represent a subset of the Tally's properties.
Parameters
----------
scores : list
A list of one or more score strings (default is []).
filters : list
A list of Filter objects (default is []).
nuclides : list
A list of Nuclide objects (default is []).
name : str
The name specified for the Tally (default is None).
id : int
The id specified for the Tally (default is None).
estimator: str
The type of estimator ('tracklength', 'analog'; default is None).
scores : list, optional
A list of one or more score strings (default is []).
filters : list, optional
A list of Filter objects (default is []).
nuclides : list, optional
A list of Nuclide objects (default is []).
name : str, optional
The name specified for the Tally (default is None).
id : int, optional
The id specified for the Tally (default is None).
estimator: str, optional
The type of estimator ('tracklength', 'analog'; default is None).
Returns
-------
A Tally object.
tally : Tally
A tally matching the specified criteria
Raises
------
LookupError : An error when a Tally meeting all of the input
parameters cannot be found in the statepoint.
LookupError
If a Tally meeting all of the input parameters cannot be found in
the statepoint.
"""
tally = None
@ -640,7 +681,7 @@ class StatePoint(object):
# Iterate over the scores requested by the user
for score in scores:
if not score in test_tally.scores:
if score not in test_tally.scores:
contains_scores = False
break
@ -653,7 +694,7 @@ class StatePoint(object):
# Iterate over the Filters requested by the user
for filter in filters:
if not filter in test_tally.filters:
if filter not in test_tally.filters:
contains_filters = False
break
@ -666,7 +707,7 @@ class StatePoint(object):
# Iterate over the Nuclides requested by the user
for nuclide in nuclides:
if not nuclide in test_tally.nuclides:
if nuclide not in test_tally.nuclides:
contains_nuclides = False
break
@ -683,14 +724,15 @@ class StatePoint(object):
return tally
def link_with_summary(self, summary):
"""Links Tallies and Filters with Summary model information.
This routine retrieves model information (materials, geometry) from a
Summary object populated with an HDF5 'summary.h5' file and inserts
it into the Tally objects. This can be helpful when viewing and
manipulating large scale Tally data.
Summary object populated with an HDF5 'summary.h5' file and inserts it
into the Tally objects. This can be helpful when viewing and
manipulating large scale Tally data. Note that it is necessary to link
against a summary to populate the Tallies with any user-specified "name"
XML tags.
Parameters
----------
@ -699,8 +741,10 @@ class StatePoint(object):
Raises
------
ValueError : An error when the argument passed to the 'summary'
parameter is not an openmc.Summary object.
ValueError
An error when the argument passed to the 'summary' parameter is not
an openmc.Summary object.
"""
if not isinstance(summary, openmc.summary.Summary):
@ -709,7 +753,6 @@ class StatePoint(object):
raise ValueError(msg)
for tally_id, tally in self.tallies.items():
# Get the Tally name from the summary file
tally.name = summary.tallies[tally_id].name
tally.with_summary = True
@ -717,7 +760,6 @@ class StatePoint(object):
nuclide_zaids = copy.deepcopy(tally.nuclides)
for nuclide_zaid in nuclide_zaids:
tally.remove_nuclide(nuclide_zaid)
if nuclide_zaid == -1:
tally.add_nuclide(openmc.Nuclide('total'))
@ -725,7 +767,6 @@ class StatePoint(object):
tally.add_nuclide(summary.nuclides[nuclide_zaid])
for filter in tally.filters:
if filter.type == 'surface':
surface_ids = []
for bin in filter.bins:
@ -749,12 +790,11 @@ class StatePoint(object):
for bin in filter.bins:
material_ids.append(summary.materials[bin].id)
filter.bins = material_ids
self._with_summary = True
def _get_data(self, n, typeCode, size):
return list(struct.unpack('={0}{1}'.format(n,typeCode),
return list(struct.unpack('={0}{1}'.format(n, typeCode),
self._f.read(n*size)))
def _get_int(self, n=1, path=None):

View file

@ -2,16 +2,18 @@ import numpy as np
import openmc
try:
import h5py
except ImportError:
msg = 'Unable to import h5py which is needed by openmc.summary'
raise ImportError(msg)
class Summary(object):
"""Information summarizing the geometry, materials, and tallies used in a
simulation.
"""
def __init__(self, filename):
# A user may not have h5py, but they can still use the rest of the
# Python API so we'll only try to import h5py if the user actually inits
# a Summary object.
import h5py
openmc.reset_auto_ids()
@ -27,9 +29,7 @@ class Summary(object):
self._read_geometry()
self._read_tallies()
def _read_metadata(self):
# Read OpenMC version
self.version = [self._f['version_major'][0],
self._f['version_minor'][0],
@ -44,9 +44,7 @@ class Summary(object):
self.gen_per_batch = self._f['gen_per_batch'][0]
self.n_procs = self._f['n_procs'][0]
def _read_geometry(self):
# Read in and initialize the Materials and Geometry
self._read_nuclides()
self._read_materials()
@ -56,9 +54,7 @@ class Summary(object):
self._read_lattices()
self._finalize_geometry()
def _read_nuclides(self):
self.n_nuclides = self._f['nuclides/n_nuclides'][0]
# Initialize dictionary for each Nuclide
@ -67,7 +63,6 @@ class Summary(object):
self.nuclides = {}
for key in self._f['nuclides'].keys():
if key == 'n_nuclides':
continue
@ -88,9 +83,7 @@ class Summary(object):
self.nuclides[zaid] = openmc.Nuclide(name=name, xs=xs)
self.nuclides[zaid].zaid = zaid
def _read_materials(self):
self.n_materials = self._f['materials/n_materials'][0]
# Initialize dictionary for each Material
@ -99,7 +92,6 @@ class Summary(object):
self.materials = {}
for key in self._f['materials'].keys():
if key == 'n_materials':
continue
@ -116,7 +108,6 @@ class Summary(object):
# Read the names of the S(a,b) tables for this Material
for i in range(1, n_sab+1):
sab_table = self._f['materials'][key]['sab_tables'][str(i)][0]
# Read the cross-section identifiers for each S(a,b) table
@ -148,9 +139,7 @@ class Summary(object):
# Add the Material to the global dictionary of all Materials
self.materials[index] = material
def _read_surfaces(self):
self.n_surfaces = self._f['geometry/n_surfaces'][0]
# Initialize dictionary for each Surface
@ -159,7 +148,6 @@ class Summary(object):
self.surfaces = {}
for key in self._f['geometry/surfaces'].keys():
if key == 'n_surfaces':
continue
@ -171,7 +159,6 @@ class Summary(object):
coeffs = self._f['geometry/surfaces'][key]['coefficients'][...]
# Create the Surface based on its type
if surf_type == 'X Plane':
x0 = coeffs[0]
surface = openmc.XPlane(surface_id, bc, x0, name)
@ -232,9 +219,7 @@ class Summary(object):
# Add Surface to global dictionary of all Surfaces
self.surfaces[index] = surface
def _read_cells(self):
self.n_cells = self._f['geometry/n_cells'][0]
# Initialize dictionary for each Cell
@ -251,7 +236,6 @@ class Summary(object):
self._cell_fills = {}
for key in self._f['geometry/cells'].keys():
if key == 'n_cells':
continue
@ -310,9 +294,7 @@ class Summary(object):
# Add the Cell to the global dictionary of all Cells
self.cells[index] = cell
def _read_universes(self):
self.n_universes = self._f['geometry/n_universes'][0]
# Initialize dictionary for each Universe
@ -321,7 +303,6 @@ class Summary(object):
self.universes = {}
for key in self._f['geometry/universes'].keys():
if key == 'n_universes':
continue
@ -340,9 +321,7 @@ class Summary(object):
# Add the Universe to the global list of Universes
self.universes[index] = universe
def _read_lattices(self):
self.n_lattices = self._f['geometry/n_lattices'][0]
# Initialize lattices for each Lattice
@ -351,7 +330,6 @@ class Summary(object):
self.lattices = {}
for key in self._f['geometry/lattices'].keys():
if key == 'n_lattices':
continue
@ -383,7 +361,6 @@ class Summary(object):
lattice.lower_left = lower_left
lattice.pitch = pitch
# If the Universe specified outer the Lattice is not void (-22)
if outer != -22:
lattice.outer = self.universes[outer]
@ -395,14 +372,14 @@ class Summary(object):
for x in range(universe_ids.shape[0]):
for y in range(universe_ids.shape[1]):
for z in range(universe_ids.shape[2]):
universes[x,y,z] = \
self.get_universe_by_id(universe_ids[x,y,z])
universes[x, y, z] = \
self.get_universe_by_id(universe_ids[x, y, z])
# Transpose, reverse y-dimension for appropriate ordering
shape = universes.shape
universes = np.transpose(universes, (1,0,2))
universes = np.transpose(universes, (1, 0, 2))
universes.shape = shape
universes = universes[:,::-1,:]
universes = universes[:, ::-1, :]
lattice.universes = universes
if offset_size > 0:
@ -420,8 +397,8 @@ class Summary(object):
pitch = self._f['geometry/lattices'][key]['pitch'][...]
outer = self._f['geometry/lattices'][key]['outer'][0]
universe_ids = \
self._f['geometry/lattices'][key]['universes'][...]
universe_ids = self._f[
'geometry/lattices'][key]['universes'][...]
# Create the Lattice
lattice = openmc.HexLattice(lattice_id=lattice_id, name=name)
@ -441,9 +418,9 @@ class Summary(object):
for i in range(universe_ids.shape[0]):
for j in range(universe_ids.shape[1]):
for k in range(universe_ids.shape[2]):
if universe_ids[i,j,k] != -1:
universes[i,j,k] = \
self.get_universe_by_id(universe_ids[i,j,k])
if universe_ids[i, j, k] != -1:
universes[i, j, k] = self.get_universe_by_id(
universe_ids[i, j, k])
if offset_size > 0:
lattice.offsets = offsets
@ -451,15 +428,12 @@ class Summary(object):
# Add the Lattice to the global dictionary of all Lattices
self.lattices[index] = lattice
def _finalize_geometry(self):
# Initialize Geometry object
self.openmc_geometry = openmc.Geometry()
# Iterate over all Cells and add fill Materials, Universes and Lattices
for cell_key in self._cell_fills.keys():
# Determine fill type ('normal', 'universe', or 'lattice') and ID
fill_type = self._cell_fills[cell_key][0]
fill_id = self._cell_fills[cell_key][1]
@ -482,16 +456,14 @@ class Summary(object):
root_universe = self.get_universe_by_id(0)
self.openmc_geometry.root_universe = root_universe
def _read_tallies(self):
# Initialize dictionaries for the Tallies
# Keys - Tally IDs
# Values - Tally objects
self.tallies = {}
# Read the number of tallies
if not 'tallies' in self._f.keys():
if 'tallies' not in self._f:
self.n_tallies = 0
return
@ -554,8 +526,11 @@ class Summary(object):
# Add Tally to the global dictionary of all Tallies
self.tallies[tally_id] = tally
def make_opencg_geometry(self):
"""Create OpenCG geometry based on the information contained in the summary
file. The geometry is stored as the 'opencg_geometry' attribute.
"""
try:
from openmc.opencg_compatible import get_opencg_geometry
@ -567,8 +542,22 @@ class Summary(object):
if self.opencg_geometry is None:
self.opencg_geometry = get_opencg_geometry(self.openmc_geometry)
def get_nuclide_by_zaid(self, zaid):
"""Return a Nuclide object given the 'zaid' identifier for the nuclide.
Parameters
----------
zaid : int
1000*Z + A, where Z is the atomic number of the nuclide and A is the
mass number. For example, the zaid for U-235 is 92235.
Returns
-------
nuclide : openmc.nuclide.Nuclide or None
Nuclide matching the specified zaid, or None if no matching object
is found.
"""
for index, nuclide in self.nuclides.items():
if nuclide._zaid == zaid:
@ -576,8 +565,20 @@ class Summary(object):
return None
def get_material_by_id(self, material_id):
"""Return a Material object given the material id
Parameters
----------
id : int
Unique identifier for the material
Returns
-------
material : openmc.material.Material
Material with given id
"""
for index, material in self.materials.items():
if material._id == material_id:
@ -585,8 +586,20 @@ class Summary(object):
return None
def get_surface_by_id(self, surface_id):
"""Return a Surface object given the surface id
Parameters
----------
id : int
Unique identifier for the surface
Returns
-------
surface : openmc.surface.Surface
Surface with given id
"""
for index, surface in self.surfaces.items():
if surface._id == surface_id:
@ -594,8 +607,20 @@ class Summary(object):
return None
def get_cell_by_id(self, cell_id):
"""Return a Cell object given the cell id
Parameters
----------
id : int
Unique identifier for the cell
Returns
-------
cell : openmc.universe.Cell
Cell with given id
"""
for index, cell in self.cells.items():
if cell._id == cell_id:
@ -603,8 +628,20 @@ class Summary(object):
return None
def get_universe_by_id(self, universe_id):
"""Return a Universe object given the universe id
Parameters
----------
id : int
Unique identifier for the universe
Returns
-------
universe : openmc.universe.Universe
Universe with given id
"""
for index, universe in self.universes.items():
if universe._id == universe_id:
@ -612,8 +649,20 @@ class Summary(object):
return None
def get_lattice_by_id(self, lattice_id):
"""Return a Lattice object given the lattice id
Parameters
----------
id : int
Unique identifier for the lattice
Returns
-------
lattice : openmc.universe.Lattice
Lattice with given id
"""
for index, lattice in self.lattices.items():
if lattice._id == lattice_id:

File diff suppressed because it is too large Load diff

File diff suppressed because it is too large Load diff

View file

@ -4,22 +4,39 @@ from openmc.checkvalue import *
class Trigger(object):
"""A criterion for when to finish a simulation based on tally uncertainties.
Parameters
----------
trigger_type : {'variance', 'std_dev', 'rel_err'}
Determine whether to trigger on the variance, standard deviation, or
relative error of scores.
threshold : float
The threshold for the trigger type.
Attributes
----------
trigger_type : {'variance', 'std_dev', 'rel_err'}
Determine whether to trigger on the variance, standard deviation, or
relative error of scores.
threshold : float
The threshold for the trigger type.
scores : list of str
Scores which should be checked against the trigger
"""
def __init__(self, trigger_type, threshold):
# Initialize Mesh class attributes
self.trigger_type = trigger_type
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
@ -36,36 +53,29 @@ class Trigger(object):
else:
return existing
@property
def trigger_type(self):
return self._trigger_type
@property
def threshold(self):
return self._threshold
@property
def scores(self):
return self._scores
@trigger_type.setter
def trigger_type(self, trigger_type):
if not trigger_type in ['variance', 'std_dev', 'rel_err']:
if trigger_type not in ['variance', 'std_dev', 'rel_err']:
msg = 'Unable to create a tally trigger with ' \
'type "{0}"'.format(trigger_type)
raise ValueError(msg)
self._trigger_type = trigger_type
@threshold.setter
def threshold(self, threshold):
if not is_float(threshold):
msg = 'Unable to set a tally trigger threshold with ' \
'threshold "{0}"'.format(threshold)
@ -73,8 +83,15 @@ class Trigger(object):
self._threshold = threshold
def add_score(self, score):
"""Add a score to the list of scores to be checked against the trigger.
Parameters
----------
score : str
Score to append
"""
if not is_string(score):
msg = 'Unable to add score "{0}" to tally trigger since ' \
@ -87,28 +104,25 @@ class Trigger(object):
else:
self._scores.append(score)
def __repr__(self):
string = 'Trigger\n'
string += '{0: <16}{1}{2}\n'.format('\tType', '=\t', self._trigger_type)
string += '{0: <16}{1}{2}\n'.format('\tThreshold', '=\t', self._threshold)
string += '{0: <16}{1}{2}\n'.format('\tScores', '=\t', self._scores)
return string
def get_trigger_xml(self, element):
"""Return XML representation of the trigger
Returns
-------
element : xml.etree.ElementTree.Element
XML element containing trigger data
"""
subelement = ET.SubElement(element, "trigger")
subelement.set("type", self._trigger_type)
subelement.set("threshold", str(self._threshold))
# Scores
if len(self._scores) != 0:
scores = ''
for score in self._scores:
scores += '{0} '.format(score)
scores.rstrip(' ')
subelement.set("scores", scores)
subelement.set("scores", ' '.join(map(str, self._scores)))

File diff suppressed because it is too large Load diff