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

merged with origin/mg-mode-delayed

Browse source
This commit is contained in:
Sam Shaner 2016-10-16 09:55:07 -04:00
commit 808787f45b
30 changed files with 1430 additions and 1348 deletions
Download patch file Download diff file Expand all files Collapse all files

2
.gitignore vendored
View file

@ -50,8 +50,6 @@ tests/ctestscript.run
# HDF5 files
*.h5
!tests/1d_mgxs.h5
!examples/xml/pincell_multigroup/mgxs.h5
# Build files
src/CMakeCache.txt

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

@ -6,7 +6,7 @@ Multi-Group Cross Section Library Format
OpenMC can be run in continuous-energy mode or multi-group mode, provided the
nuclear data is available. In continuous-energy mode, the
``cross_sections.xml`` file contains necessary meta-data for each data set,
``cross_sections.xml`` file contains necessary meta-data for each dataset,
including the name and a file system location where the complete library
can be found. In multi-group mode, the multi-group meta-data and the
nuclear data itself is contained within an ``mgxs.h5``. This portion of
@ -40,8 +40,8 @@ MGXS Library Specification
The data within <library name> contains the temperature-dependent multi-group
data for the nuclide or material that it represents.
:Attributes: - **awr** (*double*) -- The atomic weight ratio (optional, i.e. it
is not meaningful for material-wise data)
:Attributes: - **atomic_weight_ratio** (*double*) -- The atomic weight ratio (optional,
i.e. it is not meaningful for material-wise data)
- **fissionable** (*int*) -- Whether the dataset is fissionable
(1) or not (0).
- **representation** (*char[]*) -- The method used to generate and
@ -51,26 +51,27 @@ data for the nuclide or material that it represents.
data was generated with angular dependent fluxes and thus the
data is angle-dependent. Valid values are either "isotropic" or
"angle".
- **num-azimuthal** (*int*) -- Number of equal width angular bins
- **num_azimuthal** (*int*) -- Number of equal width angular bins
that the azimuthal angular domain is subdivided if the
`representation` attribute is "angle". This parameter is
ignored otherwise.
- **num-polar** (*int*) -- Number of equal width angular bins
- **num_polar** (*int*) -- Number of equal width angular bins
that the polar angular domain is subdivided if the
`representation` attribute is "angle". This parameter is
ignored otherwise.
- **scatter-format** (*char[]*) -- The representation of the
- **scatter_format** (*char[]*) -- The representation of the
scattering angular distribution. The options are either
"legendre", "histogram", or "tabular". If not provided, the
default of "legendre" will be assumed.
- **order** (*int*) -- Either the Legendre order, number of bins,
or number of points (depending on the value of `scatter-format`)
or number of points (depending on the value of `scatter_format`)
used to describe the angular distribution associated with each
group-to-group transfer probability.
**/<library name>/kTs/**
:Datasets: - **<TTT>K** (*double*) -- kT values (in MeV) for each Temperature
:Datasets:
- **<TTT>K** (*double*) -- kT values (in MeV) for each Temperature
TTT (in Kelvin), rounded to the nearest integer
**/<library name>/<TTT>K/**
@ -90,19 +91,19 @@ Temperature-dependent data, provided for temperature <TTT>K.
cross section.
This is a 1-D vector if `representation` is "isotropic", or a 3-D
vector if `representation` is "angle" with dimensions of
[groups, azimuthal, polar]. This is only required if the data set
[groups, azimuthal, polar]. This is only required if the dataset
is fissionable and fission-tallies are expected to be used.
- **kappa-fission** (*double[]* or *double[][][]*) -- Kappa-Fission
(energy-release from fission) cross section.
This is a 1-D vector if `representation` is "isotropic", or a 3-D
vector if `representation` is "angle" with dimensions of
[groups, azimuthal, polar]. This is only required if the data set
[groups, azimuthal, polar]. This is only required if the dataset
is fissionable and fission-tallies are expected to be used.
- **chi** (*double[]* or *double[][][]*) -- Fission neutron energy
spectra.
This is a 1-D vector if `representation` is "isotropic", or a 3-D
vector if `representation` is "angle" with dimensions of
[groups, azimuthal, polar]. This is only required if the data set
[groups, azimuthal, polar]. This is only required if the dataset
is fissionable and fission-tallies are expected to be used.
- **nu-fission** (*double[]* to *double[][][][]*) -- Nu-Fission
cross section.
@ -115,40 +116,40 @@ Temperature-dependent data, provided for temperature <TTT>K.
- **inverse-velocity** (*double[]*) -- Average inverse velocity
for each of the groups in the library. This dataset is optional.
**/<library name>/<TTT>K/scatter data/**
**/<library name>/<TTT>K/scatter_data/**
Data specific to neutron scattering for the temperature <TTT>K
:Datasets: - **g_min** (*int[]* or *int[][][]) --
:Datasets: - **g_min** (*int[]* or *int[][][]*) --
Minimum (most energetic) outgoing groups with non-zero values of
the scattering matrix. These group numbers use the standard
ordering where the fastest neutron energy group is group 1 while
the slowest neutron energy group is group G.
The dimensionality of `g_out bounds` is:
`g_min[g_in]`, or `g_min[num-polar][num-azimuthal][g_in]`.
The dimensionality of `g_min` is:
`g_min[g_in]`, or `g_min[num_polar][num_azimuthal][g_in]`.
The former is used when `representation` is "isotropic", and the
latter when `representation` is "angle".
- **g_max** (*int[]* or *int[][][]) --
- **g_max** (*int[]* or *int[][][]*) --
Maximum (least energetic) outgoing groups with non-zero values of
the scattering matrix. These group numbers use the standard
ordering where the fastest neutron energy group is group 1 while
the slowest neutron energy group is group G.
The dimensionality of `g_out bounds` is:
`g_max[g_in]`, or `g_max[num-polar][num-azimuthal][g_in]`.
The dimensionality of `g_max` is:
`g_max[g_in]`, or `g_max[num_polar][num_azimuthal][g_in]`.
The former is used when `representation` is "isotropic", and the
latter when `representation` is "angle".
- **scatter matrix** (*double[]*) -- Flattened representation of the
- **scatter_matrix** (*double[]*) -- Flattened representation of the
scattering moment matrices. The pre-flattened array is shaped as
follows (in row-major format):
`scatter matrix[order(+1)][g_in][g_out]`, or
`scatter matrix[num-polar][num-azimuthal][order(+1)][g_in][g_out]`
`scatter_matrix[order(+1)][g_in][g_out]`, or
`scatter_matrix[num_polar][num_azimuthal][order(+1)][g_in][g_out]`
The former is used when `representation` is "isotropic", and the
latter when `representation` is "angle". Note that if the value of
`scatter-format` is "legendre", the order dimension will be one
`scatter_format` is "legendre", the order dimension will be one
larger than the value of `order`, otherwise it will match `order`.
Finally, the g_out dimension has a dimensionality of
`g_out bounds`[0] to `g_out bounds`[1].
- **multiplicity matrix** (*double[]*) -- Flattened representation of
`g_min` to `g_max`.
- **multiplicity_matrix** (*double[]*) -- Flattened representation of
the scattering moment matrices. This dataset provides the code with
a scaling factor to account for neutrons being produced in (n,xn)
reactions. This is assumed isotropic and therefore is not repeated
@ -156,9 +157,8 @@ Data specific to neutron scattering for the temperature <TTT>K
optional, if it is not provided no multiplication (i.e., values of
1.0) will be assumed.
The pre-flattened array is shaped as follows (in row-major format):
`multiplicity matrix[g_in][g_out]`, or
`multiplicity matrix[num-polar][num-azimuthal][g_in][g_out]`
`multiplicity_matrix[g_in][g_out]`, or
`multiplicity_matrix[num_polar][num_azimuthal][g_in][g_out]`
The former is used when `representation` is "isotropic", and the
latter when `representation` is "angle". Finally, the g_out
dimension has a dimensionality of `g_out bounds`[0] to
`g_out bounds`[1].
dimension has a dimensionality of `g_min` to `g_max`.

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

@ -121,10 +121,13 @@ multi-group mode.
``<cutoff>`` Element
--------------------
The ``<cutoff>`` element indicates the weight cutoff used below which particles
undergo Russian roulette. Surviving particles are assigned a user-determined
weight. Note that weight cutoffs and Russian rouletting are not turned on by
default. This element has the following attributes/sub-elements:
The ``<cutoff>`` element indicates two kinds of cutoffs. The first is the weight
cutoff used below which particles undergo Russian roulette. Surviving particles
are assigned a user-determined weight. Note that weight cutoffs and Russian
rouletting are not turned on by default. The second is the energy cutoff which
is used to kill particles under certain energy. The energy cutoff should not be
used unless you know particles under the energy are of no importance to results
you care. This element has the following attributes/sub-elements:
:weight:
The weight below which particles undergo Russian roulette.
@ -137,6 +140,11 @@ default. This element has the following attributes/sub-elements:
*Default*: 1.0
:energy:
The energy under which particles will be killed.
*Default*: 0.0
.. _eigenvalue:
``<eigenvalue>`` Element

4
examples/xml/basic/materials.xml
View file

@ -8,8 +8,8 @@
<material id="41">
<density value="1.0" units="g/cc" />
<nuclide name="H1" ao="2.0" />
<nuclide name="O16" ao="1.0" />
<nuclide name="H1" wo="2.0" />
<nuclide name="O16" wo="1.0" />
<sab name="c_H_in_H2O"/>
</material>

BIN
examples/xml/pincell_multigroup/mgxs.h5
View file

Binary file not shown.

6
openmc/mesh.py
View file

@ -298,12 +298,12 @@ class Mesh(object):
# Would prefer to have the z ranges be the max supported float, but
# these values are apparently different between python and Fortran.
# Choosing a safe and sane default.
# Values of +/-1000 are used here as there seems to be an
# Values of +/-1e10 are used here as there seems to be an
# inconsistency between what numpy uses as the max float and what
# Fortran expects for a real(8), so this avoids code complication
# and achieves the same goal.
zplanes = [openmc.ZPlane(z0=-1000., boundary_type='reflective'),
openmc.ZPlane(z0=1000., boundary_type='reflective')]
zplanes = [openmc.ZPlane(z0=-1e10, boundary_type='reflective'),
openmc.ZPlane(z0=1e10, boundary_type='reflective')]
else:
zplanes = [openmc.ZPlane(z0=self.lower_left[2],
boundary_type=bc[4]),

99
openmc/mgxs/library.py
View file

@ -186,14 +186,14 @@ class Library(object):
@property
def domains(self):
if self._domains == 'all':
if self.domain_type == 'material':
if self._domains is 'all':
if self.domain_type is 'material':
return self.openmc_geometry.get_all_materials()
elif self.domain_type in ['cell', 'distribcell']:
return self.openmc_geometry.get_all_material_cells()
elif self.domain_type == 'universe':
elif self.domain_type is 'universe':
return self.openmc_geometry.get_all_universes()
elif self.domain_type == 'mesh':
elif self.domain_type is 'mesh':
raise ValueError('Unable to get domains for Mesh domain type')
else:
raise ValueError('Unable to get domains without a domain type')
@ -265,7 +265,7 @@ class Library(object):
@mgxs_types.setter
def mgxs_types(self, mgxs_types):
all_mgxs_types = openmc.mgxs.MGXS_TYPES + openmc.mgxs.MDGXS_TYPES
if mgxs_types == 'all':
if mgxs_types is 'all':
self._mgxs_types = all_mgxs_types
else:
cv.check_iterable_type('mgxs_types', mgxs_types, basestring)
@ -277,7 +277,7 @@ class Library(object):
def by_nuclide(self, by_nuclide):
cv.check_type('by_nuclide', by_nuclide, bool)
if by_nuclide == True and self.domain_type == 'mesh':
if by_nuclide == True and self.domain_type is 'mesh':
raise ValueError('Unable to create MGXS library by nuclide with '
'mesh domain')
@ -287,7 +287,7 @@ class Library(object):
def domain_type(self, domain_type):
cv.check_value('domain type', domain_type, openmc.mgxs.DOMAIN_TYPES)
if self.by_nuclide == True and domain_type == 'mesh':
if self.by_nuclide == True and domain_type is 'mesh':
raise ValueError('Unable to create MGXS library by nuclide with '
'mesh domain')
@ -297,21 +297,21 @@ class Library(object):
def domains(self, domains):
# Use all materials, cells or universes in the geometry as domains
if domains == 'all':
if domains is 'all':
self._domains = domains
# User specified a list of material, cell or universe domains
else:
if self.domain_type == 'material':
if self.domain_type is 'material':
cv.check_iterable_type('domain', domains, openmc.Material)
all_domains = self.openmc_geometry.get_all_materials()
elif self.domain_type in ['cell', 'distribcell']:
cv.check_iterable_type('domain', domains, openmc.Cell)
all_domains = self.openmc_geometry.get_all_material_cells()
elif self.domain_type == 'universe':
elif self.domain_type is 'universe':
cv.check_iterable_type('domain', domains, openmc.Universe)
all_domains = self.openmc_geometry.get_all_universes()
elif self.domain_type == 'mesh':
elif self.domain_type is 'mesh':
cv.check_iterable_type('domain', domains, openmc.Mesh)
# The mesh and geometry are independent, so set all_domains
@ -355,7 +355,7 @@ class Library(object):
def correction(self, correction):
cv.check_value('correction', correction, ('P0', None))
if correction == 'P0' and self.legendre_order > 0:
if correction is 'P0' and self.legendre_order > 0:
warn('The P0 correction will be ignored since the scattering '
'order "{}" is greater than zero'.format(self.legendre_order))
@ -367,7 +367,7 @@ class Library(object):
cv.check_greater_than('legendre_order', legendre_order, 0, equality=True)
cv.check_less_than('legendre_order', legendre_order, 10, equality=True)
if self.correction == 'P0' and legendre_order > 0:
if self.correction is 'P0' and legendre_order > 0:
msg = 'The P0 correction will be ignored since the scattering ' \
'order {} is greater than zero'.format(self.legendre_order)
warn(msg, RuntimeWarning)
@ -505,7 +505,7 @@ class Library(object):
self._openmc_geometry = statepoint.summary.openmc_geometry
self._nuclides = statepoint.summary.nuclides
if statepoint.run_mode == 'k-eigenvalue':
if statepoint.run_mode is 'k-eigenvalue':
self._keff = statepoint.k_combined[0]
# Load tallies for each MGXS for each domain and mgxs type
@ -543,13 +543,13 @@ class Library(object):
"""
if self.domain_type == 'material':
if self.domain_type is 'material':
cv.check_type('domain', domain, (openmc.Material, Integral))
elif self.domain_type == 'cell' or self.domain_type == 'distribcell':
elif self.domain_type is 'cell' or self.domain_type is 'distribcell':
cv.check_type('domain', domain, (openmc.Cell, Integral))
elif self.domain_type == 'universe':
elif self.domain_type is 'universe':
cv.check_type('domain', domain, (openmc.Universe, Integral))
elif self.domain_type == 'mesh':
elif self.domain_type is 'mesh':
cv.check_type('domain', domain, (openmc.Mesh, Integral))
# Check that requested domain is included in library
@ -662,7 +662,7 @@ class Library(object):
# Clone this Library to initialize the subdomain-averaged version
subdomain_avg_library = copy.deepcopy(self)
if subdomain_avg_library.domain_type == 'distribcell':
if subdomain_avg_library.domain_type is 'distribcell':
subdomain_avg_library.domain_type = 'cell'
else:
return subdomain_avg_library
@ -671,7 +671,7 @@ class Library(object):
for domain in self.domains:
for mgxs_type in self.mgxs_types:
mgxs = subdomain_avg_library.get_mgxs(domain, mgxs_type)
if mgxs.domain_type == 'distribcell':
if mgxs.domain_type is 'distribcell':
avg_mgxs = mgxs.get_subdomain_avg_xs()
subdomain_avg_library.all_mgxs[domain.id][mgxs_type] = avg_mgxs
@ -751,7 +751,7 @@ class Library(object):
for mgxs_type in self.mgxs_types:
mgxs = self.all_mgxs[domain.id][mgxs_type]
if subdomains == 'avg':
if subdomains is 'avg':
mgxs = mgxs.get_subdomain_avg_xs()
mgxs.build_hdf5_store(filename, directory, xs_type=xs_type,
@ -824,8 +824,9 @@ class Library(object):
def get_xsdata(self, domain, xsdata_name, nuclide='total', xs_type='macro',
order=None, subdomain=None):
"""Generates an openmc.XSdata object describing a multi-group cross section
data set for eventual combination in to an openmc.MGXSLibrary object
(i.e., the library). Note that this method does not build an XSdata
dataset for writing to an openmc.MGXSLibrary object.
Note that this method does not build an XSdata
object with nested temperature tables. The temperature of each
XSdata object will be left at the default value of 300K.
@ -856,7 +857,7 @@ class Library(object):
Returns
-------
xsdata : openmc.XSdata
Multi-Group Cross Section data set object.
Multi-Group Cross Section dataset object.
Raises
------
@ -913,7 +914,7 @@ class Library(object):
self.representation = 'isotropic'
if nuclide is not 'total':
xsdata.awr = self._nuclides[nuclide][1]
xsdata.atomic_weight_ratio = self._nuclides[nuclide][1]
if subdomain is None:
subdomain = 'all'
@ -921,10 +922,10 @@ class Library(object):
subdomain = [subdomain]
# Now get xs data itself
if 'nu-transport' in self.mgxs_types and self.correction == 'P0':
if 'nu-transport' in self.mgxs_types and self.correction is 'P0':
mymgxs = self.get_mgxs(domain, 'nu-transport')
xsdata.set_total_mgxs(mymgxs, xs_type=xs_type, nuclide=[nuclide],
subdomains=subdomain)
subdomain=subdomain)
elif 'total' in self.mgxs_types:
mymgxs = self.get_mgxs(domain, 'total')
@ -1044,10 +1045,11 @@ class Library(object):
nuclide=[nuclide],
subdomain=subdomain)
# Since we are not using multiplicity, then scattering
# multiplication (nu-scatter) must be accounted for
# approximately by using an adjusted absorption cross section.
if 'total' in self.mgxs_types:
# Since we are not using multiplicity, then
# scattering multiplication (nu-scatter) must be
# accounted for approximately by using an adjusted
# absorption cross section.
if 'total' in self.mgxs_types or 'transport' in self.mgxs_types:
for i in range(len(xsdata.temperatures)):
xsdata._absorption[i] \
= np.subtract(xsdata._total[i], np.sum(
@ -1057,9 +1059,11 @@ class Library(object):
def create_mg_library(self, xs_type='macro', xsdata_names=None):
"""Creates an openmc.MGXSLibrary object to contain the MGXS data for the
Multi-Group mode of OpenMC. Note that this library will not make use
of nested temperature tables. Every dataset in the library will be
treated as if it was at the same default temperature.
Multi-Group mode of OpenMC.
Note that this library will not make use of nested temperature tables.
Every dataset in the library will be treated as if it was at the same
default temperature.
Parameters
----------
@ -1106,7 +1110,7 @@ class Library(object):
mgxs_file = openmc.MGXSLibrary(self.energy_groups,
delayed_groups=self.num_delayed_groups)
if self.domain_type == 'mesh':
if self.domain_type is 'mesh':
# Create the xsdata objects and add to the mgxs_file
i = 0
for domain in self.domains:
@ -1152,14 +1156,15 @@ class Library(object):
def create_mg_mode(self, xsdata_names=None, bc=['reflective'] * 6):
"""Creates an openmc.MGXSLibrary object to contain the MGXS data for the
Multi-Group mode of OpenMC as well as the associated openmc.Materials
and openmc.Geometry objects. The created Geometry is the same as that
used to generate the MGXS data, with the only differences being
modifications to point to newly-created Materials which point to the
multi-group data. This method only creates a macroscopic
MGXS Library even if nuclidic tallies are specified in the Library.
Note that this library will not make use of nested temperature tables.
Every dataset in the library will be treated as if it was at the same
default temperature.
and openmc.Geometry objects.
The created Geometry is the same as that used to generate the MGXS
data, with the only differences being modifications to point to
newly-created Materials which point to the multi-group data. This
method only creates a macroscopic MGXS Library even if nuclidic tallies
are specified in the Library. Note that this library will not make
use of nested temperature tables. Every dataset in the library will
be treated as if it was at the same default temperature.
Parameters
----------
@ -1208,14 +1213,14 @@ class Library(object):
# the multiple meshes could be overlapping or in disparate regions
# of the continuous energy model. The next step makes sure there is
# only one before continuing.
if self.domain_type == 'mesh':
if self.domain_type is 'mesh':
cv.check_length("domains", self.domains, 1, 1)
# Get the MGXS File Data
mgxs_file = self.create_mg_library('macro', xsdata_names)
# Now move on the creating the geometry and assigning materials
if self.domain_type == 'mesh':
if self.domain_type is 'mesh':
root = openmc.Universe(name='root', universe_id=0)
# Add cells representative of the mesh with reflective BC
@ -1261,13 +1266,13 @@ class Library(object):
materials.append(material)
# Differentiate Geometry with new Material
if self.domain_type == 'material':
if self.domain_type is 'material':
# Fill all appropriate Cells with new Material
for cell in all_cells:
if cell.fill.id == domain.id:
cell.fill = material
elif self.domain_type == 'cell':
elif self.domain_type is 'cell':
for cell in all_cells:
if cell.id == domain.id:
cell.fill = material

329
openmc/mgxs_library.py
View file

@ -31,7 +31,7 @@ class XSdata(object):
representation : {'isotropic', 'angle'}, optional
Method used in generating the MGXS (isotropic or angle-dependent flux
weighting). Defaults to 'isotropic'
temperatures : numpy.ndarray
temperatures : Iterable of float
Temperatures (in units of Kelvin) of the provided datasets. Defaults
to a single temperature at 294K.
delayed_groups : int
@ -41,7 +41,7 @@ class XSdata(object):
----------
name : str
Unique identifier for the xsdata object
awr : float
aromic_weight_ratio : float
Atomic weight ratio of an isotope. That is, the ratio of the mass
of the isotope to the mass of a single neutron.
temperatures : numpy.ndarray
@ -64,7 +64,8 @@ class XSdata(object):
weighting).
num_azimuthal : int
Number of equal width angular bins that the azimuthal angular domain is
subdivided into. This only applies when ``representation`` is "angle".
subdivided into. This only applies when :attr:`XSdata.representation`
is "angle".
num_polar : int
Number of equal width angular bins that the polar angular domain is
subdivided into. This only applies when ``representation`` is "angle".
@ -84,7 +85,7 @@ class XSdata(object):
fission : dict of numpy.ndarray
Group-wise fission cross section.
kappa_fission : dict of numpy.ndarray
Group-wise kappa-fission cross section.
Group-wise kappa_fission cross section.
chi : dict of numpy.ndarray
Group-wise fission spectra ordered by increasing group index (i.e.,
fast to thermal). This attribute should be used if making the common
@ -125,7 +126,7 @@ class XSdata(object):
self.delayed_groups = delayed_groups
self.temperatures = temperatures
self.representation = representation
self._awr = None
self._atomic_weight_ratio = None
self._fissionable = False
self._scatter_format = 'legendre'
self._order = None
@ -164,8 +165,8 @@ class XSdata(object):
return self._representation
@property
def awr(self):
return self._awr
def atomic_weight_ratio(self):
return self._atomic_weight_ratio
@property
def fissionable(self):
@ -281,31 +282,18 @@ class XSdata(object):
check_value('representation', representation, _REPRESENTATIONS)
self._representation = representation
@awr.setter
def awr(self, awr):
# Check validity of type and that the awr value is > 0
check_type('awr', awr, Real)
check_greater_than('awr', awr, 0.0)
self._awr = awr
@fissionable.setter
def fissionable(self, fissionable):
check_type('fissionable', fissionable, bool)
self._fissionable = fissionable
@atomic_weight_ratio.setter
def atomic_weight_ratio(self, atomic_weight_ratio):
# Check validity of type and that the atomic_weight_ratio value is > 0
check_type('atomic_weight_ratio', atomic_weight_ratio, Real)
check_greater_than('atomic_weight_ratio', atomic_weight_ratio, 0.0)
self._atomic_weight_ratio = atomic_weight_ratio
@temperatures.setter
def temperatures(self, temperatures):
check_iterable_type('temperatures', temperatures, Real)
check_type('temperatures', temperatures, Iterable,
expected_iter_type=Real)
# Convert to a numpy array so we can easily get the shape for checking
nptemperatures = np.asarray(temperatures)
check_value('temperatures dimensionality', nptemperatures.ndim, [1])
self._temperatures = nptemperatures
self._temperatures = np.array(temperatures)
@scatter_format.setter
def scatter_format(self, scatter_format):
@ -379,8 +367,8 @@ class XSdata(object):
total: np.ndarray
Total Cross Section
temperature : float
Temperature (in units of Kelvin) of the provided dataset. Defaults
to 294K
Temperature (in Kelvin) of the data. Defaults to room temperature
(294K).
See also
--------
@ -398,13 +386,13 @@ class XSdata(object):
self.energy_groups.num_groups)]
# Convert to a numpy array so we can easily get the shape for checking
nptotal = np.asarray(total)
check_value('total shape', nptotal.shape, shapes)
total = np.asarray(total)
check_value('total shape', total.shape, shapes)
check_type('temperature', temperature, Real)
check_value('temperature', temperature, self.temperatures)
i = self.temperatures.tolist().index(temperature)
self._total[i] = nptotal
i = np.where(self.temperatures == temperature)[0][0]
self._total[i] = total
def set_absorption(self, absorption, temperature=294.):
"""This method sets the cross section for this XSdata object at the
@ -415,8 +403,8 @@ class XSdata(object):
absorption: np.ndarray
Absorption Cross Section
temperature : float
Temperature (in units of Kelvin) of the provided dataset. Defaults
to 294K
Temperature (in Kelvin) of the data. Defaults to room temperature
(294K).
See also
--------
@ -434,13 +422,13 @@ class XSdata(object):
self.energy_groups.num_groups)]
# Convert to a numpy array so we can easily get the shape for checking
npabsorption = np.asarray(absorption)
check_value('absorption shape', npabsorption.shape, shapes)
absorption = np.asarray(absorption)
check_value('absorption shape', absorption.shape, shapes)
check_type('temperature', temperature, Real)
check_value('temperature', temperature, self.temperatures)
i = self.temperatures.tolist().index(temperature)
self._absorption[i] = npabsorption
i = np.where(self.temperatures == temperature)[0][0]
self._absorption[i] = absorption
def set_fission(self, fission, temperature=294.):
"""This method sets the cross section for this XSdata object at the
@ -451,8 +439,8 @@ class XSdata(object):
fission: np.ndarray
Fission Cross Section
temperature : float
Temperature (in units of Kelvin) of the provided dataset. Defaults
to 294K
Temperature (in Kelvin) of the data. Defaults to room temperature
(294K).
See also
--------
@ -470,15 +458,15 @@ class XSdata(object):
self.energy_groups.num_groups)]
# Convert to a numpy array so we can easily get the shape for checking
npfission = np.asarray(fission)
check_value('fission shape', npfission.shape, shapes)
fission = np.asarray(fission)
check_value('fission shape', fission.shape, shapes)
check_type('temperature', temperature, Real)
check_value('temperature', temperature, self.temperatures)
i = self.temperatures.tolist().index(temperature)
self._fission[i] = npfission
i = np.where(self.temperatures == temperature)[0][0]
self._fission[i] = fission
if np.sum(npfission) > 0.0:
if np.sum(fission) > 0.0:
self._fissionable = True
def set_kappa_fission(self, kappa_fission, temperature=294.):
@ -490,8 +478,8 @@ class XSdata(object):
kappa_fission: np.ndarray
Kappa-Fission Cross Section
temperature : float
Temperature (in units of Kelvin) of the provided dataset. Defaults
to 294K
Temperature (in Kelvin) of the data. Defaults to room temperature
(294K).
See also
--------
@ -510,15 +498,15 @@ class XSdata(object):
self.energy_groups.num_groups)]
# Convert to a numpy array so we can easily get the shape for checking
npkappa_fission = np.asarray(kappa_fission)
check_value('kappa fission shape', npkappa_fission.shape, shapes)
kappa_fission = np.asarray(kappa_fission)
check_value('kappa fission shape', kappa_fission.shape, shapes)
check_type('temperature', temperature, Real)
check_value('temperature', temperature, self.temperatures)
i = self.temperatures.tolist().index(temperature)
self._kappa_fission[i] = npkappa_fission
i = np.where(self.temperatures == temperature)[0][0]
self._kappa_fission[i] = kappa_fission
if np.sum(npkappa_fission) > 0.0:
if np.sum(kappa_fission) > 0.0:
self._fissionable = True
def set_chi(self, chi, temperature=294.):
@ -530,8 +518,8 @@ class XSdata(object):
chi: np.ndarray
Fission Spectrum
temperature : float
Temperature (in units of Kelvin) of the provided dataset. Defaults
to 294K
Temperature (in Kelvin) of the data. Defaults to room temperature
(294K).
See also
--------
@ -547,13 +535,13 @@ class XSdata(object):
self.energy_groups.num_groups)]
# Convert to a numpy array so we can easily get the shape for checking
npchi = np.asarray(chi)
check_value('chi shape', npchi.shape, shapes)
chi = np.asarray(chi)
check_value('chi shape', chi.shape, shapes)
check_type('temperature', temperature, Real)
check_value('temperature', temperature, self.temperatures)
i = self.temperatures.tolist().index(temperature)
self._chi[i] = npchi
i = np.where(self.temperatures == temperature)[0][0]
self._chi[i] = chi
def set_chi_prompt(self, chi_prompt, temperature=294.):
"""This method sets the cross section for this XSdata object at the
@ -581,13 +569,13 @@ class XSdata(object):
self.energy_groups.num_groups)]
# Convert to a numpy array so we can easily get the shape for checking
npchi_prompt = np.asarray(chi_prompt)
check_value('chi prompt shape', npchi_prompt.shape, shapes)
chi_prompt = np.asarray(chi_prompt)
check_value('chi prompt shape', chi_prompt.shape, shapes)
check_type('temperature', temperature, Real)
check_value('temperature', temperature, self.temperatures)
i = self.temperatures.tolist().index(temperature)
self._chi_prompt[i] = npchi_prompt
i = np.where(self.temperatures == temperature)[0][0]
self._chi_prompt[i] = chi_prompt
def set_chi_delayed(self, chi_delayed, temperature=294.):
"""This method sets the cross section for this XSdata object at the
@ -618,13 +606,13 @@ class XSdata(object):
self.energy_groups.num_groups)]
# Convert to a numpy array so we can easily get the shape for checking
npchi_delayed = np.asarray(chi_delayed)
check_value('chi delayed shape', npchi_delayed.shape, shapes)
chi_delayed = np.asarray(chi_delayed)
check_value('chi delayed shape', chi_delayed.shape, shapes)
check_type('temperature', temperature, Real)
check_value('temperature', temperature, self.temperatures)
i = self.temperatures.tolist().index(temperature)
self._chi_delayed[i] = npchi_delayed
i = np.where(self.temperatures == temperature)[0][0]
self._chi_delayed[i] = chi_delayed
def set_beta(self, beta, temperature=294.):
"""This method sets the cross section for this XSdata object at the
@ -654,13 +642,13 @@ class XSdata(object):
(self.delayed_groups, self.num_polar, self.num_azimuthal)]
# Convert to a numpy array so we can easily get the shape for checking
npbeta = np.asarray(beta)
check_value('beta shape', npbeta.shape, shapes)
beta = np.asarray(beta)
check_value('beta shape', beta.shape, shapes)
check_type('temperature', temperature, Real)
check_value('temperature', temperature, self.temperatures)
i = self.temperatures.tolist().index(temperature)
self._beta[i] = npbeta
i = np.where(self.temperatures == temperature)[0][0]
self._beta[i] = beta
def set_decay_rate(self, decay_rate, temperature=294.):
"""This method sets the cross section for this XSdata object at the
@ -692,13 +680,13 @@ class XSdata(object):
(self.delayed_groups, self.num_polar, self.num_azimuthal)]
# Convert to a numpy array so we can easily get the shape for checking
npdecay_rate = np.asarray(decay_rate)
check_value('decay rate shape', npdecay_rate.shape, shapes)
decay_rate = np.asarray(decay_rate)
check_value('decay rate shape', decay_rate.shape, shapes)
check_type('temperature', temperature, Real)
check_value('temperature', temperature, self.temperatures)
i = self.temperatures.tolist().index(temperature)
self._decay_rate[i] = npdecay_rate
i = np.where(self.temperatures == temperature)[0][0]
self._decay_rate[i] = decay_rate
def set_scatter_matrix(self, scatter, temperature=294.):
"""This method sets the cross section for this XSdata object at the
@ -709,8 +697,8 @@ class XSdata(object):
scatter: np.ndarray
Scattering Matrix Cross Section
temperature : float
Temperature (in units of Kelvin) of the provided dataset. Defaults
to 294K
Temperature (in Kelvin) of the data. Defaults to room temperature
(294K).
See also
--------
@ -728,15 +716,15 @@ class XSdata(object):
self.energy_groups.num_groups)]
# Convert to a numpy array so we can easily get the shape for checking
npscatter = np.asarray(scatter)
check_iterable_type('scatter', npscatter, Real,
max_depth=len(npscatter.shape))
check_value('scatter shape', npscatter.shape, shapes)
scatter = np.asarray(scatter)
check_iterable_type('scatter', scatter, Real,
max_depth=len(scatter.shape))
check_value('scatter shape', scatter.shape, shapes)
check_type('temperature', temperature, Real)
check_value('temperature', temperature, self.temperatures)
i = self.temperatures.tolist().index(temperature)
self._scatter_matrix[i] = npscatter
i = np.where(self.temperatures == temperature)[0][0]
self._scatter_matrix[i] = scatter
def set_multiplicity_matrix(self, multiplicity, temperature=294.):
"""This method sets the cross section for this XSdata object at the
@ -747,8 +735,8 @@ class XSdata(object):
multiplicity: np.ndarray
Multiplicity Matrix Cross Section
temperature : float
Temperature (in units of Kelvin) of the provided dataset. Defaults
to 294K
Temperature (in Kelvin) of the data. Defaults to room temperature
(294K).
See also
--------
@ -766,15 +754,15 @@ class XSdata(object):
self.energy_groups.num_groups)]
# Convert to a numpy array so we can easily get the shape for checking
npmultiplicity = np.asarray(multiplicity)
check_iterable_type('multiplicity', npmultiplicity, Real,
max_depth=len(npmultiplicity.shape))
check_value('multiplicity shape', npmultiplicity.shape, shapes)
multiplicity = np.asarray(multiplicity)
check_iterable_type('multiplicity', multiplicity, Real,
max_depth=len(multiplicity.shape))
check_value('multiplicity shape', multiplicity.shape, shapes)
check_type('temperature', temperature, Real)
check_value('temperature', temperature, self.temperatures)
i = self.temperatures.tolist().index(temperature)
self._multiplicity_matrix[i] = npmultiplicity
i = np.where(self.temperatures == temperature)[0][0]
self._multiplicity_matrix[i] = multiplicity
def set_nu_fission(self, nu_fission, temperature=294.):
"""This method sets the cross section for this XSdata object at the
@ -785,8 +773,8 @@ class XSdata(object):
nu_fission: np.ndarray
Nu-fission Cross Section
temperature : float
Temperature (in units of Kelvin) of the provided dataset. Defaults
to 294K
Temperature (in Kelvin) of the data. Defaults to room temperature
(294K).
See also
--------
@ -807,16 +795,16 @@ class XSdata(object):
self.energy_groups.num_groups)]
# Convert to a numpy array so we can easily get the shape for checking
npnu_fission = np.asarray(nu_fission)
check_value('nu_fission shape', npnu_fission.shape, shapes)
check_iterable_type('nu_fission', npnu_fission, Real,
max_depth=len(npnu_fission.shape))
nu_fission = np.asarray(nu_fission)
check_value('nu_fission shape', nu_fission.shape, shapes)
check_iterable_type('nu_fission', nu_fission, Real,
max_depth=len(nu_fission.shape))
check_type('temperature', temperature, Real)
check_value('temperature', temperature, self.temperatures)
i = self.temperatures.tolist().index(temperature)
self._nu_fission[i] = npnu_fission
if np.sum(npnu_fission) > 0.0:
i = np.where(self.temperatures == temperature)[0][0]
self._nu_fission[i] = nu_fission
if np.sum(nu_fission) > 0.0:
self._fissionable = True
def set_prompt_nu_fission(self, prompt_nu_fission, temperature=294.):
@ -850,16 +838,16 @@ class XSdata(object):
self.energy_groups.num_groups)]
# Convert to a numpy array so we can easily get the shape for checking
npprompt_nu_fission = np.asarray(prompt_nu_fission)
check_value('prompt_nu_fission shape', npprompt_nu_fission.shape, shapes)
check_iterable_type('prompt_nu_fission', npprompt_nu_fission, Real,
max_depth=len(npprompt_nu_fission.shape))
prompt_nu_fission = np.asarray(prompt_nu_fission)
check_value('prompt_nu_fission shape', prompt_nu_fission.shape, shapes)
check_iterable_type('prompt_nu_fission', prompt_nu_fission, Real,
max_depth=len(prompt_nu_fission.shape))
check_type('temperature', temperature, Real)
check_value('temperature', temperature, self.temperatures)
i = self.temperatures.tolist().index(temperature)
self._prompt_nu_fission[i] = npprompt_nu_fission
if np.sum(npprompt_nu_fission) > 0.0:
i = np.where(self.temperatures == temperature)[0][0]
self._prompt_nu_fission[i] = prompt_nu_fission
if np.sum(prompt_nu_fission) > 0.0:
self._fissionable = True
def set_delayed_nu_fission(self, delayed_nu_fission, temperature=294.):
@ -893,17 +881,17 @@ class XSdata(object):
self.energy_groups.num_groups)]
# Convert to a numpy array so we can easily get the shape for checking
npdelayed_nu_fission = np.asarray(delayed_nu_fission)
check_value('delayed_nu_fission shape', npdelayed_nu_fission.shape,
delayed_nu_fission = np.asarray(delayed_nu_fission)
check_value('delayed_nu_fission shape', delayed_nu_fission.shape,
shapes)
check_iterable_type('delayed_nu_fission', npdelayed_nu_fission, Real,
max_depth=len(npdelayed_nu_fission.shape))
check_iterable_type('delayed_nu_fission', delayed_nu_fission, Real,
max_depth=len(delayed_nu_fission.shape))
check_type('temperature', temperature, Real)
check_value('temperature', temperature, self.temperatures)
i = self.temperatures.tolist().index(temperature)
self._delayed_nu_fission[i] = npdelayed_nu_fission
if np.sum(npdelayed_nu_fission) > 0.0:
i = np.where(self.temperatures == temperature)[0][0]
self._delayed_nu_fission[i] = delayed_nu_fission
if np.sum(delayed_nu_fission) > 0.0:
self._fissionable = True
def set_inverse_velocity(self, inv_vel, temperature=294.):
@ -915,8 +903,8 @@ class XSdata(object):
inv_vel: np.ndarray
Inverse velocity in units of sec/cm.
temperature : float
Temperature (in units of Kelvin) of the provided dataset. Defaults
to 294K
Temperature (in Kelvin) of the data. Defaults to room temperature
(294K).
"""
@ -931,13 +919,13 @@ class XSdata(object):
self.energy_groups.num_groups)]
# Convert to a numpy array so we can easily get the shape for checking
npinv_vel = np.asarray(inv_vel)
check_value('inverse_velocity shape', npinv_vel.shape, shapes)
inv_vel = np.asarray(inv_vel)
check_value('inverse_velocity shape', inv_vel.shape, shapes)
check_type('temperature', temperature, Real)
check_value('temperature', temperature, self.temperatures)
i = self.temperatures.tolist().index(temperature)
self._inverse_velocity[i] = npinv_vel
i = np.where(self.temperatures == temperature)[0][0]
self._inverse_velocity[i] = inv_vel
def set_total_mgxs(self, total, temperature=294., nuclide='total',
xs_type='macro', subdomain=None):
@ -951,8 +939,8 @@ class XSdata(object):
MGXS Object containing the total, transport or nu-transport cross
section for the domain of interest.
temperature : float
Temperature (in units of Kelvin) of the provided dataset. Defaults
to 294K
Temperature (in Kelvin) of the data. Defaults to room temperature
(294K).
nuclide : str
Individual nuclide (or 'total' if obtaining material-wise data)
to gather data for. Defaults to 'total'.
@ -978,7 +966,7 @@ class XSdata(object):
check_type('temperature', temperature, Real)
check_value('temperature', temperature, self.temperatures)
i = self.temperatures.tolist().index(temperature)
i = np.where(self.temperatures == temperature)[0][0]
if self.representation is 'isotropic':
self._total[i] = total.get_xs(nuclides=nuclide, xs_type=xs_type,
subdomains=subdomain)
@ -997,8 +985,8 @@ class XSdata(object):
MGXS Object containing the absorption cross section
for the domain of interest.
temperature : float
Temperature (in units of Kelvin) of the provided dataset. Defaults
to 294K
Temperature (in Kelvin) of the data. Defaults to room temperature
(294K).
nuclide : str
Individual nuclide (or 'total' if obtaining material-wise data)
to gather data for. Defaults to 'total'.
@ -1024,7 +1012,7 @@ class XSdata(object):
check_type('temperature', temperature, Real)
check_value('temperature', temperature, self.temperatures)
i = self.temperatures.tolist().index(temperature)
i = np.where(self.temperatures == temperature)[0][0]
if self.representation is 'isotropic':
self._absorption[i] = absorption.get_xs(nuclides=nuclide,
xs_type=xs_type,
@ -1044,8 +1032,8 @@ class XSdata(object):
MGXS Object containing the fission cross section
for the domain of interest.
temperature : float
Temperature (in units of Kelvin) of the provided dataset. Defaults
to 294K
Temperature (in Kelvin) of the data. Defaults to room temperature
(294K).
nuclide : str
Individual nuclide (or 'total' if obtaining material-wise data)
to gather data for. Defaults to 'total'.
@ -1071,7 +1059,7 @@ class XSdata(object):
check_type('temperature', temperature, Real)
check_value('temperature', temperature, self.temperatures)
i = self.temperatures.tolist().index(temperature)
i = np.where(self.temperatures == temperature)[0][0]
if self.representation is 'isotropic':
self._fission[i] = fission.get_xs(nuclides=nuclide,
xs_type=xs_type,
@ -1091,8 +1079,8 @@ class XSdata(object):
MGXS Object containing the nu-fission cross section
for the domain of interest.
temperature : float
Temperature (in units of Kelvin) of the provided dataset. Defaults
to 294K
Temperature (in Kelvin) of the data. Defaults to room temperature
(294K).
nuclide : str
Individual nuclide (or 'total' if obtaining material-wise data)
to gather data for. Defaults to 'total'.
@ -1119,7 +1107,7 @@ class XSdata(object):
check_type('temperature', temperature, Real)
check_value('temperature', temperature, self.temperatures)
i = self.temperatures.tolist().index(temperature)
i = np.where(self.temperatures == temperature)[0][0]
if self.representation is 'isotropic':
self._nu_fission[i] = nu_fission.get_xs(nuclides=nuclide,
xs_type=xs_type,
@ -1173,7 +1161,7 @@ class XSdata(object):
check_type('temperature', temperature, Real)
check_value('temperature', temperature, self.temperatures)
i = self.temperatures.tolist().index(temperature)
i = np.where(self.temperatures == temperature)[0][0]
if self.representation is 'isotropic':
self._prompt_nu_fission[i] = prompt_nu_fission.get_xs\
(nuclides=nuclide,
@ -1230,7 +1218,7 @@ class XSdata(object):
check_type('temperature', temperature, Real)
check_value('temperature', temperature, self.temperatures)
i = self.temperatures.tolist().index(temperature)
i = np.where(self.temperatures == temperature)[0][0]
if self.representation is 'isotropic':
self._delayed_nu_fission[i] = delayed_nu_fission.get_xs\
(nuclides=nuclide,
@ -1256,8 +1244,8 @@ class XSdata(object):
MGXS Object containing the kappa-fission cross section
for the domain of interest.
temperature : float
Temperature (in units of Kelvin) of the provided dataset. Defaults
to 294K
Temperature (in Kelvin) of the data. Defaults to room temperature
(294K).
nuclide : str
Individual nuclide (or 'total' if obtaining material-wise data)
to gather data for. Defaults to 'total'.
@ -1275,7 +1263,7 @@ class XSdata(object):
"""
check_type('k_fission', k_fission, openmc.mgxs.KappaFissionXS)
check_type('kappa_fission', k_fission, openmc.mgxs.KappaFissionXS)
check_value('energy_groups', k_fission.energy_groups,
[self.energy_groups])
check_value('domain_type', k_fission.domain_type,
@ -1283,7 +1271,7 @@ class XSdata(object):
check_type('temperature', temperature, Real)
check_value('temperature', temperature, self.temperatures)
i = self.temperatures.tolist().index(temperature)
i = np.where(self.temperatures == temperature)[0][0]
if self.representation is 'isotropic':
self._kappa_fission[i] = k_fission.get_xs(nuclides=nuclide,
xs_type=xs_type,
@ -1302,8 +1290,8 @@ class XSdata(object):
chi: openmc.mgxs.Chi
MGXS Object containing chi for the domain of interest.
temperature : float
Temperature (in units of Kelvin) of the provided dataset. Defaults
to 294K
Temperature (in Kelvin) of the data. Defaults to room temperature
(294K).
nuclide : str
Individual nuclide (or 'total' if obtaining material-wise data)
to gather data for. Defaults to 'total'.
@ -1328,7 +1316,7 @@ class XSdata(object):
check_type('temperature', temperature, Real)
check_value('temperature', temperature, self.temperatures)
i = self.temperatures.tolist().index(temperature)
i = np.where(self.temperatures == temperature)[0][0]
if self.representation is 'isotropic':
self._chi[i] = chi.get_xs(nuclides=nuclide,
xs_type=xs_type, subdomains=subdomain)
@ -1373,7 +1361,7 @@ class XSdata(object):
check_type('temperature', temperature, Real)
check_value('temperature', temperature, self.temperatures)
i = self.temperatures.tolist().index(temperature)
i = np.where(self.temperatures == temperature)[0][0]
if self.representation is 'isotropic':
self._chi_prompt[i] = chi_prompt.get_xs(nuclides=nuclide,
xs_type=xs_type,
@ -1421,7 +1409,7 @@ class XSdata(object):
check_type('temperature', temperature, Real)
check_value('temperature', temperature, self.temperatures)
i = self.temperatures.tolist().index(temperature)
i = np.where(self.temperatures == temperature)[0][0]
if self.representation is 'isotropic':
self._chi_delayed[i] = chi_delayed.get_xs(nuclides=nuclide,
xs_type=xs_type,
@ -1467,7 +1455,7 @@ class XSdata(object):
check_type('temperature', temperature, Real)
check_value('temperature', temperature, self.temperatures)
i = self.temperatures.tolist().index(temperature)
i = np.where(self.temperatures == temperature)[0][0]
if self.representation is 'isotropic':
self._beta[i] = beta.get_xs(nuclides=nuclide,
xs_type=xs_type,
@ -1513,7 +1501,7 @@ class XSdata(object):
check_type('temperature', temperature, Real)
check_value('temperature', temperature, self.temperatures)
i = self.temperatures.tolist().index(temperature)
i = np.where(self.temperatures == temperature)[0][0]
if self.representation is 'isotropic':
self._decay_rate[i] = decay_rate.get_xs(nuclides=nuclide,
xs_type=xs_type,
@ -1536,8 +1524,8 @@ class XSdata(object):
MGXS Object containing the scatter matrix cross section
for the domain of interest.
temperature : float
Temperature (in units of Kelvin) of the provided dataset. Defaults
to 294K
Temperature (in Kelvin) of the data. Defaults to room temperature
(294K).
nuclide : str
Individual nuclide (or 'total' if obtaining material-wise data)
to gather data for. Defaults to 'total'.
@ -1563,7 +1551,7 @@ class XSdata(object):
check_type('temperature', temperature, Real)
check_value('temperature', temperature, self.temperatures)
if (self.scatter_format != 'legendre'):
if self.scatter_format is not 'legendre':
msg = 'Anisotropic scattering representations other than ' \
'Legendre expansions have not yet been implemented in ' \
'openmc.mgxs.'
@ -1579,7 +1567,7 @@ class XSdata(object):
check_value('legendre_order', scatter.legendre_order,
[self.order])
i = self.temperatures.tolist().index(temperature)
i = np.where(self.temperatures == temperature)[0][0]
if self.representation is 'isotropic':
# Get the scattering orders in the outermost dimension
self._scatter_matrix[i] = np.zeros((self.num_orders,
@ -1617,8 +1605,8 @@ class XSdata(object):
MGXS Object containing the scattering matrix cross section
for the domain of interest.
temperature : float
Temperature (in units of Kelvin) of the provided dataset. Defaults
to 294K
Temperature (in Kelvin) of the data. Defaults to room temperature
(294K).
nuclide : str
Individual nuclide (or 'total' if obtaining material-wise data)
to gather data for. Defaults to 'total'.
@ -1657,7 +1645,7 @@ class XSdata(object):
check_value('domain_type', scatter.domain_type,
['universe', 'cell', 'material', 'mesh'])
i = self.temperatures.tolist().index(temperature)
i = np.where(self.temperatures == temperature)[0][0]
if self.representation is 'isotropic':
nuscatt = nuscatter.get_xs(nuclides=nuclide,
xs_type=xs_type, moment=0,
@ -1675,20 +1663,25 @@ class XSdata(object):
self._multiplicity_matrix[i] = \
np.nan_to_num(self._multiplicity_matrix[i])
def _get_xsdata_group(self, file):
def to_hdf5(self, file):
"""Write XSdata to an HDF5 file
Parameters
----------
file : h5py.File
HDF5 File (a root Group) to write to
"""
grp = file.create_group(self.name)
if self.awr is not None:
grp.attrs['awr'] = self.awr
if self.atomic_weight_ratio is not None:
grp.attrs['atomic_weight_ratio'] = self.atomic_weight_ratio
if self.fissionable is not None:
grp.attrs['fissionable'] = self.fissionable
if self.representation is not None:
grp.attrs['representation'] = np.array(self.representation,
dtype='S')
if self.representation == 'angle':
if self.representation is 'angle':
if self.num_azimuthal is not None:
grp.attrs['num-azimuthal'] = self.num_azimuthal
@ -1696,7 +1689,7 @@ class XSdata(object):
grp.attrs['num-polar'] = self.num_polar
if self.scatter_format is not None:
grp.attrs['scatter-format'] = np.array(self.scatter_format,
grp.attrs['scatter_format'] = np.array(self.scatter_format,
dtype='S')
if self.order is not None:
@ -1766,7 +1759,7 @@ class XSdata(object):
xs_grp.create_dataset("beta", data=self._beta[i])
if self._decay_rate[i] is not None:
xs_grp.create_dataset("decay-rate", data=self._decay_rate[i])
xs_grp.create_dataset("decay rate", data=self._decay_rate[i])
if self._scatter_matrix[i] is None:
raise ValueError('Scatter matrix must be provided when '
@ -1793,8 +1786,8 @@ class XSdata(object):
flat_scatt.append(matrix[l, g_in, g_out])
# And write it.
scatt_grp = xs_grp.create_group('scatter data')
scatt_grp.create_dataset("scatter matrix",
scatt_grp = xs_grp.create_group('scatter_data')
scatt_grp.create_dataset("scatter_matrix",
data=np.array(flat_scatt))
# Repeat for multiplicity
@ -1817,7 +1810,7 @@ class XSdata(object):
scatt_grp.create_dataset("g_min", data=g_out_bounds[:, 0])
scatt_grp.create_dataset("g_max", data=g_out_bounds[:, 1])
elif self.representation == 'angle':
elif self.representation is 'angle':
Np = self.num_polar
Na = self.num_azimuthal
g_out_bounds = np.zeros((Np, Na, G, 2), dtype=np.int)
@ -1842,8 +1835,8 @@ class XSdata(object):
flat_scatt.append(matrix[l, g_in, g_out])
# And write it.
scatt_grp = xs_grp.create_group('scatter data')
scatt_grp.create_dataset("scatter matrix",
scatt_grp = xs_grp.create_group('scatter_data')
scatt_grp.create_dataset("scatter_matrix",
data=np.array(flat_scatt))
# Repeat for multiplicity
@ -1860,7 +1853,7 @@ class XSdata(object):
flat_mult.append(matrix[g_in, g_out])
# And write it.
scatt_grp.create_dataset("multiplicity matrix",
scatt_grp.create_dataset("multiplicity_matrix",
data=np.array(flat_mult))
# And finally, adjust g_out_bounds for 1-based group counting
@ -2000,6 +1993,6 @@ class MGXSLibrary(object):
file.attrs['group structure'] = self.energy_groups.group_edges
for xsdata in self._xsdatas:
xsdata._get_xsdata_group(file)
xsdata.to_hdf5(file)
file.close()

70
openmc/settings.py
View file

@ -89,10 +89,14 @@ class Settings(object):
Seed for the linear congruential pseudorandom number generator
survival_biasing : bool
Indicate whether survival biasing is to be used
weight : float
Weight cutoff below which particle undergo Russian roulette
weight_avg : float
Weight assigned to particles that are not killed after Russian roulette
cutoff : dict
Dictionary defining weight cutoff and energy cutoff. The dictionary may
have three keys, 'weight', 'weight_avg' and 'energy'. Value for 'weight'
should be a float indicating weight cutoff below which particle undergo
Russian roulette. Value for 'weight_avg' should be a float indicating
weight assigned to particles that are not killed after Russian
roulette. Value of energy should be a float indicating energy in MeV
below which particle will be killed.
entropy_dimension : tuple or list
Number of Shannon entropy mesh cells in the x, y, and z directions,
respectively
@ -214,8 +218,7 @@ class Settings(object):
self._temperature = {}
# Cutoff subelement
self._weight = None
self._weight_avg = None
self._cutoff = None
# Uniform fission source subelement
self._ufs_dimension = 1
@ -393,12 +396,8 @@ class Settings(object):
return self._track
@property
def weight(self):
return self._weight
@property
def weight_avg(self):
return self._weight_avg
def cutoff(self):
return self._cutoff
@property
def ufs_dimension(self):
@ -635,17 +634,29 @@ class Settings(object):
cv.check_type('survival biasing', survival_biasing, bool)
self._survival_biasing = survival_biasing
@weight.setter
def weight(self, weight):
cv.check_type('weight cutoff', weight, Real)
cv.check_greater_than('weight cutoff', weight, 0.0)
self._weight = weight
@cutoff.setter
def cutoff(self, cutoff):
if not isinstance(cutoff, Mapping):
msg = 'Unable to set cutoff from "{0}" which is not a '\
' Python dictionary'.format(cutoff)
raise ValueError(msg)
for key in cutoff:
if key == 'weight':
cv.check_type('weight cutoff', cutoff['weight'], Real)
cv.check_greater_than('weight cutoff', cutoff['weight'], 0.0)
elif key == 'weight_avg':
cv.check_type('average survival weight', cutoff['weight_avg'],
Real)
cv.check_greater_than('average survival weight',
cutoff['weight_avg'], 0.0)
elif key == 'energy':
cv.check_type('energy cutoff', cutoff['energy'], Real)
cv.check_greater_than('energy cutoff', cutoff['energy'], 0.0)
else:
msg = 'Unable to set cutoff to "{0}" which is unsupported by '\
'OpenMC'.format(key)
@weight_avg.setter
def weight_avg(self, weight_avg):
cv.check_type('average survival weight', weight_avg, Real)
cv.check_greater_than('average survival weight', weight_avg, 0.0)
self._weight_avg = weight_avg
self._cutoff = cutoff
@entropy_dimension.setter
def entropy_dimension(self, dimension):
@ -1030,14 +1041,19 @@ class Settings(object):
element.text = str(self._survival_biasing).lower()
def _create_cutoff_subelement(self):
if self._weight is not None:
if self._cutoff is not None:
element = ET.SubElement(self._settings_file, "cutoff")
if 'weight' in self._cutoff:
subelement = ET.SubElement(element, "weight")
subelement.text = str(self._cutoff['weight'])
subelement = ET.SubElement(element, "weight")
subelement.text = str(self._weight)
if 'weight_avg' in self._cutoff:
subelement = ET.SubElement(element, "weight_avg")
subelement.text = str(self._cutoff['weight_avg'])
subelement = ET.SubElement(element, "weight_avg")
subelement.text = str(self._weight_avg)
if 'energy' in self._cutoff:
subelement = ET.SubElement(element, "energy")
subelement.text = str(self._cutoff['energy'])
def _create_entropy_subelement(self):
if self._entropy_lower_left is not None and \

57
scripts/openmc-update-mgxs.py → scripts/openmc-update-mgxs
View file

@ -22,6 +22,7 @@ optional arguments:
"""
from __future__ import print_function
import os
from shutil import move
import warnings
import xml.etree.ElementTree as ET
@ -50,10 +51,10 @@ def parse_args():
help='HDF5 Compression Level')
args = vars(parser.parse_args())
if args['output'] == '':
if args['output'] is '':
filename = args['input'].name
extension = filename[filename.rfind('.'):]
if extension == '.xml':
extension = filenameos.path.splitext()
if extension is '.xml':
filename = filename[:filename.rfind('.')] + '.h5'
args['output'] = filename
@ -62,21 +63,19 @@ def parse_args():
def get_data(element, entry):
try:
value = element.find(entry).text
except:
value = element.find(entry)
if value is not None:
value = value.text.strip()
else:
if entry in element.attrib:
value = element.attrib[entry]
value = element.attrib[entry].strip()
else:
value = None
if value is not None:
value = value.strip()
return value
if __name__ == '__main__':
if __name__ is '__main__':
args = parse_args()
# Parse the XML data.
@ -118,7 +117,7 @@ if __name__ == '__main__':
representation = get_data(xsdata_elem, 'representation')
if representation is None:
representation = 'isotropic'
if representation == 'angle':
if representation is 'angle':
n_azi = int(get_data(xsdata_elem, 'num_azimuthal'))
n_pol = int(get_data(xsdata_elem, 'num_polar'))
@ -139,22 +138,22 @@ if __name__ == '__main__':
try:
# It is in our list, so store that entry
i = names.index(name)
except:
except ValueError:
# It is not in our list, so add it
i = -1
xsd.append(openmc.XSdata(name, energy_groups,
temperatures=temperatures,
representation=representation))
if awr is not None:
xsd[-1].awr = awr
if representation == 'angle':
xsd[-1].atomic_weight_ratio = awr
if representation is 'angle':
xsd[-1].num_azimuthal = n_azi
xsd[-1].num_polar = n_pol
xsd[-1].scatter_format = scatter_format
xsd[-1].order = order
names.append(name)
if scatter_format == 'legendre':
if scatter_format is 'legendre':
order_dim = order + 1
else:
order_dim = order
@ -164,9 +163,9 @@ if __name__ == '__main__':
temp = get_data(xsdata_elem, 'total')
if temp is not None:
temp = np.array(temp.split())
temp = np.array(temp.split(), dtype=float)
total = temp.astype(np.float)
total = np.reshape(total, xsd[i].vector_shape)
total.shape = xsd[i].vector_shape
xsd[i].set_total(total, temperature)
if inverse_velocities is not None:
@ -175,41 +174,41 @@ if __name__ == '__main__':
temp = get_data(xsdata_elem, 'absorption')
temp = np.array(temp.split())
absorption = temp.astype(np.float)
absorption = np.reshape(absorption, xsd[i].vector_shape)
absorption.shape = xsd[i].vector_shape
xsd[i].set_absorption(absorption, temperature)
temp = get_data(xsdata_elem, 'scatter')
temp = np.array(temp.split())
temp = temp.astype(np.float)
scatter = np.reshape(temp, xsd[i].pn_matrix_shape)
scatter = temp.astype(np.float)
scatter.shape = xsd[i].pn_matrix_shape
xsd[i].set_scatter_matrix(scatter, temperature)
temp = get_data(xsdata_elem, 'multiplicity')
if temp is not None:
temp = np.array(temp.split())
temp = temp.astype(np.float)
multiplicity = np.reshape(temp, xsd[i].matrix_shape)
multiplicity = temp.astype(np.float)
multiplicity.shape = xsd[i].matrix_shape
xsd[i].set_multiplicity_matrix(multiplicity, temperature)
temp = get_data(xsdata_elem, 'fission')
if temp is not None:
temp = np.array(temp.split())
fission = temp.astype(np.float)
fission = np.reshape(fission, xsd[i].vector_shape)
fission.shape = xsd[i].vector_shape
xsd[i].set_fission(fission, temperature)
temp = get_data(xsdata_elem, 'kappa_fission')
if temp is not None:
temp = np.array(temp.split())
kappa_fission = temp.astype(np.float)
kappa_fission = np.reshape(kappa_fission, xsd[i].vector_shape)
kappa_fission.shape = xsd[i].vector_shape
xsd[i].set_kappa_fission(kappa_fission, temperature)
temp = get_data(xsdata_elem, 'chi')
if temp is not None:
temp = np.array(temp.split())
chi = temp.astype(np.float)
chi = np.reshape(chi, xsd[i].vector_shape)
chi.shape = xsd[i].vector_shape
xsd[i].set_chi(chi, temperature)
else:
chi = None
@ -217,11 +216,11 @@ if __name__ == '__main__':
temp = get_data(xsdata_elem, 'nu_fission')
if temp is not None:
temp = np.array(temp.split())
temp = temp.astype(np.float)
nu_fission = temp.astype(np.float)
if chi is not None:
nu_fission = np.reshape(temp, xsd[i].vector_shape)
nu_fission.shape = xsd[i].vector_shape
else:
nu_fission = np.reshape(temp, xsd[i].matrix_shape)
nu_fission.shape = xsd[i].matrix_shape
xsd[i].set_nu_fission(nu_fission, temperature)
# Build library as we go, but first we have enough to initialize it

4
setup.py
View file

@ -59,5 +59,9 @@ if have_cython:
kwargs.update({
'ext_modules': cythonize('openmc/data/reconstruct.pyx')
})
import numpy as np
kwargs.update({
'include_dirs': [np.get_include()]
})
setup(**kwargs)

78
src/geometry_header.F90
View file

@ -1,6 +1,12 @@
module geometry_header
use constants, only: HALF, TWO, THREE, INFINITY
use algorithm, only: find
use constants, only: HALF, TWO, THREE, INFINITY, K_BOLTZMANN, &
MATERIAL_VOID, NONE
use dict_header, only: DictCharInt, DictIntInt
use material_header, only: Material
use stl_vector, only: VectorReal
use string, only: to_lower
implicit none
@ -319,4 +325,74 @@ contains
end if
end function get_local_hex
!===============================================================================
! GET_TEMPERATURES returns a list of temperatures that each nuclide/S(a,b) table
! appears at in the model. Later, this list is used to determine the actual
! temperatures to read (which may be different if interpolation is used)
!===============================================================================
subroutine get_temperatures(cells, materials, material_dict, nuclide_dict, &
n_nucs, nuc_temps, sab_dict, n_sabs, sab_temps)
type(Cell), allocatable, intent(in) :: cells(:)
type(Material), allocatable, intent(in) :: materials(:)
type(DictIntInt), intent(in) :: material_dict
type(DictCharInt), intent(in) :: nuclide_dict
integer, intent(in) :: n_nucs
type(VectorReal), allocatable, intent(out) :: nuc_temps(:)
type(DictCharInt), optional, intent(in) :: sab_dict
integer, optional, intent(in) :: n_sabs
type(VectorReal), optional, allocatable, intent(out) :: sab_temps(:)
integer :: i, j, k
integer :: i_nuclide ! index in nuclides array
integer :: i_sab ! index in S(a,b) array
integer :: i_material
real(8) :: temperature ! temperature in Kelvin
allocate(nuc_temps(n_nucs))
if (present(n_sabs) .and. present(sab_temps)) allocate(sab_temps(n_sabs))
do i = 1, size(cells)
do j = 1, size(cells(i) % material)
! Skip any non-material cells and void materials
if (cells(i) % material(j) == NONE .or. &
cells(i) % material(j) == MATERIAL_VOID) cycle
! Get temperature of cell (rounding to nearest integer)
if (size(cells(i) % sqrtkT) > 1) then
temperature = cells(i) % sqrtkT(j)**2 / K_BOLTZMANN
else
temperature = cells(i) % sqrtkT(1)**2 / K_BOLTZMANN
end if
i_material = material_dict % get_key(cells(i) % material(j))
associate (mat => materials(i_material))
NUC_NAMES_LOOP: do k = 1, size(mat % names)
! Get index in nuc_temps array
i_nuclide = nuclide_dict % get_key(to_lower(mat % names(k)))
! Add temperature if it hasn't already been added
if (find(nuc_temps(i_nuclide), temperature) == -1) then
call nuc_temps(i_nuclide) % push_back(temperature)
end if
end do NUC_NAMES_LOOP
if (present(sab_temps) .and. present(sab_dict) .and. &
mat % n_sab > 0) then
SAB_NAMES_LOOP: do k = 1, size(mat % sab_names)
! Get index in nuc_temps array
i_sab = sab_dict % get_key(to_lower(mat % sab_names(k)))
! Add temperature if it hasn't already been added
if (find(sab_temps(i_sab), temperature) == -1) then
call sab_temps(i_sab) % push_back(temperature)
end if
end do SAB_NAMES_LOOP
end if
end associate
end do
end do
end subroutine get_temperatures
end module geometry_header

1
src/global.F90
View file

@ -310,6 +310,7 @@ module global
logical :: survival_biasing = .false.
real(8) :: weight_cutoff = 0.25_8
real(8) :: energy_cutoff = ZERO
real(8) :: weight_survive = ONE
! ============================================================================

33
src/hdf5_interface.F90
View file

@ -85,17 +85,16 @@ module hdf5_interface
public :: write_dataset
public :: read_dataset
public :: check_attribute
public :: attribute_exists
public :: write_attribute
public :: read_attribute
public :: file_create
public :: file_open
public :: file_close
public :: create_group
public :: check_group
public :: object_exists
public :: open_group
public :: close_group
public :: check_dataset
public :: open_dataset
public :: close_dataset
public :: get_shape
@ -254,7 +253,7 @@ contains
! CHECK_ATTRIBUTE Checks to see if an attribute exists in the object
!===============================================================================
function check_attribute(object_id, name) result(exists)
function attribute_exists(object_id, name) result(exists)
integer(HID_T), intent(in) :: object_id
character(*), intent(in) :: name ! name of group
logical :: exists
@ -264,13 +263,13 @@ contains
! Check if attribute exists
call h5aexists_by_name_f(object_id, '.', trim(name), exists, hdf5_err)
end function check_attribute
end function attribute_exists
!===============================================================================
! CHECK_GROUP Checks to see if a group exists in the object
!===============================================================================
function check_group(object_id, name) result(exists)
function object_exists(object_id, name) result(exists)
integer(HID_T), intent(in) :: object_id
character(*), intent(in) :: name ! name of group
logical :: exists
@ -280,21 +279,7 @@ contains
! Check if group exists
call h5ltpath_valid_f(object_id, trim(name), .true., exists, hdf5_err)
end function check_group
!===============================================================================
! CHECK_DATASET Checks to see if a dataset exists in the object
!===============================================================================
function check_dataset(object_id, name) result(exists)
integer(HID_T), intent(in) :: object_id
character(*), intent(in) :: name ! name of group
logical :: exists
! Wrap check_group since the method used there will work here too
exists = check_group(object_id, name)
end function check_dataset
end function object_exists
!===============================================================================
! GET_DATASETS Gets a list of all the datasets in a given location.
@ -368,7 +353,7 @@ contains
integer :: hdf5_err ! HDF5 error code
! Check if group exists
exists = check_group(group_id, name)
exists = object_exists(group_id, name)
! open group if it exists
if (exists) then
@ -391,7 +376,7 @@ contains
logical :: exists ! does the group exist
! Check if group exists
exists = check_group(group_id, name)
exists = object_exists(group_id, name)
! create group
if (exists) then
@ -429,7 +414,7 @@ contains
integer :: hdf5_err ! HDF5 error code
! Check if group exists
exists = check_group(group_id, name)
exists = object_exists(group_id, name)
! open group if it exists
if (exists) then

102
src/input_xml.F90
View file

@ -1,17 +1,15 @@
module input_xml
use hdf5
use algorithm, only: find
use cmfd_input, only: configure_cmfd
use constants
use dict_header, only: DictIntInt, ElemKeyValueCI
use dict_header, only: DictIntInt, DictCharInt, ElemKeyValueCI
use distribution_multivariate
use distribution_univariate
use endf, only: reaction_name
use energy_grid, only: grid_method, n_log_bins
use error, only: fatal_error, warning
use geometry_header, only: Cell, Lattice, RectLattice, HexLattice
use geometry_header, only: Cell, Lattice, RectLattice, HexLattice, &
get_temperatures
use global
use hdf5_interface
use list_header, only: ListChar, ListInt, ListReal
@ -61,6 +59,9 @@ contains
call time_read_xs % stop()
end if
! Normalize atom/weight percents
if (run_mode /= MODE_PLOTTING) call normalize_ao()
end subroutine read_input_xml
!===============================================================================
@ -669,8 +670,15 @@ contains
! Cutoffs
if (check_for_node(doc, "cutoff")) then
call get_node_ptr(doc, "cutoff", node_cutoff)
call get_node_value(node_cutoff, "weight", weight_cutoff)
call get_node_value(node_cutoff, "weight_avg", weight_survive)
if (check_for_node(node_cutoff, "weight")) then
call get_node_value(node_cutoff, "weight", weight_cutoff)
end if
if (check_for_node(node_cutoff, "weight_avg")) then
call get_node_value(node_cutoff, "weight_avg", weight_survive)
end if
if (check_for_node(node_cutoff, "energy")) then
call get_node_value(node_cutoff, "energy", energy_cutoff)
end if
end if
! Particle trace
@ -2120,7 +2128,9 @@ contains
call assign_temperatures(material_temps)
! Determine desired temperatures for each nuclide and S(a,b) table
call get_temperatures(nuc_temps, sab_temps)
call get_temperatures(cells, materials, material_dict, nuclide_dict, &
n_nuclides_total, nuc_temps, sab_dict, &
n_sab_tables, sab_temps)
! Read continuous-energy cross sections
if (run_CE .and. run_mode /= MODE_PLOTTING) then
@ -2129,9 +2139,6 @@ contains
call time_read_xs % stop()
end if
! Normalize atom/weight percents
if (run_mode /= MODE_PLOTTING) call normalize_ao()
! Clear dictionary
call library_dict % clear()
end subroutine read_materials
@ -4646,14 +4653,14 @@ contains
! Open file for reading
file_id = file_open(path_cross_sections, 'r', parallel=.true.)
if (check_attribute(file_id, "energy_groups")) then
if (attribute_exists(file_id, "energy_groups")) then
! Get neutron energy group count
call read_attribute(num_energy_groups, file_id, "energy_groups")
else
call fatal_error("'energy_groups' attribute must exist!")
end if
if (check_attribute(file_id, "delayed_groups")) then
if (attribute_exists(file_id, "delayed_groups")) then
! Get neutron delayed group count
call read_attribute(num_delayed_groups, file_id, "delayed_groups")
else
@ -4665,7 +4672,7 @@ contains
allocate(rev_energy_bins(num_energy_groups + 1))
allocate(energy_bins(num_energy_groups + 1))
if (check_attribute(file_id, "group structure")) then
if (attribute_exists(file_id, "group structure")) then
! Get neutron group structure
call read_attribute(energy_bins, file_id, "group structure")
else
@ -4698,7 +4705,7 @@ contains
libraries(i) % materials(1) = names(i)
end do
! Close MGXS HDF file
! Close MGXS HDF5 file
call file_close(file_id)
end subroutine read_mg_cross_sections_header
@ -5646,7 +5653,7 @@ contains
if (run_CE) then
awr = nuclides(mat % nuclide(j)) % awr
else
awr = ONE
awr = nuclides_MG(mat % nuclide(j)) % obj % awr
end if
! if given weight percent, convert all values so that they are divided
@ -5671,7 +5678,7 @@ contains
if (run_CE) then
awr = nuclides(mat % nuclide(j)) % awr
else
awr = ONE
awr = nuclides_MG(mat % nuclide(j)) % obj % awr
end if
x = mat % atom_density(j)
sum_percent = sum_percent + x*awr
@ -5934,67 +5941,6 @@ contains
end do
end subroutine assign_temperatures
!===============================================================================
! GET_TEMPERATURES returns a list of temperatures that each nuclide/S(a,b) table
! appears at in the model. Later, this list is used to determine the actual
! temperatures to read (which may be different if interpolation is used)
!===============================================================================
subroutine get_temperatures(nuc_temps, sab_temps)
type(VectorReal), allocatable, intent(out) :: nuc_temps(:)
type(VectorReal), allocatable, intent(out) :: sab_temps(:)
integer :: i, j, k
integer :: i_nuclide ! index in nuclides array
integer :: i_sab ! index in S(a,b) array
integer :: i_material
real(8) :: temperature ! temperature in Kelvin
allocate(nuc_temps(n_nuclides_total))
allocate(sab_temps(n_sab_tables))
do i = 1, size(cells)
do j = 1, size(cells(i) % material)
! Skip any non-material cells and void materials
if (cells(i) % material(j) == NONE .or. &
cells(i) % material(j) == MATERIAL_VOID) cycle
! Get temperature of cell (rounding to nearest integer)
if (size(cells(i) % sqrtkT) > 1) then
temperature = cells(i) % sqrtkT(j)**2 / K_BOLTZMANN
else
temperature = cells(i) % sqrtkT(1)**2 / K_BOLTZMANN
end if
i_material = material_dict % get_key(cells(i) % material(j))
associate (mat => materials(i_material))
NUC_NAMES_LOOP: do k = 1, size(mat % names)
! Get index in nuc_temps array
i_nuclide = nuclide_dict % get_key(to_lower(mat % names(k)))
! Add temperature if it hasn't already been added
if (find(nuc_temps(i_nuclide), temperature) == -1) then
call nuc_temps(i_nuclide) % push_back(temperature)
end if
end do NUC_NAMES_LOOP
if (mat % n_sab > 0) then
SAB_NAMES_LOOP: do k = 1, size(mat % sab_names)
! Get index in nuc_temps array
i_sab = sab_dict % get_key(to_lower(mat % sab_names(k)))
! Add temperature if it hasn't already been added
if (find(sab_temps(i_sab), temperature) == -1) then
call sab_temps(i_sab) % push_back(temperature)
end if
end do SAB_NAMES_LOOP
end if
end associate
end do
end do
end subroutine get_temperatures
!===============================================================================
! READ_0K_ELASTIC_SCATTERING
!===============================================================================

60
src/mgxs_data.F90
View file

@ -3,8 +3,9 @@ module mgxs_data
use constants
use algorithm, only: find
use error, only: fatal_error
use geometry_header, only: get_temperatures
use global
use hdf5
use hdf5_interface
use material_header, only: Material
use mgxs_header
use output, only: write_message
@ -48,7 +49,8 @@ contains
call write_message("Loading Cross Section Data...", 5)
! Get temperatures
call get_temperatures(temps)
call get_temperatures(cells, materials, material_dict, nuclide_dict, &
n_nuclides_total, temps)
! Open file for reading
file_id = file_open(path_cross_sections, 'r', parallel=.true.)
@ -76,7 +78,7 @@ contains
call write_message("Loading " // trim(name) // " Data...", 5)
! Check to make sure cross section set exists in the library
if (check_group(file_id, trim(name))) then
if (object_exists(file_id, trim(name))) then
xsdata_group = open_group(file_id, trim(name))
else
call fatal_error("Data for '" // trim(name) // "' does not exist in "&
@ -84,7 +86,7 @@ contains
end if
! First find out the data representation
if (check_attribute(xsdata_group, "representation")) then
if (attribute_exists(xsdata_group, "representation")) then
call read_attribute(temp_str, xsdata_group, "representation")
@ -119,6 +121,7 @@ contains
! Add name to dictionary
call already_read % add(name)
end if
end do NUCLIDE_LOOP
end do MATERIAL_LOOP
@ -230,53 +233,4 @@ contains
end subroutine get_mat_kTs
!===============================================================================
! GET_TEMPERATURES returns a list of temperatures that each MGXS table
! appears at in the model. Later, this list is used to determine the actual
! temperatures to read (which may be different if interpolation is used)
!===============================================================================
subroutine get_temperatures(temps)
type(VectorReal), allocatable, intent(out) :: temps(:)
integer :: i, j, k
integer :: i_nuclide ! index in nuclides array
integer :: i_material ! Index in materials array
real(8) :: temperature ! temperature in Kelvin
allocate(temps(n_nuclides_total))
do i = 1, size(cells)
do j = 1, size(cells(i) % material)
! Skip any non-material cells and void materials
if (cells(i) % material(j) == NONE .or. &
cells(i) % material(j) == MATERIAL_VOID) cycle
! Get temperature of cell (rounding to nearest integer)
if (size(cells(i) % sqrtkT) > 1) then
temperature = cells(i) % sqrtkT(j)**2 / K_BOLTZMANN
else
temperature = cells(i) % sqrtkT(1)**2 / K_BOLTZMANN
end if
i_material = material_dict % get_key(cells(i) % material(j))
associate (mat => materials(i_material))
NUC_NAMES_LOOP: do k = 1, size(mat % names)
! Get index in temps array
i_nuclide = nuclide_dict % get_key(to_lower(mat % names(k)))
! Add temperature if it hasn't already been added
if (find(temps(i_nuclide), temperature) == -1) then
call temps(i_nuclide) % push_back(temperature)
end if
end do NUC_NAMES_LOOP
end associate
end do
end do
end subroutine get_temperatures
end module mgxs_data

349
src/mgxs_header.F90
View file

File diff suppressed because it is too large Load diff

8
src/nuclide_header.F90
View file

@ -13,7 +13,7 @@ module nuclide_header
use error, only: fatal_error, warning
use hdf5_interface, only: read_attribute, open_group, close_group, &
open_dataset, read_dataset, close_dataset, get_shape, get_datasets, &
check_group, get_name, get_groups
object_exists, get_name, get_groups
use list_header, only: ListInt
use math, only: evaluate_legendre
use multipole_header, only: MultipoleArray
@ -350,7 +350,7 @@ module nuclide_header
call close_group(rxs_group)
! Read unresolved resonance probability tables if present
if (check_group(group_id, 'urr')) then
if (object_exists(group_id, 'urr')) then
this % urr_present = .true.
allocate(this % urr_data(n_temperature))
@ -391,7 +391,7 @@ module nuclide_header
end if
! Check for nu-total
if (check_group(group_id, 'total_nu')) then
if (object_exists(group_id, 'total_nu')) then
nu_group = open_group(group_id, 'total_nu')
! Read total nu data
@ -410,7 +410,7 @@ module nuclide_header
end if
! Read fission energy release data if present
if (check_group(group_id, 'fission_energy_release')) then
if (object_exists(group_id, 'fission_energy_release')) then
fer_group = open_group(group_id, 'fission_energy_release')
! Check to see if this is polynomial or tabulated data

7
src/physics.F90
View file

@ -127,6 +127,13 @@ contains
if (.not. p % alive) return
end if
! Kill neutron under certain energy
if (p % E < energy_cutoff) then
p % alive = .false.
p % wgt = ZERO
p % last_wgt = ZERO
end if
end subroutine sample_reaction
!===============================================================================

6
src/sab_header.F90
View file

@ -9,7 +9,7 @@ module sab_header
use error, only: warning, fatal_error
use hdf5, only: HID_T, HSIZE_T, SIZE_T
use hdf5_interface, only: read_attribute, get_shape, open_group, close_group, &
open_dataset, read_dataset, close_dataset, get_datasets, check_group, &
open_dataset, read_dataset, close_dataset, get_datasets, object_exists, &
get_name
use secondary_correlated, only: CorrelatedAngleEnergy
use stl_vector, only: VectorInt, VectorReal
@ -205,7 +205,7 @@ contains
T_group = open_group(group_id, temp_str)
! Coherent elastic data
if (check_group(T_group, 'elastic')) then
if (object_exists(T_group, 'elastic')) then
! Read cross section data
elastic_group = open_group(T_group, 'elastic')
dset_id = open_dataset(elastic_group, 'xs')
@ -247,7 +247,7 @@ contains
end if
! Inelastic data
if (check_group(T_group, 'inelastic')) then
if (object_exists(T_group, 'inelastic')) then
! Read type of inelastic data
inelastic_group = open_group(T_group, 'inelastic')

7
src/secondary_uncorrelated.F90
View file

@ -8,7 +8,8 @@ module secondary_uncorrelated
use energy_distribution, only: EnergyDistribution, LevelInelastic, &
ContinuousTabular, MaxwellEnergy, Evaporation, WattEnergy, DiscretePhoton
use error, only: warning
use hdf5_interface, only: read_attribute, open_group, close_group, check_group
use hdf5_interface, only: read_attribute, open_group, close_group, &
object_exists
use random_lcg, only: prn
!===============================================================================
@ -60,14 +61,14 @@ contains
character(MAX_WORD_LEN) :: type
! Check if angle group is present & read
if (check_group(group_id, 'angle')) then
if (object_exists(group_id, 'angle')) then
angle_group = open_group(group_id, 'angle')
call this%angle%from_hdf5(angle_group)
call close_group(angle_group)
end if
! Check if energy group is present & read
if (check_group(group_id, 'energy')) then
if (object_exists(group_id, 'energy')) then
energy_group = open_group(group_id, 'energy')
call read_attribute(type, energy_group, 'type')
select case (type)

84
src/tally.F90
View file

@ -1247,19 +1247,19 @@ contains
end if
if (i_nuclide > 0) then
score = score * nucxs % get_xs('inverse_velocity', p_g, UVW=p_uvw) &
score = score * nucxs % get_xs('inverse-velocity', p_g, UVW=p_uvw) &
/ nucxs % get_xs('total', p_g, UVW=p_uvw) * flux
else
score = score * matxs % get_xs('inverse_velocity', p_g, UVW=p_uvw) &
score = score * matxs % get_xs('inverse-velocity', p_g, UVW=p_uvw) &
/ matxs % get_xs('total', p_g, UVW=p_uvw) * flux
end if
else
if (i_nuclide > 0) then
score = flux * nucxs % get_xs('inverse_velocity', p_g, UVW=p_uvw)
score = flux * nucxs % get_xs('inverse-velocity', p_g, UVW=p_uvw)
else
score = flux * matxs % get_xs('inverse_velocity', p_g, UVW=p_uvw)
score = flux * matxs % get_xs('inverse-velocity', p_g, UVW=p_uvw)
end if
end if
@ -1445,11 +1445,11 @@ contains
score = p % absorb_wgt * flux
if (i_nuclide > 0) then
score = score * atom_density * &
nucxs % get_xs('nu_fission', p_g, UVW=p_uvw) / &
nucxs % get_xs('nu-fission', p_g, UVW=p_uvw) / &
matxs % get_xs('absorption', p_g, UVW=p_uvw)
else
score = score * &
matxs % get_xs('nu_fission', p_g, UVW=p_uvw) / &
matxs % get_xs('nu-fission', p_g, UVW=p_uvw) / &
matxs % get_xs('absorption', p_g, UVW=p_uvw)
end if
else
@ -1470,10 +1470,10 @@ contains
else
if (i_nuclide > 0) then
score = nucxs % get_xs('nu_fission', p_g, UVW=p_uvw) * &
score = nucxs % get_xs('nu-fission', p_g, UVW=p_uvw) * &
atom_density * flux
else
score = matxs % get_xs('nu_fission', p_g, UVW=p_uvw) * flux
score = matxs % get_xs('nu-fission', p_g, UVW=p_uvw) * flux
end if
end if
@ -1506,11 +1506,11 @@ contains
score = p % absorb_wgt * flux
if (i_nuclide > 0) then
score = score * atom_density * &
nucxs % get_xs('prompt_nu_fission', p_g, UVW=p_uvw) / &
nucxs % get_xs('prompt-nu-fission', p_g, UVW=p_uvw) / &
matxs % get_xs('absorption', p_g, UVW=p_uvw)
else
score = score * &
matxs % get_xs('prompt_nu_fission', p_g, UVW=p_uvw) / &
matxs % get_xs('prompt-nu-fission', p_g, UVW=p_uvw) / &
matxs % get_xs('absorption', p_g, UVW=p_uvw)
end if
else
@ -1532,10 +1532,10 @@ contains
else
if (i_nuclide > 0) then
score = nucxs % get_xs('prompt_nu_fission', p_g, UVW=p_uvw) * &
score = nucxs % get_xs('prompt-nu-fission', p_g, UVW=p_uvw) * &
atom_density * flux
else
score = matxs % get_xs('prompt_nu_fission', p_g, UVW=p_uvw) * flux
score = matxs % get_xs('prompt-nu-fission', p_g, UVW=p_uvw) * flux
end if
end if
@ -1583,11 +1583,11 @@ contains
score = p % absorb_wgt * flux
if (i_nuclide > 0) then
score = score * nucxs % get_xs('delayed_nu_fission', &
score = score * nucxs % get_xs('delayed-nu-fission', &
p_g, UVW=p_uvw, dg=d) / &
nucxs % get_xs('absorption', p_g, UVW=p_uvw)
else
score = score * matxs % get_xs('delayed_nu_fission', &
score = score * matxs % get_xs('delayed-nu-fission', &
p_g, UVW=p_uvw, dg=d) / &
matxs % get_xs('absorption', p_g, UVW=p_uvw)
end if
@ -1599,10 +1599,10 @@ contains
else
score = p % absorb_wgt * flux
if (i_nuclide > 0) then
score = score * nucxs % get_xs('delayed_nu_fission', p_g, &
score = score * nucxs % get_xs('delayed-nu-fission', p_g, &
UVW=p_uvw) / nucxs % get_xs('absorption', p_g, UVW=p_uvw)
else
score = score * matxs % get_xs('delayed_nu_fission', p_g, &
score = score * matxs % get_xs('delayed-nu-fission', p_g, &
UVW=p_uvw) / matxs % get_xs('absorption', p_g, UVW=p_uvw)
end if
end if
@ -1665,10 +1665,10 @@ contains
d = filt % groups(d_bin)
if (i_nuclide > 0) then
score = nucxs % get_xs('delayed_nu_fission', p_g, &
score = nucxs % get_xs('delayed-nu-fission', p_g, &
UVW=p_uvw, dg=d) * atom_density * flux
else
score = matxs % get_xs('delayed_nu_fission', p_g, &
score = matxs % get_xs('delayed-nu-fission', p_g, &
UVW=p_uvw, dg=d) * flux
end if
@ -1678,10 +1678,10 @@ contains
end select
else
if (i_nuclide > 0) then
score = nucxs % get_xs('delayed_nu_fission', p_g, UVW=p_uvw) &
score = nucxs % get_xs('delayed-nu-fission', p_g, UVW=p_uvw) &
* atom_density * flux
else
score = matxs % get_xs('delayed_nu_fission', p_g, UVW=p_uvw) &
score = matxs % get_xs('delayed-nu-fission', p_g, UVW=p_uvw) &
* flux
end if
end if
@ -1719,15 +1719,15 @@ contains
score = p % absorb_wgt * flux
if (i_nuclide > 0) then
score = score * nucxs % get_xs('decay_rate', p_g, &
score = score * nucxs % get_xs('decay rate', p_g, &
UVW=p_uvw, dg=d) * &
nucxs % get_xs('delayed_nu_fission', p_g, &
nucxs % get_xs('delayed-nu-fission', p_g, &
UVW=p_uvw, dg=d) / nucxs % get_xs('absorption', &
p_g, UVW=p_uvw)
else
score = score * matxs % get_xs('decay_rate', p_g, &
score = score * matxs % get_xs('decay rate', p_g, &
UVW=p_uvw, dg=d) * &
matxs % get_xs('delayed_nu_fission', p_g, &
matxs % get_xs('delayed-nu-fission', p_g, &
UVW=p_uvw, dg=d) / matxs % get_xs('absorption', &
p_g, UVW=p_uvw)
end if
@ -1747,13 +1747,13 @@ contains
do d = 1, num_delayed_groups
if (i_nuclide > 0) then
score = score + p % absorb_wgt * &
nucxs % get_xs('decay_rate', p_g, UVW=p_uvw, dg=d) * &
nucxs % get_xs('delayed_nu_fission', p_g, UVW=p_uvw, &
nucxs % get_xs('decay rate', p_g, UVW=p_uvw, dg=d) * &
nucxs % get_xs('delayed-nu-fission', p_g, UVW=p_uvw, &
dg=d) / nucxs % get_xs('absorption', p_g, UVW=p_uvw)
else
score = score + p % absorb_wgt * &
matxs % get_xs('decay_rate', p_g, UVW=p_uvw, dg=d) * &
matxs % get_xs('delayed_nu_fission', p_g, UVW=p_uvw, &
matxs % get_xs('decay rate', p_g, UVW=p_uvw, dg=d) * &
matxs % get_xs('delayed-nu-fission', p_g, UVW=p_uvw, &
dg=d) / matxs % get_xs('absorption', p_g, UVW=p_uvw)
end if
end do
@ -1783,13 +1783,13 @@ contains
if (i_nuclide > 0) then
score = score + keff * atom_density * &
fission_bank(n_bank - p % n_bank + k) % wgt * &
nucxs % get_xs('decay_rate', p_g, UVW=p_uvw, dg=g) * &
nucxs % get_xs('decay rate', p_g, UVW=p_uvw, dg=g) * &
nucxs % get_xs('fission', p_g, UVW=p_uvw) / &
matxs % get_xs('fission', p_g, UVW=p_uvw)
else
score = score + keff * &
fission_bank(n_bank - p % n_bank + k) % wgt * &
matxs % get_xs('decay_rate', p_g, UVW=p_uvw, dg=g)
matxs % get_xs('decay rate', p_g, UVW=p_uvw, dg=g)
end if
! if the delayed group filter is present, tally to corresponding
@ -1841,12 +1841,12 @@ contains
d = filt % groups(d_bin)
if (i_nuclide > 0) then
score = nucxs % get_xs('decay_rate', p_g, UVW=p_uvw, dg=d) * &
nucxs % get_xs('delayed_nu_fission', p_g, UVW=p_uvw, &
score = nucxs % get_xs('decay rate', p_g, UVW=p_uvw, dg=d) * &
nucxs % get_xs('delayed-nu-fission', p_g, UVW=p_uvw, &
dg=d) * atom_density * flux
else
score = matxs % get_xs('decay_rate', p_g, UVW=p_uvw, dg=d) * &
matxs % get_xs('delayed_nu_fission', p_g, UVW=p_uvw, &
score = matxs % get_xs('decay rate', p_g, UVW=p_uvw, dg=d) * &
matxs % get_xs('delayed-nu-fission', p_g, UVW=p_uvw, &
dg=d) * flux
end if
@ -1864,12 +1864,12 @@ contains
do d = 1, num_delayed_groups
if (i_nuclide > 0) then
score = score + atom_density * flux * &
nucxs % get_xs('decay_rate', p_g, UVW=p_uvw, dg=d) * &
nucxs % get_xs('delayed_nu_fission', p_g, UVW=p_uvw, dg=d)
nucxs % get_xs('decay rate', p_g, UVW=p_uvw, dg=d) * &
nucxs % get_xs('delayed-nu-fission', p_g, UVW=p_uvw, dg=d)
else
score = score + flux * &
matxs % get_xs('decay_rate', p_g, UVW=p_uvw, dg=d) * &
matxs % get_xs('delayed_nu_fission', p_g, UVW=p_uvw, dg=d)
matxs % get_xs('decay rate', p_g, UVW=p_uvw, dg=d) * &
matxs % get_xs('delayed-nu-fission', p_g, UVW=p_uvw, dg=d)
end if
end do
end if
@ -1901,19 +1901,19 @@ contains
end if
if (i_nuclide > 0) then
score = score * atom_density * &
nucxs % get_xs('kappa_fission', p_g, UVW=p_uvw) / &
nucxs % get_xs('kappa-fission', p_g, UVW=p_uvw) / &
matxs % get_xs('absorption', p_g, UVW=p_uvw)
else
score = score * &
matxs % get_xs('kappa_fission', p_g, UVW=p_uvw) / &
matxs % get_xs('kappa-fission', p_g, UVW=p_uvw) / &
matxs % get_xs('absorption', p_g, UVW=p_uvw)
end if
else
if (i_nuclide > 0) then
score = nucxs % get_xs('kappa_fission', p_g, UVW=p_uvw) * &
score = nucxs % get_xs('kappa-fission', p_g, UVW=p_uvw) * &
atom_density * flux
else
score = matxs % get_xs('kappa_fission', p_g, UVW=p_uvw) * flux
score = matxs % get_xs('kappa-fission', p_g, UVW=p_uvw) * flux
end if
end if

BIN
tests/1d_mgxs.h5
View file

Binary file not shown.

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

@ -0,0 +1 @@
4dc6a7b131f6757ecc9b06e93d425712f0813c546ed9eaa0aafa4990383b2f3a5b742a5a39ef5f27300d2c861b344e78de2d2c569b9439e7893dc06b03bf3b02

3
tests/test_energy_cutoff/results_true.dat Normal file
View file

@ -0,0 +1,3 @@
tally 1:
sum = 0.000000E+00
sum_sq = 0.000000E+00

84
tests/test_energy_cutoff/test_energy_cutoff.py Executable file
View file

@ -0,0 +1,84 @@
#!/usr/bin/env python
import os
import sys
sys.path.insert(0, os.pardir)
from testing_harness import PyAPITestHarness
import openmc
class EnergyCutoffTestHarness(PyAPITestHarness):
def _build_inputs(self):
# Set energy cutoff
energy_cutoff = 4e-6
# Material is composed of H-1
mat = openmc.Material(material_id=1, name='mat')
mat.set_density('atom/b-cm', 0.069335)
mat.add_nuclide('H1', 40.0)
materials_file = openmc.Materials([mat])
materials_file.export_to_xml()
# Cell is box with reflective boundary
x1 = openmc.XPlane(surface_id=1, x0=-1)
x2 = openmc.XPlane(surface_id=2, x0=1)
y1 = openmc.YPlane(surface_id=3, y0=-1)
y2 = openmc.YPlane(surface_id=4, y0=1)
z1 = openmc.ZPlane(surface_id=5, z0=-1)
z2 = openmc.ZPlane(surface_id=6, z0=1)
for surface in [x1, x2, y1, y2, z1, z2]:
surface.boundary_type = 'reflective'
box = openmc.Cell(cell_id=1, name='box')
box.region = +x1 & -x2 & +y1 & -y2 & +z1 & -z2
box.fill = mat
root = openmc.Universe(universe_id=0, name='root universe')
root.add_cell(box)
geometry = openmc.Geometry(root)
geometry.export_to_xml()
# Set the running parameters
settings_file = openmc.Settings()
settings_file.run_mode = 'fixed source'
settings_file.batches = 10
settings_file.particles = 100
settings_file.cutoff = {'energy': energy_cutoff}
bounds = [-1, -1, -1, 1, 1, 1]
uniform_dist = openmc.stats.Box(bounds[:3], bounds[3:])
watt_dist = openmc.stats.Watt()
settings_file.source = openmc.source.Source(space=uniform_dist,
energy=watt_dist)
settings_file.export_to_xml()
# Tally flux under energy cutoff
tallies = openmc.Tallies()
tally = openmc.Tally(1)
tally.scores = ['flux']
energy_filter = openmc.filter.EnergyFilter((0.0, energy_cutoff))
tally.filters = [energy_filter]
tallies.append(tally)
tallies.export_to_xml()
def _get_results(self):
"""Digest info in the statepoint and return as a string."""
# Read the statepoint file.
sp = openmc.StatePoint(self._sp_name)
# Write out tally data.
outstr = ''
t = sp.get_tally()
outstr += 'tally {0}:\n'.format(t.id)
outstr += 'sum = {0:12.6E}\n'.format(t.sum[0, 0, 0])
outstr += 'sum_sq = {0:12.6E}\n'.format(t.sum_sq[0, 0, 0])
return outstr
def _cleanup(self):
super(EnergyCutoffTestHarness, self)._cleanup()
f = os.path.join(os.getcwd(), 'tallies.xml')
if os.path.exists(f):
os.remove(f)
if __name__ == '__main__':
harness = EnergyCutoffTestHarness('statepoint.10.h5', True)
harness.main()

2
tests/test_mg_basic/results_true.dat
View file

@ -1,2 +1,2 @@
k-combined:
1.091835E+00 2.710745E-03
1.070804E+00 1.141438E-02

2
tests/test_mg_max_order/results_true.dat
View file

@ -1,2 +1,2 @@
k-combined:
1.066827E+00 2.533280E-02
1.093516E+00 1.593236E-02

2
tests/test_mg_nuclide/results_true.dat
View file

@ -1,2 +1,2 @@
k-combined:
1.091835E+00 2.710745E-03
1.070804E+00 1.141438E-02

1306
tests/test_mg_tallies/results_true.dat
View file

File diff suppressed because it is too large Load diff