openmc.Tally

class openmc.Tally(tally_id=None, name='')[source]

A tally defined by a set of scores that are accumulated for a list of nuclides given a set of filters.

Parameters:
  • tally_id (int, optional) – Unique identifier for the tally. If none is specified, an identifier will automatically be assigned
  • name (str, optional) – Name of the tally. If not specified, the name is the empty string.
Variables:
  • id (int) – Unique identifier for the tally
  • name (str) – Name of the tally
  • filters (list of openmc.Filter) – List of specified filters for the tally
  • nuclides (list of openmc.Nuclide) – List of nuclides to score results for
  • scores (list of str) – List of defined scores, e.g. ‘flux’, ‘fission’, etc.
  • estimator ({'analog', 'tracklength', 'collision'}) – Type of estimator for the tally
  • triggers (list of openmc.Trigger) – List of tally triggers
  • num_scores (int) – Total number of scores
  • num_filter_bins (int) – Total number of filter bins accounting for all filters
  • num_bins (int) – Total number of bins for the tally
  • shape (3-tuple of int) – The shape of the tally data array ordered as the number of filter bins, nuclide bins and score bins
  • filter_strides (list of int) – Stride in memory for each filter
  • num_realizations (int) – Total number of realizations
  • with_summary (bool) – Whether or not a Summary has been linked
  • sum (numpy.ndarray) – An array containing the sum of each independent realization for each bin
  • sum_sq (numpy.ndarray) – An array containing the sum of each independent realization squared for each bin
  • mean (numpy.ndarray) – An array containing the sample mean for each bin
  • std_dev (numpy.ndarray) – An array containing the sample standard deviation for each bin
  • derived (bool) – Whether or not the tally is derived from one or more other tallies
  • sparse (bool) – Whether or not the tally uses SciPy’s LIL sparse matrix format for compressed data storage
  • derivative (openmc.TallyDerivative) – A material perturbation derivative to apply to all scores in the tally.
average(scores=[], filter_type=None, filter_bins=[], nuclides=[], remove_filter=False)[source]

Vectorized average of tally data across scores, filter bins and/or nuclides using tally aggregation.

This method constructs a new tally to encapsulate the average of the data represented by the average of the data in this tally. The tally data average is determined by the scores, filter bins and nuclides specified in the input parameters.

Parameters:
  • scores (list of str) – A list of one or more score strings to average across (e.g., [‘absorption’, ‘nu-fission’]; default is [])
  • filter_type (openmc.FilterMeta) – Type of the filter, e.g. MeshFilter
  • filter_bins (Iterable of int or tuple) – A list of the filter bins corresponding to the filter_type parameter Each bin in the list is the integer ID for ‘material’, ‘surface’, ‘cell’, ‘cellborn’, and ‘universe’ Filters. Each bin is an integer for the cell instance ID for ‘distribcell’ Filters. Each bin is a 2-tuple of floats for ‘energy’ and ‘energyout’ filters corresponding to the energy boundaries of the bin of interest. Each bin is an (x,y,z) 3-tuple for ‘mesh’ filters corresponding to the mesh cell of interest.
  • nuclides (list of str) – A list of nuclide name strings to average across (e.g., [‘U235’, ‘U238’]; default is [])
  • remove_filter (bool) – If a filter is being averaged over, this bool indicates whether to remove that filter in the returned tally. Default is False.
Returns:

A new tally which encapsulates the average of data requested.

Return type:

openmc.Tally

can_merge(other)[source]

Determine if another tally can be merged with this one

If results have been loaded from a statepoint, then tallies are only mergeable along one and only one of filter bins, nuclides or scores.

Parameters:other (openmc.Tally) – Tally to check for merging
contains_filter(filter_type)[source]

Looks for a filter in the tally that matches a specified type

Parameters:filter_type (openmc.FilterMeta) – Type of the filter, e.g. MeshFilter
Returns:filter_found – True if the tally contains a filter of the requested type; otherwise false
Return type:bool
diagonalize_filter(new_filter, filter_position=-1)[source]

Diagonalize the tally data array along a new axis of filter bins.

This is a helper method for the tally arithmetic methods. This method adds the new filter to a derived tally constructed copied from this one. The data in the derived tally arrays is “diagonalized” along the bins in the new filter. This functionality is used by the openmc.mgxs module; to transport-correct scattering matrices by subtracting a ‘scatter-P1’ reaction rate tally with an energy filter from a ‘scatter’ reaction rate tally with both energy and energyout filters.

Parameters:
  • new_filter (Filter) – The filter along which to diagonalize the data in the new
  • filter_position (int) – Where to place the new filter in the Tally.filters list. Defaults to last position.
Returns:

A new derived Tally with data diagaonalized along the new filter.

Return type:

openmc.Tally

find_filter(filter_type)[source]

Return a filter in the tally that matches a specified type

Parameters:filter_type (openmc.FilterMeta) – Type of the filter, e.g. MeshFilter
Returns:filter_found – Filter from this tally with matching type, or None if no matching Filter is found
Return type:openmc.Filter
Raises:ValueError – If no matching Filter is found
get_filter_indices(filters=[], filter_bins=[])[source]

Get indices into the filter axis of this tally’s data arrays.

This is a helper method for the Tally.get_values(…) method to extract tally data. This method returns the indices into the filter axis of the tally’s data array (axis=0) for particular combinations of filters and their corresponding bins.

Parameters:
  • filters (Iterable of openmc.FilterMeta) – An iterable of filter types (e.g., [MeshFilter, EnergyFilter]; default is [])
  • filter_bins (Iterable of tuple) – A list of tuples of filter bins corresponding to the filter_types parameter (e.g., [(1,), ((0., 0.625e-6),)]; default is []). Each tuple contains bins for the corresponding filter type in the filters parameter. Each bin is an integer ID for Material-, Surface-, Cell-, Cellborn-, and Universe- Filters. Each bin is an integer for the cell instance ID for DistribcellFilters. Each bin is a 2-tuple of floats for Energy- and Energyout- Filters corresponding to the energy boundaries of the bin of interest. The bin is an (x,y,z) 3-tuple for MeshFilters corresponding to the mesh cell of interest. The order of the bins in the list must correspond to the filter_types parameter.
Returns:

A NumPy array of the filter indices

Return type:

numpy.ndarray

get_nuclide_index(nuclide)[source]

Returns the index in the Tally’s results array for a Nuclide bin

Parameters:nuclide (str) – The name of the Nuclide (e.g., ‘H1’, ‘U238’)
Returns:nuclide_index – The index in the Tally data array for this nuclide.
Return type:int
Raises:KeyError – When the argument passed to the ‘nuclide’ parameter cannot be found in the Tally.
get_nuclide_indices(nuclides)[source]

Get indices into the nuclide axis of this tally’s data arrays.

This is a helper method for the Tally.get_values(…) method to extract tally data. This method returns the indices into the nuclide axis of the tally’s data array (axis=1) for one or more nuclides.

Parameters:nuclides (list of str) – A list of nuclide name strings (e.g., [‘U235’, ‘U238’]; default is [])
Returns:A NumPy array of the nuclide indices
Return type:numpy.ndarray
get_pandas_dataframe(filters=True, nuclides=True, scores=True, derivative=True, paths=True, float_format='{:.2e}')[source]

Build a Pandas DataFrame for the Tally data.

This method constructs a Pandas DataFrame object for the Tally data with columns annotated by filter, nuclide and score bin information.

This capability has been tested for Pandas >=0.13.1. However, it is recommended to use v0.16 or newer versions of Pandas since this method uses the Multi-index Pandas feature.

Parameters:
  • filters (bool) – Include columns with filter bin information (default is True).
  • nuclides (bool) – Include columns with nuclide bin information (default is True).
  • scores (bool) – Include columns with score bin information (default is True).
  • derivative (bool) – Include columns with differential tally info (default is True).
  • paths (bool, optional) – Construct columns for distribcell tally filters (default is True). The geometric information in the Summary object is embedded into a Multi-index column with a geometric “path” to each distribcell instance.
  • float_format (str) – All floats in the DataFrame will be formatted using the given format string before printing.
Returns:

A Pandas DataFrame with each column annotated by filter, nuclide and score bin information (if these parameters are True), and the mean and standard deviation of the Tally’s data.

Return type:

pandas.DataFrame

Raises:

KeyError – When this method is called before the Tally is populated with data

get_reshaped_data(value='mean')[source]

Returns an array of tally data with one dimension per filter.

The tally data in OpenMC is stored as a 3D array with the dimensions corresponding to filters, nuclides and scores. As a result, tally data can be opaque for a user to directly index (i.e., without use of openmc.Tally.get_values()) since one must know how to properly use the number of bins and strides for each filter to index into the first (filter) dimension.

This builds and returns a reshaped version of the tally data array with unique dimensions corresponding to each tally filter. For example, suppose this tally has arrays of data with shape (8,5,5) corresponding to two filters (2 and 4 bins, respectively), five nuclides and five scores. This method will return a version of the data array with the with a new shape of (2,4,5,5) such that the first two dimensions correspond directly to the two filters with two and four bins.

Parameters:value (str) – A string for the type of value to return - ‘mean’ (default), ‘std_dev’, ‘rel_err’, ‘sum’, or ‘sum_sq’ are accepted
Returns:The tally data array indexed by filters, nuclides and scores.
Return type:numpy.ndarray
get_score_index(score)[source]

Returns the index in the Tally’s results array for a score bin

Parameters:score (str) – The score string (e.g., ‘absorption’, ‘nu-fission’)
Returns:score_index – The index in the Tally data array for this score.
Return type:int
Raises:ValueError – When the argument passed to the ‘score’ parameter cannot be found in the Tally.
get_score_indices(scores)[source]

Get indices into the score axis of this tally’s data arrays.

This is a helper method for the Tally.get_values(…) method to extract tally data. This method returns the indices into the score axis of the tally’s data array (axis=2) for one or more scores.

Parameters:scores (list of str or openmc.CrossScore) – A list of one or more score strings (e.g., [‘absorption’, ‘nu-fission’]; default is [])
Returns:A NumPy array of the score indices
Return type:numpy.ndarray
get_slice(scores=[], filters=[], filter_bins=[], nuclides=[], squeeze=False)[source]

Build a sliced tally for the specified filters, scores and nuclides.

This method constructs a new tally to encapsulate a subset of the data represented by this tally. The subset of data to include in the tally slice is determined by the scores, filters and nuclides specified in the input parameters.

Parameters:
  • scores (list of str) – A list of one or more score strings (e.g., [‘absorption’, ‘nu-fission’]
  • filters (Iterable of openmc.FilterMeta) – An iterable of filter types (e.g., [MeshFilter, EnergyFilter])
  • filter_bins (list of Iterables) – A list of iterables of filter bins corresponding to the specified filter types (e.g., [(1,), ((0., 0.625e-6),)]). Each iterable contains bins to slice for the corresponding filter type in the filters parameter. Each bin is the integer ID for ‘material’, ‘surface’, ‘cell’, ‘cellborn’, and ‘universe’ Filters. Each bin is an integer for the cell instance ID for ‘distribcell’ Filters. Each bin is a 2-tuple of floats for ‘energy’ and ‘energyout’ filters corresponding to the energy boundaries of the bin of interest. The bin is an (x,y,z) 3-tuple for ‘mesh’ filters corresponding to the mesh cell of interest. The order of the bins in the list must correspond to the filters argument.
  • nuclides (list of str) – A list of nuclide name strings (e.g., [‘U235’, ‘U238’])
  • squeeze (bool) – Whether to remove filters with only a single bin in the sliced tally
Returns:

A new tally which encapsulates the subset of data requested in the order each filter, nuclide and score is listed in the parameters.

Return type:

openmc.Tally

Raises:

ValueError – When this method is called before the Tally is populated with data.

get_values(scores=[], filters=[], filter_bins=[], nuclides=[], value='mean')[source]

Returns one or more tallied values given a list of scores, filters, filter bins and nuclides.

This method constructs a 3D NumPy array for the requested Tally data indexed by filter bin, nuclide bin, and score index. The method will order the data in the array as specified in the parameter lists.

Parameters:
  • scores (list of str) – A list of one or more score strings (e.g., [‘absorption’, ‘nu-fission’]; default is [])
  • filters (Iterable of openmc.FilterMeta) – An iterable of filter types (e.g., [MeshFilter, EnergyFilter]; default is [])
  • filter_bins (list of Iterables) – A list of tuples of filter bins corresponding to the filter_types parameter (e.g., [(1,), ((0., 0.625e-6),)]; default is []). Each tuple contains bins for the corresponding filter type in the filters parameter. Each bins is the integer ID for ‘material’, ‘surface’, ‘cell’, ‘cellborn’, and ‘universe’ Filters. Each bin is an integer for the cell instance ID for ‘distribcell’ Filters. Each bin is a 2-tuple of floats for ‘energy’ and ‘energyout’ filters corresponding to the energy boundaries of the bin of interest. The bin is an (x,y,z) 3-tuple for ‘mesh’ filters corresponding to the mesh cell of interest. The order of the bins in the list must correspond to the filter_types parameter.
  • nuclides (list of str) – A list of nuclide name strings (e.g., [‘U235’, ‘U238’]; default is [])
  • value (str) – A string for the type of value to return - ‘mean’ (default), ‘std_dev’, ‘rel_err’, ‘sum’, or ‘sum_sq’ are accepted
Returns:

A scalar or NumPy array of the Tally data indexed in the order each filter, nuclide and score is listed in the parameters.

Return type:

float or numpy.ndarray

Raises:

ValueError – When this method is called before the Tally is populated with data, or the input parameters do not correspond to the Tally’s attributes, e.g., if the score(s) do not match those in the Tally.

hybrid_product(other, binary_op, filter_product=None, nuclide_product=None, score_product=None)[source]

Combines filters, scores and nuclides with another tally.

This is a helper method for the tally arithmetic operator overloaded methods. It is called a “hybrid product” because it performs a combination of tensor (or Kronecker) and entrywise (or Hadamard) products. The filters from both tallies are combined using an entrywise (or Hadamard) product on matching filters. By default, if all nuclides are identical in the two tallies, the entrywise product is performed across nuclides; else the tensor product is performed. By default, if all scores are identical in the two tallies, the entrywise product is performed across scores; else the tensor product is performed. Users can also call the method explicitly and specify the desired product.

Parameters:
  • other (openmc.Tally) – The tally on the right hand side of the hybrid product
  • binary_op ({'+', '-', '*', '/', '^'}) – The binary operation in the hybrid product
  • filter_product ({'tensor', 'entrywise' or None}) – The type of product (tensor or entrywise) to be performed between filter data. The default is the entrywise product. Currently only the entrywise product is supported since a tally cannot contain two of the same filter.
  • nuclide_product ({'tensor', 'entrywise' or None}) – The type of product (tensor or entrywise) to be performed between nuclide data. The default is the entrywise product if all nuclides between the two tallies are the same; otherwise the default is the tensor product.
  • score_product ({'tensor', 'entrywise' or None}) – The type of product (tensor or entrywise) to be performed between score data. The default is the entrywise product if all scores between the two tallies are the same; otherwise the default is the tensor product.
Returns:

A new Tally that is the hybrid product with this one.

Return type:

openmc.Tally

Raises:

ValueError – When this method is called before the other tally is populated with data.

merge(other)[source]

Merge another tally with this one

If results have been loaded from a statepoint, then tallies are only mergeable along one and only one of filter bins, nuclides or scores.

Parameters:other (openmc.Tally) – Tally to merge with this one
Returns:merged_tally – Merged tallies
Return type:openmc.Tally
remove_filter(old_filter)[source]

Remove a filter from the tally

Parameters:old_filter (openmc.Filter) – Filter to remove
remove_nuclide(nuclide)[source]

Remove a nuclide from the tally

Parameters:nuclide (openmc.Nuclide) – Nuclide to remove
remove_score(score)[source]

Remove a score from the tally

Parameters:score (str) – Score to remove
summation(scores=[], filter_type=None, filter_bins=[], nuclides=[], remove_filter=False)[source]

Vectorized sum of tally data across scores, filter bins and/or nuclides using tally aggregation.

This method constructs a new tally to encapsulate the sum of the data represented by the summation of the data in this tally. The tally data sum is determined by the scores, filter bins and nuclides specified in the input parameters.

Parameters:
  • scores (list of str) – A list of one or more score strings to sum across (e.g., [‘absorption’, ‘nu-fission’]; default is [])
  • filter_type (openmc.FilterMeta) – Type of the filter, e.g. MeshFilter
  • filter_bins (Iterable of int or tuple) – A list of the filter bins corresponding to the filter_type parameter Each bin in the list is the integer ID for ‘material’, ‘surface’, ‘cell’, ‘cellborn’, and ‘universe’ Filters. Each bin is an integer for the cell instance ID for ‘distribcell’ Filters. Each bin is a 2-tuple of floats for ‘energy’ and ‘energyout’ filters corresponding to the energy boundaries of the bin of interest. Each bin is an (x,y,z) 3-tuple for ‘mesh’ filters corresponding to the mesh cell of interest.
  • nuclides (list of str) – A list of nuclide name strings to sum across (e.g., [‘U235’, ‘U238’]; default is [])
  • remove_filter (bool) – If a filter is being summed over, this bool indicates whether to remove that filter in the returned tally. Default is False.
Returns:

A new tally which encapsulates the sum of data requested.

Return type:

openmc.Tally

to_xml_element()[source]

Return XML representation of the tally

Returns:element – XML element containing tally data
Return type:xml.etree.ElementTree.Element