openmc.deplete.Results¶
- class openmc.deplete.Results(filename=None)[source]¶
Results from a depletion simulation
The
Results
class acts as a list that stores the results from each depletion step and provides extra methods for interrogating these results.Changed in version 0.13.1: Name changed from
ResultsList
toResults
- Parameters
filename (str) – Path to depletion result file
- export_to_materials(burnup_index: int, nuc_with_data: Optional[Iterable[str]] = None, path: Union[str, PathLike] = 'materials.xml') Materials [source]¶
Return openmc.Materials object based on results at a given step
New in version 0.12.1.
- Parameters
burn_index (int) – Index of burnup step to evaluate. See also: get_step_where for obtaining burnup step indices from other data such as the time.
nuc_with_data (Iterable of str, optional) – Nuclides to include in resulting materials. This can be specified if not all nuclides appearing in depletion results have associated neutron cross sections, and as such cannot be used in subsequent transport calculations. If not provided, nuclides from the cross_sections element of materials.xml will be used. If that element is not present, nuclides from openmc.config[‘cross_sections’] will be used.
path (PathLike) – Path to materials XML file to read. Defaults to ‘materials.xml’.
- Returns
mat_file – A modified Materials instance containing depleted material data and original isotopic compositions of non-depletable materials
- Return type
- classmethod from_hdf5(filename: Union[str, PathLike])[source]¶
Load in depletion results from a previous file
- get_atoms(mat: Union[Material, str], nuc: str, nuc_units: str = 'atoms', time_units: str = 's') Tuple[ndarray, ndarray] [source]¶
Get number of nuclides over time from a single material
Note
Initial values for some isotopes that do not appear in initial concentrations may be non-zero, depending on the value of the
openmc.deplete.CoupledOperator.dilute_initial
attribute. Theopenmc.deplete.CoupledOperator
class adds isotopes according to this setting, which can be set to zero.- Parameters
mat (openmc.Material, str) – Material object or material id to evaluate
nuc (str) – Nuclide name to evaluate
nuc_units ({"atoms", "atom/b-cm", "atom/cm3"}, optional) –
Units for the returned concentration. Default is
"atoms"
New in version 0.12.
time_units ({"s", "min", "h", "d", "a"}, optional) –
Units for the returned time array. Default is
"s"
to return the value in seconds. Other options are minutes"min"
, hours"h"
, days"d"
, and Julian years"a"
.New in version 0.12.
- Returns
times (numpy.ndarray) – Array of times in units of
time_units
concentrations (numpy.ndarray) – Concentration of specified nuclide in units of
nuc_units
- get_depletion_time() ndarray [source]¶
Return an array of the average time to deplete a material
Note
The return value will have one fewer values than several other methods, such as
get_keff()
, because no depletion is performed at the final transport stage.- Returns
times – Vector of average time to deplete a single material across all processes and materials.
- Return type
- get_keff(time_units: str = 's') Tuple[ndarray, ndarray] [source]¶
Evaluates the eigenvalue from a results list.
New in version 0.13.1.
- Parameters
time_units ({"s", "d", "min", "h", "a"}, optional) – Desired units for the times array. Options are seconds
"s"
, minutes"min"
, hours"h"
, days"d"
, and Julian years"a"
.- Returns
times (numpy.ndarray) – Array of times in specified units
eigenvalues (numpy.ndarray) – k-eigenvalue at each time. Column 0 contains the eigenvalue, while column 1 contains the associated uncertainty
- get_reaction_rate(mat: Union[Material, str], nuc: str, rx: str) Tuple[ndarray, ndarray] [source]¶
Get reaction rate in a single material/nuclide over time
Note
Initial values for some isotopes that do not appear in initial concentrations may be non-zero, depending on the value of
openmc.deplete.CoupledOperator
dilute_initial
Theopenmc.deplete.CoupledOperator
adds isotopes according to this setting, which can be set to zero.- Parameters
mat (openmc.Material, str) – Material object or material id to evaluate
nuc (str) – Nuclide name to evaluate
rx (str) – Reaction rate to evaluate
- Returns
times (numpy.ndarray) – Array of times in [s]
rates (numpy.ndarray) – Array of reaction rates
- get_step_where(time, time_units: str = 'd', atol: float = 1e-06, rtol: float = 0.001) int [source]¶
Return the index closest to a given point in time
In the event
time
lies exactly between two points, the lower index will be returned. It is possible that the index will be at most one past the point in time requested, but only according to tolerances requested.Passing
atol=math.inf
andrtol=math.inf
will return the closest index to the requested point.New in version 0.12.1.
- Parameters
time (float) – Desired point in time
time_units ({"s", "d", "min", "h", "a"}, optional) – Units on
time
. Default: days"d"
. Other options are seconds"s"
, minutes"min"
, hours"h"
and Julian years"a"
.atol (float, optional) – Absolute tolerance (in
time_units
) iftime
is not found.rtol (float, optional) – Relative tolerance if
time
is not found.
- Return type
- get_times(time_units: str = 'd') ndarray [source]¶
Return the points in time that define the depletion schedule
New in version 0.12.1.
- Parameters
time_units ({"s", "d", "min", "h", "a"}, optional) – Return the vector in these units. Default is to convert to days
"d"
. Other options are seconds"s"
, minutes"min"
, hours"h"
, days"d"
, and Julian years"a"
.- Returns
1-D vector of time points
- Return type