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

Merge branch 'mg-mode-delayed' into mg-mode-delayed-matrix

Browse source
This commit is contained in:
Sam Shaner 2016-10-18 09:08:36 -04:00
commit ad19b2092b
25 changed files with 293 additions and 3227 deletions
Download patch file Download diff file Expand all files Collapse all files

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

@ -9,7 +9,7 @@ nuclear data is available. In continuous-energy mode, the
``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
nuclear data itself is contained within an ``mgxs.h5`` file. This portion of
the manual describes the format of the multi-group data library required
to be used in the ``mgxs.h5`` file.
@ -40,10 +40,10 @@ MGXS Library Specification
The data within <library name> contains the temperature-dependent multi-group
data for the nuclide or material that it represents.
: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).
:Attributes: - **atomic_weight_ratio** (*double*) -- The atomic weight ratio
(optional, i.e. it is not meaningful for material-wise data).
- **fissionable** (*bool*) -- Whether the dataset is fissionable
(True) or not (False).
- **representation** (*char[]*) -- The method used to generate and
represent the multi-group cross sections. That is, whether they
were generated with scalar flux weighting (or reduced to a
@ -67,6 +67,15 @@ data for the nuclide or material that it represents.
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.
- **scatter_shape** (*char[]*) -- The shape of the provided
scatter and multiplicity matrix. The values provided are strings
describing the ordering the scattering array is provided in
row-major (i.e., C/C++ and Python) indexing. Valid values are
"[Order][G][G']" or "[Order][G'][G]" where "G'" denotes the
secondary/outgoing energy groups, "G" denotes the incoming
energy groups, and "Order" is the angular distribution index.
This value is not required; if not the default value of
"[Order][G][G']" will be assumed.
**/<library name>/kTs/**
@ -81,29 +90,29 @@ Temperature-dependent data, provided for temperature <TTT>K.
:Datasets: - **total** (*double[]* or *double[][][]*) -- Total 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].
[groups][azimuthal][polar].
- **absorption** (*double[]* or *double[][][]*) -- Absorption
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].
[groups][azimuthal][polar].
- **fission** (*double[]* or *double[][][]*) -- 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 dataset
[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 dataset
[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 dataset
[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.
@ -121,34 +130,30 @@ Temperature-dependent data, provided for temperature <TTT>K.
Data specific to neutron scattering for the temperature <TTT>K
: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
Minimum (most energetic) groups with non-zero values of
the scattering matrix provided. If `scatter_shape` is
"[Order][G][G']" then `g_min` will describe the minimum values
of "G'" for each "G"; if `scatter_shape` is "[Order][G'][G]"
then `g_min` will describe the minimum values of "G" for each "G'".
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_min` is:
`g_min[g_in]`, or `g_min[num_polar][num_azimuthal][g_in]`.
`g_min[g]`, or `g_min[num_polar][num_azimuthal][g]`.
The former is used when `representation` is "isotropic", and the
latter when `representation` is "angle".
- **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_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".
Similar to `g_min`, except this dataset describes the maximum
(least energetic) groups with non-zero values of
the scattering matrix.
- **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]`
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
larger than the value of `order`, otherwise it will match `order`.
Finally, the g_out dimension has a dimensionality of
`g_min` to `g_max`.
scattering moment matrices. The pre-flattened array corresponds to
the shape provied in `scatter_shape`, but if `representation` is
"angle" the dimensionality in `scatter_shape` is prepended by
"[num_polar][num_azimuthal]" dimensions. The right-most energy
group dimension will only include the entries between `g_min` and
`g_max`.
dimension has a dimensionality 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)
@ -156,9 +161,6 @@ Data specific to neutron scattering for the temperature <TTT>K
for every Legendre moment or histogram/tabular bin. This dataset is
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]`
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_min` to `g_max`.
The pre-flattened array is shapes consistent with `scatter_matrix`
except the "[Order]" dimension in `scatter_shape` is ignored since
this data is assumed isotropic.

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

Binary file not shown.

4
openmc/mesh.py
View file

@ -286,8 +286,8 @@ class Mesh(object):
openmc.XPlane(x0=self.upper_right[0],
boundary_type=bc[1])]
if len(self.dimension) == 1:
yplanes = [openmc.YPlane(y0=1000., boundary_type='reflective'),
openmc.YPlane(y0=1000., boundary_type='reflective')]
yplanes = [openmc.YPlane(y0=-1e10, boundary_type='reflective'),
openmc.YPlane(y0=1e10, boundary_type='reflective')]
else:
yplanes = [openmc.YPlane(y0=self.lower_left[1],
boundary_type=bc[2]),

60
openmc/mgxs/library.py
View file

@ -186,14 +186,14 @@ class Library(object):
@property
def domains(self):
if self._domains is 'all':
if self.domain_type is 'material':
if self._domains == 'all':
if self.domain_type == '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 is 'universe':
elif self.domain_type == 'universe':
return self.openmc_geometry.get_all_universes()
elif self.domain_type is 'mesh':
elif self.domain_type == '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 is 'all':
if mgxs_types == '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 is 'mesh':
if by_nuclide == True and self.domain_type == '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 is 'mesh':
if self.by_nuclide == True and domain_type == '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 is 'all':
if domains == 'all':
self._domains = domains
# User specified a list of material, cell or universe domains
else:
if self.domain_type is 'material':
if self.domain_type == '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 is 'universe':
elif self.domain_type == 'universe':
cv.check_iterable_type('domain', domains, openmc.Universe)
all_domains = self.openmc_geometry.get_all_universes()
elif self.domain_type is 'mesh':
elif self.domain_type == '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 is 'P0' and self.legendre_order > 0:
if correction == '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 is 'P0' and legendre_order > 0:
if self.correction == '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 is 'k-eigenvalue':
if statepoint.run_mode == '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 is 'material':
if self.domain_type == 'material':
cv.check_type('domain', domain, (openmc.Material, Integral))
elif self.domain_type is 'cell' or self.domain_type is 'distribcell':
elif self.domain_type == 'cell' or self.domain_type == 'distribcell':
cv.check_type('domain', domain, (openmc.Cell, Integral))
elif self.domain_type is 'universe':
elif self.domain_type == 'universe':
cv.check_type('domain', domain, (openmc.Universe, Integral))
elif self.domain_type is 'mesh':
elif self.domain_type == '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 is 'distribcell':
if subdomain_avg_library.domain_type == '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 is 'distribcell':
if mgxs.domain_type == '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 is 'avg':
if subdomains == 'avg':
mgxs = mgxs.get_subdomain_avg_xs()
mgxs.build_hdf5_store(filename, directory, xs_type=xs_type,
@ -896,7 +896,7 @@ class Library(object):
# Build & add metadata to XSdata object
name = xsdata_name
if nuclide is not 'total':
if nuclide != 'total':
name += '_' + nuclide
xsdata = openmc.XSdata(name, self.energy_groups)
xsdata.delayed_groups = self.num_delayed_groups
@ -913,7 +913,7 @@ class Library(object):
self.scatter_format = 'legendre'
self.representation = 'isotropic'
if nuclide is not 'total':
if nuclide != 'total':
xsdata.atomic_weight_ratio = self._nuclides[nuclide][1]
if subdomain is None:
@ -922,7 +922,7 @@ class Library(object):
subdomain = [subdomain]
# Now get xs data itself
if 'nu-transport' in self.mgxs_types and self.correction is 'P0':
if 'nu-transport' in self.mgxs_types and self.correction == 'P0':
mymgxs = self.get_mgxs(domain, 'nu-transport')
xsdata.set_total_mgxs(mymgxs, xs_type=xs_type, nuclide=[nuclide],
subdomain=subdomain)
@ -1110,7 +1110,7 @@ class Library(object):
mgxs_file = openmc.MGXSLibrary(self.energy_groups,
delayed_groups=self.num_delayed_groups)
if self.domain_type is 'mesh':
if self.domain_type == 'mesh':
# Create the xsdata objects and add to the mgxs_file
i = 0
for domain in self.domains:
@ -1143,7 +1143,7 @@ class Library(object):
xsdata_name = 'set' + str(i + 1)
else:
xsdata_name = xsdata_names[i]
if nuclide is not 'total':
if nuclide != 'total':
xsdata_name += '_' + nuclide
xsdata = self.get_xsdata(domain, xsdata_name,
@ -1213,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 is 'mesh':
if self.domain_type == '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 is 'mesh':
if self.domain_type == 'mesh':
root = openmc.Universe(name='root', universe_id=0)
# Add cells representative of the mesh with reflective BC
@ -1266,13 +1266,13 @@ class Library(object):
materials.append(material)
# Differentiate Geometry with new Material
if self.domain_type is 'material':
if self.domain_type == 'material':
# Fill all appropriate Cells with new Material
for cell in all_cells:
if cell.fill.id == domain.id:
cell.fill = material
elif self.domain_type is 'cell':
elif self.domain_type == 'cell':
for cell in all_cells:
if cell.id == domain.id:
cell.fill = material

63
openmc/mgxs_library.py
View file

@ -8,7 +8,7 @@ import h5py
import openmc
import openmc.mgxs
from openmc.checkvalue import check_type, check_value, check_greater_than, \
check_less_than, check_iterable_type
check_iterable_type
if sys.version_info[0] >= 3:
basestring = str
@ -16,6 +16,7 @@ if sys.version_info[0] >= 3:
# Supported incoming particle MGXS angular treatment representations
_REPRESENTATIONS = ['isotropic', 'angle']
_SCATTER_TYPES = ['tabular', 'legendre', 'histogram']
_SCATTER_SHAPES = ["[Order][G][G']"]
class XSdata(object):
@ -55,6 +56,8 @@ class XSdata(object):
Whether or not this is a fissionable data set.
scatter_format : {'legendre', 'histogram', or 'tabular'}
Angular distribution representation (legendre, histogram, or tabular)
scatter_shapes : {"[Order][G][G']"}
Dimensionality of the scattering and multiplicity matrices
order : int
Either the Legendre order, number of bins, or number of points used to
describe the angular distribution associated with each group-to-group
@ -129,6 +132,7 @@ class XSdata(object):
self._atomic_weight_ratio = None
self._fissionable = False
self._scatter_format = 'legendre'
self._scatter_shape = "[Order][G][G']"
self._order = None
self._num_polar = None
self._num_azimuthal = None
@ -180,6 +184,10 @@ class XSdata(object):
def scatter_format(self):
return self._scatter_format
@property
def scatter_shape(self):
return self._scatter_shape
@property
def order(self):
return self._order
@ -302,6 +310,12 @@ class XSdata(object):
check_value('scatter_format', scatter_format, _SCATTER_TYPES)
self._scatter_format = scatter_format
@scatter_shape.setter
def scatter_shape(self, scatter_shape):
# check to see it is of a valid type and value
check_value('scatter_shape', scatter_shape, _SCATTER_SHAPES)
self._scatter_shape = scatter_shape
@order.setter
def order(self, order):
@ -967,10 +981,10 @@ class XSdata(object):
check_value('temperature', temperature, self.temperatures)
i = np.where(self.temperatures == temperature)[0][0]
if self.representation is 'isotropic':
if self.representation == 'isotropic':
self._total[i] = total.get_xs(nuclides=nuclide, xs_type=xs_type,
subdomains=subdomain)
elif self.representation is 'angle':
elif self.representation == 'angle':
msg = 'Angular-Dependent MGXS have not yet been implemented'
raise ValueError(msg)
@ -1013,11 +1027,11 @@ class XSdata(object):
check_value('temperature', temperature, self.temperatures)
i = np.where(self.temperatures == temperature)[0][0]
if self.representation is 'isotropic':
if self.representation == 'isotropic':
self._absorption[i] = absorption.get_xs(nuclides=nuclide,
xs_type=xs_type,
subdomains=subdomain)
elif self.representation is 'angle':
elif self.representation == 'angle':
msg = 'Angular-Dependent MGXS have not yet been implemented'
raise ValueError(msg)
@ -1060,11 +1074,11 @@ class XSdata(object):
check_value('temperature', temperature, self.temperatures)
i = np.where(self.temperatures == temperature)[0][0]
if self.representation is 'isotropic':
if self.representation == 'isotropic':
self._fission[i] = fission.get_xs(nuclides=nuclide,
xs_type=xs_type,
subdomains=subdomain)
elif self.representation is 'angle':
elif self.representation == 'angle':
msg = 'Angular-Dependent MGXS have not yet been implemented'
raise ValueError(msg)
@ -1108,11 +1122,11 @@ class XSdata(object):
check_value('temperature', temperature, self.temperatures)
i = np.where(self.temperatures == temperature)[0][0]
if self.representation is 'isotropic':
if self.representation == 'isotropic':
self._nu_fission[i] = nu_fission.get_xs(nuclides=nuclide,
xs_type=xs_type,
subdomains=subdomain)
elif self.representation is 'angle':
elif self.representation == 'angle':
msg = 'Angular-Dependent MGXS have not yet been implemented'
raise ValueError(msg)
@ -1272,11 +1286,11 @@ class XSdata(object):
check_value('temperature', temperature, self.temperatures)
i = np.where(self.temperatures == temperature)[0][0]
if self.representation is 'isotropic':
if self.representation == 'isotropic':
self._kappa_fission[i] = k_fission.get_xs(nuclides=nuclide,
xs_type=xs_type,
subdomains=subdomain)
elif self.representation is 'angle':
elif self.representation == 'angle':
msg = 'Angular-Dependent MGXS have not yet been implemented'
raise ValueError(msg)
@ -1317,10 +1331,10 @@ class XSdata(object):
check_value('temperature', temperature, self.temperatures)
i = np.where(self.temperatures == temperature)[0][0]
if self.representation is 'isotropic':
if self.representation == 'isotropic':
self._chi[i] = chi.get_xs(nuclides=nuclide,
xs_type=xs_type, subdomains=subdomain)
elif self.representation is 'angle':
elif self.representation == 'angle':
msg = 'Angular-Dependent MGXS have not yet been implemented'
raise ValueError(msg)
@ -1551,7 +1565,7 @@ class XSdata(object):
check_type('temperature', temperature, Real)
check_value('temperature', temperature, self.temperatures)
if self.scatter_format is not 'legendre':
if self.scatter_format != 'legendre':
msg = 'Anisotropic scattering representations other than ' \
'Legendre expansions have not yet been implemented in ' \
'openmc.mgxs.'
@ -1568,7 +1582,7 @@ class XSdata(object):
[self.order])
i = np.where(self.temperatures == temperature)[0][0]
if self.representation is 'isotropic':
if self.representation == 'isotropic':
# Get the scattering orders in the outermost dimension
self._scatter_matrix[i] = np.zeros((self.num_orders,
self.energy_groups.num_groups,
@ -1578,7 +1592,7 @@ class XSdata(object):
scatter.get_xs(nuclides=nuclide, xs_type=xs_type,
moment=moment, subdomains=subdomain)
elif self.representation is 'angle':
elif self.representation == 'angle':
msg = 'Angular-Dependent MGXS have not yet been implemented'
raise ValueError(msg)
@ -1646,7 +1660,7 @@ class XSdata(object):
['universe', 'cell', 'material', 'mesh'])
i = np.where(self.temperatures == temperature)[0][0]
if self.representation is 'isotropic':
if self.representation == 'isotropic':
nuscatt = nuscatter.get_xs(nuclides=nuclide,
xs_type=xs_type, moment=0,
subdomains=subdomain)
@ -1657,7 +1671,7 @@ class XSdata(object):
xs_type=xs_type, moment=0,
subdomains=subdomain)
self._multiplicity_matrix[i] = np.divide(nuscatt, scatt)
elif self.representation is 'angle':
elif self.representation == 'angle':
msg = 'Angular-Dependent MGXS have not yet been implemented'
raise ValueError(msg)
self._multiplicity_matrix[i] = \
@ -1679,9 +1693,8 @@ class XSdata(object):
grp.attrs['fissionable'] = self.fissionable
if self.representation is not None:
grp.attrs['representation'] = np.array(self.representation,
dtype='S')
if self.representation is 'angle':
grp.attrs['representation'] = np.string_(self.representation)
if self.representation == 'angle':
if self.num_azimuthal is not None:
grp.attrs['num-azimuthal'] = self.num_azimuthal
@ -1689,9 +1702,9 @@ 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,
dtype='S')
grp.attrs['scatter_format'] = np.string_(self.scatter_format)
if self.scatter_shape is not None:
grp.attrs['scatter_shape'] = np.string_(self.scatter_shape)
if self.order is not None:
grp.attrs['order'] = self.order
@ -1810,7 +1823,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 is 'angle':
elif self.representation == 'angle':
Np = self.num_polar
Na = self.num_azimuthal
g_out_bounds = np.zeros((Np, Na, G, 2), dtype=np.int)

14
scripts/openmc-update-mgxs
View file

@ -47,14 +47,12 @@ def parse_args():
help='input XML file')
parser.add_argument('-o', '--output', nargs='?', default='',
help='output file, in HDF5 format')
parser.add_argument('-c', '--compression', type=int,
help='HDF5 Compression Level')
args = vars(parser.parse_args())
if args['output'] is '':
if args['output'] == '':
filename = args['input'].name
extension = filenameos.path.splitext()
if extension is '.xml':
if extension == '.xml':
filename = filename[:filename.rfind('.')] + '.h5'
args['output'] = filename
@ -75,7 +73,7 @@ def get_data(element, entry):
return value
if __name__ is '__main__':
if __name__ == '__main__':
args = parse_args()
# Parse the XML data.
@ -117,7 +115,7 @@ if __name__ is '__main__':
representation = get_data(xsdata_elem, 'representation')
if representation is None:
representation = 'isotropic'
if representation is 'angle':
if representation == 'angle':
n_azi = int(get_data(xsdata_elem, 'num_azimuthal'))
n_pol = int(get_data(xsdata_elem, 'num_polar'))
@ -146,14 +144,14 @@ if __name__ is '__main__':
representation=representation))
if awr is not None:
xsd[-1].atomic_weight_ratio = awr
if representation is 'angle':
if representation == '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 is 'legendre':
if scatter_format == 'legendre':
order_dim = order + 1
else:
order_dim = order

6
setup.py Normal file → Executable file
View file

@ -1,6 +1,7 @@
#!/usr/bin/env python
import glob
import numpy as np
try:
from setuptools import setup
have_setuptools = True
@ -57,10 +58,7 @@ if have_setuptools:
# If Cython is present, add resonance reconstruction capability
if have_cython:
kwargs.update({
'ext_modules': cythonize('openmc/data/reconstruct.pyx')
})
import numpy as np
kwargs.update({
'ext_modules': cythonize('openmc/data/reconstruct.pyx'),
'include_dirs': [np.get_include()]
})

2
src/input_xml.F90
View file

@ -3293,7 +3293,7 @@ contains
! Search through nuclides
pair_list => nuclide_dict % keys()
do while (associated(pair_list))
if (starts_with(pair_list % key, word)) then
if (trim(pair_list % key) == trim(word)) then
word = pair_list % key(1:150)
exit
end if

7
src/mgxs_header.F90
View file

@ -351,6 +351,13 @@ module mgxs_header
else
this % scatter_format = ANGLE_LEGENDRE
end if
if (attribute_exists(xs_id, "scatter_shape")) then
call read_attribute(temp_str, xs_id, "scatter_shape")
temp_str = trim(temp_str)
if (to_lower(temp_str) /= "[order][g][g']") then
call fatal_error("Invalid scatter_shape option!")
end if
end if
if (attribute_exists(xs_id, "fissionable")) then
call read_attribute(this % fissionable, xs_id, "fissionable")
else

9
src/output.F90
View file

@ -856,8 +856,13 @@ contains
write(UNIT=unit_tally, FMT='(1X,2A,1X,A)') repeat(" ", indent), &
"Total Material"
else
write(UNIT=unit_tally, FMT='(1X,2A,1X,A)') repeat(" ", indent), &
trim(nuclides(i_nuclide) % name)
if (run_CE) then
write(UNIT=unit_tally, FMT='(1X,2A,1X,A)') repeat(" ", indent), &
trim(nuclides(i_nuclide) % name)
else
write(UNIT=unit_tally, FMT='(1X,2A,1X,A)') repeat(" ", indent), &
trim(nuclides_MG(i_nuclide) % obj % name)
end if
end if
indent = indent + 2

17
src/state_point.F90
View file

@ -251,11 +251,20 @@ contains
allocate(str_array(tally % n_nuclide_bins))
NUCLIDE_LOOP: do j = 1, tally % n_nuclide_bins
if (tally % nuclide_bins(j) > 0) then
i_xs = index(nuclides(tally % nuclide_bins(j)) % name, '.')
if (i_xs > 0) then
str_array(j) = nuclides(tally % nuclide_bins(j)) % name(1 : i_xs-1)
if (run_CE) then
i_xs = index(nuclides(tally % nuclide_bins(j)) % name, '.')
if (i_xs > 0) then
str_array(j) = nuclides(tally % nuclide_bins(j)) % name(1 : i_xs-1)
else
str_array(j) = nuclides(tally % nuclide_bins(j)) % name
end if
else
str_array(j) = nuclides(tally % nuclide_bins(j)) % name
i_xs = index(nuclides_MG(tally % nuclide_bins(j)) % obj % name, '.')
if (i_xs > 0) then
str_array(j) = nuclides_MG(tally % nuclide_bins(j)) % obj % name(1 : i_xs-1)
else
str_array(j) = nuclides_MG(tally % nuclide_bins(j)) % obj % name
end if
end if
else
str_array(j) = 'total'

29
src/tally.F90
View file

@ -1164,6 +1164,8 @@ contains
! Do same for nucxs, point it to the microscopic nuclide data of interest
if (i_nuclide > 0) then
nucxs => nuclides_MG(i_nuclide) % obj
! And since we haven't calculated this temperature index yet, do so now
call nucxs % find_temperature(p % sqrtkT)
end if
i = 0
@ -2301,17 +2303,24 @@ contains
i_nuclide = t % nuclide_bins(k)
! Check to see if this nuclide was in the material of our collision.
do m = 1, mat % n_nuclides
if (mat % nuclide(m) == i_nuclide) then
atom_density = mat % atom_density(m)
exit
end if
end do
if (i_nuclide > 0) then
atom_density = -ONE
! Check to see if this nuclide was in the material of our collision
do m = 1, mat % n_nuclides
if (mat % nuclide(m) == i_nuclide) then
atom_density = mat % atom_density(m)
exit
end if
end do
else
atom_density = ZERO
end if
! Determine score for each bin
call score_general(p, t, (k-1)*t % n_score_bins, filter_index, &
i_nuclide, atom_density, filter_weight)
! If we found the nuclide, determine the score for each bin
if (atom_density >= ZERO) then
call score_general(p, t, (k-1)*t % n_score_bins, filter_index, &
i_nuclide, atom_density, filter_weight)
end if
end do NUCLIDE_LOOP

BIN
tests/1d_mgxs.h5
View file

Binary file not shown.

90
tests/input_set.py
View file

@ -1,8 +1,9 @@
import numpy as np
import openmc
from openmc.source import Source
from openmc.stats import Box
import numpy as np
class InputSet(object):
def __init__(self):
@ -728,60 +729,59 @@ class AssemblyInputSet(object):
class MGInputSet(InputSet):
def build_default_materials_and_geometry(self):
def build_default_materials_and_geometry(self, reps=None, as_macro=True):
# Define materials needed for 1D/1G slab problem
uo2_data = openmc.Macroscopic('uo2_iso')
uo2 = openmc.Material(name='UO2', material_id=1)
uo2.set_density('macro', 1.0)
uo2.add_macroscopic(uo2_data)
mat_names = ['uo2', 'clad', 'lwtr']
mgxs_reps = ['ang', 'ang_mu', 'iso', 'iso_mu']
clad_data = openmc.Macroscopic('clad_ang_mu')
clad = openmc.Material(name='Clad', material_id=2)
clad.set_density('macro', 1.0)
clad.add_macroscopic(clad_data)
if reps is None:
reps = mgxs_reps
water_data = openmc.Macroscopic('lwtr_iso_mu')
water = openmc.Material(name='LWTR', material_id=3)
water.set_density('macro', 1.0)
water.add_macroscopic(water_data)
xs = []
mats = []
i = 0
for mat in mat_names:
for rep in reps:
i += 1
if as_macro:
xs.append(openmc.Macroscopic(mat + '_' + rep))
mats.append(openmc.Material(name=str(i)))
mats[-1].set_density('macro', 1.)
mats[-1].add_macroscopic(xs[-1])
else:
xs.append(openmc.Nuclide(mat + '_' + rep))
mats.append(openmc.Material(name=str(i)))
mats[-1].set_density('atom/b-cm', 1.)
mats[-1].add_nuclide(xs[-1].name, 1.0, 'ao')
# Define the materials file.
self.materials += (uo2, clad, water)
# Define the materials file
self.xs_data = xs
self.materials += mats
# Define surfaces.
# Assembly/Problem Boundary
left = openmc.XPlane(x0=0.0, surface_id=200,
boundary_type='reflective')
right = openmc.XPlane(x0=10.0, surface_id=201,
boundary_type='reflective')
bottom = openmc.YPlane(y0=0.0, surface_id=300,
boundary_type='reflective')
top = openmc.YPlane(y0=10.0, surface_id=301,
boundary_type='reflective')
left = openmc.XPlane(x0=0.0, boundary_type='reflective')
right = openmc.XPlane(x0=10.0, boundary_type='reflective')
bottom = openmc.YPlane(y0=0.0, boundary_type='reflective')
top = openmc.YPlane(y0=10.0, boundary_type='reflective')
# for each material add a plane
planes = [openmc.ZPlane(z0=0.0, boundary_type='reflective')]
dz = round(5. / float(len(mats)), 4)
for i in range(len(mats) - 1):
planes.append(openmc.ZPlane(z0=dz * float(i + 1)))
planes.append(openmc.ZPlane(z0=5.0, boundary_type='reflective'))
down = openmc.ZPlane(z0=0.0, surface_id=0,
boundary_type='reflective')
fuel_clad_intfc = openmc.ZPlane(z0=2.0, surface_id=1)
clad_lwtr_intfc = openmc.ZPlane(z0=2.4, surface_id=2)
up = openmc.ZPlane(z0=5.0, surface_id=3,
boundary_type='reflective')
# Define cells
c1 = openmc.Cell(cell_id=1)
c1.region = +left & -right & +bottom & -top & +down & -fuel_clad_intfc
c1.fill = uo2
c2 = openmc.Cell(cell_id=2)
c2.region = +left & -right & +bottom & -top & +fuel_clad_intfc & -clad_lwtr_intfc
c2.fill = clad
c3 = openmc.Cell(cell_id=3)
c3.region = +left & -right & +bottom & -top & +clad_lwtr_intfc & -up
c3.fill = water
# Define cells for each material
cells = []
xy = +left & -right & +bottom & -top
for i, mat in enumerate(mats):
cells.append(openmc.Cell())
cells[-1].region = xy & +planes[i] & -planes[i + 1]
cells[-1].fill = mat
# Define root universe.
root = openmc.Universe(universe_id=0, name='root universe')
root.add_cells((c1, c2, c3))
root.add_cells(cells)
# Assign root universe to geometry
self.geometry.root_universe = root
@ -791,7 +791,7 @@ class MGInputSet(InputSet):
self.settings.inactive = 5
self.settings.particles = 100
self.settings.source = Source(space=Box([0.0, 0.0, 0.0],
[10.0, 10.0, 2.0]))
[10.0, 10.0, 5.]))
self.settings.energy_mode = "multi-group"
self.settings.cross_sections = "../1d_mgxs.h5"

2
tests/test_mg_basic/inputs_true.dat
View file

@ -1 +1 @@
002fff0d8ba33340c9039f7650ece41bd26324d2e5c2816acc232adb298ae9dbd0873c5e28520df33ebe5af4c58e2ad338547ca337b2b42f53d255bff1ae314b
115218221500e60ea43145875324a9a6800366fcc97357ef4cd8182f9253e114cbe1c7481dd3bdd8d7dce5983e2cb757ae0a2e05ede236d94cba7909703ba7f7

2
tests/test_mg_basic/results_true.dat
View file

@ -1,2 +1,2 @@
k-combined:
1.070804E+00 1.141438E-02
5.788739E-01 1.945182E-02

2
tests/test_mg_max_order/inputs_true.dat
View file

@ -1 +1 @@
8752d3167989d876c6f29b4b1ca7eed4b91793d8579c92de85633fd90235c67d63ee76d49081006047e0078f6f86eb6e3420c1f8fe44564e7a1fb62e828587d0
139d760cd83eab001ed3fd52d7146fb5f30b7021b26ada6a1e425f5446185ac6fa61ac2ac58aa9e895981e10540160b7ae7428c833024679ef1ca762715e2d51

2
tests/test_mg_max_order/results_true.dat
View file

@ -1,2 +1,2 @@
k-combined:
1.093516E+00 1.593236E-02
1.077247E+00 2.412834E-02

75
tests/test_mg_max_order/test_mg_max_order.py
View file

@ -2,83 +2,26 @@
import os
import sys
sys.path.insert(0, os.pardir)
from testing_harness import TestHarness, PyAPITestHarness
from testing_harness import PyAPITestHarness
from input_set import MGInputSet
import openmc
class MGNuclideInputSet(MGInputSet):
def build_default_materials_and_geometry(self):
# Define materials needed for 1D/1G slab problem
uo2_data = openmc.Macroscopic('uo2_iso')
uo2 = openmc.Material(name='UO2', material_id=1)
uo2.set_density('macro', 1.0)
uo2.add_macroscopic(uo2_data)
clad_data = openmc.Macroscopic('clad_iso')
clad = openmc.Material(name='Clad', material_id=2)
clad.set_density('macro', 1.0)
clad.add_macroscopic(clad_data)
water_data = openmc.Macroscopic('lwtr_iso')
water = openmc.Material(name='LWTR', material_id=3)
water.set_density('macro', 1.0)
water.add_macroscopic(water_data)
# Define the materials file.
self.materials += (uo2, clad, water)
# Define surfaces.
# Assembly/Problem Boundary
left = openmc.XPlane(x0=0.0, surface_id=200,
boundary_type='reflective')
right = openmc.XPlane(x0=10.0, surface_id=201,
boundary_type='reflective')
bottom = openmc.YPlane(y0=0.0, surface_id=300,
boundary_type='reflective')
top = openmc.YPlane(y0=10.0, surface_id=301,
boundary_type='reflective')
down = openmc.ZPlane(z0=0.0, surface_id=0,
boundary_type='reflective')
fuel_clad_intfc = openmc.ZPlane(z0=2.0, surface_id=1)
clad_lwtr_intfc = openmc.ZPlane(z0=2.4, surface_id=2)
up = openmc.ZPlane(z0=5.0, surface_id=3,
boundary_type='reflective')
# Define cells
c1 = openmc.Cell(cell_id=1)
c1.region = +left & -right & +bottom & -top & +down & -fuel_clad_intfc
c1.fill = uo2
c2 = openmc.Cell(cell_id=2)
c2.region = +left & -right & +bottom & -top & +fuel_clad_intfc & -clad_lwtr_intfc
c2.fill = clad
c3 = openmc.Cell(cell_id=3)
c3.region = +left & -right & +bottom & -top & +clad_lwtr_intfc & -up
c3.fill = water
# Define root universe.
root = openmc.Universe(universe_id=0, name='root universe')
root.add_cells((c1,c2,c3))
# Define the geometry file.
geometry = openmc.Geometry()
geometry.root_universe = root
self.geometry = geometry
class MGMaxOrderTestHarness(PyAPITestHarness):
def __init__(self, statepoint_name, tallies_present, mg=False):
PyAPITestHarness.__init__(self, statepoint_name, tallies_present)
self._input_set = MGNuclideInputSet()
self._input_set = MGInputSet()
def _build_inputs(self):
"""Write input XML files."""
reps = ['iso']
self._input_set.build_default_materials_and_geometry(reps=reps)
self._input_set.build_default_settings()
# Set P1 scattering
self._input_set.settings.max_order = 1
# Call standard input build
super(MGMaxOrderTestHarness, self)._build_inputs()
self._input_set.export()
if __name__ == '__main__':
harness = MGMaxOrderTestHarness('statepoint.10.*', False, mg=True)

2
tests/test_mg_nuclide/inputs_true.dat
View file

@ -1 +1 @@
10b64fa97cc4feb1b1b8ba401b25b68ab682aba55ad4303a76a9a293ae79cdfff2f65d74310eb917da0aac243c24713d91eafae7ec5f0c4ae535d49fc928a869
6d03b988671fe21dfe2a908685536edb3059e5bdda7974ddfa65c7f63cecb7f6b3fce49bc1533ef2787e284121779653cb14a66cfd709258e4b07d8bd7571843

2
tests/test_mg_nuclide/results_true.dat
View file

@ -1,2 +1,2 @@
k-combined:
1.070804E+00 1.141438E-02
5.788739E-01 1.945182E-02

73
tests/test_mg_nuclide/test_mg_nuclide.py
View file

@ -2,81 +2,22 @@
import os
import sys
sys.path.insert(0, os.pardir)
from testing_harness import TestHarness, PyAPITestHarness
from testing_harness import PyAPITestHarness
from input_set import MGInputSet
import openmc
class MGNuclideInputSet(MGInputSet):
def build_default_materials_and_geometry(self):
# Define materials needed for 1D/1G slab problem
# This time do using nuclide, not macroscopic
uo2 = openmc.Material(name='UO2', material_id=1)
uo2.set_density('sum', 1.0)
uo2.add_nuclide("uo2_iso", 1.0)
clad = openmc.Material(name='Clad', material_id=2)
clad.set_density('sum', 1.0)
clad.add_nuclide("clad_ang_mu", 1.0)
# water_data = openmc.Nuclide('lwtr_iso_mu', '71c')
water = openmc.Material(name='LWTR', material_id=3)
water.set_density('sum', 1.0)
water.add_nuclide("lwtr_iso_mu", 1.0)
# Define the materials file.
self.materials.default_xs = '71c'
self.materials += (uo2, clad, water)
# Define surfaces.
# Assembly/Problem Boundary
left = openmc.XPlane(x0=0.0, surface_id=200,
boundary_type='reflective')
right = openmc.XPlane(x0=10.0, surface_id=201,
boundary_type='reflective')
bottom = openmc.YPlane(y0=0.0, surface_id=300,
boundary_type='reflective')
top = openmc.YPlane(y0=10.0, surface_id=301,
boundary_type='reflective')
down = openmc.ZPlane(z0=0.0, surface_id=0,
boundary_type='reflective')
fuel_clad_intfc = openmc.ZPlane(z0=2.0, surface_id=1)
clad_lwtr_intfc = openmc.ZPlane(z0=2.4, surface_id=2)
up = openmc.ZPlane(z0=5.0, surface_id=3,
boundary_type='reflective')
# Define cells
c1 = openmc.Cell(cell_id=1)
c1.region = +left & -right & +bottom & -top & +down & -fuel_clad_intfc
c1.fill = uo2
c2 = openmc.Cell(cell_id=2)
c2.region = +left & -right & +bottom & -top & +fuel_clad_intfc & -clad_lwtr_intfc
c2.fill = clad
c3 = openmc.Cell(cell_id=3)
c3.region = +left & -right & +bottom & -top & +clad_lwtr_intfc & -up
c3.fill = water
# Define root universe.
root = openmc.Universe(universe_id=0, name='root universe')
root.add_cells((c1,c2,c3))
# Define the geometry file.
geometry = openmc.Geometry()
geometry.root_universe = root
self.geometry = geometry
class MGNuclideTestHarness(PyAPITestHarness):
def __init__(self, statepoint_name, tallies_present, mg=False):
PyAPITestHarness.__init__(self, statepoint_name, tallies_present)
self._input_set = MGNuclideInputSet()
self._input_set = MGInputSet()
def _build_inputs(self):
super(MGNuclideTestHarness, self)._build_inputs()
"""Write input XML files."""
self._input_set.build_default_materials_and_geometry(as_macro=False)
self._input_set.build_default_settings()
self._input_set.export()
if __name__ == '__main__':

2
tests/test_mg_tallies/inputs_true.dat
View file

@ -1 +1 @@
615fe0688980afc02bc11f41ac065b9479d2f6618acf6d20fcb91cd78209aed1d8fbaf73b06f153368a16aa90249d6d3ce06f840d9071a40b2118d36178fd84f
8c496e13f9456bedcc0789ffdca81913596cca3c7867717ced38c3063b6cdc8949bd204c98723b72c1b456e3fc85fc90e6710482b0c74e5e29283902a5d9ed50

2907
tests/test_mg_tallies/results_true.dat
View file

File diff suppressed because it is too large Load diff

74
tests/test_mg_tallies/test_mg_tallies.py
View file

@ -3,38 +3,84 @@
import os
import sys
sys.path.insert(0, os.pardir)
from testing_harness import TestHarness, PyAPITestHarness
from testing_harness import HashedPyAPITestHarness
import openmc
class MGTalliesTestHarness(PyAPITestHarness):
class MGTalliesTestHarness(HashedPyAPITestHarness):
def _build_inputs(self):
"""Write input XML files."""
self._input_set.build_default_materials_and_geometry(as_macro=False)
self._input_set.build_default_settings()
# Instantiate a tally mesh
mesh = openmc.Mesh(mesh_id=1)
mesh.type = 'regular'
mesh.dimension = [17, 17, 1]
mesh.dimension = [1, 1, 10]
mesh.lower_left = [0.0, 0.0, 0.0]
mesh.upper_right = [21.42, 21.42, 100.0]
mesh.upper_right = [10, 10, 5]
# Instantiate some tally filters
energy_filter = openmc.EnergyFilter([0.0, 20.0])
energyout_filter = openmc.EnergyoutFilter([0.0, 20.0])
mesh_filter = openmc.MeshFilter(mesh)
mat_filter = openmc.MaterialFilter([1,2,3])
mat_ids = [mat.id for mat in self._input_set.materials]
mat_filter = openmc.MaterialFilter(mat_ids)
tally1 = openmc.Tally(tally_id=1)
tally1.filters = [mesh_filter]
tally1.scores = ['total', 'absorption', 'flux',
'fission', 'nu-fission']
nuclides = [xs.name for xs in self._input_set.xs_data]
tally2 = openmc.Tally(tally_id=2)
tally2.filters = [mat_filter, energy_filter, energyout_filter]
tally2.scores = ['scatter', 'nu-scatter']
scores= {False: ['total', 'absorption', 'flux', 'fission', 'nu-fission'],
True: ['total', 'absorption', 'fission', 'nu-fission']}
self._input_set.tallies = openmc.Tallies([tally1, tally2])
tallies = []
for do_nuclides in [False, True]:
tallies.append(openmc.Tally())
tallies[-1].filters = [mesh_filter]
tallies[-1].estimator = 'analog'
tallies[-1].scores = scores[do_nuclides]
if do_nuclides:
tallies[-1].nuclides = nuclides
super(MGTalliesTestHarness, self)._build_inputs()
tallies.append(openmc.Tally())
tallies[-1].filters = [mesh_filter]
tallies[-1].estimator = 'tracklength'
tallies[-1].scores = scores[do_nuclides]
if do_nuclides:
tallies[-1].nuclides = nuclides
tallies.append(openmc.Tally())
tallies[-1].filters = [mat_filter, energy_filter]
tallies[-1].estimator = 'analog'
tallies[-1].scores = scores[do_nuclides] + ['scatter',
'nu-scatter']
if do_nuclides:
tallies[-1].nuclides = nuclides
tallies.append(openmc.Tally())
tallies[-1].filters = [mat_filter, energy_filter]
tallies[-1].estimator = 'collision'
tallies[-1].scores = scores[do_nuclides]
if do_nuclides:
tallies[-1].nuclides = nuclides
tallies.append(openmc.Tally())
tallies[-1].filters = [mat_filter, energy_filter]
tallies[-1].estimator = 'tracklength'
tallies[-1].scores = scores[do_nuclides]
if do_nuclides:
tallies[-1].nuclides = nuclides
tallies.append(openmc.Tally())
tallies[-1].filters = [mat_filter, energy_filter,
energyout_filter]
tallies[-1].scores = ['scatter', 'nu-scatter', 'nu-fission']
if do_nuclides:
tallies[-1].nuclides = nuclides
self._input_set.tallies = openmc.Tallies(tallies)
self._input_set.export()
def _cleanup(self):
super(MGTalliesTestHarness, self)._cleanup()