openmc.mgxs.Library

class openmc.mgxs.Library(geometry, by_nuclide=False, mgxs_types=None, name='')[source]

A multi-energy-group and multi-delayed-group cross section library for some energy group structure.

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.

This class helps automate the generation of MGXS and MDGXS objects for some energy group structure and domain type. The Library serves as a collection for MGXS and MDGXS objects with routines to automate the initialization of tallies for input files, the loading of tally data from statepoint files, data storage, energy group condensation and more.

Parameters:
  • geometry (openmc.Geometry) – A geometry which has been initialized with a root universe
  • by_nuclide (bool) – If true, computes cross sections for each nuclide in each domain
  • mgxs_types (Iterable of str) – The types of cross sections in the library (e.g., [‘total’, ‘scatter’])
  • name (str, optional) – Name of the multi-group cross section library. Used as a label to identify tallies in OpenMC ‘tallies.xml’ file.
Variables:
  • geometry (openmc.Geometry) – An geometry which has been initialized with a root universe
  • by_nuclide (bool) – If true, computes cross sections for each nuclide in each domain
  • mgxs_types (Iterable of str) – The types of cross sections in the library (e.g., [‘total’, ‘scatter’])
  • domain_type ({'material', 'cell', 'distribcell', 'universe', 'mesh'}) – Domain type for spatial homogenization
  • domains (Iterable of openmc.Material, openmc.Cell, openmc.Universe or openmc.RegularMesh) – The spatial domain(s) for which MGXS in the Library are computed
  • correction ({'P0', None}) – Apply the P0 correction to scattering matrices if set to ‘P0’
  • scatter_format ({'legendre', 'histogram'}) – Representation of the angular scattering distribution (default is ‘legendre’)
  • legendre_order (int) – The highest Legendre moment in the scattering matrix; this is used if ScatterMatrixXS.scatter_format is ‘legendre’. (default is 0)
  • histogram_bins (int) – The number of equally-spaced bins for the histogram representation of the angular scattering distribution; this is used if ScatterMatrixXS.scatter_format is ‘histogram’. (default is 16)
  • energy_groups (openmc.mgxs.EnergyGroups) – 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.
  • tally_trigger (openmc.Trigger) – An (optional) tally precision trigger given to each tally used to compute the cross section
  • all_mgxs (collections.OrderedDict) – MGXS objects keyed by domain ID and cross section type
  • sp_filename (str) – The filename of the statepoint with tally data used to the compute cross sections
  • keff (Real or None) – The combined keff from the statepoint file with tally data used to compute cross sections (for eigenvalue calculations only)
  • name (str, optional) – Name of the multi-group cross section library. Used as a label to identify tallies in OpenMC ‘tallies.xml’ file.
  • sparse (bool) – Whether or not the Library’s tallies use SciPy’s LIL sparse matrix format for compressed data storage
add_to_tallies_file(tallies_file, merge=True)[source]

Add all tallies from all MGXS objects to a tallies file.

NOTE: This assumes that Library.build_library() has been called

Parameters:
  • tallies_file (openmc.Tallies) – A Tallies collection to add each MGXS’ tallies to generate a ‘tallies.xml’ input file for OpenMC
  • merge (bool) – Indicate whether tallies should be merged when possible. Defaults to True.
build_hdf5_store(filename='mgxs.h5', directory='mgxs', subdomains='all', nuclides='all', xs_type='macro', row_column='inout', libver='earliest')[source]

Export the multi-group cross section library to an HDF5 binary file.

This method constructs an HDF5 file which stores the library’s multi-group cross section data. The data is stored in a hierarchy of HDF5 groups from the domain type, domain id, subdomain id (for distribcell domains), nuclides and cross section types. Two datasets for the mean and standard deviation are stored for each subdomain entry in the HDF5 file. The number of groups is stored as a file attribute.

NOTE: This requires the h5py Python package.

Parameters:
  • filename (str) – Filename for the HDF5 file. Defaults to ‘mgxs.h5’.
  • directory (str) – Directory for the HDF5 file. Defaults to ‘mgxs’.
  • subdomains ({'all', 'avg'}) – Report all subdomains or the average of all subdomain cross sections in the report. Defaults to ‘all’.
  • nuclides ({'all', 'sum'}) – The nuclides of the cross-sections to include in the report. This may be a list of nuclide name strings (e.g., [‘U235’, ‘U238’]). 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'}) – Store the macro or micro cross section in units of cm^-1 or barns. Defaults to ‘macro’.
  • row_column ({'inout', 'outin'}) – Store scattering matrices indexed first by incoming group and second by outgoing group (‘inout’), or vice versa (‘outin’). Defaults to ‘inout’.
  • libver ({'earliest', 'latest'}) – Compatibility mode for the HDF5 file. ‘latest’ will produce files that are less backwards compatible but have performance benefits.
Raises:

ValueError – When this method is called before a statepoint has been loaded

See also

MGXS.build_hdf5_store(), directory(), xs_type()

build_library()[source]

Initialize MGXS objects in each domain and for each reaction type in the library.

This routine will populate the all_mgxs instance attribute dictionary with MGXS subclass objects keyed by each domain ID (e.g., Material IDs) and cross section type (e.g., ‘nu-fission’, ‘total’, etc.).

check_library_for_openmc_mgxs()[source]

This routine will check the MGXS Types within a Library to ensure the MGXS types provided can be used to create a MGXS Library for OpenMC’s Multi-Group mode.

The rules to check include:

  • Either total or transport must be present.
    • Both can be available if one wants, but we should use whatever corresponds to Library.correction (if P0: transport)
  • Absorption is required.
  • A nu-fission cross section and chi values are not required as a fixed source problem could be the target.
  • Fission and kappa-fission are not required as they are only needed to support tallies the user may wish to request.
  • Scattering multiplicity should have been tallied for increased model accuracy, either using a multiplicity or scatter and nu-scatter matrix tally.
create_mg_library(xs_type='macro', xsdata_names=None, apply_domain_chi=False)[source]

Creates an openmc.MGXSLibrary object to contain the MGXS data for the Multi-Group mode of OpenMC.

Note that this library will not make use of nested temperature tables. Every dataset in the library will be treated as if it was at the same default temperature.

Parameters:
  • xs_type ({'macro', 'micro'}) – Provide the macro or micro cross section in units of cm^-1 or barns. Defaults to ‘macro’. If the Library object is not tallied by nuclide this will be set to ‘macro’ regardless.
  • xsdata_names (Iterable of str) – List of names to apply to the “xsdata” entries in the resultant mgxs data file. Defaults to ‘set1’, ‘set2’, …
  • apply_domain_chi (bool) – This parameter sets whether (True) or not (False) the domain-averaged values of chi, chi-prompt, and chi-delayed are to be applied to each of the nuclide-dependent fission energy spectra of a domain. In effect, if this is True, then every nuclide in the domain receives the same flux-weighted Chi. This is useful for downstream multigroup solvers that precompute a material-specific chi before the transport solve provides group-wise fluxes. Defaults to False.
Returns:

mgxs_file – Multi-Group Cross Section File that is ready to be printed to the file of choice by the user.

Return type:

openmc.MGXSLibrary

Raises:

ValueError – When the Library object is initialized with insufficient types of cross sections for the Library.

create_mg_mode(xsdata_names=None, bc=['reflective', 'reflective', 'reflective', 'reflective', 'reflective', 'reflective'], apply_domain_chi=False)[source]

Creates an openmc.MGXSLibrary object to contain the MGXS data for the Multi-Group mode of OpenMC as well as the associated openmc.Materials and openmc.Geometry objects.

The created Geometry is the same as that used to generate the MGXS data, with the only differences being modifications to point to newly-created Materials which point to the multi-group data. This method only creates a macroscopic MGXS Library even if nuclidic tallies are specified in the Library. Note that this library will not make use of nested temperature tables. Every dataset in the library will be treated as if it was at the same default temperature.

Parameters:
  • xsdata_names (Iterable of str) – List of names to apply to the “xsdata” entries in the resultant mgxs data file. Defaults to ‘set1’, ‘set2’, …
  • bc (iterable of {'reflective', 'periodic', 'transmission', or 'vacuum'}) – Boundary conditions for each of the four faces of a rectangle (if applying to a 2D mesh) or six faces of a parallelepiped (if applying to a 3D mesh) provided in the following order: [x min, x max, y min, y max, z min, z max]. 2-D cells do not contain the z min and z max entries.
  • apply_domain_chi (bool) – This parameter sets whether (True) or not (False) the domain-averaged values of chi, chi-prompt, and chi-delayed are to be applied to each of the nuclide-dependent fission energy spectra of a domain. In effect, if this is True, then every nuclide in the domain receives the same flux-weighted Chi. This is useful for downstream multigroup solvers that precompute a material-specific chi before the transport solve provides group-wise fluxes. Defaults to False.
Returns:

  • mgxs_file (openmc.MGXSLibrary) – Multi-Group Cross Section File that is ready to be printed to the file of choice by the user.
  • materials (openmc.Materials) – Materials file ready to be printed with all the macroscopic data present within this Library.
  • geometry (openmc.Geometry) – Materials file ready to be printed with all the macroscopic data present within this Library.

Raises:

ValueError – When the Library object is initialized with insufficient types of cross sections for the Library.

dump_to_file(filename='mgxs', directory='mgxs')[source]

Store this Library object in a pickle binary file.

Parameters:
  • filename (str) – Filename for the pickle file. Defaults to ‘mgxs’.
  • directory (str) – Directory for the pickle file. Defaults to ‘mgxs’.

See also

Library.load_from_file(), directory()

get_condensed_library(coarse_groups)[source]

Construct an energy-condensed version of this library.

This routine condenses each of the multi-group cross sections in the library to a coarse energy group structure. NOTE: This routine must be called after the load_from_statepoint(…) routine loads the tallies from the statepoint into each of the cross sections.

Parameters:coarse_groups (openmc.mgxs.EnergyGroups) – The coarse energy group structure of interest
Returns:A new multi-group cross section library condensed to the group structure of interest
Return type:openmc.mgxs.Library
Raises:ValueError – When this method is called before a statepoint has been loaded
get_mgxs(domain, mgxs_type)[source]

Return the MGXS object for some domain and reaction rate type.

This routine searches the library for an MGXS object for the spatial domain and reaction rate type requested by the user.

NOTE: This routine must be called after the build_library() routine.

Parameters:
  • domain (openmc.Material or openmc.Cell or openmc.Universe or openmc.RegularMesh or Integral) – The material, cell, or universe object of interest (or its ID)
  • mgxs_type (str) – The type of multi-group cross section object to return; allowable values are those MGXS to the Library and present in the mgxs_types attribute.
Returns:

The MGXS object for the requested domain and reaction rate type

Return type:

openmc.mgxs.MGXS

Raises:

ValueError – If no MGXS object can be found for the requested domain or multi-group cross section type

get_subdomain_avg_library()[source]

Construct a subdomain-averaged version of this library.

This routine averages each multi-group cross section across distribcell instances. The method performs spatial homogenization to compute the scalar flux-weighted average cross section across the subdomains.

NOTE: This method is only relevant for distribcell domain types and simplys returns a deep copy of the library for all other domains types.

Returns:A new multi-group cross section library averaged across subdomains
Return type:openmc.mgxs.Library
Raises:ValueError – When this method is called before a statepoint has been loaded
get_xsdata(domain, xsdata_name, nuclide='total', xs_type='macro', subdomain=None, apply_domain_chi=False)[source]

Generates an openmc.XSdata object describing a multi-group cross section dataset for writing to an openmc.MGXSLibrary object.

Note that this method does not build an XSdata object with nested temperature tables. The temperature of each XSdata object will be left at the default value of 294K.

Parameters:
  • domain (openmc.Material or openmc.Cell or openmc.Universe or openmc.RegularMesh) – The domain for spatial homogenization
  • xsdata_name (str) – Name to apply to the “xsdata” entry produced by this method
  • nuclide (str) – A nuclide name string (e.g., ‘U235’). Defaults to ‘total’ to obtain a material-wise macroscopic cross section.
  • xs_type ({'macro', 'micro'}) – Provide the macro or micro cross section in units of cm^-1 or barns. Defaults to ‘macro’. If the Library object is not tallied by nuclide this will be set to ‘macro’ regardless.
  • subdomain (iterable of int) – This parameter is not used unless using a mesh domain. In that case, the subdomain is an [i,j,k] index (1-based indexing) of the mesh cell of interest in the openmc.RegularMesh object. Note: this parameter currently only supports subdomains within a mesh, and not the subdomains of a distribcell.
  • apply_domain_chi (bool) – This parameter sets whether (True) or not (False) the domain-averaged values of chi, chi-prompt, and chi-delayed are to be applied to each of the nuclide-dependent fission energy spectra of a domain. In effect, if this is True, then every nuclide in the domain receives the same flux-weighted Chi. This is useful for downstream multigroup solvers that precompute a material-specific chi before the transport solve provides group-wise fluxes. Defaults to False.
Returns:

xsdata – Multi-Group Cross Section dataset object.

Return type:

openmc.XSdata

Raises:

ValueError – When the Library object is initialized with insufficient types of cross sections for the Library.

static load_from_file(filename='mgxs', directory='mgxs')[source]

Load a Library object from a pickle binary file.

Parameters:
  • filename (str) – Filename for the pickle file. Defaults to ‘mgxs’.
  • directory (str) – Directory for the pickle file. Defaults to ‘mgxs’.
Returns:

A Library object loaded from the pickle binary file

Return type:

openmc.mgxs.Library

See also

Library.dump_to_file(), filename(), directory()

load_from_statepoint(statepoint)[source]

Extracts tallies in an OpenMC StatePoint with the data needed to compute multi-group cross sections.

This method is needed to compute cross section data from tallies in an OpenMC StatePoint object.

NOTE: The statepoint must first be linked with an OpenMC Summary object.

Parameters:statepoint (openmc.StatePoint) – An OpenMC StatePoint object with tally data
Raises:ValueError – When this method is called with a statepoint that has not been linked with a summary object.