From a04f1e22f4d04f1d923c495c1dd5c987f0577e71 Mon Sep 17 00:00:00 2001
From: Marc Berliner <34451391+MarcBerliner@users.noreply.github.com>
Date: Thu, 3 Sep 2026 10:37:06 -0400
Subject: [PATCH 1/6] feat: input parameters for experiment step durations and
 terminations

Two things a step can now take as symbols, resolved at solve time rather
than at construction:

- `duration=pybamm.InputParameter("...")` (or any expression), evaluated
  against the inputs passed to `solve` to give the step's final time.
  Drive cycles support this too, because the interpolant no longer
  depends on the duration: instead of tiling the cycle out to a fixed
  duration, it wraps the step time with a modulo, which also means one
  built model now serves any duration.
- `termination=pybamm.CoupledVariable("Voltage [V]") > pybamm.InputParameter("V hold")`,
  or any inequality over model variables, input parameters and numbers.
  A heaviside is "left < right", so `left - right` is already in the
  event convention; the CoupledVariables in it are resolved against the
  model variables and the whole thing becomes a CustomTermination.

A symbolic duration round-trips through `to_config`/`from_config` via the
serialisation kernel. DiffSL export bakes the schedule into generated
code, so it raises for a symbolic duration.
---
 CHANGELOG.md                                  |   2 +
 .../pybamm/discretisations/discretisation.py  |  18 +--
 .../src/pybamm/experiment/step/base_step.py   | 135 +++++++++++-------
 .../experiment/step/step_termination.py       |  55 +++++--
 .../expression_tree/operations/diffsl.py      |   5 +
 .../expression_tree/operations/serialise.py   |  25 +++-
 .../pybamm/src/pybamm/models/base_model.py    |  22 +--
 .../src/pybamm/models/symbol_processor.py     |  25 +++-
 .../src/pybamm/simulation/simulation.py       |  24 ++--
 .../unit/test_experiments/test_base_step.py   |  97 ++++++++++++-
 .../test_experiment_step_termination.py       |  65 +++++++++
 .../test_simulation_with_experiment.py        |  63 ++++++++
 .../unit/test_models/test_symbol_processor.py |  19 +++
 .../test_serialisation/test_serialisation.py  |  30 ++++
 14 files changed, 467 insertions(+), 118 deletions(-)

diff --git a/CHANGELOG.md b/CHANGELOG.md
index 7a158f7d3f..7a1fdcefe8 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -8,6 +8,8 @@
 
 ## Features
 
+- Experiment steps accept a symbolic `duration` (e.g. a `pybamm.InputParameter`), resolved against the inputs passed to `Simulation.solve`, so one built model serves any duration. Drive cycles now repeat by wrapping the step time, rather than being tiled out to a fixed duration.
+- Experiment steps accept a symbolic inequality as a `termination`, e.g. `pybamm.CoupledVariable("Voltage [V]") > pybamm.InputParameter("Voltage hold [V]")`, so a cut-off can be set on any model variable and its threshold supplied at solve time.
 - Added unstructured mesh support (`UnstructuredSubMesh`, generators, and interface coupling) for arbitrary 2D/3D domains. Hexahedra must have planar faces (warped hexes raise a `GeometryError`), and `UserSuppliedUnstructuredMesh` accepts tetrahedral, triangular, and quadrilateral cells only. ([#5687](https://github.com/pybamm-team/PyBaMM/pull/5687))
 - Generalised `VectorField` to N components and added `Component`/`Norm` operators for multi-dimensional vector fields. ([#5686](https://github.com/pybamm-team/PyBaMM/pull/5686))
 - Removed the left sidebar from the documentation home page for a cleaner landing experience. ([#5699](https://github.com/pybamm-team/PyBaMM/pull/5699))
diff --git a/packages/pybamm/src/pybamm/discretisations/discretisation.py b/packages/pybamm/src/pybamm/discretisations/discretisation.py
index cdb404505d..4e3e525f7e 100644
--- a/packages/pybamm/src/pybamm/discretisations/discretisation.py
+++ b/packages/pybamm/src/pybamm/discretisations/discretisation.py
@@ -300,23 +300,7 @@ def _resolve_coupled_variables_in_model(self, model):
         """Resolve CoupledVariables in rhs, algebraic, initial_conditions, and boundary_conditions."""
 
         def resolve_symbol(symbol):
-            if isinstance(symbol, pybamm.CoupledVariable):
-                if symbol.name not in model.variables:
-                    raise pybamm.DiscretisationError(
-                        f"CoupledVariable '{symbol.name}' not found in model.variables."
-                    )
-                return resolve_symbol(model.variables[symbol.name])
-            elif hasattr(symbol, "children") and symbol.children:
-                new_children = []
-                changed = False
-                for child in symbol.children:
-                    new_child = resolve_symbol(child)
-                    new_children.append(new_child)
-                    if new_child is not child:
-                        changed = True
-                if changed:
-                    return symbol.create_copy(new_children=new_children)
-            return symbol
+            return pybamm.SymbolProcessor.resolve(symbol, model.variables)
 
         for var, expr in list(model.rhs.items()):
             resolved = resolve_symbol(expr)
diff --git a/packages/pybamm/src/pybamm/experiment/step/base_step.py b/packages/pybamm/src/pybamm/experiment/step/base_step.py
index 0cf88b05a9..d674b5e55c 100644
--- a/packages/pybamm/src/pybamm/experiment/step/base_step.py
+++ b/packages/pybamm/src/pybamm/experiment/step/base_step.py
@@ -51,11 +51,16 @@ class BaseStep:
     value : float
         The value of the step, corresponding to the type of step. Can be a number, a
         2-tuple (for cccv_ode), a 2-column array (for drive cycles), or a 1-argument function of t
-    duration : float, optional
-        The duration of the step in seconds.
-    termination : str or list, optional
-        A string or list of strings indicating the condition(s) that will terminate the
-        step. If a list, the step will terminate when any of the conditions are met.
+    duration : float or str or :class:`pybamm.Symbol`, optional
+        The duration of the step in seconds. A symbolic duration (e.g. a
+        :class:`pybamm.InputParameter`) is evaluated at solve time against the inputs
+        passed to :meth:`pybamm.Simulation.solve`.
+    termination : str or list or :class:`pybamm.Symbol`, optional
+        A condition, or list of conditions, that will terminate the step; the step ends
+        when any of them is met. Each condition is a string (e.g. ``"4.2V"``), a
+        :class:`pybamm.step.BaseTermination`, or a symbolic inequality such as
+        ``pybamm.CoupledVariable("Voltage [V]") > pybamm.InputParameter("V hold")``,
+        whose threshold may be an input parameter resolved at solve time.
     period : float or string, optional
         The period of the step. If a float, the value is in seconds. If a string, the
         value should be a valid time string, e.g. "1 hour".
@@ -135,28 +140,19 @@ def __init__(
             duration = self.default_duration(value)
         self.duration = _convert_time_to_seconds(duration)
 
-        # If drive cycle, repeat the drive cycle until the end of the experiment,
-        # and create an interpolant
+        # If drive cycle, repeat the drive cycle until the end of the experiment by
+        # wrapping the step time, so that the duration (which may be symbolic) does
+        # not enter the interpolant, and create an interpolant
         if self.is_drive_cycle:
-            t_max = self.duration
-            if t_max > value[-1, 0]:
-                # duration longer than drive cycle values so loop
-                nloop = np.ceil(t_max / value[-1, 0]).astype(int)
-                tstep = np.diff(value[:, 0])[0]
-                t = []
-                y = []
-                for i in range(nloop):
-                    t.append(value[:, 0] + ((value[-1, 0] + tstep) * i))
-                    y.append(value[:, 1])
-                t = np.asarray(t).flatten()
-                y = np.asarray(y).flatten()
-            else:
-                t, y = value[:, 0], value[:, 1]
-
+            t, y = value[:, 0], value[:, 1]
+            # Each repeat starts one sample after the last one ended, so the cycle's
+            # period is one sample longer than its last time point
+            cycle_period = t[-1] + np.diff(t)[0]
+            step_time = pybamm.t - pybamm.InputParameter("start time")
             self.value = pybamm.Interpolant(
-                t,
-                y,
-                pybamm.t - pybamm.InputParameter("start time"),
+                np.append(t, cycle_period),
+                np.append(y, y[0]),
+                step_time % cycle_period,
                 name="Drive Cycle",
             )
 
@@ -185,15 +181,13 @@ def __init__(
             termination = []
         elif not isinstance(termination, list):
             termination = [termination]
-        self.termination = []
-        for term in termination:
-            term_obj = None
+
+        def _build_termination(term):
             if isinstance(term, str):
-                operator, typ, val = _parse_termination(term, self.value)
-                term_obj = _read_termination((operator, typ, val))
-            else:
-                term_obj = _read_termination(term)
-            self.termination.append(term_obj)
+                term = _parse_termination(term, self.value)
+            return _read_termination(term)
+
+        self.termination = [_build_termination(term) for term in termination]
 
         if (
             hasattr(self, "calculate_charge_or_discharge")
@@ -217,6 +211,18 @@ def __init__(
         self.next_start_time = None
         self.end_time = None
 
+    def duration_seconds(self, parameter_values, inputs=None):
+        """The duration in seconds, resolving a symbolic duration at solve time."""
+        duration = parameter_values.process_symbol(
+            pybamm.convert_to_symbol(self.duration)
+        )
+        return float(duration.evaluate(inputs=inputs))
+
+    def period_seconds(self, parameter_values, inputs=None):
+        """The period in seconds, resolving a symbolic period at solve time."""
+        period = parameter_values.process_symbol(pybamm.convert_to_symbol(self.period))
+        return float(period.evaluate(inputs=inputs))
+
     @staticmethod
     def is_implicit() -> bool:
         return False
@@ -352,9 +358,9 @@ def default_duration(self, value):
     def default_period():
         return 60.0  # seconds
 
-    def default_time_vector(self, solver, tf, t0=0):
+    def default_time_vector(self, solver, tf, t0=0, parameter_values=None, inputs=None):
         if self.period is not None:
-            period = self.period
+            period = self.period_seconds(parameter_values, inputs)
         elif self.is_drive_cycle and solver.supports_interp:
             # Infer the period from the drive cycle
             period = np.diff(self.value.x[0]).min()
@@ -364,7 +370,9 @@ def default_time_vector(self, solver, tf, t0=0):
 
         return np.linspace(t0, tf, npts)
 
-    def setup_timestepping(self, solver, tf, t_interp=None):
+    def setup_timestepping(
+        self, solver, tf, t_interp=None, parameter_values=None, inputs=None
+    ):
         """
         Setup timestepping for the model.
 
@@ -378,11 +386,17 @@ def setup_timestepping(self, solver, tf, t_interp=None):
             The time points at which to interpolate the solution
         """
         if solver.supports_interp:
-            return self._setup_timestepping(solver, tf, t_interp)
+            return self._setup_timestepping(
+                solver, tf, t_interp, parameter_values, inputs
+            )
         else:
-            return self._setup_timestepping_dense_t_eval(solver, tf, t_interp)
+            return self._setup_timestepping_dense_t_eval(
+                solver, tf, t_interp, parameter_values, inputs
+            )
 
-    def _setup_timestepping(self, solver, tf, t_interp):
+    def _setup_timestepping(
+        self, solver, tf, t_interp, parameter_values=None, inputs=None
+    ):
         """
         Setup timestepping for the model. This returns a t_eval vector that stops
         only at the first and last time points. If t_interp and the period are
@@ -400,13 +414,7 @@ def _setup_timestepping(self, solver, tf, t_interp):
             The time points at which to interpolate the solution
         """
         if self.is_drive_cycle:
-            t_eval = self.value.x[0]
-            # If the drive cycle is longer than the final time,
-            # then truncate the drive cycle
-            if t_eval[-1] > tf:
-                t_eval = t_eval[t_eval <= tf]
-            if t_eval[-1] != tf:
-                t_eval = np.append(t_eval, tf)
+            t_eval = self._drive_cycle_time_points(tf)
         else:
             t_eval = np.array([0, tf])
 
@@ -415,13 +423,28 @@ def _setup_timestepping(self, solver, tf, t_interp):
 
         if t_interp is None:
             if self.period is not None:
-                t_interp = self.default_time_vector(solver, tf)
+                t_interp = self.default_time_vector(
+                    solver, tf, parameter_values=parameter_values, inputs=inputs
+                )
             else:
                 t_interp = solver.process_t_interp(t_interp)
 
         return t_eval, t_interp
 
-    def _setup_timestepping_dense_t_eval(self, solver, tf, t_interp):
+    def _drive_cycle_time_points(self, tf):
+        """The drive cycle's sample times, repeated to cover ``[0, tf]``."""
+        # The last point of ``x`` is the wrap point, which repeats the first sample
+        samples, period = self.value.x[0][:-1], self.value.x[0][-1]
+        nloop = max(int(np.ceil(tf / period)), 1)
+        t_eval = np.concatenate([samples + period * i for i in range(nloop)])
+        t_eval = t_eval[t_eval <= tf]
+        if t_eval[-1] != tf:
+            t_eval = np.append(t_eval, tf)
+        return t_eval
+
+    def _setup_timestepping_dense_t_eval(
+        self, solver, tf, t_interp, parameter_values=None, inputs=None
+    ):
         """
         Setup timestepping for the model. By default, this returns a dense t_eval which
         stops the solver at each point in the t_eval vector. This method is for solvers
@@ -436,7 +459,9 @@ def _setup_timestepping_dense_t_eval(self, solver, tf, t_interp):
         t_interp: np.array | None
             The time points at which to interpolate the solution
         """
-        t_eval = self.default_time_vector(solver, tf)
+        t_eval = self.default_time_vector(
+            solver, tf, parameter_values=parameter_values, inputs=inputs
+        )
 
         t_interp = solver.process_t_interp(t_interp)
 
@@ -543,6 +568,14 @@ def record_tags(
         direction,
     ):
         """Record all the args for repr and hash"""
+        # A Symbol has no truth value and no stable repr, but its id identifies the
+        # expression it stands for
+        if isinstance(duration, pybamm.Symbol):
+            duration = str(duration)
+        if isinstance(termination, pybamm.Symbol):
+            termination = str(termination)
+        if isinstance(period, pybamm.Symbol):
+            period = str(period)
         repr_args = f"{value}, duration={duration}"
         hash_args = f"{value}"
         if termination:
@@ -664,6 +697,10 @@ def _convert_time_to_seconds(time_and_units):
     if time_and_units is None:
         return time_and_units
 
+    # A symbolic time is resolved at solve time, by `_evaluate_time`
+    if isinstance(time_and_units, pybamm.Symbol):
+        return time_and_units
+
     # If the time is a number, assume it is in seconds
     if isinstance(time_and_units, numbers.Number):
         if time_and_units <= 0:
diff --git a/packages/pybamm/src/pybamm/experiment/step/step_termination.py b/packages/pybamm/src/pybamm/experiment/step/step_termination.py
index 4f1d84a5b2..1077f3010c 100644
--- a/packages/pybamm/src/pybamm/experiment/step/step_termination.py
+++ b/packages/pybamm/src/pybamm/experiment/step/step_termination.py
@@ -1,3 +1,4 @@
+from dataclasses import dataclass
 from warnings import warn
 
 import pybamm
@@ -220,15 +221,49 @@ def __hash__(self) -> int:
         return hash((type(self).__name__, self.name, id(self.event_function)))
 
 
-def _read_termination(termination, operator=None):
-    if isinstance(termination, tuple):
-        op, typ, value = termination
-    else:
+_TERMINATION_MAP: dict[str, type[BaseTermination]] = {
+    "current": CurrentTermination,
+    "voltage": VoltageTermination,
+    "C-rate": CRateTermination,
+}
+
+
+def _parse_termination_class(termination_tuple) -> BaseTermination:
+    if len(termination_tuple) != 3:
+        raise ValueError(
+            f"Termination tuple must be of the form (operator, type, value), "
+            f"but got {termination_tuple}"
+        )
+    op, typ, value = termination_tuple
+    cls = _TERMINATION_MAP[typ]
+    return cls(value, operator=op)
+
+
+@dataclass(frozen=True, slots=True)
+class _HeavsideResidual:
+    symbol: pybamm.Symbol
+
+    def __call__(self, variables: dict[str, pybamm.Symbol]) -> pybamm.Symbol:
+        # A heaviside is "left < right", so left - right is the event residual
+        residual = self.symbol.left - self.symbol.right
+        return pybamm.SymbolProcessor.resolve(residual, variables)
+
+
+def _parse_termination_symbol(termination: pybamm.Symbol) -> CustomTermination:
+    if not isinstance(termination, (pybamm.EqualHeaviside, pybamm.NotEqualHeaviside)):
+        raise TypeError(
+            "A symbolic termination must be an inequality between symbols, e.g. "
+            'pybamm.CoupledVariable("Voltage [V]") > 3.0, but got '
+            f"'{termination!s}'. Note that inequalities are only exact when "
+            'pybamm.settings.heaviside_smoothing is "exact".'
+        )
+    return CustomTermination(str(termination), _HeavsideResidual(termination))
+
+
+def _read_termination(termination):
+    if isinstance(termination, pybamm.Symbol):
+        return _parse_termination_symbol(termination)
+    if not isinstance(termination, tuple):
         return termination
 
-    termination_class = {
-        "current": CurrentTermination,
-        "voltage": VoltageTermination,
-        "C-rate": CRateTermination,
-    }[typ]
-    return termination_class(value, operator=op)
+    return _parse_termination_class(termination)
diff --git a/packages/pybamm/src/pybamm/expression_tree/operations/diffsl.py b/packages/pybamm/src/pybamm/expression_tree/operations/diffsl.py
index debc6e6802..a5d7a6e50e 100644
--- a/packages/pybamm/src/pybamm/expression_tree/operations/diffsl.py
+++ b/packages/pybamm/src/pybamm/expression_tree/operations/diffsl.py
@@ -407,6 +407,11 @@ def _normalise_schedule_value(value):
     @staticmethod
     def _effective_step_duration(step: pybamm.step.BaseStep, initial_start_time):
         effective_duration = step.duration
+        if isinstance(effective_duration, pybamm.Symbol):
+            raise NotImplementedError(
+                "Steps with a symbolic duration cannot be exported to DiffSL: the "
+                "schedule is baked into the generated code, so it needs a number."
+            )
         if step.end_time is not None and initial_start_time is not None:
             start_dt = (step.start_time - initial_start_time).total_seconds()
             end_dt = (step.end_time - initial_start_time).total_seconds()
diff --git a/packages/pybamm/src/pybamm/expression_tree/operations/serialise.py b/packages/pybamm/src/pybamm/expression_tree/operations/serialise.py
index 5c25087b1f..a7e7f930d2 100644
--- a/packages/pybamm/src/pybamm/expression_tree/operations/serialise.py
+++ b/packages/pybamm/src/pybamm/expression_tree/operations/serialise.py
@@ -1789,6 +1789,8 @@ def serialise_experiment(experiment) -> dict:
             ``CustomStepExplicit``, ``CustomStepImplicit``); these have no
             JSON representation.
         """
+        from pybamm.expression_tree.operations.serialise_kernel import encode
+
         step_type_map = {
             "Current": "current",
             "Rest": "rest",
@@ -1831,7 +1833,10 @@ def _serialise_step(step):
 
             step_config: dict = {"type": step_type}
             # Use ``input_duration`` so ``uses_default_duration`` round-trips.
-            if step.input_duration is not None:
+            if isinstance(step.input_duration, pybamm.Symbol):
+                # A symbolic duration is an expression tree, not a time string
+                step_config["duration"] = encode(step.input_duration)
+            elif step.input_duration is not None:
                 step_config["duration"] = step.input_duration
 
             if step_type != "rest":
@@ -1882,7 +1887,9 @@ def _serialise_step(step):
             for field, default in field_defaults.items():
                 value = getattr(step, field, None)
                 if value is not None and value != default:
-                    step_config[field] = value
+                    step_config[field] = (
+                        encode(value) if isinstance(value, pybamm.Symbol) else value
+                    )
             tags = getattr(step, "tags", None)
             if tags:
                 step_config["tags"] = tags
@@ -1938,6 +1945,8 @@ def deserialise_experiment(data: dict):
         -------
         :class:`pybamm.Experiment`
         """
+        from pybamm.expression_tree.operations.serialise_kernel import decode
+
         step_func_map = _experiment_step_factories()
         term_class_map = {
             "voltage": pybamm.step.VoltageTermination,
@@ -1986,7 +1995,10 @@ def _parse_step(step_dict):
             # ``uses_default_duration`` (used by infeasibility handling).
             duration_kwargs = {}
             if "duration" in step_dict and step_dict["duration"] is not None:
-                duration_kwargs["duration"] = step_dict["duration"]
+                duration = step_dict["duration"]
+                if isinstance(duration, dict):
+                    duration = decode(duration)
+                duration_kwargs["duration"] = duration
             terminations = None
             if step_dict.get("terminations"):
                 terminations = [
@@ -2002,7 +2014,12 @@ def _parse_step(step_dict):
                 "direction",
             ):
                 if step_dict.get(field) is not None:
-                    extra_kwargs[field] = step_dict[field]
+                    field_value = step_dict[field]
+                    extra_kwargs[field] = (
+                        decode(field_value)
+                        if isinstance(field_value, dict)
+                        else field_value
+                    )
             if step_dict.get("start_time") is not None:
                 extra_kwargs["start_time"] = datetime.fromisoformat(
                     step_dict["start_time"]
diff --git a/packages/pybamm/src/pybamm/models/base_model.py b/packages/pybamm/src/pybamm/models/base_model.py
index 21ff2f0f75..26184fc1a3 100644
--- a/packages/pybamm/src/pybamm/models/base_model.py
+++ b/packages/pybamm/src/pybamm/models/base_model.py
@@ -419,30 +419,10 @@ def process_and_register_variable(self, name: str, symbol: pybamm.Symbol):
                 f"Cannot process variable '{name}' without a `symbol_processor`."
             )
 
-        symbol = self._resolve_coupled_variables(symbol)
+        symbol = SymbolProcessor.resolve(symbol, self.variables)
         value = self.symbol_processor(name=name, symbol=symbol)
         self._variables_processed[name] = value
 
-    def _resolve_coupled_variables(self, symbol: pybamm.Symbol) -> pybamm.Symbol:
-        """Resolve CoupledVariables by looking up their targets in self._variables."""
-        if isinstance(symbol, pybamm.CoupledVariable):
-            if symbol.name not in self._variables:
-                raise ValueError(
-                    f"CoupledVariable '{symbol.name}' not found in model.variables"
-                )
-            return self._resolve_coupled_variables(self._variables[symbol.name])
-        elif hasattr(symbol, "children") and symbol.children:
-            new_children = []
-            changed = False
-            for child in symbol.children:
-                new_child = self._resolve_coupled_variables(child)
-                new_children.append(new_child)
-                if new_child is not child:
-                    changed = True
-            if changed:
-                return symbol.create_copy(new_children=new_children)
-        return symbol
-
     def update_processed_variables(self, processed_vars: dict[str, pybamm.Symbol]):
         """
         Update the _variables_processed dict with new processed variables.
diff --git a/packages/pybamm/src/pybamm/models/symbol_processor.py b/packages/pybamm/src/pybamm/models/symbol_processor.py
index 27ec085cce..b9354f55a7 100644
--- a/packages/pybamm/src/pybamm/models/symbol_processor.py
+++ b/packages/pybamm/src/pybamm/models/symbol_processor.py
@@ -15,7 +15,9 @@ class SymbolProcessor:
     This class provides a convenient way to process symbols using both
     :class:`pybamm.ParameterValues` and :class:`pybamm.Discretisation` objects.
     Once both are set, calling the processor on a symbol will first substitute
-    parameters, then discretise the result.
+    parameters, then discretise the result. :meth:`resolve` turns a symbol written
+    against variable names (:class:`pybamm.CoupledVariable`) into one over the model's
+    own symbols, ready for either.
 
     Attributes
     ----------
@@ -40,6 +42,27 @@ def __init__(self):
         self._discretisation = None
         self._parameter_values = None
 
+    @staticmethod
+    def resolve(symbol: pybamm.Symbol, variables: dict) -> pybamm.Symbol:
+        """Replace each CoupledVariable in ``symbol`` with the entry of ``variables``
+        of that name, so that a symbol written against variable names can be
+        processed like any other model equation."""
+        resolve = pybamm.SymbolProcessor.resolve
+        if isinstance(symbol, pybamm.CoupledVariable):
+            if symbol.name not in variables:
+                raise ValueError(
+                    f"CoupledVariable '{symbol.name}' not found in model.variables"
+                )
+            return resolve(variables[symbol.name], variables)
+        if symbol.children:
+            children = [resolve(child, variables) for child in symbol.children]
+            if any(
+                new is not old
+                for new, old in zip(children, symbol.children, strict=True)
+            ):
+                return symbol.create_copy(new_children=children)
+        return symbol
+
     def __call__(self, name: str, symbol: pybamm.Symbol) -> pybamm.Symbol:
         """
         Process a symbol by applying parameter values and discretisation.
diff --git a/packages/pybamm/src/pybamm/simulation/simulation.py b/packages/pybamm/src/pybamm/simulation/simulation.py
index de112cb39b..e1d6857bc7 100644
--- a/packages/pybamm/src/pybamm/simulation/simulation.py
+++ b/packages/pybamm/src/pybamm/simulation/simulation.py
@@ -415,19 +415,22 @@ def _build_experiment_step_inputs(
                 )
         return inputs
 
+    @staticmethod
+    def _get_experiment_step_index(step_or_key):
+        if isinstance(step_or_key, str):
+            return step_or_key
+        return step_or_key.basic_repr()
+
     def _get_built_experiment_model(self, step_or_key):
         if self._experiment_uses_unified_model:
             return self._built_experiment_model
-        if isinstance(step_or_key, str):
-            return self.steps_to_built_models[step_or_key]
-        return self.steps_to_built_models[step_or_key.basic_repr()]
+        return self.steps_to_built_models[self._get_experiment_step_index(step_or_key)]
 
     def _get_built_experiment_solver(self, step_or_key):
         if self._experiment_uses_unified_model:
             return self._built_experiment_solver
-        if isinstance(step_or_key, str):
-            return self.steps_to_built_solvers[step_or_key]
-        return self.steps_to_built_solvers[step_or_key.basic_repr()]
+        key = self._get_experiment_step_index(step_or_key)
+        return self.steps_to_built_solvers[key]
 
     def _evaluate_step_termination_expression_from_solution(
         self, term, step_solution, step
@@ -478,7 +481,8 @@ def _decode_combined_step_termination(self, step_solution, step, model, inputs):
             )
 
             try:
-                value = event.expression.evaluate(t=t, y=y, inputs=inputs)
+                processed = model.process_symbol(event.expression)
+                value = processed.evaluate(t=t, y=y, inputs=inputs)
             except NotImplementedError:  # pragma: no cover
                 # If the raw expression still contains unevaluated symbols, fall back to
                 # the processed variables on the solved step. This is slower, but it works
@@ -1008,7 +1012,7 @@ def solve(
                 step = experiment_steps[idx]
                 start_time = current_solution.t[-1]
 
-                dt = step.duration
+                dt = step.duration_seconds(self._parameter_values, user_inputs)
                 if step.end_time is not None:
                     remaining = (
                         step.end_time
@@ -1029,7 +1033,7 @@ def solve(
 
                 logs["step number"] = (step_num, cycle_length)
                 logs["step operating conditions"] = step_str
-                logs["step duration"] = step.duration
+                logs["step duration"] = dt
                 callbacks.on_step_start(logs)
 
                 active_step_index = step_indices[idx] if uses_unified else None
@@ -1042,7 +1046,7 @@ def solve(
                 )
 
                 t_eval, t_interp_processed = step.setup_timestepping(
-                    solver, dt, t_interp
+                    solver, dt, t_interp, self._parameter_values, user_inputs
                 )
 
                 state_mapper = self._get_state_mapper_for_solution(
diff --git a/packages/pybamm/tests/unit/test_experiments/test_base_step.py b/packages/pybamm/tests/unit/test_experiments/test_base_step.py
index dec3cc7b0a..ba3bd156e2 100644
--- a/packages/pybamm/tests/unit/test_experiments/test_base_step.py
+++ b/packages/pybamm/tests/unit/test_experiments/test_base_step.py
@@ -67,12 +67,18 @@ def test_drive_cycle_validation_and_looping():
 
     step = pybamm.step.current(np.array([[0.0, 0.0], [2.0, 1.0]]), duration=5)
 
-    np.testing.assert_array_equal(
-        step.value.x[0], np.array([0.0, 2.0, 4.0, 6.0, 8.0, 10.0])
-    )
-    np.testing.assert_array_equal(
-        step.value.y, np.array([0.0, 1.0, 0.0, 1.0, 0.0, 1.0])
-    )
+    # One cycle, plus the wrap point that repeats the first sample: the step time is
+    # wrapped into this window rather than the cycle being tiled out to the duration
+    np.testing.assert_array_equal(step.value.x[0], np.array([0.0, 2.0, 4.0]))
+    np.testing.assert_array_equal(step.value.y, np.array([0.0, 1.0, 0.0]))
+
+    # The drive cycle repeats every 4 s, for as long as the step runs
+    def current_at(t):
+        return step.value.evaluate(t=t, inputs={"start time": 0.0})
+
+    for t in (0.0, 1.0, 2.0, 3.0):
+        assert current_at(t) == current_at(t + 4.0) == current_at(t + 8.0)
+    assert current_at(1.0) == 0.5
 
 
 def test_default_time_vector_uses_drive_cycle_period_only_when_supported():
@@ -133,3 +139,82 @@ def test_parse_termination_requires_operator_for_input_parameter():
         pybamm.experiment.step.base_step._parse_termination(
             "2A", pybamm.InputParameter("I_app")
         )
+
+
+@pytest.mark.parametrize(
+    "duration, inputs",
+    [
+        # Resolved from the solver inputs directly...
+        (pybamm.InputParameter("d"), {"d": 1200}),
+        # ...or through a Parameter whose value in parameter_values is "[input]"
+        (pybamm.Parameter("d"), {"d": 1200}),
+        # ...or through a Parameter that simply has a value
+        (pybamm.Parameter("fixed d"), {}),
+        # A number needs neither
+        (1200, {}),
+    ],
+)
+def test_duration_is_resolved_at_solve_time(duration, inputs):
+    parameter_values = pybamm.ParameterValues({"d": "[input]", "fixed d": 1200})
+    step = pybamm.step.c_rate(1, duration=duration)
+
+    assert not step.uses_default_duration
+    assert step.duration_seconds(parameter_values, inputs) == 1200.0
+
+
+def test_symbolic_duration_needs_its_input():
+    step = pybamm.step.c_rate(1, duration=pybamm.InputParameter("d"))
+
+    with pytest.raises(KeyError, match="Input parameter 'd' not found"):
+        step.duration_seconds(pybamm.ParameterValues({}))
+
+
+def test_drive_cycle_repeats_independently_of_its_duration():
+    drive_cycle = np.array([[0.0, 1.0], [1.0, -1.0], [2.0, 0.5]])
+    symbolic = pybamm.step.current(drive_cycle, duration=pybamm.InputParameter("d"))
+    numeric = pybamm.step.current(drive_cycle, duration=10)
+
+    # The interpolant wraps the step time, so it does not depend on the duration
+    assert symbolic.value == numeric.value
+    assert symbolic.duration_seconds(pybamm.ParameterValues({}), {"d": 10}) == 10
+    # Sample times are repeated up to the resolved final time, one cycle being 3 s
+    np.testing.assert_array_equal(
+        symbolic._drive_cycle_time_points(8),
+        np.array([0.0, 1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0, 8.0]),
+    )
+
+
+def test_symbolic_duration_and_termination_are_recorded():
+    termination = pybamm.CoupledVariable("Voltage [V]") < pybamm.InputParameter("Vmin")
+    duration = pybamm.InputParameter("d")
+    step = pybamm.step.c_rate(1, duration=duration, termination=termination)
+
+    # A Symbol has no truth value, so its string form is what gets recorded
+    assert f"duration={duration}" in repr(step)
+    assert f"termination={termination}" in step.basic_repr()
+    # Steps differing only in duration still share a model
+    assert (
+        step.basic_repr()
+        == pybamm.step.c_rate(1, duration=1800, termination=termination).basic_repr()
+    )
+
+
+@pytest.mark.parametrize(
+    "period, inputs",
+    [
+        (pybamm.InputParameter("sample"), {"sample": 300}),
+        (pybamm.Parameter("sample"), {"sample": 300}),
+        ("5 minutes", {}),
+        (300, {}),
+    ],
+)
+def test_period_is_resolved_at_solve_time(period, inputs):
+    parameter_values = pybamm.ParameterValues({"sample": "[input]"})
+    step = pybamm.step.c_rate(1, duration=1800, period=period)
+
+    assert step.period_seconds(parameter_values, inputs) == 300.0
+    # The period sets the output sampling, so it must be resolved before t_interp
+    _, t_interp = step.setup_timestepping(
+        DummySolver(), 1800, parameter_values=parameter_values, inputs=inputs
+    )
+    np.testing.assert_array_equal(t_interp, np.arange(0, 1801, 300.0))
diff --git a/packages/pybamm/tests/unit/test_experiments/test_experiment_step_termination.py b/packages/pybamm/tests/unit/test_experiments/test_experiment_step_termination.py
index 1939df7b91..5ef54cff7e 100644
--- a/packages/pybamm/tests/unit/test_experiments/test_experiment_step_termination.py
+++ b/packages/pybamm/tests/unit/test_experiments/test_experiment_step_termination.py
@@ -205,3 +205,68 @@ def test_unique_steps_independent_of_cycle_count(self):
                 f"unique_steps={len(exp.unique_steps)} for n_cycles={n_cycles}, "
                 f"expected {n_unique} (steps must scale with template, not cycles)"
             )
+
+
+class TestInequalityTermination:
+    @pytest.mark.parametrize(
+        "expression, residual",
+        [
+            # `a > b` is the heaviside `b < a`, so the sides come back swapped
+            (
+                pybamm.CoupledVariable("V") > pybamm.InputParameter("V hold"),
+                lambda v: pybamm.InputParameter("V hold") - v,
+            ),
+            (
+                pybamm.CoupledVariable("V") <= pybamm.InputParameter("V hold"),
+                lambda v: v - pybamm.InputParameter("V hold"),
+            ),
+            (
+                pybamm.CoupledVariable("V") * 2 < 7.0,
+                lambda v: v * 2 - pybamm.Scalar(7.0),
+            ),
+        ],
+    )
+    def test_inequality_becomes_a_custom_termination(self, expression, residual):
+        term = pybamm.step.base_step._read_termination(expression)
+        v = pybamm.Variable("V")
+
+        assert isinstance(term, pybamm.step.CustomTermination)
+        assert term.name == f"{expression} [experiment]"
+        # A heaviside is "left < right", so left - right is positive before the
+        # inequality holds and negative once it does: the event convention. The
+        # CoupledVariable is looked up in the variables the termination is handed,
+        # exactly as every other termination finds its variable.
+        assert term.get_event({"V": v}, None).expression == residual(v)
+
+    def test_inequality_termination_over_a_custom_variable(self):
+        # A variable that is not a standard model output, referenced by name
+        model = pybamm.lithium_ion.SPM()
+        model.variables["Headroom [V]"] = model.variables["Voltage [V]"] - 3.0
+        step = pybamm.step.c_rate(
+            1, duration=3600, termination=pybamm.CoupledVariable("Headroom [V]") < 0.6
+        )
+        sim = pybamm.Simulation(
+            model,
+            experiment=pybamm.Experiment([step]),
+            solver=pybamm.IDAKLUSolver(),
+        )
+
+        sol = sim.solve(calc_esoh=False)
+
+        assert sol.termination == f"event: {step.termination[0].name}"
+        assert sol["Headroom [V]"].data[-1] == pytest.approx(0.6, abs=1e-3)
+
+    def test_inequality_termination_rejects_unknown_variable(self):
+        step = pybamm.step.c_rate(
+            1, duration=3600, termination=pybamm.CoupledVariable("Not a variable") < 1
+        )
+        sim = pybamm.Simulation(
+            pybamm.lithium_ion.SPM(), experiment=pybamm.Experiment([step])
+        )
+
+        with pytest.raises(ValueError, match="'Not a variable' not found"):
+            sim.solve(calc_esoh=False)
+
+    def test_symbolic_termination_must_be_an_inequality(self):
+        with pytest.raises(TypeError, match="must be an inequality between symbols"):
+            pybamm.step.c_rate(1, duration=1, termination=pybamm.InputParameter("Vmin"))
diff --git a/packages/pybamm/tests/unit/test_experiments/test_simulation_with_experiment.py b/packages/pybamm/tests/unit/test_experiments/test_simulation_with_experiment.py
index 13f6c719e8..ea45711ac0 100644
--- a/packages/pybamm/tests/unit/test_experiments/test_simulation_with_experiment.py
+++ b/packages/pybamm/tests/unit/test_experiments/test_simulation_with_experiment.py
@@ -2598,3 +2598,66 @@ def test_repeated_solves_refresh_initial_soc(self, experiment_model_mode):
 
         # Reusing the same Simulation must refresh experiment ICs when SOC changes.
         assert ic1 != ic2
+
+    @pytest.mark.parametrize("experiment_model_mode", ["unified", "legacy"])
+    def test_run_experiment_with_symbolic_duration(self, experiment_model_mode):
+        experiment = pybamm.Experiment(
+            [pybamm.step.c_rate(1, duration=pybamm.InputParameter("Step duration [s]"))]
+        )
+        sim = pybamm.Simulation(
+            pybamm.lithium_ion.SPM(),
+            experiment=experiment,
+            solver=pybamm.IDAKLUSolver(),
+            experiment_model_mode=experiment_model_mode,
+        )
+
+        for duration in (600, 1200):
+            sol = sim.solve(inputs={"Step duration [s]": duration}, calc_esoh=False)
+            assert sol.termination == "final time"
+            assert sol.t[-1] == pytest.approx(duration)
+
+    @pytest.mark.parametrize("experiment_model_mode", ["unified", "legacy"])
+    def test_run_experiment_with_symbolic_termination(self, experiment_model_mode):
+        step = pybamm.step.c_rate(
+            1,
+            duration=3600,
+            termination=pybamm.CoupledVariable("Voltage [V]")
+            < pybamm.InputParameter("Voltage cut-off [V]"),
+        )
+        experiment = pybamm.Experiment([step])
+        sim = pybamm.Simulation(
+            pybamm.lithium_ion.SPM(),
+            experiment=experiment,
+            solver=pybamm.IDAKLUSolver(),
+            experiment_model_mode=experiment_model_mode,
+        )
+
+        # The threshold is only read at solve time, so one built model serves both
+        for cut_off in (3.6, 3.5):
+            sol = sim.solve(inputs={"Voltage cut-off [V]": cut_off}, calc_esoh=False)
+            assert sol.termination == f"event: {step.termination[0].name}"
+            assert sol["Voltage [V]"].data[-1] == pytest.approx(cut_off, abs=1e-3)
+            assert sol.t[-1] < 3600
+
+    def test_run_experiment_drive_cycle_with_symbolic_duration(self):
+        drive_cycle = np.column_stack([np.arange(0, 501, 100.0), [1, 2, -1, 0, 1, -2]])
+        experiment = pybamm.Experiment(
+            [pybamm.step.current(drive_cycle, duration=pybamm.InputParameter("dur"))]
+        )
+        sim = pybamm.Simulation(
+            pybamm.lithium_ion.SPM(),
+            experiment=experiment,
+            solver=pybamm.IDAKLUSolver(),
+        )
+
+        # The interpolant wraps the step time, so the same built model serves any
+        # duration: the drive cycle repeats for as long as the step runs
+        sol = sim.solve(inputs={"dur": 1700}, calc_esoh=False)
+        assert sol.termination == "final time"
+        assert sol.t[-1] == pytest.approx(1700)
+        current = sol["Current [A]"]
+        for t in (0.0, 150.0, 400.0):
+            assert current(t) == pytest.approx(current(t + 600.0))
+
+        sol = sim.solve(inputs={"dur": 300}, calc_esoh=False)
+        assert sol.t[-1] == pytest.approx(300)
diff --git a/packages/pybamm/tests/unit/test_models/test_symbol_processor.py b/packages/pybamm/tests/unit/test_models/test_symbol_processor.py
index 5238b60e2d..b0c39beecb 100644
--- a/packages/pybamm/tests/unit/test_models/test_symbol_processor.py
+++ b/packages/pybamm/tests/unit/test_models/test_symbol_processor.py
@@ -144,3 +144,22 @@ def test_call_without_setup(self):
 
         with pytest.raises(ValueError, match=r"Cannot process a symbol"):
             processor("test", symbol)
+
+
+class TestSymbolProcessorResolve:
+    def test_resolve_replaces_coupled_variables_by_name(self):
+        a = pybamm.Variable("a")
+        variables = {"a": a, "twice a": 2 * pybamm.CoupledVariable("a")}
+
+        # Nested references resolve all the way down to the model's own symbols
+        resolved = pybamm.SymbolProcessor.resolve(
+            pybamm.CoupledVariable("twice a") + pybamm.CoupledVariable("a"), variables
+        )
+
+        assert resolved == 2 * a + a
+        # A symbol with nothing to resolve is returned as-is, not copied
+        assert pybamm.SymbolProcessor.resolve(a, variables) is a
+
+    def test_resolve_raises_for_an_unknown_name(self):
+        with pytest.raises(ValueError, match="CoupledVariable 'b' not found"):
+            pybamm.SymbolProcessor.resolve(pybamm.CoupledVariable("b"), {})
diff --git a/packages/pybamm/tests/unit/test_serialisation/test_serialisation.py b/packages/pybamm/tests/unit/test_serialisation/test_serialisation.py
index a4a95e6cd3..804e468bcb 100644
--- a/packages/pybamm/tests/unit/test_serialisation/test_serialisation.py
+++ b/packages/pybamm/tests/unit/test_serialisation/test_serialisation.py
@@ -2932,6 +2932,36 @@ def test_input_parameter_value_round_trip(self):
         exp2 = pybamm.Experiment.from_config(config)
         assert isinstance(exp2.steps[0].value, pybamm.InputParameter)
 
+    @pytest.mark.parametrize(
+        "duration",
+        [
+            pybamm.InputParameter("Step duration [s]"),
+            2 * pybamm.InputParameter("d"),
+        ],
+    )
+    def test_symbolic_duration_round_trip(self, duration):
+        exp = pybamm.Experiment([pybamm.step.c_rate(1, duration=duration)])
+
+        # The config must survive a trip through JSON, not just through Python
+        config = json.loads(json.dumps(exp.to_config()))
+        exp2 = pybamm.Experiment.from_config(config)
+
+        assert exp2.steps[0].duration == duration
+        assert not exp2.steps[0].uses_default_duration
+
+    @pytest.mark.parametrize(
+        "period", [pybamm.InputParameter("sample"), 2 * pybamm.Parameter("sample")]
+    )
+    def test_symbolic_period_round_trip(self, period):
+        exp = pybamm.Experiment([pybamm.step.c_rate(1, duration=1800, period=period)])
+
+        config = json.loads(json.dumps(exp.to_config()))
+        exp2 = pybamm.Experiment.from_config(config)
+
+        assert exp2.steps[0].period == period
+        # The step's own value must survive alongside it
+        assert exp2.steps[0].value == 1.0
+
     def test_legacy_steps_format(self):
         """from_config also accepts flat {'steps': [...]} format."""
         config = {

From 0f23dca49228bf55a41b49ad3e18e8db7d356b19 Mon Sep 17 00:00:00 2001
From: Marc Berliner <34451391+MarcBerliner@users.noreply.github.com>
Date: Fri, 4 Sep 2026 12:50:31 -0400
Subject: [PATCH 2/6] feat: symbolic step settings and terminations

Step duration, period, temperature and control value may be a
`pybamm.Parameter` or `pybamm.InputParameter`. The parameter pass runs
once at setup (`BaseStep.process_parameters`); solving evaluates each
setting with the solver inputs.

`SymbolicTermination` terminates a step on any inequality over model
variables. The voltage, current and C-rate terminations build on it,
and every termination serialises, symbolic thresholds included.

`SymbolProcessor.resolve` is the one place `CoupledVariable`s are looked
up, replacing the copies in `BaseModel` and `Discretisation`.

Co-Authored-By: Claude Fable 5.1 <noreply@anthropic.com>
---
 CHANGELOG.md                                  |   6 +
 .../api/experiment/experiment_steps.rst       |   5 +
 .../custom-experiments.ipynb                  |  13 +-
 .../src/pybamm/experiment/step/__init__.py    |   2 +-
 .../src/pybamm/experiment/step/base_step.py   | 194 ++++++++++--------
 .../experiment/step/step_termination.py       | 151 +++++++-------
 .../expression_tree/operations/diffsl.py      |  17 +-
 .../expression_tree/operations/serialise.py   |  33 +--
 .../src/pybamm/models/symbol_processor.py     |   2 +-
 .../src/pybamm/simulation/simulation.py       |  28 +--
 .../unit/test_experiments/test_base_step.py   |  57 ++++-
 .../test_experiment_step_termination.py       |  45 ++--
 .../test_simulation_with_experiment.py        |   4 +-
 .../test_serialisation/test_serialisation.py  |  20 ++
 14 files changed, 330 insertions(+), 247 deletions(-)

diff --git a/CHANGELOG.md b/CHANGELOG.md
index 7a1fdcefe8..b3298893dd 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -1,5 +1,11 @@
 # [Unreleased](https://github.com/pybamm-team/PyBaMM/)
 
+## Features
+
+- Step `duration`, `period` and `temperature` accept a `pybamm.Parameter` or `pybamm.InputParameter`.
+- Added `SymbolicTermination`: terminate a step on any inequality over model variables.
+- Voltage, current and C-rate terminations build on `SymbolicTermination` and serialise symbolic thresholds.
+
 ## Bug fixes
 
 - The `integration` nox session no longer installs the `pydiffsol` extra on macOS Intel CI runners, where it has no working build. ([#5726](https://github.com/pybamm-team/PyBaMM/pull/5726))
diff --git a/docs/source/api/experiment/experiment_steps.rst b/docs/source/api/experiment/experiment_steps.rst
index 3a1d2d78d3..b73e5cf8b1 100644
--- a/docs/source/api/experiment/experiment_steps.rst
+++ b/docs/source/api/experiment/experiment_steps.rst
@@ -33,6 +33,11 @@ Custom steps can be defined using either explicit or implicit control:
 Step terminations
 -----------------
 
+Any inequality over model variables can terminate a step:
+
+.. autoclass:: pybamm.step.SymbolicTermination
+    :members:
+
 Standard step termination events are implemented by the following classes, which are
 called when the termination is specified by a specific string. These classes can be
 either be called directly or via the string format specified in the class docstring
diff --git a/docs/source/examples/notebooks/simulations_and_experiments/custom-experiments.ipynb b/docs/source/examples/notebooks/simulations_and_experiments/custom-experiments.ipynb
index 37f0b94fb4..a6fdead898 100644
--- a/docs/source/examples/notebooks/simulations_and_experiments/custom-experiments.ipynb
+++ b/docs/source/examples/notebooks/simulations_and_experiments/custom-experiments.ipynb
@@ -31,7 +31,7 @@
    "source": [
     "## Custom termination\n",
     "\n",
-    "Termination of a step can be specified using a few standard strings (e.g. \"4.2V\" for voltage, \"1 A\" for current, \"C/2\" for C-rate), or via a custom termination step. The custom termination step can be specified based on any variable in the model.\n",
+    "Termination of a step can be specified using a few standard strings (e.g. \"4.2V\" for voltage, \"1 A\" for current, \"C/2\" for C-rate), or as an inequality on any model variable, referenced by name with `pybamm.CoupledVariable`. The threshold can be a number, a `pybamm.Parameter` or a `pybamm.InputParameter`.\n",
     "Below, we show an example where we specify a custom termination step based on keeping the anode potential above 0V, which is a common limit used to avoid lithium plating,"
    ]
   },
@@ -63,16 +63,9 @@
     "parameter_values = pybamm.ParameterValues(\"Chen2020\")\n",
     "\n",
     "\n",
-    "# Create a custom termination event for the anode potential cut-off at 0.02V\n",
+    "# Terminate when the anode potential drops to 0.02V\n",
     "# We use 0.02V instead of 0V to give a safety factor\n",
-    "def anode_potential_cutoff(variables):\n",
-    "    return variables[\"Anode potential [V]\"] - 0.02\n",
-    "\n",
-    "\n",
-    "# The CustomTermination class takes a name and function\n",
-    "anode_potential_termination = pybamm.step.CustomTermination(\n",
-    "    name=\"Anode potential cut-off [V]\", event_function=anode_potential_cutoff\n",
-    ")\n",
+    "anode_potential_termination = pybamm.CoupledVariable(\"Anode potential [V]\") < 0.02\n",
     "\n",
     "# Provide a list of termination events, each step will stop whenever the first\n",
     "# termination event is reached\n",
diff --git a/packages/pybamm/src/pybamm/experiment/step/__init__.py b/packages/pybamm/src/pybamm/experiment/step/__init__.py
index 0ee7f4a1ee..9cf7b834da 100644
--- a/packages/pybamm/src/pybamm/experiment/step/__init__.py
+++ b/packages/pybamm/src/pybamm/experiment/step/__init__.py
@@ -1,5 +1,5 @@
 from .steps import *
 from .base_step import BaseStep, BaseStepExplicit, BaseStepImplicit
-from .step_termination import BaseTermination, CurrentTermination, VoltageTermination, CustomTermination, CRateTermination, CrateTermination, _read_termination
+from .step_termination import BaseTermination, CurrentTermination, VoltageTermination, SymbolicTermination, CustomTermination, CRateTermination, CrateTermination, _read_termination
 
 __all__ = ['base_step', 'step_termination', 'steps']
diff --git a/packages/pybamm/src/pybamm/experiment/step/base_step.py b/packages/pybamm/src/pybamm/experiment/step/base_step.py
index d674b5e55c..8117ae246f 100644
--- a/packages/pybamm/src/pybamm/experiment/step/base_step.py
+++ b/packages/pybamm/src/pybamm/experiment/step/base_step.py
@@ -39,6 +39,23 @@ class ControlKind(str, Enum):
     """
 
 
+class _Direction(str, Enum):
+    """The potential directions of a step."""
+
+    CHARGE = "charge"
+    DISCHARGE = "discharge"
+    REST = "rest"
+
+
+class _SymbolicInput(str, Enum):
+    """Step inputs which may be a symbolic function of parameters or real numbers."""
+
+    CONTROL_TARGET = "_control_target"
+    TEMPERATURE = "temperature"
+    PERIOD = "period"
+    DURATION = "duration"
+
+
 class BaseStep:
     """
     Class representing one step in an experiment.
@@ -74,7 +91,7 @@ class BaseStep:
     description : str, optional
         A description of the step.
     direction : str, optional
-        The direction of the step, e.g. "Charge" or "Discharge" or "Rest".
+        The direction of the step, e.g. "charge" or "discharge" or "rest".
     skip_ok : bool, optional
         If True, the step will be skipped if it is infeasible at the initial conditions.
         Default is True.
@@ -93,10 +110,12 @@ def __init__(
         direction: str | None = None,
         skip_ok: bool = True,
     ):
-        potential_directions = ["charge", "discharge", "rest", None]
-        if direction not in potential_directions:
+        # Filled by `process_parameters` at setup, see `evaluate`
+        self._processed_variables = {}
+        if direction not in _Direction and direction is not None:
+            _DIRECTIONS = [d.value for d in _Direction] + [None]
             raise ValueError(
-                f"Invalid direction: {direction}. Must be one of {potential_directions}"
+                f"Invalid direction: {direction}. Must be one of {_DIRECTIONS}"
             )
         self.input_duration = duration
         self.input_value = value
@@ -211,17 +230,49 @@ def _build_termination(term):
         self.next_start_time = None
         self.end_time = None
 
-    def duration_seconds(self, parameter_values, inputs=None):
-        """The duration in seconds, resolving a symbolic duration at solve time."""
-        duration = parameter_values.process_symbol(
-            pybamm.convert_to_symbol(self.duration)
-        )
-        return float(duration.evaluate(inputs=inputs))
+    def process_parameters(self, parameter_values: pybamm.ParameterValues):
+        """Substitute parameter values into the step's symbolic settings. Done once at
+        setup, so that at solve time each needs only the solver inputs."""
+        for parameter in _SymbolicInput:
+            name = parameter.value
+            value = getattr(self, name)
+            self._processed_variables[name] = parameter_values.process_symbol(value)
+
+    def _evaluate(self, name: str, inputs: dict[str, float]) -> float | None:
+        # A step that has not been set up still answers for its numeric settings
+        value = self._processed_variables.get(name, getattr(self, name))
+        if isinstance(value, pybamm.Symbol):
+            return float(value.evaluate(inputs=inputs))
+        return value
+
+    def period_seconds(self, inputs: dict[str, float]) -> float | None:
+        """Evaluate the step's period in seconds, using the solver inputs if it is
+        symbolic."""
+        return self._evaluate(_SymbolicInput.PERIOD, inputs)
+
+    def duration_seconds(self, inputs: dict[str, float]) -> float | None:
+        """Evaluate the step's duration in seconds, using the solver inputs if it is
+        symbolic."""
+        return self._evaluate(_SymbolicInput.DURATION, inputs)
+
+    def temperature_kelvin(self, inputs: dict[str, float]) -> float | None:
+        """Evaluate the step's temperature in Kelvin, using the solver inputs if it is
+        symbolic."""
+        return self._evaluate(_SymbolicInput.TEMPERATURE, inputs)
+
+    def control_target_value(self, inputs: dict[str, float]) -> float | None:
+        """The control target as a number, or ``None`` if it is not a constant."""
+        return self._evaluate(_SymbolicInput.CONTROL_TARGET, inputs)
 
-    def period_seconds(self, parameter_values, inputs=None):
-        """The period in seconds, resolving a symbolic period at solve time."""
-        period = parameter_values.process_symbol(pybamm.convert_to_symbol(self.period))
-        return float(period.evaluate(inputs=inputs))
+    @property
+    def has_symbolic_control_target(self) -> bool:
+        """Whether the control target is not a constant once parameters are substituted
+        (an input parameter, a drive cycle, or a custom step with no target)."""
+        name = _SymbolicInput.CONTROL_TARGET
+        target = self._processed_variables.get(name, getattr(self, name))
+        if isinstance(target, pybamm.Symbol):
+            return not target.is_constant()
+        return target is None
 
     @staticmethod
     def is_implicit() -> bool:
@@ -322,10 +373,6 @@ def _control_target(self):
         """The control target as a symbol (the prescribed current, voltage, ...)."""
         return self.value
 
-    def control_target_value(self, parameter_values):
-        """The control target as a number, or ``None`` if it is not a constant."""
-        return _constant_value(parameter_values, self._control_target)
-
     def unified_branch_repr(self):
         """Branch identity in unified mode: control kind and everything but the value."""
         parts = [self.control_kind or type(self).__name__]
@@ -358,9 +405,9 @@ def default_duration(self, value):
     def default_period():
         return 60.0  # seconds
 
-    def default_time_vector(self, solver, tf, t0=0, parameter_values=None, inputs=None):
+    def default_time_vector(self, solver, tf, t0=0, inputs=None):
         if self.period is not None:
-            period = self.period_seconds(parameter_values, inputs)
+            period = self.period_seconds(inputs)
         elif self.is_drive_cycle and solver.supports_interp:
             # Infer the period from the drive cycle
             period = np.diff(self.value.x[0]).min()
@@ -370,9 +417,7 @@ def default_time_vector(self, solver, tf, t0=0, parameter_values=None, inputs=No
 
         return np.linspace(t0, tf, npts)
 
-    def setup_timestepping(
-        self, solver, tf, t_interp=None, parameter_values=None, inputs=None
-    ):
+    def setup_timestepping(self, solver, tf, t_interp=None, inputs=None):
         """
         Setup timestepping for the model.
 
@@ -386,17 +431,11 @@ def setup_timestepping(
             The time points at which to interpolate the solution
         """
         if solver.supports_interp:
-            return self._setup_timestepping(
-                solver, tf, t_interp, parameter_values, inputs
-            )
+            return self._setup_timestepping(solver, tf, t_interp, inputs)
         else:
-            return self._setup_timestepping_dense_t_eval(
-                solver, tf, t_interp, parameter_values, inputs
-            )
+            return self._setup_timestepping_dense_t_eval(solver, tf, t_interp, inputs)
 
-    def _setup_timestepping(
-        self, solver, tf, t_interp, parameter_values=None, inputs=None
-    ):
+    def _setup_timestepping(self, solver, tf, t_interp, inputs=None):
         """
         Setup timestepping for the model. This returns a t_eval vector that stops
         only at the first and last time points. If t_interp and the period are
@@ -423,9 +462,7 @@ def _setup_timestepping(
 
         if t_interp is None:
             if self.period is not None:
-                t_interp = self.default_time_vector(
-                    solver, tf, parameter_values=parameter_values, inputs=inputs
-                )
+                t_interp = self.default_time_vector(solver, tf, inputs=inputs)
             else:
                 t_interp = solver.process_t_interp(t_interp)
 
@@ -442,9 +479,7 @@ def _drive_cycle_time_points(self, tf):
             t_eval = np.append(t_eval, tf)
         return t_eval
 
-    def _setup_timestepping_dense_t_eval(
-        self, solver, tf, t_interp, parameter_values=None, inputs=None
-    ):
+    def _setup_timestepping_dense_t_eval(self, solver, tf, t_interp, inputs=None):
         """
         Setup timestepping for the model. By default, this returns a dense t_eval which
         stops the solver at each point in the t_eval vector. This method is for solvers
@@ -459,9 +494,7 @@ def _setup_timestepping_dense_t_eval(
         t_interp: np.array | None
             The time points at which to interpolate the solution
         """
-        t_eval = self.default_time_vector(
-            solver, tf, parameter_values=parameter_values, inputs=inputs
-        )
+        t_eval = self.default_time_vector(solver, tf, inputs=inputs)
 
         t_interp = solver.process_t_interp(t_interp)
 
@@ -549,11 +582,11 @@ def value_based_charge_or_discharge(self):
             init_curr = self.value
         sign = np.sign(init_curr)
         if sign == 0:
-            return "rest"
+            return _Direction.REST.value
         elif sign > 0:
-            return "discharge"
+            return _Direction.DISCHARGE.value
         else:
-            return "charge"
+            return _Direction.CHARGE.value
 
     def record_tags(
         self,
@@ -570,32 +603,28 @@ def record_tags(
         """Record all the args for repr and hash"""
         # A Symbol has no truth value and no stable repr, but its id identifies the
         # expression it stands for
-        if isinstance(duration, pybamm.Symbol):
-            duration = str(duration)
-        if isinstance(termination, pybamm.Symbol):
-            termination = str(termination)
-        if isinstance(period, pybamm.Symbol):
-            period = str(period)
-        repr_args = f"{value}, duration={duration}"
-        hash_args = f"{value}"
-        if termination:
-            repr_args += f", termination={termination}"
-            hash_args += f", termination={termination}"
-        if period:
-            repr_args += f", period={period}"
-        if temperature:
-            repr_args += f", temperature={temperature}"
-            hash_args += f", temperature={temperature}"
-        if tags:
-            repr_args += f", tags={tags}"
-        if start_time:
-            repr_args += f", start_time={start_time}"
-        if description:
-            repr_args += f", description={description}"
-        if direction:
-            repr_args += f", direction={direction}"
-            hash_args += f", direction={direction}"
-        return repr_args, hash_args
+
+        reprs = [str(value)]
+        hashes = [str(value)]
+
+        def record(name, item, *, hashed):
+            if isinstance(item, pybamm.Symbol):
+                item = str(item)
+            if not item:
+                return
+            reprs.append(f"{name}={item}")
+            if hashed:
+                hashes.append(f"{name}={item}")
+
+        record("duration", duration, hashed=False)
+        record("termination", termination, hashed=True)
+        record("period", period, hashed=False)
+        record("temperature", temperature, hashed=True)
+        record("tags", tags, hashed=False)
+        record("start_time", start_time, hashed=False)
+        record("description", description, hashed=False)
+        record("direction", direction, hashed=True)
+        return ", ".join(reprs), ", ".join(hashes)
 
 
 class BaseStepExplicit(BaseStep):
@@ -678,16 +707,6 @@ def set_up(self, new_model, new_parameter_values):
 }
 
 
-def _constant_value(parameter_values, target):
-    """Numeric value of ``target`` if it is a state/time-independent constant (after
-    parameter substitution), else ``None``. ``target`` may be ``None`` (custom steps),
-    a number, or a symbol (e.g. a drive-cycle interpolant, which is not constant)."""
-    if target is None:
-        return None
-    processed = parameter_values.process_symbol(pybamm.convert_to_symbol(target))
-    return float(processed.evaluate()) if processed.is_constant() else None
-
-
 def get_unit_from(a_string: str) -> str:
     return a_string.lstrip("0123456789.-eE ")
 
@@ -726,8 +745,12 @@ def _convert_time_to_seconds(time_and_units):
 
 def _convert_temperature_to_kelvin(temperature_and_units):
     """Convert a temperature in Celsius or Kelvin to a temperature in Kelvin"""
-    # If the temperature is a number, assume it is in Kelvin
-    if isinstance(temperature_and_units, int | float) or temperature_and_units is None:
+    # If the temperature is a number, assume it is in Kelvin; a symbol is resolved at
+    # solve time, see `BaseStep.process_parameters`
+    if (
+        isinstance(temperature_and_units, int | float | pybamm.Symbol)
+        or temperature_and_units is None
+    ):
         return temperature_and_units
 
     # Split number and units
@@ -795,9 +818,8 @@ def _parse_termination(term_str, value):
 
 
 def _check_input_params(value):
-    """Check if self.value is a function of input parameters"""
+    """Check if a step's value depends on parameters that are only known at solve time"""
     leaves = value.post_order(filter=lambda node: len(node.children) == 0)
-    contains_input_parameter = any(
-        isinstance(leaf, pybamm.InputParameter) for leaf in leaves
+    return any(
+        isinstance(leaf, pybamm.InputParameter | pybamm.Parameter) for leaf in leaves
     )
-    return contains_input_parameter
diff --git a/packages/pybamm/src/pybamm/experiment/step/step_termination.py b/packages/pybamm/src/pybamm/experiment/step/step_termination.py
index 1077f3010c..1e845a3865 100644
--- a/packages/pybamm/src/pybamm/experiment/step/step_termination.py
+++ b/packages/pybamm/src/pybamm/experiment/step/step_termination.py
@@ -1,4 +1,3 @@
-from dataclasses import dataclass
 from warnings import warn
 
 import pybamm
@@ -8,8 +7,8 @@ class BaseTermination:
     """
     Base class for a termination event for an experiment step. To create a custom
     termination, a class must implement `get_event_expression` to return the symbolic
-    expression for the event. In most cases the class
-    :class:`pybamm.step.CustomTermination` can be used to assist with this.
+    expression for the event. In most cases :class:`pybamm.step.SymbolicTermination`
+    can be used instead.
 
     Parameters
     ----------
@@ -60,18 +59,52 @@ def __hash__(self) -> int:
         return hash((type(self).__name__, self.value, self.operator))
 
 
-class CRateTermination(BaseTermination):
+class SymbolicTermination(BaseTermination):
+    """
+    Termination defined by an inequality over model variables, input parameters and
+    parameters, e.g. ``pybamm.CoupledVariable("Voltage [V]") > pybamm.InputParameter(
+    "V hold")``: the step terminates once it holds. Model variables are referenced by
+    name as :class:`pybamm.CoupledVariable`. The voltage, current and C-rate
+    terminations are special cases that build their inequality from a threshold.
+
+    Parameters
+    ----------
+    value : :class:`pybamm.EqualHeaviside` or :class:`pybamm.NotEqualHeaviside`
+        The inequality.
+    """
+
+    def inequality(self, step):
+        """The inequality, or ``None`` if the step has no such event."""
+        return self.value
+
+    def get_event_name(self, step):
+        return f"{self.inequality(step)} [experiment]"
+
+    def get_event_expression(self, variables, step):
+        inequality = self.inequality(step)
+        if inequality is None:
+            return None
+        # An inequality is "left < right", so left - right is positive before it holds
+        # and zero when the event triggers
+        residual = pybamm.SymbolProcessor.resolve(
+            inequality.left - inequality.right, variables
+        )
+        # A termination over known numbers is a number
+        return residual.evaluate() if residual.is_constant() else residual
+
+
+class CRateTermination(SymbolicTermination):
     """
     Termination based on C-rate, created when a string termination of the C-rate type
     (e.g. "C/10") is provided
     """
 
+    def inequality(self, step):
+        return abs(pybamm.CoupledVariable("C-rate")) < self.value
+
     def get_event_name(self, step):
         return "C-rate cut-off [experiment]"
 
-    def get_event_expression(self, variables, step):
-        return abs(variables["C-rate"]) - self.value
-
 
 class CrateTermination(CRateTermination):
     """
@@ -87,71 +120,44 @@ def __init__(self, value, operator=None):
         warn(warning, stacklevel=2)
 
 
-class CurrentTermination(BaseTermination):
+class CurrentTermination(SymbolicTermination):
     """
     Termination based on current, created when a string termination of the current type
     (e.g. "1A") is provided
     """
 
-    def get_event_name(self, step):
-        operator = self.operator
-        if operator == ">":
-            return f"Current [A] > {self.value} [A] [experiment]"
-        elif operator == "<":
-            return f"Current [A] < {self.value} [A] [experiment]"
-        else:
-            return f"abs(Current [A]) < {self.value} [A] [experiment]"
+    def inequality(self, step):
+        current = pybamm.CoupledVariable("Current [A]")
+        if self.operator == ">":
+            return current > self.value
+        if self.operator == "<":
+            return current < self.value
+        return abs(current) < self.value
 
-    def get_event_expression(self, variables, step):
-        operator = self.operator
-        if operator == ">":
-            return self.value - variables["Current [A]"]
-        elif operator == "<":
-            return variables["Current [A]"] - self.value
-        else:
-            return abs(variables["Current [A]"]) - self.value
+    def get_event_name(self, step):
+        return str(self.inequality(step)) + " [experiment]"
 
 
-class VoltageTermination(BaseTermination):
+class VoltageTermination(SymbolicTermination):
     """
     Termination based on voltage, created when a string termination of the voltage type
-    (e.g. "4.2V") is provided
+    (e.g. "4.2V") is provided. Without an operator, the step's direction decides
+    whether the cut-off is above or below the voltage.
     """
 
-    def _get_operator(self, step):
-        operator = self.operator
-        if operator is None:
-            direction = step.direction
-            if direction == "charge":
-                operator = ">"
-            elif direction == "discharge":
-                operator = "<"
-            else:
-                return None
-        return operator
-
-    def get_event_name(self, step):
-        operator = self._get_operator(step)
-        if operator is None:
-            return None
-        return f"Voltage {operator} {self.value} [V] [experiment]"
+    def _operator(self, step):
+        return self.operator or {"charge": ">", "discharge": "<"}.get(step.direction)
 
-    def get_event_expression(self, variables, step):
-        # The voltage event should be positive at the start of charge/
-        # discharge. We use the sign of the current or power input to
-        # figure out whether the voltage event is greater than the starting
-        # voltage (charge) or less (discharge) and set the sign of the
-        # event accordingly
-        operator = self._get_operator(step)
+    def inequality(self, step):
+        operator = self._operator(step)
         if operator is None:
             return None
+        voltage = pybamm.CoupledVariable("Battery voltage [V]")
+        return voltage > self.value if operator == ">" else voltage < self.value
 
-        if operator == ">":
-            sign = -1
-        else:
-            sign = 1
-
-        return sign * (variables["Battery voltage [V]"] - self.value)
+    def get_event_name(self, step):
+        operator = self._operator(step)
+        return operator and f"Voltage {operator} {self.value} [V] [experiment]"
 
 
 class Voltage:
@@ -239,31 +245,14 @@ def _parse_termination_class(termination_tuple) -> BaseTermination:
     return cls(value, operator=op)
 
 
-@dataclass(frozen=True, slots=True)
-class _HeavsideResidual:
-    symbol: pybamm.Symbol
-
-    def __call__(self, variables: dict[str, pybamm.Symbol]) -> pybamm.Symbol:
-        # A heaviside is "left < right", so left - right is the event residual
-        residual = self.symbol.left - self.symbol.right
-        return pybamm.SymbolProcessor.resolve(residual, variables)
-
-
-def _parse_termination_symbol(termination: pybamm.Symbol) -> CustomTermination:
-    if not isinstance(termination, (pybamm.EqualHeaviside, pybamm.NotEqualHeaviside)):
+def _read_termination(termination):
+    if isinstance(termination, pybamm.EqualHeaviside | pybamm.NotEqualHeaviside):
+        return SymbolicTermination(termination)
+    if isinstance(termination, pybamm.Symbol):
         raise TypeError(
             "A symbolic termination must be an inequality between symbols, e.g. "
-            'pybamm.CoupledVariable("Voltage [V]") > 3.0, but got '
-            f"'{termination!s}'. Note that inequalities are only exact when "
-            'pybamm.settings.heaviside_smoothing is "exact".'
+            f'pybamm.CoupledVariable("Voltage [V]") > 3.0, but got "{termination!s}".'
         )
-    return CustomTermination(str(termination), _HeavsideResidual(termination))
-
-
-def _read_termination(termination):
-    if isinstance(termination, pybamm.Symbol):
-        return _parse_termination_symbol(termination)
-    if not isinstance(termination, tuple):
-        return termination
-
-    return _parse_termination_class(termination)
+    if isinstance(termination, tuple):
+        return _parse_termination_class(termination)
+    return termination
diff --git a/packages/pybamm/src/pybamm/expression_tree/operations/diffsl.py b/packages/pybamm/src/pybamm/expression_tree/operations/diffsl.py
index a5d7a6e50e..466347887a 100644
--- a/packages/pybamm/src/pybamm/expression_tree/operations/diffsl.py
+++ b/packages/pybamm/src/pybamm/expression_tree/operations/diffsl.py
@@ -450,10 +450,14 @@ def _experiment_schedule_key(
         padding_duration = self._padding_step_duration(
             step, effective_duration, initial_start_time
         )
-        target = step.control_target_value(sim._parameter_values)
-        ambient = (
-            step.temperature or sim._parameter_values[sim._AMBIENT_TEMPERATURE_INPUT]
-        )
+        if step.has_symbolic_control_target:
+            target = None
+        else:
+            target = step.control_target_value(inputs=None)
+
+        ambient = step.temperature
+        if ambient is None:
+            ambient = sim._parameter_values[sim._AMBIENT_TEMPERATURE_INPUT]
 
         return (
             branch_index,
@@ -496,7 +500,10 @@ def _get_unified_experiment_schedule_states(
                 stop_expr = duration_stop
             else:
                 stop_expr = pybamm.minimum(duration_stop, branch)
-            target = step.control_target_value(sim._parameter_values)
+            if step.has_symbolic_control_target:
+                target = None
+            else:
+                target = step.control_target_value(inputs=None)
             schedule_states.append(
                 _ExperimentScheduleState(
                     len(schedule_states),
diff --git a/packages/pybamm/src/pybamm/expression_tree/operations/serialise.py b/packages/pybamm/src/pybamm/expression_tree/operations/serialise.py
index a7e7f930d2..456caf38c2 100644
--- a/packages/pybamm/src/pybamm/expression_tree/operations/serialise.py
+++ b/packages/pybamm/src/pybamm/expression_tree/operations/serialise.py
@@ -1804,6 +1804,7 @@ def serialise_experiment(experiment) -> dict:
             "CurrentTermination": "current",
             "CrateTermination": "c-rate",
             "CRateTermination": "c-rate",
+            "SymbolicTermination": "symbolic",
         }
 
         # Top-level defaults; per-step values are emitted only when they differ.
@@ -1858,24 +1859,24 @@ def _serialise_step(step):
                 terminations = []
                 for term in step.termination:
                     term_class_name = term.__class__.__name__
-                    if term_class_name == "CustomTermination":
-                        raise NotImplementedError(
-                            "CustomTermination cannot be serialised: it "
-                            "carries a user-supplied Python callable "
-                            "(``event_function``) that has no JSON "
-                            "representation."
-                        )
                     if term_class_name not in termination_type_map:
                         raise NotImplementedError(
                             f"Cannot serialise termination of type "
                             f"{term_class_name!r}: only the built-in "
                             f"termination classes "
                             f"({sorted(termination_type_map)!r}) are "
-                            f"supported."
+                            f"supported. A CustomTermination carries a Python "
+                            "callable; use a SymbolicTermination instead."
                         )
-                    term_type = termination_type_map[term_class_name]
-                    term_config = {"type": term_type, "value": term.value}
-                    if hasattr(term, "operator") and term.operator:
+                    # A symbolic threshold or event is an expression tree
+                    value = term.value
+                    if isinstance(value, pybamm.Symbol):
+                        value = encode(value)
+                    term_config = {
+                        "type": termination_type_map[term_class_name],
+                        "value": value,
+                    }
+                    if term.operator:
                         term_config["operator"] = term.operator
                     terminations.append(term_config)
                 step_config["terminations"] = terminations
@@ -1956,14 +1957,16 @@ def deserialise_experiment(data: dict):
 
         def _parse_termination(term_dict):
             term_type = term_dict.get("type")
+            value = term_dict["value"]
+            value = decode(value) if isinstance(value, dict) else float(value)
+            if term_type == "symbolic":
+                return pybamm.step.SymbolicTermination(value)
             if term_type not in term_class_map:
                 raise ValueError(
                     f"Unknown termination type: {term_type!r}. "
-                    f"Expected one of {list(term_class_map)!r}."
+                    f"Expected one of {['symbolic', *term_class_map]!r}."
                 )
-            value = float(term_dict["value"])
-            operator = term_dict.get("operator")
-            return term_class_map[term_type](value, operator=operator)
+            return term_class_map[term_type](value, operator=term_dict.get("operator"))
 
         def _parse_step(step_dict):
             step_type = step_dict.get("type")
diff --git a/packages/pybamm/src/pybamm/models/symbol_processor.py b/packages/pybamm/src/pybamm/models/symbol_processor.py
index b9354f55a7..4b56e76239 100644
--- a/packages/pybamm/src/pybamm/models/symbol_processor.py
+++ b/packages/pybamm/src/pybamm/models/symbol_processor.py
@@ -53,7 +53,7 @@ def resolve(symbol: pybamm.Symbol, variables: dict) -> pybamm.Symbol:
                 raise ValueError(
                     f"CoupledVariable '{symbol.name}' not found in model.variables"
                 )
-            return resolve(variables[symbol.name], variables)
+            return resolve(pybamm.convert_to_symbol(variables[symbol.name]), variables)
         if symbol.children:
             children = [resolve(child, variables) for child in symbol.children]
             if any(
diff --git a/packages/pybamm/src/pybamm/simulation/simulation.py b/packages/pybamm/src/pybamm/simulation/simulation.py
index e1d6857bc7..d734a3fd4e 100644
--- a/packages/pybamm/src/pybamm/simulation/simulation.py
+++ b/packages/pybamm/src/pybamm/simulation/simulation.py
@@ -281,7 +281,7 @@ def _set_up_unified_experiment_model(self, parameter_values):
         value_input = pybamm.InputParameter(self._STEP_VALUE_INPUT)
 
         def is_collapsible(step):
-            return step.control_target_value(self._parameter_values) is not None
+            return not step.has_symbolic_control_target
 
         def branch_key(step):
             return (
@@ -394,14 +394,13 @@ def _build_experiment_step_inputs(
 
         inputs[self._START_TIME_INPUT] = start_time
         if include_temperature:
-            ambient = (
-                step.temperature
-                or self._parameter_values[self._AMBIENT_TEMPERATURE_INPUT]
-            )
-            # The unified model reads ambient temperature as a solver input, so it
-            # must be numeric; the parameter value can be a pybamm.Scalar.
-            if isinstance(ambient, pybamm.Scalar):
-                ambient = ambient.value
+            ambient = step.temperature_kelvin(user_inputs)
+            if ambient is None:
+                ambient = self._parameter_values[self._AMBIENT_TEMPERATURE_INPUT]
+                # The unified model reads ambient temperature as a solver input, so
+                # it must be numeric; the parameter value can be a pybamm.Scalar.
+                if isinstance(ambient, pybamm.Scalar):
+                    ambient = ambient.value
             inputs[self._AMBIENT_TEMPERATURE_INPUT] = ambient
 
         if self._experiment_uses_unified_model:
@@ -409,9 +408,10 @@ def _build_experiment_step_inputs(
             if self._experiment_uses_value_input:
                 # Collapsible steps read their control target from this shared input;
                 # a non-collapsible active step doesn't read it, so pass the dummy.
-                target = step.control_target_value(self._parameter_values)
                 inputs[self._STEP_VALUE_INPUT] = (
-                    target if target is not None else DUMMY_INPUT_PARAMETER_VALUE
+                    DUMMY_INPUT_PARAMETER_VALUE
+                    if step.has_symbolic_control_target
+                    else step.control_target_value(user_inputs)
                 )
         return inputs
 
@@ -540,6 +540,8 @@ def _set_up_and_parameterise_experiment(self, solve_kwargs=None):
         parameter_values = self._parameter_values.copy()
 
         self._check_experiment_input_parameters(parameter_values)
+        for step in self.experiment.steps:
+            step.process_parameters(parameter_values)
 
         if (
             solve_kwargs is not None
@@ -1012,7 +1014,7 @@ def solve(
                 step = experiment_steps[idx]
                 start_time = current_solution.t[-1]
 
-                dt = step.duration_seconds(self._parameter_values, user_inputs)
+                dt = step.duration_seconds(user_inputs)
                 if step.end_time is not None:
                     remaining = (
                         step.end_time
@@ -1046,7 +1048,7 @@ def solve(
                 )
 
                 t_eval, t_interp_processed = step.setup_timestepping(
-                    solver, dt, t_interp, self._parameter_values, user_inputs
+                    solver, dt, t_interp, inputs=user_inputs
                 )
 
                 state_mapper = self._get_state_mapper_for_solution(
diff --git a/packages/pybamm/tests/unit/test_experiments/test_base_step.py b/packages/pybamm/tests/unit/test_experiments/test_base_step.py
index ba3bd156e2..86015d70ee 100644
--- a/packages/pybamm/tests/unit/test_experiments/test_base_step.py
+++ b/packages/pybamm/tests/unit/test_experiments/test_base_step.py
@@ -155,18 +155,25 @@ def test_parse_termination_requires_operator_for_input_parameter():
     ],
 )
 def test_duration_is_resolved_at_solve_time(duration, inputs):
-    parameter_values = pybamm.ParameterValues({"d": "[input]", "fixed d": 1200})
+    parameter_values = pybamm.ParameterValues(
+        {"d": "[input]", "fixed d": 1200, "Nominal cell capacity [A.h]": 1}
+    )
     step = pybamm.step.c_rate(1, duration=duration)
-
     assert not step.uses_default_duration
-    assert step.duration_seconds(parameter_values, inputs) == 1200.0
+
+    # The parameter pass happens once at setup; solving then needs only the inputs
+    step.process_parameters(parameter_values)
+    assert step.duration_seconds(inputs) == 1200.0
+    # The user-facing attribute stays as it was written
+    assert step.duration == duration
 
 
 def test_symbolic_duration_needs_its_input():
     step = pybamm.step.c_rate(1, duration=pybamm.InputParameter("d"))
+    step.process_parameters(pybamm.ParameterValues({"Nominal cell capacity [A.h]": 1}))
 
     with pytest.raises(KeyError, match="Input parameter 'd' not found"):
-        step.duration_seconds(pybamm.ParameterValues({}))
+        step.duration_seconds({})
 
 
 def test_drive_cycle_repeats_independently_of_its_duration():
@@ -176,7 +183,8 @@ def test_drive_cycle_repeats_independently_of_its_duration():
 
     # The interpolant wraps the step time, so it does not depend on the duration
     assert symbolic.value == numeric.value
-    assert symbolic.duration_seconds(pybamm.ParameterValues({}), {"d": 10}) == 10
+    symbolic.process_parameters(pybamm.ParameterValues({}))
+    assert symbolic.duration_seconds({"d": 10}) == 10
     # Sample times are repeated up to the resolved final time, one cycle being 3 s
     np.testing.assert_array_equal(
         symbolic._drive_cycle_time_points(8),
@@ -209,12 +217,41 @@ def test_symbolic_duration_and_termination_are_recorded():
     ],
 )
 def test_period_is_resolved_at_solve_time(period, inputs):
-    parameter_values = pybamm.ParameterValues({"sample": "[input]"})
+    parameter_values = pybamm.ParameterValues(
+        {"sample": "[input]", "Nominal cell capacity [A.h]": 1}
+    )
     step = pybamm.step.c_rate(1, duration=1800, period=period)
 
-    assert step.period_seconds(parameter_values, inputs) == 300.0
+    step.process_parameters(parameter_values)
+    assert step.period_seconds(inputs) == 300.0
     # The period sets the output sampling, so it must be resolved before t_interp
-    _, t_interp = step.setup_timestepping(
-        DummySolver(), 1800, parameter_values=parameter_values, inputs=inputs
-    )
+    _, t_interp = step.setup_timestepping(DummySolver(), 1800, inputs=inputs)
     np.testing.assert_array_equal(t_interp, np.arange(0, 1801, 300.0))
+
+
+def test_every_step_setting_is_processed_the_same_way():
+    parameter_values = pybamm.ParameterValues(
+        {"I": "[input]", "T": "[input]", "d": "[input]", "p": "[input]"}
+    )
+    step = pybamm.step.current(
+        pybamm.Parameter("I"),
+        temperature=pybamm.Parameter("T"),
+        duration=pybamm.Parameter("d"),
+        period=pybamm.Parameter("p"),
+    )
+    inputs = {"I": 2.0, "T": 300.0, "d": 1800.0, "p": 60.0}
+
+    step.process_parameters(parameter_values)
+
+    assert step.temperature_kelvin(inputs) == 300.0
+    assert step.duration_seconds(inputs) == 1800.0
+    assert step.period_seconds(inputs) == 60.0
+    # An input-dependent control target is not a constant, so it is not collapsible
+    assert step.has_symbolic_control_target is True
+    # A constant one is, once its parameters are substituted
+    c_rate = pybamm.step.c_rate(0.5, duration=1)
+    c_rate.process_parameters(pybamm.ParameterValues("Marquis2019"))
+    assert c_rate.control_target_value(inputs) == pytest.approx(0.5 * 0.680616)
+    # The user-facing attributes are left exactly as written
+    assert step.temperature == pybamm.Parameter("T")
+    assert step.duration == pybamm.Parameter("d")
diff --git a/packages/pybamm/tests/unit/test_experiments/test_experiment_step_termination.py b/packages/pybamm/tests/unit/test_experiments/test_experiment_step_termination.py
index 5ef54cff7e..2c9090e220 100644
--- a/packages/pybamm/tests/unit/test_experiments/test_experiment_step_termination.py
+++ b/packages/pybamm/tests/unit/test_experiments/test_experiment_step_termination.py
@@ -29,19 +29,19 @@ def test_c_rate_termination(self):
         term = pybamm.step.CRateTermination(0.02)
         assert term.value == 0.02
         assert term.operator is None
-        variables = {"C-rate": pybamm.Scalar(0.02)}
-        assert term.get_event(variables, None).evaluate() == 0
+        variables = {"C-rate": 0.02}
+        assert term.get_event(variables, None).expression == 0
         with pytest.warns(DeprecationWarning):
             term_old = pybamm.step.CrateTermination(0.02)
         assert (
-            term.get_event(variables, None).evaluate()
-            == term_old.get_event(variables, None).evaluate()
+            term.get_event(variables, None).expression
+            == term_old.get_event(variables, None).expression
         )
 
     def test_current_and_voltage_termination_operator_branches(self):
         variables = {
-            "Current [A]": pybamm.Scalar(0.5),
-            "Battery voltage [V]": pybamm.Scalar(4.2),
+            "Current [A]": 0.5,
+            "Battery voltage [V]": 4.2,
         }
         rest_step = SimpleNamespace(direction="rest")
 
@@ -49,23 +49,19 @@ def test_current_and_voltage_termination_operator_branches(self):
         current_lt = pybamm.step.CurrentTermination(0.6, operator="<")
         assert current_gt.get_event_name(None) == "Current [A] > 0.4 [A] [experiment]"
         assert current_lt.get_event_name(None) == "Current [A] < 0.6 [A] [experiment]"
-        assert current_gt.get_event_expression(
-            variables, None
-        ).evaluate() == pytest.approx(-0.1)
-        assert current_lt.get_event_expression(
-            variables, None
-        ).evaluate() == pytest.approx(-0.1)
+        assert current_gt.get_event_expression(variables, None) == pytest.approx(-0.1)
+        assert current_lt.get_event_expression(variables, None) == pytest.approx(-0.1)
 
         voltage_gt = pybamm.step.VoltageTermination(4.1, operator=">")
         voltage_lt = pybamm.step.VoltageTermination(4.3, operator="<")
         assert voltage_gt.get_event_name(rest_step) == "Voltage > 4.1 [V] [experiment]"
         assert voltage_lt.get_event_name(rest_step) == "Voltage < 4.3 [V] [experiment]"
-        assert voltage_gt.get_event_expression(
-            variables, rest_step
-        ).evaluate() == pytest.approx(-0.1)
-        assert voltage_lt.get_event_expression(
-            variables, rest_step
-        ).evaluate() == pytest.approx(-0.1)
+        assert voltage_gt.get_event_expression(variables, rest_step) == pytest.approx(
+            -0.1
+        )
+        assert voltage_lt.get_event_expression(variables, rest_step) == pytest.approx(
+            -0.1
+        )
 
         assert (pybamm.step.step_termination.Current() > 0.4) == current_gt
         assert (pybamm.step.step_termination.Current() < 0.6) == current_lt
@@ -75,7 +71,7 @@ def test_current_and_voltage_termination_operator_branches(self):
     def test_voltage_termination_returns_none_without_charge_or_discharge(self):
         term = pybamm.step.VoltageTermination(4.2)
         step = SimpleNamespace(direction="rest")
-        variables = {"Battery voltage [V]": pybamm.Scalar(4.2)}
+        variables = {"Battery voltage [V]": 4.2}
 
         assert term.get_event_name(step) is None
         assert term.get_event_expression(variables, step) is None
@@ -226,12 +222,13 @@ class TestInequalityTermination:
             ),
         ],
     )
-    def test_inequality_becomes_a_custom_termination(self, expression, residual):
+    def test_inequality_becomes_a_symbolic_termination(self, expression, residual):
         term = pybamm.step.base_step._read_termination(expression)
         v = pybamm.Variable("V")
 
-        assert isinstance(term, pybamm.step.CustomTermination)
-        assert term.name == f"{expression} [experiment]"
+        assert isinstance(term, pybamm.step.SymbolicTermination)
+        # The name is the inequality as written
+        assert term.get_event_name(None) == f"{expression} [experiment]"
         # A heaviside is "left < right", so left - right is positive before the
         # inequality holds and negative once it does: the event convention. The
         # CoupledVariable is looked up in the variables the termination is handed,
@@ -253,7 +250,7 @@ def test_inequality_termination_over_a_custom_variable(self):
 
         sol = sim.solve(calc_esoh=False)
 
-        assert sol.termination == f"event: {step.termination[0].name}"
+        assert sol.termination == f"event: {step.termination[0].get_event_name(step)}"
         assert sol["Headroom [V]"].data[-1] == pytest.approx(0.6, abs=1e-3)
 
     def test_inequality_termination_rejects_unknown_variable(self):
@@ -268,5 +265,5 @@ def test_inequality_termination_rejects_unknown_variable(self):
             sim.solve(calc_esoh=False)
 
     def test_symbolic_termination_must_be_an_inequality(self):
-        with pytest.raises(TypeError, match="must be an inequality between symbols"):
+        with pytest.raises(TypeError, match="must be an inequality"):
             pybamm.step.c_rate(1, duration=1, termination=pybamm.InputParameter("Vmin"))
diff --git a/packages/pybamm/tests/unit/test_experiments/test_simulation_with_experiment.py b/packages/pybamm/tests/unit/test_experiments/test_simulation_with_experiment.py
index ea45711ac0..a5e8e090d9 100644
--- a/packages/pybamm/tests/unit/test_experiments/test_simulation_with_experiment.py
+++ b/packages/pybamm/tests/unit/test_experiments/test_simulation_with_experiment.py
@@ -2635,7 +2635,9 @@ def test_run_experiment_with_symbolic_termination(self, experiment_model_mode):
         # The threshold is only read at solve time, so one built model serves both
         for cut_off in (3.6, 3.5):
             sol = sim.solve(inputs={"Voltage cut-off [V]": cut_off}, calc_esoh=False)
-            assert sol.termination == f"event: {step.termination[0].name}"
+            assert (
+                sol.termination == f"event: {step.termination[0].get_event_name(step)}"
+            )
             assert sol["Voltage [V]"].data[-1] == pytest.approx(cut_off, abs=1e-3)
             assert sol.t[-1] < 3600
 
diff --git a/packages/pybamm/tests/unit/test_serialisation/test_serialisation.py b/packages/pybamm/tests/unit/test_serialisation/test_serialisation.py
index 804e468bcb..475a383122 100644
--- a/packages/pybamm/tests/unit/test_serialisation/test_serialisation.py
+++ b/packages/pybamm/tests/unit/test_serialisation/test_serialisation.py
@@ -2962,6 +2962,26 @@ def test_symbolic_period_round_trip(self, period):
         # The step's own value must survive alongside it
         assert exp2.steps[0].value == 1.0
 
+    @pytest.mark.parametrize(
+        "termination",
+        [
+            # An inequality over a model variable and an input parameter
+            pybamm.CoupledVariable("Voltage [V]") > pybamm.InputParameter("V hold"),
+            # A named termination whose threshold is symbolic
+            pybamm.step.VoltageTermination(pybamm.Parameter("V cut"), operator="<"),
+            pybamm.step.CurrentTermination(0.05),
+        ],
+    )
+    def test_symbolic_termination_round_trip(self, termination):
+        exp = pybamm.Experiment(
+            [pybamm.step.current(1, duration=3600, termination=termination)]
+        )
+
+        config = json.loads(json.dumps(exp.to_config()))
+        exp2 = pybamm.Experiment.from_config(config)
+
+        assert exp2.steps[0].termination == exp.steps[0].termination
+
     def test_legacy_steps_format(self):
         """from_config also accepts flat {'steps': [...]} format."""
         config = {

From f490940131ec4c9258767ccd7ceac2b3559fcd19 Mon Sep 17 00:00:00 2001
From: Marc Berliner <34451391+MarcBerliner@users.noreply.github.com>
Date: Fri, 4 Sep 2026 12:57:35 -0400
Subject: [PATCH 3/6] docs: trim changelog entry

---
 CHANGELOG.md | 1 -
 1 file changed, 1 deletion(-)

diff --git a/CHANGELOG.md b/CHANGELOG.md
index b3298893dd..e143baa4d1 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -4,7 +4,6 @@
 
 - Step `duration`, `period` and `temperature` accept a `pybamm.Parameter` or `pybamm.InputParameter`.
 - Added `SymbolicTermination`: terminate a step on any inequality over model variables.
-- Voltage, current and C-rate terminations build on `SymbolicTermination` and serialise symbolic thresholds.
 
 ## Bug fixes
 

From 4c1acdac04e7aeb25512620b77516d894f2bd451 Mon Sep 17 00:00:00 2001
From: Marc Berliner <34451391+MarcBerliner@users.noreply.github.com>
Date: Fri, 4 Sep 2026 13:01:14 -0400
Subject: [PATCH 4/6] docs: changelog PR link

---
 CHANGELOG.md | 4 ++--
 1 file changed, 2 insertions(+), 2 deletions(-)

diff --git a/CHANGELOG.md b/CHANGELOG.md
index e143baa4d1..491d93ceb2 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -2,8 +2,8 @@
 
 ## Features
 
-- Step `duration`, `period` and `temperature` accept a `pybamm.Parameter` or `pybamm.InputParameter`.
-- Added `SymbolicTermination`: terminate a step on any inequality over model variables.
+- Step `duration`, `period` and `temperature` accept a `pybamm.Parameter` or `pybamm.InputParameter` ([#5744](https://github.com/pybamm-team/PyBaMM/pull/5744))
+- Added `SymbolicTermination`: terminate a step on any inequality over model variables ([#5744](https://github.com/pybamm-team/PyBaMM/pull/5744))
 
 ## Bug fixes
 

From 20eaa5b2048b235305a01c49548bc4f72dbb7038 Mon Sep 17 00:00:00 2001
From: Marc Berliner <34451391+MarcBerliner@users.noreply.github.com>
Date: Fri, 4 Sep 2026 13:18:50 -0400
Subject: [PATCH 5/6] fix: Python 3.10 enum membership; drop redundant event
 name override

`None in _Direction` raises TypeError on 3.10, so compare against the
values. CurrentTermination's get_event_name matched the inherited one.

Co-Authored-By: Claude Fable 5.1 <noreply@anthropic.com>
---
 packages/pybamm/src/pybamm/experiment/step/base_step.py     | 6 +++---
 .../pybamm/src/pybamm/experiment/step/step_termination.py   | 3 ---
 .../test_experiments/test_experiment_step_termination.py    | 4 ++--
 .../tests/unit/test_experiments/test_experiment_steps.py    | 2 +-
 .../test_experiments/test_simulation_with_experiment.py     | 6 ++----
 5 files changed, 8 insertions(+), 13 deletions(-)

diff --git a/packages/pybamm/src/pybamm/experiment/step/base_step.py b/packages/pybamm/src/pybamm/experiment/step/base_step.py
index 8117ae246f..f2b4f232c7 100644
--- a/packages/pybamm/src/pybamm/experiment/step/base_step.py
+++ b/packages/pybamm/src/pybamm/experiment/step/base_step.py
@@ -112,10 +112,10 @@ def __init__(
     ):
         # Filled by `process_parameters` at setup, see `evaluate`
         self._processed_variables = {}
-        if direction not in _Direction and direction is not None:
-            _DIRECTIONS = [d.value for d in _Direction] + [None]
+        directions = [d.value for d in _Direction] + [None]
+        if direction not in directions:
             raise ValueError(
-                f"Invalid direction: {direction}. Must be one of {_DIRECTIONS}"
+                f"Invalid direction: {direction}. Must be one of {directions}"
             )
         self.input_duration = duration
         self.input_value = value
diff --git a/packages/pybamm/src/pybamm/experiment/step/step_termination.py b/packages/pybamm/src/pybamm/experiment/step/step_termination.py
index 1e845a3865..6ff397ef74 100644
--- a/packages/pybamm/src/pybamm/experiment/step/step_termination.py
+++ b/packages/pybamm/src/pybamm/experiment/step/step_termination.py
@@ -134,9 +134,6 @@ def inequality(self, step):
             return current < self.value
         return abs(current) < self.value
 
-    def get_event_name(self, step):
-        return str(self.inequality(step)) + " [experiment]"
-
 
 class VoltageTermination(SymbolicTermination):
     """
diff --git a/packages/pybamm/tests/unit/test_experiments/test_experiment_step_termination.py b/packages/pybamm/tests/unit/test_experiments/test_experiment_step_termination.py
index 2c9090e220..8c6de592ed 100644
--- a/packages/pybamm/tests/unit/test_experiments/test_experiment_step_termination.py
+++ b/packages/pybamm/tests/unit/test_experiments/test_experiment_step_termination.py
@@ -47,8 +47,8 @@ def test_current_and_voltage_termination_operator_branches(self):
 
         current_gt = pybamm.step.CurrentTermination(0.4, operator=">")
         current_lt = pybamm.step.CurrentTermination(0.6, operator="<")
-        assert current_gt.get_event_name(None) == "Current [A] > 0.4 [A] [experiment]"
-        assert current_lt.get_event_name(None) == "Current [A] < 0.6 [A] [experiment]"
+        assert current_gt.get_event_name(None) == "0.4 < Current [A] [experiment]"
+        assert current_lt.get_event_name(None) == "Current [A] < 0.6 [experiment]"
         assert current_gt.get_event_expression(variables, None) == pytest.approx(-0.1)
         assert current_lt.get_event_expression(variables, None) == pytest.approx(-0.1)
 
diff --git a/packages/pybamm/tests/unit/test_experiments/test_experiment_steps.py b/packages/pybamm/tests/unit/test_experiments/test_experiment_steps.py
index 7a62b83827..720bfba220 100644
--- a/packages/pybamm/tests/unit/test_experiments/test_experiment_steps.py
+++ b/packages/pybamm/tests/unit/test_experiments/test_experiment_steps.py
@@ -395,7 +395,7 @@ def test_symbolic_termination_expression_helpers(self):
         events = step.get_termination_events(variables)
         assert [event.name for event in events] == [
             "Voltage < 2.5 [V] [experiment]",
-            "abs(Current [A]) < 0.05 [A] [experiment]",
+            "abs(Current [A]) < 0.05 [experiment]",
         ]
         np.testing.assert_allclose(
             step.get_combined_termination_expression(variables).evaluate(),
diff --git a/packages/pybamm/tests/unit/test_experiments/test_simulation_with_experiment.py b/packages/pybamm/tests/unit/test_experiments/test_simulation_with_experiment.py
index a5e8e090d9..a04bf92192 100644
--- a/packages/pybamm/tests/unit/test_experiments/test_simulation_with_experiment.py
+++ b/packages/pybamm/tests/unit/test_experiments/test_simulation_with_experiment.py
@@ -1123,9 +1123,7 @@ def test_run_multi_termination_step_unified_matches_legacy(self):
         legacy_hold = legacy_sol.cycles[0].steps[1]
         unified_hold = unified_sol.cycles[0].steps[1]
 
-        assert (
-            legacy_hold.termination == "event: abs(Current [A]) < 0.5 [A] [experiment]"
-        )
+        assert legacy_hold.termination == "event: abs(Current [A]) < 0.5 [experiment]"
         assert legacy_hold.termination == unified_hold.termination
         np.testing.assert_allclose(
             legacy_hold.t[-1], unified_hold.t[-1], rtol=5e-5, atol=5e-4
@@ -1242,7 +1240,7 @@ def test_skip_ok_with_multiple_infeasible_terminations_in_unified_model(self):
         assert len(sol.cycles[0].steps) == 1
         assert (
             sol.cycles[0].steps[0].termination
-            == "event: abs(Current [A]) < 0.01 [A] [experiment]"
+            == "event: abs(Current [A]) < 0.01 [experiment]"
         )
 
     def test_all_empty_solution_errors(self):

From 3e3304c64b3f3c76d7d8128c20b7a98761eb72b1 Mon Sep 17 00:00:00 2001
From: Marc Berliner <34451391+MarcBerliner@users.noreply.github.com>
Date: Fri, 4 Sep 2026 13:22:22 -0400
Subject: [PATCH 6/6] Update base_step.py

---
 packages/pybamm/src/pybamm/experiment/step/base_step.py | 8 +++++---
 1 file changed, 5 insertions(+), 3 deletions(-)

diff --git a/packages/pybamm/src/pybamm/experiment/step/base_step.py b/packages/pybamm/src/pybamm/experiment/step/base_step.py
index f2b4f232c7..de078fd85a 100644
--- a/packages/pybamm/src/pybamm/experiment/step/base_step.py
+++ b/packages/pybamm/src/pybamm/experiment/step/base_step.py
@@ -47,6 +47,9 @@ class _Direction(str, Enum):
     REST = "rest"
 
 
+_DIRECTIONS = frozenset([d.value for d in _Direction] + [None])
+
+
 class _SymbolicInput(str, Enum):
     """Step inputs which may be a symbolic function of parameters or real numbers."""
 
@@ -112,10 +115,9 @@ def __init__(
     ):
         # Filled by `process_parameters` at setup, see `evaluate`
         self._processed_variables = {}
-        directions = [d.value for d in _Direction] + [None]
-        if direction not in directions:
+        if direction not in _DIRECTIONS:
             raise ValueError(
-                f"Invalid direction: {direction}. Must be one of {directions}"
+                f"Invalid direction: {direction}. Must be one of {list(_DIRECTIONS)}"
             )
         self.input_duration = duration
         self.input_value = value