add continue feature for depletion (#3272)

Co-authored-by: Connor Moreno <camoreno@wisc.edu>
Co-authored-by: Patrick Shriwise <pshriwise@gmail.com>
Co-authored-by: Lewis Gross <ligross@cnerg-docker-04.neep.wisc.edu>
Co-authored-by: Paul Romano <paul.k.romano@gmail.com>
This commit is contained in:
Lewis Gross 2025-03-05 19:07:50 -06:00 committed by GitHub
parent ced8929128
commit 557b714d87
No known key found for this signature in database
GPG key ID: B5690EEEBB952194
3 changed files with 169 additions and 16 deletions

View file

@ -1,6 +1,7 @@
"""abc module.
This module contains Abstract Base Classes for implementing operator, integrator, depletion system solver, and operator helper classes
This module contains Abstract Base Classes for implementing operator,
integrator, depletion system solver, and operator helper classes
"""
from __future__ import annotations
@ -24,7 +25,8 @@ from openmc.utility_funcs import change_directory
from openmc import Material
from .stepresult import StepResult
from .chain import Chain
from .results import Results
from .results import Results, _SECONDS_PER_MINUTE, _SECONDS_PER_HOUR, \
_SECONDS_PER_DAY, _SECONDS_PER_JULIAN_YEAR
from .pool import deplete
from .reaction_rates import ReactionRates
from .transfer_rates import TransferRates
@ -36,12 +38,6 @@ __all__ = [
"Integrator", "SIIntegrator", "DepSystemSolver", "add_params"]
_SECONDS_PER_MINUTE = 60
_SECONDS_PER_HOUR = 60*60
_SECONDS_PER_DAY = 24*60*60
_SECONDS_PER_JULIAN_YEAR = 365.25*24*60*60
def _normalize_timesteps(
timesteps: Sequence[float] | Sequence[tuple[float, str]],
source_rates: float | Sequence[float],
@ -572,6 +568,18 @@ class Integrator(ABC):
:attr:`solver`.
.. versionadded:: 0.12
continue_timesteps : bool, optional
Whether or not to treat the current solve as a continuation of a
previous simulation. Defaults to `False`. When `False`, the depletion
steps provided are appended to any previous steps. If `True`, the
timesteps provided to the `Integrator` must exacly match any that
exist in the `prev_results` passed to the `Operator`. The `power`,
`power_density`, or `source_rates` must match as well. The
method of specifying `power`, `power_density`, or
`source_rates` should be the same as the initial run.
.. versionadded:: 0.15.1
Attributes
----------
operator : openmc.deplete.abc.TransportOperator
@ -613,7 +621,8 @@ class Integrator(ABC):
power_density: Optional[Union[float, Sequence[float]]] = None,
source_rates: Optional[Union[float, Sequence[float]]] = None,
timestep_units: str = 's',
solver: str = "cram48"
solver: str = "cram48",
continue_timesteps: bool = False,
):
# Check number of stages previously used
if operator.prev_res is not None:
@ -625,6 +634,8 @@ class Integrator(ABC):
"this uses {}".format(
self.__class__.__name__, res.data.shape[0],
self._num_stages))
elif continue_timesteps:
raise ValueError("Continuation run requires passing prev_results.")
self.operator = operator
self.chain = operator.chain
@ -642,6 +653,35 @@ class Integrator(ABC):
# Normalize timesteps and source rates
seconds, source_rates = _normalize_timesteps(
timesteps, source_rates, timestep_units, operator)
if continue_timesteps:
# Get timesteps and source rates from previous results
prev_times = operator.prev_res.get_times(timestep_units)
prev_source_rates = operator.prev_res.get_source_rates()
prev_timesteps = np.diff(prev_times)
# Make sure parameters from the previous results are consistent with
# those passed to operator
num_prev = len(prev_timesteps)
if not np.array_equal(prev_timesteps, timesteps[:num_prev]):
raise ValueError(
"You are attempting to continue a run in which the previous timesteps "
"do not have the same initial timesteps as those provided to the "
"Integrator. Please make sure you are using the correct timesteps."
)
if not np.array_equal(prev_source_rates, source_rates[:num_prev]):
raise ValueError(
"You are attempting to continue a run in which the previous results "
"do not have the same initial source rates, powers, or power densities "
"as those provided to the Integrator. Please make sure you are using "
"the correct powers, power densities, or source rates and previous "
"results file."
)
# Run with only the new time steps and source rates provided
seconds = seconds[num_prev:]
source_rates = source_rates[num_prev:]
self.timesteps = np.asarray(seconds)
self.source_rates = np.asarray(source_rates)
@ -932,6 +972,18 @@ class SIIntegrator(Integrator):
:attr:`solver`.
.. versionadded:: 0.12
continue_timesteps : bool, optional
Whether or not to treat the current solve as a continuation of a
previous simulation. Defaults to `False`. If `False`, all time
steps and source rates will be run in an append fashion and will run
after whatever time steps exist, if any. If `True`, the timesteps
provided to the `Integrator` must match exactly those that exist
in the `prev_results` passed to the `Opereator`. The `power`,
`power_density`, or `source_rates` must match as well. The
method of specifying `power`, `power_density`, or
`source_rates` should be the same as the initial run.
.. versionadded:: 0.15.1
Attributes
----------
@ -973,13 +1025,14 @@ class SIIntegrator(Integrator):
source_rates: Optional[Sequence[float]] = None,
timestep_units: str = 's',
n_steps: int = 10,
solver: str = "cram48"
solver: str = "cram48",
continue_timesteps: bool = False,
):
check_type("n_steps", n_steps, Integral)
check_greater_than("n_steps", n_steps, 0)
super().__init__(
operator, timesteps, power, power_density, source_rates,
timestep_units=timestep_units, solver=solver)
timestep_units=timestep_units, solver=solver, continue_timesteps=continue_timesteps)
self.n_steps = n_steps
def _get_bos_data_from_operator(self, step_index, step_power, n_bos):

View file

@ -17,6 +17,11 @@ from openmc.checkvalue import PathLike
__all__ = ["Results", "ResultsList"]
_SECONDS_PER_MINUTE = 60
_SECONDS_PER_HOUR = 60*60
_SECONDS_PER_DAY = 24*60*60
_SECONDS_PER_JULIAN_YEAR = 365.25*24*60*60 # 365.25 due to the leap year
def _get_time_as(seconds: float, units: str) -> float:
"""Converts the time in seconds to time in different units
@ -31,13 +36,13 @@ def _get_time_as(seconds: float, units: str) -> float:
"""
if units == "a":
return seconds / (60 * 60 * 24 * 365.25) # 365.25 due to the leap year
return seconds / _SECONDS_PER_JULIAN_YEAR
if units == "d":
return seconds / (60 * 60 * 24)
return seconds / _SECONDS_PER_DAY
elif units == "h":
return seconds / (60 * 60)
return seconds / _SECONDS_PER_HOUR
elif units == "min":
return seconds / 60
return seconds / _SECONDS_PER_MINUTE
else:
return seconds
@ -71,7 +76,6 @@ class Results(list):
data.append(StepResult.from_hdf5(fh, i))
super().__init__(data)
@classmethod
def from_hdf5(cls, filename: PathLike):
"""Load in depletion results from a previous file
@ -460,6 +464,26 @@ class Results(list):
return _get_time_as(times, time_units)
def get_source_rates(self) -> np.ndarray:
"""
.. versionadded:: 0.15.1
Returns
-------
numpy.ndarray
1-D vector of source rates at each point in the depletion simulation
with the units originally defined by the user.
"""
# Results duplicate the final source rate at the final simulation time
source_rates = np.fromiter(
(r.source_rate for r in self),
dtype=self[0].source_rate.dtype,
count=len(self)-1,
)
return source_rates
def get_step_where(
self, time, time_units: str = "d", atol: float = 1e-6, rtol: float = 1e-3
) -> int:

View file

@ -0,0 +1,76 @@
"""Unit tests for openmc.deplete continue run capability.
These tests run in two steps: first a normal run and then a continue run using the previous results
"""
import pytest
import openmc.deplete
from tests import dummy_operator
def test_continue(run_in_tmpdir):
"""Test to ensure that a properly defined continue run works"""
# set up the problem
bundle = dummy_operator.SCHEMES['predictor']
operator = dummy_operator.DummyOperator()
# initial depletion
bundle.solver(operator, [1.0, 2.0], [1.0, 2.0]).integrate()
# set up continue run
prev_res = openmc.deplete.Results(operator.output_dir / "depletion_results.h5")
operator = dummy_operator.DummyOperator(prev_res)
# if continue run happens, test passes
bundle.solver(operator, [1.0, 2.0, 3.0, 4.0], [1.0, 2.0, 3.0, 4.0],
continue_timesteps=True).integrate()
def test_mismatched_initial_times(run_in_tmpdir):
"""Test to ensure that a continue run with different initial steps is properly caught"""
# set up the problem
bundle = dummy_operator.SCHEMES['predictor']
operator = dummy_operator.DummyOperator()
# perform initial steps
bundle.solver(operator, [0.75, 0.75], [1.0, 1.0]).integrate()
# restart
prev_res = openmc.deplete.Results(operator.output_dir / "depletion_results.h5")
operator = dummy_operator.DummyOperator(prev_res)
with pytest.raises(
ValueError,
match="You are attempting to continue a run in which the previous timesteps "
"do not have the same initial timesteps as those provided to the "
"Integrator. Please make sure you are using the correct timesteps.",
):
bundle.solver(
operator, [0.75, 0.5, 0.75], [1.0, 1.0, 1.0], continue_timesteps=True
).integrate()
def test_mismatched_initial_source_rates(run_in_tmpdir):
"""Test to ensure that a continue run with different initial steps is properly caught"""
# set up the problem
bundle = dummy_operator.SCHEMES['predictor']
operator = dummy_operator.DummyOperator()
# perform initial steps
bundle.solver(operator, [0.75, 0.75], [1.0, 1.0]).integrate()
# restart
prev_res = openmc.deplete.Results(operator.output_dir / "depletion_results.h5")
operator = dummy_operator.DummyOperator(prev_res)
with pytest.raises(
ValueError,
match="You are attempting to continue a run in which the previous results "
"do not have the same initial source rates, powers, or power densities "
"as those provided to the Integrator. Please make sure you are using "
"the correct powers, power densities, or source rates and previous results file.",
):
bundle.solver(
operator, [0.75, 0.75, 0.75], [1.0, 2.0, 1.0], continue_timesteps=True
).integrate()