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

Merge pull request #778 from nelsonag/mg_angle

Browse source 
Incorporation of angle-dependent MGXS capabilities in `openmc.MGXS`
This commit is contained in:
Will Boyd 2017-02-14 21:24:20 -05:00 • committed by GitHub
commit 2a86def754
20 changed files with 1815 additions and 973 deletions
Download patch file Download diff file Expand all files Collapse all files

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

@ -90,7 +90,7 @@ 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].
[polar][azimuthal][groups].
- **absorption** (*double[]* or *double[][][]*) -- Absorption
cross section.
This is a 1-D vector if `representation` is "isotropic", or a 3-D
@ -100,19 +100,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 dataset
[polar][azimuthal][groups]. 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
[polar][azimuthal][groups]. 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
[polar][azimuthal][groups]. 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.
@ -122,8 +122,11 @@ Temperature-dependent data, provided for temperature <TTT>K.
spectra as well and thus will have one additional dimension
for the outgoing energy group. In this case, `nu-fission` has the
same dimensionality as `multiplicity matrix`.
- **inverse-velocity** (*double[]*) -- Average inverse velocity
for each of the groups in the library. This dataset is optional.
- **inverse-velocity** (*double[]* or *double[][][]*) --
Average inverse velocity for each of the groups in the library.
This dataset is optional. This is a 1-D vector if `representation`
is "isotropic", or a 3-D vector if `representation` is "angle"
with dimensions of [polar][azimuthal][groups].
**/<library name>/<TTT>K/scatter_data/**

205
docs/source/pythonapi/examples/mdgxs-part-i.ipynb
View file

File diff suppressed because one or more lines are too long

62
docs/source/pythonapi/examples/mgxs-part-i.ipynb
View file

@ -421,17 +421,13 @@
"OrderedDict([('flux', Tally\n",
" \tID =\t10000\n",
" \tName =\t\n",
" \tFilters =\t\n",
" \t\tCellFilter\t[1]\n",
" \t\tEnergyFilter\t[ 0.00000000e+00 6.25000000e-01 2.00000000e+07]\n",
" \tFilters =\tCellFilter, EnergyFilter\n",
" \tNuclides =\ttotal \n",
" \tScores =\t['flux']\n",
" \tEstimator =\ttracklength), ('absorption', Tally\n",
" \tID =\t10001\n",
" \tName =\t\n",
" \tFilters =\t\n",
" \t\tCellFilter\t[1]\n",
" \t\tEnergyFilter\t[ 0.00000000e+00 6.25000000e-01 2.00000000e+07]\n",
" \tFilters =\tCellFilter, EnergyFilter\n",
" \tNuclides =\ttotal \n",
" \tScores =\t['absorption']\n",
" \tEstimator =\ttracklength)])"
@ -521,12 +517,12 @@
" %%%%%%%%%%%\n",
"\n",
" | The OpenMC Monte Carlo Code\n",
" Copyright | 2011-2016 Massachusetts Institute of Technology\n",
" Copyright | 2011-2017 Massachusetts Institute of Technology\n",
" License | http://openmc.readthedocs.io/en/latest/license.html\n",
" Version | 0.8.0\n",
" Git SHA1 | da5563eddb5f2c2d6b2c9839d518de40962b78f2\n",
" Date/Time | 2016-10-31 12:23:45\n",
" OpenMP Threads | 4\n",
" Git SHA1 | 54b65c8bda6af5788bd762b8cf9855d1a8008238\n",
" Date/Time | 2017-02-12 13:36:24\n",
" OpenMP Threads | 8\n",
"\n",
" ===========================================================================\n",
" ========================> INITIALIZATION <=========================\n",
@ -536,11 +532,11 @@
" Reading geometry XML file...\n",
" Reading materials XML file...\n",
" Reading cross sections XML file...\n",
" Reading H1 from /home/romano/openmc/scripts/nndc_hdf5/H1.h5\n",
" Reading O16 from /home/romano/openmc/scripts/nndc_hdf5/O16.h5\n",
" Reading U235 from /home/romano/openmc/scripts/nndc_hdf5/U235.h5\n",
" Reading U238 from /home/romano/openmc/scripts/nndc_hdf5/U238.h5\n",
" Reading Zr90 from /home/romano/openmc/scripts/nndc_hdf5/Zr90.h5\n",
" Reading H1 from /opt/xsdata/nndc/H1.h5\n",
" Reading O16 from /opt/xsdata/nndc/O16.h5\n",
" Reading U235 from /opt/xsdata/nndc/U235.h5\n",
" Reading U238 from /opt/xsdata/nndc/U238.h5\n",
" Reading Zr90 from /opt/xsdata/nndc/Zr90.h5\n",
" Maximum neutron transport energy: 2.00000E+07 eV for H1\n",
" Reading tallies XML file...\n",
" Building neighboring cells lists for each surface...\n",
@ -611,20 +607,20 @@
"\n",
" =======================> TIMING STATISTICS <=======================\n",
"\n",
" Total time for initialization = 7.3858E-01 seconds\n",
" Reading cross sections = 5.3599E-01 seconds\n",
" Total time in simulation = 2.1031E+01 seconds\n",
" Time in transport only = 1.9960E+01 seconds\n",
" Time in inactive batches = 2.6543E+00 seconds\n",
" Time in active batches = 1.8376E+01 seconds\n",
" Time synchronizing fission bank = 5.8389E-03 seconds\n",
" Sampling source sites = 4.2676E-03 seconds\n",
" SEND/RECV source sites = 1.4523E-03 seconds\n",
" Time accumulating tallies = 1.5633E-04 seconds\n",
" Total time for finalization = 7.9179E-04 seconds\n",
" Total time elapsed = 2.1784E+01 seconds\n",
" Calculation Rate (inactive) = 9418.63 neutrons/second\n",
" Calculation Rate (active) = 5441.81 neutrons/second\n",
" Total time for initialization = 4.1327E-01 seconds\n",
" Reading cross sections = 3.2638E-01 seconds\n",
" Total time in simulation = 2.2324E+00 seconds\n",
" Time in transport only = 2.1226E+00 seconds\n",
" Time in inactive batches = 3.0650E-01 seconds\n",
" Time in active batches = 1.9259E+00 seconds\n",
" Time synchronizing fission bank = 2.7640E-03 seconds\n",
" Sampling source sites = 2.0198E-03 seconds\n",
" SEND/RECV source sites = 7.0929E-04 seconds\n",
" Time accumulating tallies = 4.5355E-05 seconds\n",
" Total time for finalization = 4.1885E-04 seconds\n",
" Total time elapsed = 2.6534E+00 seconds\n",
" Calculation Rate (inactive) = 81567.1 neutrons/second\n",
" Calculation Rate (active) = 51923.2 neutrons/second\n",
"\n",
" ============================> RESULTS <============================\n",
"\n",
@ -832,7 +828,7 @@
"cell_type": "code",
"execution_count": 20,
"metadata": {
"collapsed": true
"collapsed": false
},
"outputs": [],
"source": [
@ -915,7 +911,7 @@
" <td>2.000000e+07</td>\n",
" <td>total</td>\n",
" <td>(((total / flux) - (absorption / flux)) - (sca...</td>\n",
" <td>-7.771561e-16</td>\n",
" <td>-3.330669e-16</td>\n",
" <td>0.002570</td>\n",
" </tr>\n",
" </tbody>\n",
@ -929,7 +925,7 @@
"\n",
" score mean std. dev. \n",
"0 (((total / flux) - (absorption / flux)) - (sca... -2.66e-15 1.13e-02 \n",
"1 (((total / flux) - (absorption / flux)) - (sca... -7.77e-16 2.57e-03 "
"1 (((total / flux) - (absorption / flux)) - (sca... -3.33e-16 2.57e-03 "
]
},
"execution_count": 22,
@ -1192,7 +1188,7 @@
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.5.2"
"version": "3.6.0"
}
},
"nbformat": 4,

91
docs/source/pythonapi/examples/mgxs-part-ii.ipynb
View file

File diff suppressed because one or more lines are too long

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

File diff suppressed because one or more lines are too long

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

File diff suppressed because one or more lines are too long

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

@ -13,7 +13,7 @@ inactive = 10
particles = 1000
###############################################################################
# Exporting to OpenMC mgxs.xml file
# Exporting to OpenMC mgxs.h5 file
###############################################################################
# Instantiate the energy group data

8
openmc/filter.py
View file

@ -1008,7 +1008,7 @@ class EnergyFilter(RealFilter):
# them as necessary to account for other filters.
lo_bins = np.repeat(self.bins[:-1], self.stride)
hi_bins = np.repeat(self.bins[1:], self.stride)
tile_factor = data_size / len(lo_bins)
tile_factor = int(data_size / len(lo_bins))
lo_bins = np.tile(lo_bins, tile_factor)
hi_bins = np.tile(hi_bins, tile_factor)
@ -1453,7 +1453,7 @@ class PolarFilter(RealFilter):
msg = 'Unable to add bin edge "{0}" to a "{1}" ' \
'since it is less than 0'.format(edge, type(self))
raise ValueError(msg)
elif edge > np.pi:
elif not np.isclose(edge, np.pi) and edge > np.pi:
msg = 'Unable to add bin edge "{0}" to a "{1}" ' \
'since it is greater than pi'.format(edge, type(self))
raise ValueError(msg)
@ -1552,11 +1552,11 @@ class AzimuthalFilter(RealFilter):
'since it is a non-integer or floating point ' \
'value'.format(edge, type(self))
raise ValueError(msg)
elif edge < -np.pi:
elif not np.isclose(edge, -np.pi) and edge < -np.pi:
msg = 'Unable to add bin edge "{0}" to a "{1}" ' \
'since it is less than -pi'.format(edge, type(self))
raise ValueError(msg)
elif edge > np.pi:
elif not np.isclose(edge, np.pi) and edge > np.pi:
msg = 'Unable to add bin edge "{0}" to a "{1}" ' \
'since it is greater than pi'.format(edge, type(self))
raise ValueError(msg)

78
openmc/mgxs/library.py
View file

@ -73,6 +73,10 @@ class Library(object):
Energy group structure for energy condensation
num_delayed_groups : int
Number of delayed groups
num_polar : Integral
Number of equi-width polar angle bins for angle discretization
num_azimuthal : Integral
Number of equi-width azimuthal angle bins for angle discretization
estimator : str or None
The tally estimator used to compute multi-group cross sections.
If None, the default for each MGXS type is used.
@ -107,6 +111,8 @@ class Library(object):
self._domain_type = None
self._domains = 'all'
self._energy_groups = None
self._num_polar = 1
self._num_azimuthal = 1
self._num_delayed_groups = 0
self._correction = 'P0'
self._scatter_format = 'legendre'
@ -144,6 +150,8 @@ class Library(object):
clone._legendre_order = self.legendre_order
clone._histogram_bins = self.histogram_bins
clone._energy_groups = copy.deepcopy(self.energy_groups, memo)
clone._num_polar = self.num_polar
clone._num_azimuthal = self.num_azimuthal
clone._num_delayed_groups = self.num_delayed_groups
clone._tally_trigger = copy.deepcopy(self.tally_trigger, memo)
clone._all_mgxs = copy.deepcopy(self.all_mgxs)
@ -217,6 +225,14 @@ class Library(object):
def num_delayed_groups(self):
return self._num_delayed_groups
@property
def num_polar(self):
return self._num_polar
@property
def num_azimuthal(self):
return self._num_azimuthal
@property
def correction(self):
return self._correction
@ -353,6 +369,18 @@ class Library(object):
equality=True)
self._num_delayed_groups = num_delayed_groups
@num_polar.setter
def num_polar(self, num_polar):
cv.check_type('num_polar', num_polar, Integral)
cv.check_greater_than('num_polar', num_polar, 0)
self._num_polar = num_polar
@num_azimuthal.setter
def num_azimuthal(self, num_azimuthal):
cv.check_type('num_azimuthal', num_azimuthal, Integral)
cv.check_greater_than('num_azimuthal', num_azimuthal, 0)
self._num_azimuthal = num_azimuthal
@correction.setter
def correction(self, correction):
cv.check_value('correction', correction, ('P0', None))
@ -470,9 +498,13 @@ class Library(object):
self.all_mgxs[domain.id] = OrderedDict()
for mgxs_type in self.mgxs_types:
if mgxs_type in openmc.mgxs.MDGXS_TYPES:
mgxs = openmc.mgxs.MDGXS.get_mgxs(mgxs_type, name=self.name)
mgxs = openmc.mgxs.MDGXS.get_mgxs(
mgxs_type, name=self.name, num_polar=self.num_polar,
num_azimuthal=self.num_azimuthal)
else:
mgxs = openmc.mgxs.MGXS.get_mgxs(mgxs_type, name=self.name)
mgxs = openmc.mgxs.MGXS.get_mgxs(
mgxs_type, name=self.name, num_polar=self.num_polar,
num_azimuthal=self.num_azimuthal)
mgxs.domain = domain
mgxs.domain_type = self.domain_type
@ -530,7 +562,7 @@ class Library(object):
mgxs.delayed_groups = None
else:
mgxs.delayed_groups \
= list(range(1, self.num_delayed_groups+1))
= list(range(1, self.num_delayed_groups + 1))
for tally in mgxs.tallies.values():
tallies_file.append(tally, merge=merge)
@ -953,11 +985,16 @@ class Library(object):
name = xsdata_name
if nuclide != 'total':
name += '_' + nuclide
xsdata = openmc.XSdata(name, self.energy_groups)
if self.num_polar > 1 or self.num_azimuthal > 1:
representation = 'angle'
else:
representation = 'isotropic'
xsdata = openmc.XSdata(name, self.energy_groups,
representation=representation)
xsdata.num_delayed_groups = self.num_delayed_groups
# Right now only isotropic weighting is supported
self.representation = 'isotropic'
if self.num_polar > 1 or self.num_azimuthal > 1:
xsdata.num_polar = self.num_polar
xsdata.num_azimuthal = self.num_azimuthal
if nuclide != 'total':
xsdata.atomic_weight_ratio = self._nuclides[nuclide][1]
@ -1098,15 +1135,28 @@ class Library(object):
if 'total' in self.mgxs_types or 'transport' in self.mgxs_types:
if xsdata.scatter_format == 'legendre':
for i in range(len(xsdata.temperatures)):
xsdata._absorption[i] = \
np.subtract(xsdata._total[i], np.sum(
xsdata._scatter_matrix[i][0, :, :], axis=1))
if representation == 'isotropic':
xsdata._absorption[i] = \
np.subtract(xsdata._total[i], np.sum(
xsdata._scatter_matrix[i][:, :, 0],
axis=1))
elif representation == 'angle':
xsdata._absorption[i] = \
np.subtract(xsdata._total[i], np.sum(
xsdata._scatter_matrix[i][:, :, :, :, 0],
axis=3))
elif xsdata.scatter_format == 'histogram':
for i in range(len(xsdata.temperatures)):
xsdata._absorption[i] = \
np.subtract(xsdata._total[i], np.sum(np.sum(
xsdata._scatter_matrix[i][:, :, :], axis=0),
axis=1))
if representation == 'isotropic':
xsdata._absorption[i] = \
np.subtract(xsdata._total[i], np.sum(np.sum(
xsdata._scatter_matrix[i][:, :, :],
axis=2), axis=1))
elif representation == 'angle':
xsdata._absorption[i] = \
np.subtract(xsdata._total[i], np.sum(np.sum(
xsdata._scatter_matrix[i][:, :, :, :, :],
axis=4), axis=3))
return xsdata

496
openmc/mgxs/mdgxs.py
View file

@ -1,6 +1,7 @@
from __future__ import division
from collections import Iterable, OrderedDict
import itertools
from numbers import Integral
import warnings
import os
@ -55,6 +56,12 @@ class MDGXS(MGXS):
tallies in OpenMC 'tallies.xml' file.
delayed_groups : list of int
Delayed groups to filter out the xs
num_polar : Integral, optional
Number of equi-width polar angle bins for angle discretization;
defaults to one bin
num_azimuthal : Integral, optional
Number of equi-width azimuthal angle bins for angle discretization;
defaults to one bin
Attributes
----------
@ -72,6 +79,10 @@ class MDGXS(MGXS):
Energy group structure for energy condensation
delayed_groups : list of int
Delayed groups to filter out the xs
num_polar : Integral
Number of equi-width polar angle bins for angle discretization
num_azimuthal : Integral
Number of equi-width azimuthal angle bins for angle discretization
tally_trigger : openmc.Trigger
An (optional) tally precision trigger given to each tally used to
compute the cross section
@ -118,10 +129,12 @@ class MDGXS(MGXS):
The key used to index multi-group cross sections in an HDF5 data store
"""
def __init__(self, domain=None, domain_type=None, energy_groups=None,
delayed_groups=None, by_nuclide=False, name=''):
delayed_groups=None, by_nuclide=False, name='',
num_polar=1, num_azimuthal=1):
super(MDGXS, self).__init__(domain, domain_type, energy_groups,
by_nuclide, name)
by_nuclide, name, num_polar, num_azimuthal)
self._delayed_groups = None
@ -142,6 +155,8 @@ class MDGXS(MGXS):
clone._domain_type = self.domain_type
clone._energy_groups = copy.deepcopy(self.energy_groups, memo)
clone._delayed_groups = copy.deepcopy(self.delayed_groups, memo)
clone._num_polar = self.num_polar
clone._num_azimuthal = self.num_azimuthal
clone._tally_trigger = copy.deepcopy(self.tally_trigger, memo)
clone._rxn_rate_tally = copy.deepcopy(self._rxn_rate_tally, memo)
clone._xs_tally = copy.deepcopy(self._xs_tally, memo)
@ -160,13 +175,23 @@ class MDGXS(MGXS):
else:
return existing
@property
def _dont_squeeze(self):
"""Create a tuple of axes which should not be removed during the get_xs
process
"""
if self.num_polar > 1 or self.num_azimuthal > 1:
return (0, 1, 3, 4)
else:
return (1, 2)
@property
def delayed_groups(self):
return self._delayed_groups
@property
def num_delayed_groups(self):
if self.delayed_groups == None:
if self.delayed_groups is None:
return 1
else:
return len(self.delayed_groups)
@ -174,7 +199,7 @@ class MDGXS(MGXS):
@delayed_groups.setter
def delayed_groups(self, delayed_groups):
if delayed_groups != None:
if delayed_groups is not None:
cv.check_type('delayed groups', delayed_groups, list, int)
cv.check_greater_than('num delayed groups', len(delayed_groups), 0)
@ -196,14 +221,15 @@ class MDGXS(MGXS):
if self.delayed_groups != None:
delayed_filter = openmc.DelayedGroupFilter(self.delayed_groups)
return [[energy_filter], [delayed_filter, energy_filter]]
filters = [[energy_filter], [delayed_filter, energy_filter]]
else:
return [[energy_filter], [energy_filter]]
filters = [[energy_filter], [energy_filter]]
return self._add_angle_filters(filters)
@staticmethod
def get_mgxs(mdgxs_type, domain=None, domain_type=None,
energy_groups=None, delayed_groups=None,
by_nuclide=False, name=''):
def get_mgxs(mdgxs_type, domain=None, domain_type=None, energy_groups=None,
delayed_groups=None, by_nuclide=False, name='',
num_polar=1, num_azimuthal=1):
"""Return a MDGXS subclass object for some energy group structure within
some spatial domain for some reaction type.
@ -229,7 +255,12 @@ class MDGXS(MGXS):
tallies in OpenMC 'tallies.xml' file. Defaults to the empty string.
delayed_groups : list of int
Delayed groups to filter out the xs
num_polar : Integral, optional
Number of equi-width polar angle bins for angle discretization;
defaults to one bin
num_azimuthal : Integral, optional
Number of equi-width azimuthal angle bins for angle discretization;
defaults to one bin
Returns
-------
openmc.mgxs.MDGXS
@ -256,6 +287,8 @@ class MDGXS(MGXS):
mdgxs.by_nuclide = by_nuclide
mdgxs.name = name
mdgxs.num_polar = num_polar
mdgxs.num_azimuthal = num_azimuthal
return mdgxs
def get_xs(self, groups='all', subdomains='all', nuclides='all',
@ -386,21 +419,29 @@ class MDGXS(MGXS):
else:
num_delayed_groups = len(delayed_groups)
# Reshape tally data array with separate axes for domain, energy groups,
# delayed groups, and nuclides
num_subdomains = int(xs.shape[0] / (num_groups * num_delayed_groups))
new_shape = (num_subdomains, num_delayed_groups, num_groups)
# Reshape tally data array with separate axes for domain,
# energy groups, delayed groups, and nuclides
# Accommodate the polar and azimuthal bins if needed
num_subdomains = \
int(xs.shape[0] / (num_groups * num_delayed_groups *
self.num_polar * self.num_azimuthal))
if self.num_polar > 1 or self.num_azimuthal > 1:
new_shape = (self.num_polar, self.num_azimuthal, num_subdomains,
num_delayed_groups, num_groups)
else:
new_shape = (num_subdomains, num_delayed_groups, num_groups)
new_shape += xs.shape[1:]
xs = np.reshape(xs, new_shape)
# Reverse data if user requested increasing energy groups since
# tally data is stored in order of increasing energies
if order_groups == 'increasing':
xs = xs[:, :, ::-1, :]
xs = xs[..., ::-1, :]
if squeeze:
xs = np.squeeze(xs)
xs = np.atleast_1d(xs)
# We want to squeeze out everything but the polar, azimuthal,
# delayed group, and energy group data.
xs = self._squeeze_xs(xs)
return xs
@ -539,7 +580,7 @@ class MDGXS(MGXS):
"""
if self.delayed_groups == None:
if self.delayed_groups is None:
super(MDGXS, self).print_xs(subdomains, nuclides, xs_type)
return
@ -549,7 +590,7 @@ class MDGXS(MGXS):
elif self.domain_type == 'distribcell':
subdomains = np.arange(self.num_subdomains, dtype=np.int)
elif self.domain_type == 'mesh':
xyz = [range(1, x+1) for x in self.domain.dimension]
xyz = [range(1, x + 1) for x in self.domain.dimension]
subdomains = list(itertools.product(*xyz))
else:
subdomains = [self.domain.id]
@ -581,6 +622,14 @@ class MDGXS(MGXS):
print(string)
return
# Set polar/azimuthal bins
if self.num_polar > 1 or self.num_azimuthal > 1:
polar_bins = np.linspace(0., np.pi, num=self.num_polar + 1,
endpoint=True)
azimuthal_bins = np.linspace(-np.pi, np.pi,
num=self.num_azimuthal + 1,
endpoint=True)
# Loop over all subdomains
for subdomain in subdomains:
@ -605,20 +654,45 @@ class MDGXS(MGXS):
template = '{0: <12}Group {1} [{2: <10} - {3: <10}eV]:\t'
# Loop over energy groups ranges
for group in range(1, self.num_groups+1):
bounds = self.energy_groups.get_group_bounds(group)
string += template.format('', group, bounds[0], bounds[1])
average = self.get_xs([group], [subdomain], [nuclide],
xs_type=xs_type, value='mean',
delayed_groups=[delayed_group])
rel_err = self.get_xs([group], [subdomain], [nuclide],
xs_type=xs_type, value='rel_err',
delayed_groups=[delayed_group])
average = average.flatten()[0]
rel_err = rel_err.flatten()[0] * 100.
string += '{:.2e} +/- {:1.2e}%'.format(average, rel_err)
string += '\n'
average_xs = self.get_xs(nuclides=[nuclide],
subdomains=[subdomain],
xs_type=xs_type, value='mean',
delayed_groups=[delayed_group])
rel_err_xs = self.get_xs(nuclides=[nuclide],
subdomains=[subdomain],
xs_type=xs_type, value='rel_err',
delayed_groups=[delayed_group])
rel_err_xs = rel_err_xs * 100.
if self.num_polar > 1 or self.num_azimuthal > 1:
# Loop over polar, azimuthal, and energy group ranges
for pol in range(len(polar_bins) - 1):
pol_low, pol_high = polar_bins[pol: pol + 2]
for azi in range(len(azimuthal_bins) - 1):
azi_low, azi_high = azimuthal_bins[azi: azi + 2]
string += '\t\tPolar Angle: [{0:5f} - {1:5f}]'.format(
pol_low, pol_high) + \
'\tAzimuthal Angle: [{0:5f} - {1:5f}]'.format(
azi_low, azi_high) + '\n'
for group in range(1, self.num_groups + 1):
bounds = \
self.energy_groups.get_group_bounds(group)
string += '\t' + template.format('', group,
bounds[0],
bounds[1])
string += '{0:.2e} +/- {1:.2e}%'.format(
average_xs[pol, azi, group - 1],
rel_err_xs[pol, azi, group - 1])
string += '\n'
string += '\n'
else:
# Loop over energy groups ranges
for group in range(1, self.num_groups+1):
bounds = self.energy_groups.get_group_bounds(group)
string += template.format('', group, bounds[0], bounds[1])
string += '{0:.2e} +/- {1:.2e}%'.format(
average_xs[group - 1], rel_err_xs[group - 1])
string += '\n'
string += '\n'
string += '\n'
@ -755,7 +829,6 @@ class MDGXS(MGXS):
if self.by_nuclide and nuclides == 'sum':
# Use tally summation to sum across all nuclides
query_nuclides = [nuclides]
xs_tally = self.xs_tally.summation(nuclides=self.get_nuclides())
df = xs_tally.get_pandas_dataframe(
distribcell_paths=distribcell_paths)
@ -768,14 +841,12 @@ class MDGXS(MGXS):
# If the user requested a specific set of nuclides
elif self.by_nuclide and nuclides != 'all':
query_nuclides = nuclides
xs_tally = self.xs_tally.get_slice(nuclides=nuclides)
df = xs_tally.get_pandas_dataframe(
distribcell_paths=distribcell_paths)
# If the user requested all nuclides, keep nuclide column in dataframe
else:
query_nuclides = self.nuclides
df = self.xs_tally.get_pandas_dataframe(
distribcell_paths=distribcell_paths)
@ -785,40 +856,9 @@ class MDGXS(MGXS):
else:
df = df.drop('score', axis=1)
# Override energy groups bounds with indices
all_groups = np.arange(self.num_groups, 0, -1, dtype=np.int)
all_groups = np.repeat(all_groups, len(query_nuclides))
if 'energy low [eV]' in df and 'energyout low [eV]' in df:
df.rename(columns={'energy low [eV]': 'group in'},
inplace=True)
in_groups = np.tile(all_groups, int(self.num_subdomains *
self.num_delayed_groups))
in_groups = np.repeat(in_groups, int(df.shape[0] / in_groups.size))
df['group in'] = in_groups
del df['energy high [eV]']
df.rename(columns={'energyout low [eV]': 'group out'},
inplace=True)
out_groups = np.repeat(all_groups, self.xs_tally.num_scores)
out_groups = np.tile(out_groups, int(df.shape[0] / out_groups.size))
df['group out'] = out_groups
del df['energyout high [eV]']
columns = ['group in', 'group out']
elif 'energyout low [eV]' in df:
df.rename(columns={'energyout low [eV]': 'group out'},
inplace=True)
in_groups = np.tile(all_groups, int(df.shape[0] / all_groups.size))
df['group out'] = in_groups
del df['energyout high [eV]']
columns = ['group out']
elif 'energy low [eV]' in df:
df.rename(columns={'energy low [eV]': 'group in'}, inplace=True)
in_groups = np.tile(all_groups, int(df.shape[0] / all_groups.size))
df['group in'] = in_groups
del df['energy high [eV]']
columns = ['group in']
# Convert azimuthal, polar, energy in and energy out bin values in to
# bin indices
columns = self._df_convert_columns_to_bins(df)
# Select out those groups the user requested
if not isinstance(groups, string_types):
@ -834,7 +874,7 @@ class MDGXS(MGXS):
else:
densities = self.get_nuclide_densities('sum')
densities = np.repeat(densities, len(self.rxn_rate_tally.scores))
tile_factor = df.shape[0] / len(densities)
tile_factor = int(df.shape[0] / len(densities))
df['mean'] /= np.tile(densities, tile_factor)
df['std. dev.'] /= np.tile(densities, tile_factor)
@ -842,7 +882,7 @@ class MDGXS(MGXS):
# energy groups such that data is from fast to thermal
if self.domain_type == 'mesh':
mesh_str = 'mesh {0}'.format(self.domain.id)
df.sort_values(by=[(mesh_str, 'x'), (mesh_str, 'y'), \
df.sort_values(by=[(mesh_str, 'x'), (mesh_str, 'y'),
(mesh_str, 'z')] + columns, inplace=True)
else:
df.sort_values(by=[self.domain_type] + columns, inplace=True)
@ -896,6 +936,12 @@ class ChiDelayed(MDGXS):
tallies in OpenMC 'tallies.xml' file.
delayed_groups : list of int
Delayed groups to filter out the xs
num_polar : Integral, optional
Number of equi-width polar angle bins for angle discretization;
defaults to one bin
num_azimuthal : Integral, optional
Number of equi-width azimuthal angle bins for angle discretization;
defaults to one bin
Attributes
----------
@ -913,6 +959,10 @@ class ChiDelayed(MDGXS):
Energy group structure for energy condensation
delayed_groups : list of int
Delayed groups to filter out the xs
num_polar : Integral
Number of equi-width polar angle bins for angle discretization
num_azimuthal : Integral
Number of equi-width azimuthal angle bins for angle discretization
tally_trigger : openmc.Trigger
An (optional) tally precision trigger given to each tally used to
compute the cross section
@ -962,9 +1012,11 @@ class ChiDelayed(MDGXS):
"""
def __init__(self, domain=None, domain_type=None, energy_groups=None,
delayed_groups=None, by_nuclide=False, name=''):
delayed_groups=None, by_nuclide=False, name='',
num_polar=1, num_azimuthal=1):
super(ChiDelayed, self).__init__(domain, domain_type, energy_groups,
delayed_groups, by_nuclide, name)
delayed_groups, by_nuclide, name,
num_polar, num_azimuthal)
self._rxn_type = 'chi-delayed'
self._estimator = 'analog'
@ -978,11 +1030,13 @@ class ChiDelayed(MDGXS):
group_edges = self.energy_groups.group_edges
energyout = openmc.EnergyoutFilter(group_edges)
energyin = openmc.EnergyFilter([group_edges[0], group_edges[-1]])
if self.delayed_groups != None:
if self.delayed_groups is not None:
delayed_filter = openmc.DelayedGroupFilter(self.delayed_groups)
return [[delayed_filter, energyin], [delayed_filter, energyout]]
filters = [[delayed_filter, energyin], [delayed_filter, energyout]]
else:
return [[energyin], [energyout]]
filters = [[energyin], [energyout]]
return self._add_angle_filters(filters)
@property
def tally_keys(self):
@ -1326,21 +1380,28 @@ class ChiDelayed(MDGXS):
else:
num_delayed_groups = len(delayed_groups)
# Reshape tally data array with separate axes for domain, energy groups,
# delayed groups, and nuclides
num_subdomains = int(xs.shape[0] / (num_groups * num_delayed_groups))
new_shape = (num_subdomains, num_delayed_groups, num_groups)
# Reshape tally data array with separate axes for domain, energy
# groups, and accomodate the polar and azimuthal bins if needed
num_subdomains = int(xs.shape[0] / (num_delayed_groups *
num_groups * self.num_polar *
self.num_azimuthal))
if self.num_polar > 1 or self.num_azimuthal > 1:
new_shape = (self.num_polar, self.num_azimuthal, num_subdomains,
num_delayed_groups, num_groups)
else:
new_shape = (num_subdomains, num_delayed_groups, num_groups)
new_shape += xs.shape[1:]
xs = np.reshape(xs, new_shape)
# Reverse data if user requested increasing energy groups since
# tally data is stored in order of increasing energies
if order_groups == 'increasing':
xs = xs[:, :, ::-1, :]
xs = xs[..., ::-1, :]
if squeeze:
xs = np.squeeze(xs)
xs = np.atleast_1d(xs)
# We want to squeeze out everything but the polar, azimuthal,
# and energy group data.
xs = self._squeeze_xs(xs)
return xs
@ -1389,6 +1450,12 @@ class DelayedNuFissionXS(MDGXS):
tallies in OpenMC 'tallies.xml' file.
delayed_groups : list of int
Delayed groups to filter out the xs
num_polar : Integral, optional
Number of equi-width polar angle bins for angle discretization;
defaults to one bin
num_azimuthal : Integral, optional
Number of equi-width azimuthal angle bins for angle discretization;
defaults to one bin
Attributes
----------
@ -1406,6 +1473,10 @@ class DelayedNuFissionXS(MDGXS):
Energy group structure for energy condensation
delayed_groups : list of int
Delayed groups to filter out the xs
num_polar : Integral
Number of equi-width polar angle bins for angle discretization
num_azimuthal : Integral
Number of equi-width azimuthal angle bins for angle discretization
tally_trigger : openmc.Trigger
An (optional) tally precision trigger given to each tally used to
compute the cross section
@ -1455,10 +1526,12 @@ class DelayedNuFissionXS(MDGXS):
"""
def __init__(self, domain=None, domain_type=None, energy_groups=None,
delayed_groups=None, by_nuclide=False, name=''):
delayed_groups=None, by_nuclide=False, name='',
num_polar=1, num_azimuthal=1):
super(DelayedNuFissionXS, self).__init__(domain, domain_type,
energy_groups, delayed_groups,
by_nuclide, name)
by_nuclide, name, num_polar,
num_azimuthal)
self._rxn_type = 'delayed-nu-fission'
@ -1513,6 +1586,12 @@ class Beta(MDGXS):
tallies in OpenMC 'tallies.xml' file.
delayed_groups : list of int
Delayed groups to filter out the xs
num_polar : Integral, optional
Number of equi-width polar angle bins for angle discretization;
defaults to one bin
num_azimuthal : Integral, optional
Number of equi-width azimuthal angle bins for angle discretization;
defaults to one bin
Attributes
----------
@ -1530,6 +1609,10 @@ class Beta(MDGXS):
Energy group structure for energy condensation
delayed_groups : list of int
Delayed groups to filter out the xs
num_polar : Integral
Number of equi-width polar angle bins for angle discretization
num_azimuthal : Integral
Number of equi-width azimuthal angle bins for angle discretization
tally_trigger : openmc.Trigger
An (optional) tally precision trigger given to each tally used to
compute the cross section
@ -1579,9 +1662,11 @@ class Beta(MDGXS):
"""
def __init__(self, domain=None, domain_type=None, energy_groups=None,
delayed_groups=None, by_nuclide=False, name=''):
delayed_groups=None, by_nuclide=False, name='',
num_polar=1, num_azimuthal=1):
super(Beta, self).__init__(domain, domain_type, energy_groups,
delayed_groups, by_nuclide, name)
delayed_groups, by_nuclide, name, num_polar,
num_azimuthal)
self._rxn_type = 'beta'
@property
@ -1685,6 +1770,12 @@ class DecayRate(MDGXS):
tallies in OpenMC 'tallies.xml' file.
delayed_groups : list of int
Delayed groups to filter out the xs
num_polar : Integral, optional
Number of equi-width polar angle bins for angle discretization;
defaults to one bin
num_azimuthal : Integral, optional
Number of equi-width azimuthal angle bins for angle discretization;
defaults to one bin
Attributes
----------
@ -1702,6 +1793,10 @@ class DecayRate(MDGXS):
Energy group structure for energy condensation
delayed_groups : list of int
Delayed groups to filter out the xs
num_polar : Integral
Number of equi-width polar angle bins for angle discretization
num_azimuthal : Integral
Number of equi-width azimuthal angle bins for angle discretization
tally_trigger : openmc.Trigger
An (optional) tally precision trigger given to each tally used to
compute the cross section
@ -1751,9 +1846,11 @@ class DecayRate(MDGXS):
"""
def __init__(self, domain=None, domain_type=None, energy_groups=None,
delayed_groups=None, by_nuclide=False, name=''):
delayed_groups=None, by_nuclide=False, name='',
num_polar=1, num_azimuthal=1):
super(DecayRate, self).__init__(domain, domain_type, energy_groups,
delayed_groups, by_nuclide, name)
delayed_groups, by_nuclide, name,
num_polar, num_azimuthal)
self._rxn_type = 'decay-rate'
@property
@ -1771,11 +1868,14 @@ class DecayRate(MDGXS):
group_edges = self.energy_groups.group_edges
energy_filter = openmc.EnergyFilter(group_edges)
if self.delayed_groups != None:
if self.delayed_groups is not None:
delayed_filter = openmc.DelayedGroupFilter(self.delayed_groups)
return [[delayed_filter, energy_filter], [delayed_filter, energy_filter]]
filters = [[delayed_filter, energy_filter], [delayed_filter,
energy_filter]]
else:
return [[energy_filter], [energy_filter]]
filters = [[energy_filter], [energy_filter]]
return self._add_angle_filters(filters)
@property
def xs_tally(self):
@ -1847,6 +1947,12 @@ class MatrixMDGXS(MDGXS):
tallies in OpenMC 'tallies.xml' file.
delayed_groups : list of int
Delayed groups to filter out the xs
num_polar : Integral, optional
Number of equi-width polar angle bins for angle discretization;
defaults to one bin
num_azimuthal : Integral, optional
Number of equi-width azimuthal angle bins for angle discretization;
defaults to one bin
Attributes
----------
@ -1864,6 +1970,10 @@ class MatrixMDGXS(MDGXS):
Energy group structure for energy condensation
delayed_groups : list of int
Delayed groups to filter out the xs
num_polar : Integral
Number of equi-width polar angle bins for angle discretization
num_azimuthal : Integral
Number of equi-width azimuthal angle bins for angle discretization
tally_trigger : openmc.Trigger
An (optional) tally precision trigger given to each tally used to
compute the cross section
@ -1911,6 +2021,16 @@ class MatrixMDGXS(MDGXS):
"""
@property
def _dont_squeeze(self):
"""Create a tuple of axes which should not be removed during the get_xs
process
"""
if self.num_polar > 1 or self.num_azimuthal > 1:
return (0, 1, 3, 4, 5)
else:
return (1, 2, 3)
@property
def filters(self):
# Create the non-domain specific Filters for the Tallies
@ -1920,9 +2040,11 @@ class MatrixMDGXS(MDGXS):
if self.delayed_groups is not None:
delayed = openmc.DelayedGroupFilter(self.delayed_groups)
return [[energy], [delayed, energy, energyout]]
filters = [[energy], [delayed, energy, energyout]]
else:
return [[energy], [energy, energyout]]
filters = [[energy], [energy, energyout]]
return self._add_angle_filters(filters)
def get_xs(self, in_groups='all', out_groups='all',
subdomains='all', nuclides='all',
@ -2076,25 +2198,39 @@ class MatrixMDGXS(MDGXS):
num_delayed_groups = len(delayed_groups)
# Reshape tally data array with separate axes for domain and energy
num_subdomains = int(xs.shape[0] / (num_in_groups * num_out_groups *
num_delayed_groups))
new_shape = (num_subdomains, num_delayed_groups, num_in_groups,
num_out_groups)
new_shape += xs.shape[1:]
xs = np.reshape(xs, new_shape)
# Accomodate the polar and azimuthal bins if needed
num_subdomains = int(xs.shape[0] / (num_delayed_groups *
num_in_groups * num_out_groups *
self.num_polar *
self.num_azimuthal))
if self.num_polar > 1 or self.num_azimuthal > 1:
new_shape = (self.num_polar, self.num_azimuthal, num_subdomains,
num_delayed_groups, num_in_groups, num_out_groups)
new_shape += xs.shape[1:]
xs = np.reshape(xs, new_shape)
# Transpose the matrix if requested by user
if row_column == 'outin':
xs = np.swapaxes(xs, 2, 3)
# Transpose the matrix if requested by user
if row_column == 'outin':
xs = np.swapaxes(xs, 4, 5)
else:
new_shape = (num_subdomains, num_delayed_groups, num_in_groups,
num_out_groups)
new_shape += xs.shape[1:]
xs = np.reshape(xs, new_shape)
# Transpose the matrix if requested by user
if row_column == 'outin':
xs = np.swapaxes(xs, 2, 3)
# Reverse data if user requested increasing energy groups since
# tally data is stored in order of increasing energies
if order_groups == 'increasing':
xs = xs[:, :, ::-1, ::-1, :]
xs = xs[..., ::-1, ::-1, :]
if squeeze:
xs = np.squeeze(xs)
xs = np.atleast_2d(xs)
# We want to squeeze out everything but the polar, azimuthal,
# and in/out energy group data.
xs = self._squeeze_xs(xs)
return xs
@ -2184,7 +2320,7 @@ class MatrixMDGXS(MDGXS):
elif self.domain_type == 'distribcell':
subdomains = np.arange(self.num_subdomains, dtype=np.int)
elif self.domain_type == 'mesh':
xyz = [range(1, x+1) for x in self.domain.dimension]
xyz = [range(1, x + 1) for x in self.domain.dimension]
subdomains = list(itertools.product(*xyz))
else:
subdomains = [self.domain.id]
@ -2217,13 +2353,20 @@ class MatrixMDGXS(MDGXS):
return
string += '{0: <16}\n'.format('\tEnergy Groups:')
template = '{0: <12}Group {1} [{2: <10} - {3: <10}MeV]\n'
template = '{0: <12}Group {1} [{2: <10} - {3: <10}eV]\n'
# Loop over energy groups ranges
for group in range(1, self.num_groups + 1):
bounds = self.energy_groups.get_group_bounds(group)
string += template.format('', group, bounds[0], bounds[1])
# Set polar and azimuthal bins if necessary
if self.num_polar > 1 or self.num_azimuthal > 1:
pol_bins = np.linspace(0., np.pi, num=self.num_polar + 1,
endpoint=True)
azi_bins = np.linspace(-np.pi, np.pi, num=self.num_azimuthal + 1,
endpoint=True)
# Loop over all subdomains
for subdomain in subdomains:
@ -2250,47 +2393,97 @@ class MatrixMDGXS(MDGXS):
template = '{0: <12}Group {1} -> Group {2}:\t\t'
# Loop over incoming/outgoing energy groups ranges
for in_group in range(1, self.num_groups + 1):
for out_group in range(1, self.num_groups + 1):
string += template.format('', in_group, out_group)
average = self.get_xs([in_group], [out_group],
[subdomain], [nuclide],
xs_type=xs_type,
value='mean',
delayed_groups=[delayed_group])
rel_err = self.get_xs([in_group], [out_group],
[subdomain], [nuclide],
xs_type=xs_type,
value='rel_err',
delayed_groups=[delayed_group])
average = average.flatten()[0]
rel_err = rel_err.flatten()[0] * 100.
string += '{:.2e} +/- {:.2e}%'.format(average,
rel_err)
average_xs = self.get_xs(nuclides=[nuclide],
subdomains=[subdomain],
xs_type=xs_type, value='mean',
delayed_groups=[delayed_group])
rel_err_xs = self.get_xs(nuclides=[nuclide],
subdomains=[subdomain],
xs_type=xs_type,
value='rel_err',
delayed_groups=[delayed_group])
rel_err_xs = rel_err_xs * 100.
if self.num_polar > 1 or self.num_azimuthal > 1:
# Loop over polar, azi, and in/out group ranges
for pol in range(len(pol_bins) - 1):
pol_low, pol_high = pol_bins[pol: pol + 2]
for azi in range(len(azi_bins) - 1):
azi_low, azi_high = azi_bins[azi: azi + 2]
string += '\t\tPolar Angle: [{0:5f} - {1:5f}]'.format(
pol_low, pol_high) + \
'\tAzimuthal Angle: [{0:5f} - {1:5f}]'.format(
azi_low, azi_high) + '\n'
for in_group in range(1, self.num_groups + 1):
for out_group in range(1, self.num_groups + 1):
string += '\t' + template.format(
'', in_group, out_group)
string += '{0:.2e} +/- {1:.2e}%'.format(
average_xs[pol, azi, in_group - 1,
out_group - 1],
rel_err_xs[pol, azi, in_group - 1,
out_group - 1])
string += '\n'
string += '\n'
string += '\n'
else:
# Loop over incoming/outgoing energy groups ranges
for in_group in range(1, self.num_groups + 1):
for out_group in range(1, self.num_groups + 1):
string += template.format(
'', in_group, out_group)
string += '{:.2e} +/- {:.2e}%'.format(
average_xs[in_group-1, out_group-1],
rel_err_xs[in_group-1, out_group-1])
string += '\n'
string += '\n'
string += '\n'
string += '\n'
else:
template = '{0: <12}Group {1} -> Group {2}:\t\t'
# Loop over incoming/outgoing energy groups ranges
for in_group in range(1, self.num_groups + 1):
for out_group in range(1, self.num_groups + 1):
string += template.format('', in_group, out_group)
average = self.get_xs([in_group], [out_group],
[subdomain], [nuclide],
xs_type=xs_type, value='mean')
rel_err = self.get_xs([in_group], [out_group],
[subdomain], [nuclide],
xs_type=xs_type, value='rel_err')
average = average.flatten()[0]
rel_err = rel_err.flatten()[0] * 100.
string += '{:.2e} +/- {:.2e}%'.format(average,
rel_err)
average_xs = self.get_xs(nuclides=[nuclide],
subdomains=[subdomain],
xs_type=xs_type, value='mean')
rel_err_xs = self.get_xs(nuclides=[nuclide],
subdomains=[subdomain],
xs_type=xs_type, value='rel_err')
rel_err_xs = rel_err_xs * 100.
if self.num_polar > 1 or self.num_azimuthal > 1:
# Loop over polar, azi, and in/out energy group ranges
for pol in range(len(pol_bins) - 1):
pol_low, pol_high = pol_bins[pol: pol + 2]
for azi in range(len(azi_bins) - 1):
azi_low, azi_high = azi_bins[azi: azi + 2]
string += '\t\tPolar Angle: [{0:5f} - {1:5f}]'.format(
pol_low, pol_high) + \
'\tAzimuthal Angle: [{0:5f} - {1:5f}]'.format(
azi_low, azi_high) + '\n'
for in_group in range(1, self.num_groups + 1):
for out_group in range(1, self.num_groups + 1):
string += '\t' + template.format(
'', in_group, out_group)
string += '{0:.2e} +/- {1:.2e}%'.format(
average_xs[pol, azi, in_group - 1,
out_group - 1],
rel_err_xs[pol, azi, in_group - 1,
out_group - 1])
string += '\n'
string += '\n'
string += '\n'
else:
# Loop over incoming/outgoing energy groups ranges
for in_group in range(1, self.num_groups + 1):
for out_group in range(1, self.num_groups + 1):
string += template.format('', in_group,
out_group)
string += '{0:.2e} +/- {1:.2e}%'.format(
average_xs[in_group - 1, out_group - 1],
rel_err_xs[in_group - 1, out_group - 1])
string += '\n'
string += '\n'
string += '\n'
string += '\n'
string += '\n'
string += '\n'
@ -2346,6 +2539,12 @@ class DelayedNuFissionMatrixXS(MatrixMDGXS):
tallies in OpenMC 'tallies.xml' file.
delayed_groups : list of int
Delayed groups to filter out the xs
num_polar : Integral, optional
Number of equi-width polar angle bins for angle discretization;
defaults to one bin
num_azimuthal : Integral, optional
Number of equi-width azimuthal angle bins for angle discretization;
defaults to one bin
Attributes
----------
@ -2363,6 +2562,10 @@ class DelayedNuFissionMatrixXS(MatrixMDGXS):
Energy group structure for energy condensation
delayed_groups : list of int
Delayed groups to filter out the xs
num_polar : Integral
Number of equi-width polar angle bins for angle discretization
num_azimuthal : Integral
Number of equi-width azimuthal angle bins for angle discretization
tally_trigger : openmc.Trigger
An (optional) tally precision trigger given to each tally used to
compute the cross section
@ -2412,11 +2615,14 @@ class DelayedNuFissionMatrixXS(MatrixMDGXS):
"""
def __init__(self, domain=None, domain_type=None, energy_groups=None,
delayed_groups=None, by_nuclide=False, name=''):
delayed_groups=None, by_nuclide=False, name='',
num_polar=1, num_azimuthal=1):
super(DelayedNuFissionMatrixXS, self).__init__(domain, domain_type,
energy_groups,
delayed_groups,
by_nuclide, name)
by_nuclide, name,
num_polar,
num_azimuthal)
self._rxn_type = 'delayed-nu-fission'
self._hdf5_key = 'delayed-nu-fission matrix'
self._estimator = 'analog'

1014
openmc/mgxs/mgxs.py
View file

File diff suppressed because it is too large Load diff

188
openmc/mgxs_library.py
View file

@ -937,12 +937,8 @@ class XSdata(object):
check_value('temperature', temperature, self.temperatures)
i = np.where(self.temperatures == temperature)[0][0]
if self.representation == 'isotropic':
self._total[i] = total.get_xs(nuclides=nuclide, xs_type=xs_type,
subdomains=subdomain)
elif self.representation == 'angle':
msg = 'Angular-Dependent MGXS have not yet been implemented'
raise ValueError(msg)
self._total[i] = total.get_xs(nuclides=nuclide, xs_type=xs_type,
subdomains=subdomain)
def set_absorption_mgxs(self, absorption, temperature=294.,
nuclide='total', xs_type='macro', subdomain=None):
@ -983,13 +979,9 @@ class XSdata(object):
check_value('temperature', temperature, self.temperatures)
i = np.where(self.temperatures == temperature)[0][0]
if self.representation == 'isotropic':
self._absorption[i] = absorption.get_xs(nuclides=nuclide,
xs_type=xs_type,
subdomains=subdomain)
elif self.representation == 'angle':
msg = 'Angular-Dependent MGXS have not yet been implemented'
raise ValueError(msg)
self._absorption[i] = absorption.get_xs(nuclides=nuclide,
xs_type=xs_type,
subdomains=subdomain)
def set_fission_mgxs(self, fission, temperature=294., nuclide='total',
xs_type='macro', subdomain=None):
@ -1030,13 +1022,9 @@ class XSdata(object):
check_value('temperature', temperature, self.temperatures)
i = np.where(self.temperatures == temperature)[0][0]
if self.representation == 'isotropic':
self._fission[i] = fission.get_xs(nuclides=nuclide,
xs_type=xs_type,
subdomains=subdomain)
elif self.representation == 'angle':
msg = 'Angular-Dependent MGXS have not yet been implemented'
raise ValueError(msg)
self._fission[i] = fission.get_xs(nuclides=nuclide,
xs_type=xs_type,
subdomains=subdomain)
def set_nu_fission_mgxs(self, nu_fission, temperature=294.,
nuclide='total', xs_type='macro', subdomain=None):
@ -1078,13 +1066,9 @@ class XSdata(object):
check_value('temperature', temperature, self.temperatures)
i = np.where(self.temperatures == temperature)[0][0]
if self.representation == 'isotropic':
self._nu_fission[i] = nu_fission.get_xs(nuclides=nuclide,
xs_type=xs_type,
subdomains=subdomain)
elif self.representation == 'angle':
msg = 'Angular-Dependent MGXS have not yet been implemented'
raise ValueError(msg)
self._nu_fission[i] = nu_fission.get_xs(nuclides=nuclide,
xs_type=xs_type,
subdomains=subdomain)
if np.sum(self._nu_fission) > 0.0:
self._fissionable = True
@ -1134,14 +1118,8 @@ class XSdata(object):
check_value('temperature', temperature, self.temperatures)
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,
xs_type=xs_type,
subdomains=subdomain)
elif self.representation is 'angle':
msg = 'Angular-Dependent MGXS have not yet been implemented'
raise ValueError(msg)
self._prompt_nu_fission[i] = prompt_nu_fission.get_xs(
nuclides=nuclide, xs_type=xs_type, subdomains=subdomain)
if np.sum(self._prompt_nu_fission) > 0.0:
self._fissionable = True
@ -1191,14 +1169,8 @@ class XSdata(object):
check_value('temperature', temperature, self.temperatures)
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,
xs_type=xs_type,
subdomains=subdomain)
elif self.representation is 'angle':
msg = 'Angular-Dependent MGXS have not yet been implemented'
raise ValueError(msg)
self._delayed_nu_fission[i] = delayed_nu_fission.get_xs(
nuclides=nuclide, xs_type=xs_type, subdomains=subdomain)
if np.sum(self._delayed_nu_fission) > 0.0:
self._fissionable = True
@ -1244,13 +1216,9 @@ class XSdata(object):
check_value('temperature', temperature, self.temperatures)
i = np.where(self.temperatures == temperature)[0][0]
if self.representation == 'isotropic':
self._kappa_fission[i] = k_fission.get_xs(nuclides=nuclide,
xs_type=xs_type,
subdomains=subdomain)
elif self.representation == 'angle':
msg = 'Angular-Dependent MGXS have not yet been implemented'
raise ValueError(msg)
self._kappa_fission[i] = k_fission.get_xs(nuclides=nuclide,
xs_type=xs_type,
subdomains=subdomain)
def set_chi_mgxs(self, chi, temperature=294., nuclide='total',
xs_type='macro', subdomain=None):
@ -1288,15 +1256,11 @@ class XSdata(object):
check_value('temperature', temperature, self.temperatures)
i = np.where(self.temperatures == temperature)[0][0]
if self.representation == 'isotropic':
self._chi[i] = chi.get_xs(nuclides=nuclide,
xs_type=xs_type, subdomains=subdomain)
elif self.representation == 'angle':
msg = 'Angular-Dependent MGXS have not yet been implemented'
raise ValueError(msg)
self._chi[i] = chi.get_xs(nuclides=nuclide, xs_type=xs_type,
subdomains=subdomain)
def set_chi_prompt_mgxs(self, chi_prompt, temperature=294., nuclide='total',
xs_type='macro', subdomain=None):
def set_chi_prompt_mgxs(self, chi_prompt, temperature=294.,
nuclide='total', xs_type='macro', subdomain=None):
"""This method allows for an openmc.mgxs.ChiPrompt
to be used to set chi-prompt for this XSdata object.
@ -1333,13 +1297,9 @@ class XSdata(object):
check_value('temperature', temperature, self.temperatures)
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,
subdomains=subdomain)
elif self.representation is 'angle':
msg = 'Angular-Dependent MGXS have not yet been implemented'
raise ValueError(msg)
self._chi_prompt[i] = chi_prompt.get_xs(nuclides=nuclide,
xs_type=xs_type,
subdomains=subdomain)
def set_chi_delayed_mgxs(self, chi_delayed, temperature=294.,
nuclide='total', xs_type='macro', subdomain=None):
@ -1381,13 +1341,9 @@ class XSdata(object):
check_value('temperature', temperature, self.temperatures)
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,
subdomains=subdomain)
elif self.representation is 'angle':
msg = 'Angular-Dependent MGXS have not yet been implemented'
raise ValueError(msg)
self._chi_delayed[i] = chi_delayed.get_xs(nuclides=nuclide,
xs_type=xs_type,
subdomains=subdomain)
def set_beta_mgxs(self, beta, temperature=294.,
nuclide='total', xs_type='macro', subdomain=None):
@ -1426,13 +1382,9 @@ class XSdata(object):
check_value('temperature', temperature, self.temperatures)
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,
subdomains=subdomain)
elif self.representation is 'angle':
msg = 'Angular-Dependent MGXS have not yet been implemented'
raise ValueError(msg)
self._beta[i] = beta.get_xs(nuclides=nuclide,
xs_type=xs_type,
subdomains=subdomain)
def set_decay_rate_mgxs(self, decay_rate, temperature=294.,
nuclide='total', xs_type='macro', subdomain=None):
@ -1472,13 +1424,9 @@ class XSdata(object):
check_value('temperature', temperature, self.temperatures)
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,
subdomains=subdomain)
elif self.representation is 'angle':
msg = 'Angular-Dependent MGXS have not yet been implemented'
raise ValueError(msg)
self._decay_rate[i] = decay_rate.get_xs(nuclides=nuclide,
xs_type=xs_type,
subdomains=subdomain)
def set_scatter_matrix_mgxs(self, scatter, temperature=294.,
nuclide='total', xs_type='macro',
@ -1543,22 +1491,24 @@ class XSdata(object):
[self.order])
i = np.where(self.temperatures == temperature)[0][0]
if self.representation == 'isotropic':
if self.scatter_format == 'legendre':
# Get the scattering orders in the outermost dimension
self._scatter_matrix[i] = \
np.zeros(self.xs_shapes["[G][G'][Order]"])
if self.scatter_format == 'legendre':
self._scatter_matrix[i] = \
np.zeros(self.xs_shapes["[G][G'][Order]"])
# Get the scattering orders in the outermost dimension
if self.representation == 'isotropic':
for moment in range(self.num_orders):
self._scatter_matrix[i][:, :, moment] = \
scatter.get_xs(nuclides=nuclide, xs_type=xs_type,
moment=moment, subdomains=subdomain)
else:
self._scatter_matrix[i] = \
scatter.get_xs(nuclides=nuclide, xs_type=xs_type,
subdomains=subdomain)
elif self.representation == 'angle':
msg = 'Angular-Dependent MGXS have not yet been implemented'
raise ValueError(msg)
elif self.representation == 'angle':
for moment in range(self.num_orders):
self._scatter_matrix[i][:, :, :, :, moment] = \
scatter.get_xs(nuclides=nuclide, xs_type=xs_type,
moment=moment, subdomains=subdomain)
else:
self._scatter_matrix[i] = \
scatter.get_xs(nuclides=nuclide, xs_type=xs_type,
subdomains=subdomain)
def set_multiplicity_matrix_mgxs(self, nuscatter, scatter=None,
temperature=294., nuclide='total',
@ -1623,24 +1573,21 @@ class XSdata(object):
check_value('domain_type', scatter.domain_type,
openmc.mgxs.DOMAIN_TYPES)
i = np.where(self.temperatures == temperature)[0][0]
if self.representation == 'isotropic':
nuscatt = nuscatter.get_xs(nuclides=nuclide,
xs_type=xs_type, moment=0,
subdomains=subdomain)
if isinstance(nuscatter, openmc.mgxs.MultiplicityMatrixXS):
self._multiplicity_matrix[i] = nuscatt
else:
scatt = scatter.get_xs(nuclides=nuclide,
xs_type=xs_type, moment=0,
subdomains=subdomain)
if scatter.scatter_format == 'histogram':
scatt = np.sum(scatt, axis=0)
if nuscatter.scatter_format == 'histogram':
nuscatt = np.sum(nuscatt, axis=0)
self._multiplicity_matrix[i] = np.divide(nuscatt, scatt)
elif self.representation == 'angle':
msg = 'Angular-Dependent MGXS have not yet been implemented'
raise ValueError(msg)
nuscatt = nuscatter.get_xs(nuclides=nuclide,
xs_type=xs_type, moment=0,
subdomains=subdomain)
if isinstance(nuscatter, openmc.mgxs.MultiplicityMatrixXS):
self._multiplicity_matrix[i] = nuscatt
else:
scatt = scatter.get_xs(nuclides=nuclide,
xs_type=xs_type, moment=0,
subdomains=subdomain)
if scatter.scatter_format == 'histogram':
scatt = np.sum(scatt, axis=0)
if nuscatter.scatter_format == 'histogram':
nuscatt = np.sum(nuscatt, axis=0)
self._multiplicity_matrix[i] = np.divide(nuscatt, scatt)
self._multiplicity_matrix[i] = \
np.nan_to_num(self._multiplicity_matrix[i])
@ -1684,13 +1631,8 @@ class XSdata(object):
check_value('temperature', temperature, self.temperatures)
i = np.where(self.temperatures == temperature)[0][0]
if self.representation == 'isotropic':
self._inverse_velocity[i] = inverse_velocity.get_xs\
(nuclides=nuclide, xs_type=xs_type,
subdomains=subdomain)
elif self.representation == 'angle':
msg = 'Angular-Dependent MGXS have not yet been implemented'
raise ValueError(msg)
self._inverse_velocity[i] = inverse_velocity.get_xs(
nuclides=nuclide, xs_type=xs_type, subdomains=subdomain)
def to_hdf5(self, file):
"""Write XSdata to an HDF5 file

4
src/input_xml.F90
View file

@ -3228,7 +3228,7 @@ contains
end do
filt % bins(Nangle + 1) = PI
else
call fatal_error("Number of bins for mu filter must be&
call fatal_error("Number of bins for polar filter must be&
& greater than 1 on tally " &
// trim(to_str(t % id)) // ".")
end if
@ -3262,7 +3262,7 @@ contains
end do
filt % bins(Nangle + 1) = PI
else
call fatal_error("Number of bins for mu filter must be&
call fatal_error("Number of bins for azimuthal filter must be&
& greater than 1 on tally " &
// trim(to_str(t % id)) // ".")
end if

2
src/mgxs_data.F90
View file

@ -183,7 +183,7 @@ contains
end select
! Do not read materials which we do not actually use in the problem to
! save space
! reduce storage
if (allocated(kTs(i_mat) % data)) then
call macro_xs(i_mat) % obj % combine(kTs(i_mat), mat, nuclides_MG, &
num_energy_groups, num_delayed_groups, max_order, &

36
src/mgxs_header.F90
View file

@ -438,14 +438,14 @@ module mgxs_header
real(8) :: dmu, mu, norm, chi_sum
integer :: order, order_dim, gin, gout, l, imu, length
type(VectorInt) :: temps_to_read
integer :: t, dg
integer :: t, dg, order_data
type(Jagged2D), allocatable :: input_scatt(:), scatt_coeffs(:)
type(Jagged1D), allocatable :: temp_mult(:)
integer, allocatable :: gmin(:), gmax(:)
! Call generic data gathering routine (will populate the metadata)
call mgxs_from_hdf5(this, xs_id, temperature, method, tolerance, &
temps_to_read, order_dim)
temps_to_read, order_data)
! Set the number of delayed groups
this % num_delayed_groups = delayed_groups
@ -1009,7 +1009,7 @@ module mgxs_header
length = 0
do gin = 1, energy_groups
length = length + order_dim * (gmax(gin) - gmin(gin) + 1)
length = length + order_data * (gmax(gin) - gmin(gin) + 1)
end do
! Allocate flattened array
@ -1022,8 +1022,10 @@ module mgxs_header
! Compare the number of orders given with the maximum order of the
! problem. Strip off the supefluous orders if needed.
if (this % scatter_format == ANGLE_LEGENDRE) then
order = min(order_dim - 1, max_order)
order = min(order_data - 1, max_order)
order_dim = order + 1
else
order_dim = order_data
end if
! Convert temp_arr to a jagged array ((gin) % data(l, gout)) for
@ -1038,6 +1040,8 @@ module mgxs_header
input_scatt(gin) % data(l, gout) = temp_arr(index)
index = index + 1
end do
! Adjust index for the orders we didnt take
index = index + (order_data - order_dim)
end do
end do
@ -1108,7 +1112,7 @@ module mgxs_header
deallocate(input_scatt)
! Now get the multiplication matrix
if (object_exists(scatt_grp, 'multiplicity matrix')) then
if (object_exists(scatt_grp, 'multiplicity_matrix')) then
! Now use this information to find the length of a container array
! to hold the flattened data
@ -1120,7 +1124,7 @@ module mgxs_header
! Allocate flattened array
allocate(temp_arr(length))
call read_dataset(temp_arr, scatt_grp, "multiplicity matrix")
call read_dataset(temp_arr, scatt_grp, "multiplicity_matrix")
! Convert temp_arr to a jagged array ((gin) % data(gout)) for
! passing to ScattData
@ -1215,14 +1219,14 @@ module mgxs_header
real(8) :: dmu, mu, norm, chi_sum
integer :: order, order_dim, gin, gout, l, imu, dg
type(VectorInt) :: temps_to_read
integer :: t, length, ipol, iazi
integer :: t, length, ipol, iazi, order_data
type(Jagged2D), allocatable :: input_scatt(:, :, :), scatt_coeffs(:, :, :)
type(Jagged1D), allocatable :: temp_mult(:, :, :)
integer, allocatable :: gmin(:, :, :), gmax(:, :, :)
! Call generic data gathering routine (will populate the metadata)
call mgxs_from_hdf5(this, xs_id, temperature, method, tolerance, &
temps_to_read, order_dim)
temps_to_read, order_data)
! Set the number of delayed groups
this % num_delayed_groups = delayed_groups
@ -1975,7 +1979,7 @@ module mgxs_header
do ipol = 1, this % n_pol
do iazi = 1, this % n_azi
do gin = 1, energy_groups
length = length + order_dim * (gmax(gin, iazi, ipol) - &
length = length + order_data * (gmax(gin, iazi, ipol) - &
gmin(gin, iazi, ipol) + 1)
end do
end do
@ -1991,8 +1995,10 @@ module mgxs_header
! Compare the number of orders given with the maximum order of the
! problem. Strip off the superfluous orders if needed.
if (this % scatter_format == ANGLE_LEGENDRE) then
order = min(order_dim - 1, max_order)
order = min(order_data - 1, max_order)
order_dim = order + 1
else
order_dim = order_data
end if
! Convert temp_1d to a jagged array ((gin) % data(l, gout)) for
@ -2011,6 +2017,8 @@ module mgxs_header
temp_1d(index)
index = index + 1
end do ! gout
! Adjust index for the orders we didnt take
index = index + (order_data - order_dim)
end do ! order
end do ! gin
end do ! iazi
@ -2092,7 +2100,7 @@ module mgxs_header
deallocate(input_scatt)
! Now get the multiplication matrix
if (object_exists(scatt_grp, 'multiplicity matrix')) then
if (object_exists(scatt_grp, 'multiplicity_matrix')) then
! Now use this information to find the length of a container array
! to hold the flattened data
@ -2108,7 +2116,7 @@ module mgxs_header
! Allocate flattened array
allocate(temp_1d(length))
call read_dataset(temp_1d, scatt_grp, "multiplicity matrix")
call read_dataset(temp_1d, scatt_grp, "multiplicity_matrix")
! Convert temp_1d to a jagged array ((gin) % data(gout)) for passing
! to ScattData
@ -2959,9 +2967,9 @@ module mgxs_header
! Now create our jagged data from the dense data
call jagged_from_dense_2D(scatt_coeffs(:, :, :, iazi, ipol), &
jagged_scatt)
jagged_scatt, gmin, gmax)
call jagged_from_dense_1D(temp_mult(:, :, iazi, ipol), &
jagged_mult, gmin, gmax)
jagged_mult)
! Initialize the ScattData Object
call this % xs(t) % scatter(iazi, ipol) % obj % init(gmin, &

8
src/scattdata_header.F90
View file

@ -266,12 +266,12 @@ contains
allocate(matrix(groups))
do gin = 1, groups
allocate(matrix(gin) % data(order, gmin(gin):gmax(gin)))
matrix(gin) % data = coeffs(gin) % data
matrix(gin) % data(:, :) = coeffs(gin) % data(:, :)
end do
! Get scattxs value
allocate(this % scattxs(groups))
! Get this by summing the un-normalized P0 coefficient in matrix
! Get this by summing the un-normalized angular distribution in matrix
! over all outgoing groups
do gin = 1, groups
this % scattxs(gin) = sum(matrix(gin) % data(:, :))
@ -317,7 +317,7 @@ contains
this % dist(gin) % data(imu - 1, gout)
end do
! Now make sure integral norms to zero
! Normalize the integral to unity
norm = this % dist(gin) % data(order, gout)
if (norm > ZERO) then
this % fmu(gin) % data(:, gout) = &
@ -578,7 +578,7 @@ contains
imu = 1
else
imu = binary_search(this % dist(gin) % data(:, gout), &
size(this % dist(gin) % data(:, gout)), xi)
size(this % dist(gin) % data(:, gout)), xi) + 1
end if
! Randomly select a mu in this bin.

2
tests/test_mg_basic/results_true.dat
View file

@ -1,2 +1,2 @@
k-combined:
1.086852E+00 2.677280E-02
1.073147E+00 1.602384E-02

2
tests/test_mg_max_order/results_true.dat
View file

@ -1,2 +1,2 @@
k-combined:
1.077247E+00 2.412834E-02
1.074551E+00 1.871525E-02

2
tests/test_mg_nuclide/results_true.dat
View file

@ -1,2 +1,2 @@
k-combined:
1.086852E+00 2.677280E-02
1.073147E+00 1.602384E-02

2
tests/test_mg_tallies/results_true.dat
View file

@ -1 +1 @@
864328b2c4f3c4bfa9c80c756acbedac6fbdf3d502c03c2f2a3e7461b774729c47a55341c9515eff246e1bd445485f0d0a05a5712663ee499046afdbe0c77aef
1817642a0d20d8ef437c7970cae2e77c265a566ef9843ac1b7cd33f76af09365dc8bd5c045f4297360b38d509864bd629b06dfc4d043b4b1ae7d8a5e2e93fae3