From 9b8dbe6e942f32095b9594f63dfcd82500025306 Mon Sep 17 00:00:00 2001 From: Adam Nelson Date: Sun, 22 May 2016 05:56:29 -0400 Subject: [PATCH] added NuFissionMatrixXS class to MGXS types and incorporated in to Library and XsData --- openmc/mgxs/library.py | 22 ++- openmc/mgxs/mgxs.py | 430 ++++++++++++++++++++++++++++++++++++++++- openmc/mgxs_library.py | 46 +++-- 3 files changed, 474 insertions(+), 24 deletions(-) diff --git a/openmc/mgxs/library.py b/openmc/mgxs/library.py index 4eae2cb067..816b78e3d3 100644 --- a/openmc/mgxs/library.py +++ b/openmc/mgxs/library.py @@ -460,7 +460,10 @@ class Library(object): ---------- domain : Material or Cell or Universe or Integral The material, cell, or universe object of interest (or its ID) - mgxs_type : {'total', 'transport', 'nu-transport', 'absorption', 'capture', 'fission', 'nu-fission', 'kappa-fission', 'scatter', 'nu-scatter', 'scatter matrix', 'nu-scatter matrix', 'chi'} + mgxs_type : {'total', 'transport', 'nu-transport', 'absorption', + 'capture', 'fission', 'nu-fission', 'kappa-fission', + 'scatter', 'nu-scatter', 'scatter matrix', + 'nu-scatter matrix', 'nu-fission matrix', chi'} The type of multi-group cross section object to return Returns @@ -853,13 +856,20 @@ class Library(object): mymgxs = self.get_mgxs(domain, 'kappa-fission') xsdata.set_kappa_fission_mgxs(mymgxs, xs_type=xs_type, nuclide=[nuclide]) - if 'chi' in self.mgxs_types: - mymgxs = self.get_mgxs(domain, 'chi') - xsdata.set_chi_mgxs(mymgxs, xs_type=xs_type, nuclide=[nuclide]) - if 'nu-fission' in self.mgxs_types: - mymgxs = self.get_mgxs(domain, 'nu-fission') + # For chi and nu-fission we can either have only a nu-fission matrix + # provided, or vectors of chi and nu-fission provided + if 'nu-fission matrix' in self.mgxs_types: + mymgxs = self.get_mgxs(domain, 'nu-fission matrix') xsdata.set_nu_fission_mgxs(mymgxs, xs_type=xs_type, nuclide=[nuclide]) + else: + if 'chi' in self.mgxs_types: + mymgxs = self.get_mgxs(domain, 'chi') + xsdata.set_chi_mgxs(mymgxs, xs_type=xs_type, nuclide=[nuclide]) + if 'nu-fission' in self.mgxs_types: + mymgxs = self.get_mgxs(domain, 'nu-fission') + xsdata.set_nu_fission_mgxs(mymgxs, xs_type=xs_type, + nuclide=[nuclide]) # multiplicity requires scatter and nu-scatter if ((('scatter matrix' in self.mgxs_types) and ('nu-scatter matrix' in self.mgxs_types))): diff --git a/openmc/mgxs/mgxs.py b/openmc/mgxs/mgxs.py index 5be84bb2ce..7173349376 100644 --- a/openmc/mgxs/mgxs.py +++ b/openmc/mgxs/mgxs.py @@ -32,6 +32,7 @@ MGXS_TYPES = ['total', 'nu-scatter', 'scatter matrix', 'nu-scatter matrix', + 'nu-fission matrix', 'chi'] @@ -425,7 +426,10 @@ class MGXS(object): Parameters ---------- - mgxs_type : {'total', 'transport', 'nu-transport', 'absorption', 'capture', 'fission', 'nu-fission', 'kappa-fission', 'scatter', 'nu-scatter', 'scatter matrix', 'nu-scatter matrix', 'chi'} + mgxs_type : {'total', 'transport', 'nu-transport', 'absorption', + 'capture', 'fission', 'nu-fission', 'kappa-fission', + 'scatter', 'nu-scatter', 'scatter matrix', + 'nu-scatter matrix', 'nu-fission matrix', 'chi'} The type of multi-group cross section object to return domain : openmc.Material or openmc.Cell or openmc.Universe The domain for spatial homogenization @@ -474,6 +478,8 @@ class MGXS(object): mgxs = ScatterMatrixXS(domain, domain_type, energy_groups) elif mgxs_type == 'nu-scatter matrix': mgxs = NuScatterMatrixXS(domain, domain_type, energy_groups) + elif mgxs_type == 'nu-fission matrix': + mgxs = NuFissionMatrixXS(domain, domain_type, energy_groups) elif mgxs_type == 'chi': mgxs = Chi(domain, domain_type, energy_groups) @@ -3182,6 +3188,428 @@ class NuScatterMatrixXS(ScatterMatrixXS): self._hdf5_key = 'nu-scatter matrix' +class NuFissionMatrixXS(MGXS): + """A fission production matrix multi-group cross section. + + This class can be used for both OpenMC input generation and tally data + post-processing to compute spatially-homogenized and energy-integrated + multi-group cross sections for deterministic neutronics calculations. + + Parameters + ---------- + domain : openmc.Material or openmc.Cell or openmc.Universe + The domain for spatial homogenization + domain_type : {'material', 'cell', 'distribcell', 'universe'} + The domain type for spatial homogenization + energy_groups : openmc.mgxs.EnergyGroups + The energy group structure for energy condensation + by_nuclide : bool + If true, computes cross sections for each nuclide in domain + name : str, optional + Name of the multi-group cross section. Used as a label to identify + tallies in OpenMC 'tallies.xml' file. + + Attributes + ---------- + name : str, optional + Name of the multi-group cross section + rxn_type : str + Reaction type (e.g., 'total', 'nu-fission', etc.) + by_nuclide : bool + If true, computes cross sections for each nuclide in domain + domain : Material or Cell or Universe + Domain for spatial homogenization + domain_type : {'material', 'cell', 'distribcell', 'universe'} + Domain type for spatial homogenization + energy_groups : openmc.mgxs.EnergyGroups + Energy group structure for energy condensation + tally_trigger : openmc.Trigger + An (optional) tally precision trigger given to each tally used to + compute the cross section + scores : list of str + The scores in each tally used to compute the multi-group cross section + filters : list of openmc.Filter + The filters in each tally used to compute the multi-group cross section + tally_keys : list of str + The keys into the tallies dictionary for each tally used to compute + the multi-group cross section + estimator : {'tracklength', 'analog'} + The tally estimator used to compute the multi-group cross section + tallies : collections.OrderedDict + OpenMC tallies needed to compute the multi-group cross section + rxn_rate_tally : openmc.Tally + Derived tally for the reaction rate tally used in the numerator to + compute the multi-group cross section. This attribute is None + unless the multi-group cross section has been computed. + xs_tally : openmc.Tally + Derived tally for the multi-group cross section. This attribute + is None unless the multi-group cross section has been computed. + num_subdomains : int + The number of subdomains is unity for 'material', 'cell' and 'universe' + domain types. When the This is equal to the number of cell instances + for 'distribcell' domain types (it is equal to unity prior to loading + tally data from a statepoint file). + num_nuclides : int + The number of nuclides for which the multi-group cross section is + being tracked. This is unity if the by_nuclide attribute is False. + nuclides : Iterable of str or 'sum' + The optional user-specified nuclides for which to compute cross + sections (e.g., 'U-238', 'O-16'). If by_nuclide is True but nuclides + are not specified by the user, all nuclides in the spatial domain + are included. This attribute is 'sum' if by_nuclide is false. + sparse : bool + Whether or not the MGXS' tallies use SciPy's LIL sparse matrix format + for compressed data storage + loaded_sp : bool + Whether or not a statepoint file has been loaded with tally data + derived : bool + Whether or not the MGXS is merged from one or more other MGXS + hdf5_key : str + The key used to index multi-group cross sections in an HDF5 data store + + """ + + def __init__(self, domain=None, domain_type=None, + groups=None, by_nuclide=False, name=''): + super(NuFissionMatrixXS, self).__init__(domain, domain_type, + groups, by_nuclide, name) + self._rxn_type = 'nu-fission matrix' + self._hdf5_key = 'nu-fission matrix' + + def __deepcopy__(self, memo): + clone = super(NuFissionMatrixXS, self).__deepcopy__(memo) + return clone + + @property + def scores(self): + scores = ['flux', 'nu-fission'] + + return scores + + @property + def filters(self): + group_edges = self.energy_groups.group_edges + energy = openmc.Filter('energy', group_edges) + energyout = openmc.Filter('energyout', group_edges) + + filters = [[energy], [energy, energyout]] + + return filters + + @property + def estimator(self): + return 'analog' + + @property + def rxn_rate_tally(self): + if self._rxn_rate_tally is None: + self._rxn_rate_tally = self.tallies['nu-fission'] + self._rxn_rate_tally.sparse = self.sparse + return self._rxn_rate_tally + + def get_slice(self, nuclides=[], in_groups=[], out_groups=[]): + """Build a sliced NuFissionMatrix for the specified nuclides and + energy groups. + + This method constructs a new MGXS to encapsulate a subset of the data + represented by this MGXS. The subset of data to include in the tally + slice is determined by the nuclides and energy groups specified in + the input parameters. + + Parameters + ---------- + nuclides : list of str + A list of nuclide name strings + (e.g., ['U-235', 'U-238']; default is []) + in_groups : list of int + A list of incoming energy group indices starting at 1 for the high + energies (e.g., [1, 2, 3]; default is []) + out_groups : list of int + A list of outgoing energy group indices starting at 1 for the high + energies (e.g., [1, 2, 3]; default is []) + + Returns + ------- + openmc.mgxs.MGXS + A new tally which encapsulates the subset of data requested for the + nuclide(s) and/or energy group(s) requested in the parameters. + + """ + + # Call super class method and null out derived tallies + slice_xs = super(NuFissionMatrixXS, self).get_slice(nuclides, + in_groups) + slice_xs._rxn_rate_tally = None + slice_xs._xs_tally = None + + # Slice outgoing energy groups if needed + if len(out_groups) != 0: + filter_bins = [] + for group in out_groups: + group_bounds = self.energy_groups.get_group_bounds(group) + filter_bins.append(group_bounds) + filter_bins = [tuple(filter_bins)] + + # Slice each of the tallies across energyout groups + for tally_type, tally in slice_xs.tallies.items(): + if tally.contains_filter('energyout'): + tally_slice = tally.get_slice(filters=['energyout'], + filter_bins=filter_bins) + slice_xs.tallies[tally_type] = tally_slice + + slice_xs.sparse = self.sparse + return slice_xs + + def get_xs(self, in_groups='all', out_groups='all', + subdomains='all', nuclides='all', + xs_type='macro', order_groups='increasing', + row_column='inout', value='mean', **kwargs): + r"""Returns an array of multi-group cross sections. + + This method constructs a 2D NumPy array for the requested scattering + matrix data data for one or more energy groups and subdomains. + + NOTE: The scattering moments are not multiplied by the :math:`(2l+1)/2` + prefactor in the expansion of the scattering source into Legendre + moments in the neutron transport equation. + + Parameters + ---------- + in_groups : Iterable of Integral or 'all' + Incoming energy groups of interest. Defaults to 'all'. + out_groups : Iterable of Integral or 'all' + Outgoing energy groups of interest. Defaults to 'all'. + subdomains : Iterable of Integral or 'all' + Subdomain IDs of interest. Defaults to 'all'. + nuclides : Iterable of str or 'all' or 'sum' + A list of nuclide name strings (e.g., ['U-235', 'U-238']). The + special string 'all' will return the cross sections for all nuclides + in the spatial domain. The special string 'sum' will return the + cross section summed over all nuclides. Defaults to 'all'. + xs_type: {'macro', 'micro'} + Return the macro or micro cross section in units of cm^-1 or barns. + Defaults to 'macro'. + order_groups: {'increasing', 'decreasing'} + Return the cross section indexed according to increasing or + decreasing energy groups (decreasing or increasing energies). + Defaults to 'increasing'. + row_column: {'inout', 'outin'} + Return the cross section indexed first by incoming group and + second by outgoing group ('inout'), or vice versa ('outin'). + Defaults to 'inout'. + value : str + A string for the type of value to return - 'mean', 'std_dev', or + 'rel_err' are accepted. Defaults to the empty string. + + Returns + ------- + ndarray + A NumPy array of the multi-group cross section indexed in the order + each group and subdomain is listed in the parameters. + + Raises + ------ + ValueError + When this method is called before the multi-group cross section is + computed from tally data. + + """ + + cv.check_value('value', value, ['mean', 'std_dev', 'rel_err']) + cv.check_value('xs_type', xs_type, ['macro', 'micro']) + + filters = [] + filter_bins = [] + + # Construct a collection of the domain filter bins + if not isinstance(subdomains, basestring): + cv.check_iterable_type('subdomains', subdomains, Integral, max_depth=2) + for subdomain in subdomains: + filters.append(self.domain_type) + filter_bins.append((subdomain,)) + + # Construct list of energy group bounds tuples for all requested groups + if not isinstance(in_groups, basestring): + cv.check_iterable_type('groups', in_groups, Integral) + for group in in_groups: + filters.append('energy') + filter_bins.append((self.energy_groups.get_group_bounds(group),)) + + # Construct list of energy group bounds tuples for all requested groups + if not isinstance(out_groups, basestring): + cv.check_iterable_type('groups', out_groups, Integral) + for group in out_groups: + filters.append('energyout') + filter_bins.append((self.energy_groups.get_group_bounds(group),)) + + # Construct a collection of the nuclides to retrieve from the xs tally + if self.by_nuclide: + if nuclides == 'all' or nuclides == 'sum' or nuclides == ['sum']: + query_nuclides = self.get_all_nuclides() + else: + query_nuclides = nuclides + else: + query_nuclides = ['total'] + + # Use tally summation if user requested the sum for all nuclides + if nuclides == 'sum' or nuclides == ['sum']: + xs_tally = self.xs_tally.summation(nuclides=query_nuclides) + xs = xs_tally.get_values(filters=filters, filter_bins=filter_bins, + value=value) + else: + xs = self.xs_tally.get_values(filters=filters, + filter_bins=filter_bins, + nuclides=query_nuclides, value=value) + + xs = np.nan_to_num(xs) + + # Divide by atom number densities for microscopic cross sections + if xs_type == 'micro': + if self.by_nuclide: + densities = self.get_nuclide_densities(nuclides) + else: + densities = self.get_nuclide_densities('sum') + if value == 'mean' or value == 'std_dev': + xs /= densities[np.newaxis, :, np.newaxis] + + # Reverse data if user requested increasing energy groups since + # tally data is stored in order of increasing energies + if order_groups == 'increasing': + if in_groups == 'all': + num_in_groups = self.num_groups + else: + num_in_groups = len(in_groups) + if out_groups == 'all': + num_out_groups = self.num_groups + else: + num_out_groups = len(out_groups) + + # Reshape tally data array with separate axes for domain and energy + num_subdomains = int(xs.shape[0] / + (num_in_groups * num_out_groups)) + new_shape = (num_subdomains, num_in_groups, num_out_groups) + new_shape += xs.shape[1:] + xs = np.reshape(xs, new_shape) + + # Transpose the scattering matrix if requested by user + if row_column == 'outin': + xs = np.swapaxes(xs, 1, 2) + + # Reverse energies to align with increasing energy groups + xs = xs[:, ::-1, ::-1, :] + + # Eliminate trivial dimensions + xs = np.squeeze(xs) + xs = np.atleast_2d(xs) + + return xs + + def print_xs(self, subdomains='all', nuclides='all', + xs_type='macro'): + """Prints a string representation for the multi-group cross section. + + Parameters + ---------- + subdomains : Iterable of Integral or 'all' + The subdomain IDs of the cross sections to include in the report. + Defaults to 'all'. + nuclides : Iterable of str or 'all' or 'sum' + The nuclides of the cross-sections to include in the report. This + may be a list of nuclide name strings (e.g., ['U-235', 'U-238']). + The special string 'all' will report the cross sections for all + nuclides in the spatial domain. The special string 'sum' will report + the cross sections summed over all nuclides. Defaults to 'all'. + xs_type: {'macro', 'micro'} + Return the macro or micro cross section in units of cm^-1 or barns. + Defaults to 'macro'. + + """ + + # Construct a collection of the subdomains to report + if not isinstance(subdomains, basestring): + cv.check_iterable_type('subdomains', subdomains, Integral) + elif self.domain_type == 'distribcell': + subdomains = np.arange(self.num_subdomains, dtype=np.int) + else: + subdomains = [self.domain.id] + + # Construct a collection of the nuclides to report + if self.by_nuclide: + if nuclides == 'all': + nuclides = self.get_all_nuclides() + if nuclides == 'sum': + nuclides = ['sum'] + else: + cv.check_iterable_type('nuclides', nuclides, basestring) + else: + nuclides = ['sum'] + + cv.check_value('xs_type', xs_type, ['macro', 'micro']) + + # Build header for string with type and domain info + string = 'Multi-Group XS\n' + string += '{0: <16}=\t{1}\n'.format('\tReaction Type', self.rxn_type) + string += '{0: <16}=\t{1}\n'.format('\tDomain Type', self.domain_type) + string += '{0: <16}=\t{1}\n'.format('\tDomain ID', self.domain.id) + + # If cross section data has not been computed, only print string header + if self.tallies is None: + print(string) + return + + string += '{0: <16}\n'.format('\tEnergy Groups:') + template = '{0: <12}Group {1} [{2: <10} - {3: <10}MeV]\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]) + + # Loop over all subdomains + for subdomain in subdomains: + + if self.domain_type == 'distribcell': + string += \ + '{0: <16}=\t{1}\n'.format('\tSubdomain', subdomain) + + # Loop over all Nuclides + for nuclide in nuclides: + + # Build header for nuclide type + if xs_type != 'sum': + string += '{0: <16}=\t{1}\n'.format('\tNuclide', nuclide) + + # Build header for cross section type + if xs_type == 'macro': + string += '{0: <16}\n'.format('\tCross Sections [cm^-1]:') + else: + string += '{0: <16}\n'.format('\tCross Sections [barns]:') + + 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 += '{:1.2e} +/- {:1.2e}%'.format(average, rel_err) + string += '\n' + string += '\n' + string += '\n' + string += '\n' + + print(string) + + class Chi(MGXS): """The fission spectrum. diff --git a/openmc/mgxs_library.py b/openmc/mgxs_library.py index 88ae05808a..7559b427a1 100644 --- a/openmc/mgxs_library.py +++ b/openmc/mgxs_library.py @@ -131,6 +131,8 @@ class XSdata(object): num_polar : int Number of equal width angular bins that the polar angular domain is subdivided into. This only applies when ``representation`` is "angle". + use_chi : bool + Whether or not a chi vector or nu-fission matrix was used. vector_shape : iterable of int Dimensionality of vector multi-group cross sections (e.g., the total cross section). The return result depends on the value of @@ -292,6 +294,10 @@ class XSdata(object): def num_azimuthal(self): return self._num_azimuthal + @property + def use_chi(self): + return self._use_chi + @property def total(self): return self._total @@ -461,6 +467,11 @@ class XSdata(object): check_greater_than('num_azimuthal', num_azimuthal, 0) self._num_azimuthal = num_azimuthal + @use_chi.setter + def use_chi(self, use_chi): + check_type('use_chi', use_chi, bool) + self._use_chi = use_chi + @total.setter def total(self, total): check_type('total', total, Iterable, expected_iter_type=Real) @@ -512,8 +523,8 @@ class XSdata(object): @chi.setter def chi(self, chi): - if self._use_chi is not None: - if not self._use_chi: + if self.use_chi is not None: + if not self.use_chi: msg = 'Providing chi when nu_fission already provided as a' \ 'matrix' raise ValueError(msg) @@ -529,8 +540,8 @@ class XSdata(object): self._chi = npchi - if self._use_chi is not None: - self._use_chi = True + if self.use_chi is not None: + self.use_chi = True @scatter.setter def scatter(self, scatter): @@ -574,8 +585,8 @@ class XSdata(object): check_iterable_type('nu_fission', npnu_fission, Real, max_depth=len(npnu_fission.shape)) - if self._use_chi is not None: - if self._use_chi: + if self.use_chi is not None: + if self.use_chi: check_value('nu_fission shape', npnu_fission.shape, [self.vector_shape]) else: @@ -587,9 +598,9 @@ class XSdata(object): # Find out if we have a nu-fission matrix or vector # and set a flag to allow other methods to check this later. if npnu_fission.shape == self.vector_shape: - self._use_chi = True + self.use_chi = True else: - self._use_chi = False + self.use_chi = False self._nu_fission = npnu_fission if np.sum(self._nu_fission) > 0.0: @@ -728,10 +739,8 @@ class XSdata(object): """ - # The NuFissionXS class does not have the capability to produce - # a fission matrix and therefore if this path is pursued, we know - # chi must be used. - check_type('nu_fission', nu_fission, openmc.mgxs.NuFissionXS) + check_type('nu_fission', nu_fission, (openmc.mgxs.NuFissionXS, + openmc.mgxs.NuFissionMatrixXS)) check_value('energy_groups', nu_fission.energy_groups, [self.energy_groups]) check_value('domain_type', nu_fission.domain_type, @@ -744,7 +753,10 @@ class XSdata(object): msg = 'Angular-Dependent MGXS have not yet been implemented' raise ValueError(msg) - self._use_chi = True + if type(nu_fission) is openmc.mgxs.NuFissionMatrixXS: + self.use_chi = False + else: + self.use_chi = True if np.sum(self._nu_fission) > 0.0: self._fissionable = True @@ -809,8 +821,8 @@ class XSdata(object): """ - if self._use_chi is not None: - if not self._use_chi: + if self.use_chi is not None: + if not self.use_chi: msg = 'Providing chi when nu_fission already provided as a ' \ 'matrix!' raise ValueError(msg) @@ -827,8 +839,8 @@ class XSdata(object): msg = 'Angular-Dependent MGXS have not yet been implemented' raise ValueError(msg) - if self._use_chi is not None: - self._use_chi = True + if self.use_chi is not None: + self.use_chi = True def set_scatter_mgxs(self, scatter, nuclide='total', xs_type='macro'): """This method allows for an openmc.mgxs.ScatterMatrixXS