Repository: pybamm-team/PyBaMM
Pull request #4993 (2025-04-24 20:21:48)

Add a solver option to change `on_extrapolation` behavior
# Description

Let's you enable `"error"`, `"warn"`, or `"ignore"` on extrapolation events.

## Type of change

Please add a line in the relevant section of [CHANGELOG.md](https://github.com/pybamm-team/PyBaMM/blob/develop/CHANGELOG.md) to document the change (include PR #)

# Important checks:

Please confirm the following before marking the PR as ready for review:
- No style issues: `nox -s pre-commit`
- All tests pass: `nox -s tests`
- The documentation builds: `nox -s doctests`
- Code is commented for hard-to-understand areas
- Tests added that prove fix is effective or that feature works


codecov[bot]: ## [Codecov](https://app.codecov.io/gh/pybamm-team/PyBaMM/pull/4993?dropdown=coverage&src=pr&el=h1&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=pybamm-team) Report
All modified and coverable lines are covered by tests :white_check_mark:
> Project coverage is 98.57%. Comparing base [(`b4c9df8`)](https://app.codecov.io/gh/pybamm-team/PyBaMM/commit/b4c9df8bc210d6d68bf53d38939864d5f856539d?dropdown=coverage&el=desc&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=pybamm-team) to head [(`88e0168`)](https://app.codecov.io/gh/pybamm-team/PyBaMM/commit/88e01683563d3629810f128ce0eb33944adc9735?dropdown=coverage&el=desc&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=pybamm-team).
> Report is 2 commits behind head on develop.

<details><summary>Additional details and impacted files</summary>


```diff
@@             Coverage Diff             @@
##           develop    #4993      +/-   ##
===========================================
- Coverage    98.60%   98.57%   -0.04%     
===========================================
  Files          304      304              
  Lines        23809    23652     -157     
===========================================
- Hits         23478    23316     -162     
- Misses         331      336       +5     
```

</details>

[:umbrella: View full report in Codecov by Sentry](https://app.codecov.io/gh/pybamm-team/PyBaMM/pull/4993?dropdown=coverage&src=pr&el=continue&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=pybamm-team).   
:loudspeaker: Have feedback on the report? [Share it here](https://about.codecov.io/codecov-pr-comment-feedback/?utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=pybamm-team).

<details><summary> :rocket: New features to boost your workflow: </summary>

- :snowflake: [Test Analytics](https://docs.codecov.com/docs/test-analytics): Detect flaky tests, report on failures, and find test suite problems.
</details>

Patch:
diff --git a/CHANGELOG.md b/CHANGELOG.md
index 406cfc4ed..57f32cc96 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -2,6 +2,7 @@
 
 ## Optimizations
 
+- Add a solver option to change `on_extrapolation` behavior to `"error"`, `"warn"`, or `"ignore"` on extrapolation events. ([#4993](https://github.com/pybamm-team/PyBaMM/pull/4993))
 - Improve reliability of `CasadiAlgebraicSolver` and added an option for the `step_tol` of the Newton iteration. ([#4985](https://github.com/pybamm-team/PyBaMM/pull/4985))
 
 ## Bug fixes
diff --git a/src/pybamm/solvers/base_solver.py b/src/pybamm/solvers/base_solver.py
index 6f13db2b6..1539252c0 100644
--- a/src/pybamm/solvers/base_solver.py
+++ b/src/pybamm/solvers/base_solver.py
@@ -39,6 +39,9 @@ class BaseSolver:
     output_variables : list[str], optional
         List of variables to calculate and return. If none are specified then
         the complete state vector is returned (can be very large) (default is [])
+    on_extrapolation : str, optional
+        What to do if the solver is extrapolating. Options are "warn", "error", or "ignore".
+        Default is "warn".
     """
 
     def __init__(
@@ -49,6 +52,7 @@ class BaseSolver:
         root_method=None,
         root_tol=1e-6,
         extrap_tol=None,
+        on_extrapolation=None,
         output_variables=None,
     ):
         self.method = method
@@ -58,6 +62,7 @@ class BaseSolver:
         self.root_method = root_method
         self.extrap_tol = extrap_tol or -1e-10
         self.output_variables = [] if output_variables is None else output_variables
+        self.on_extrapolation = on_extrapolation or "warn"
         self._model_set_up = {}
 
         # Defaults, can be overwritten by specific solver
@@ -65,7 +70,6 @@ class BaseSolver:
         self._ode_solver = False
         self._algebraic_solver = False
         self._supports_interp = False
-        self._on_extrapolation = "warn"
         self.computed_var_fcns = {}
         self._mp_context = self.get_platform_context(platform.system())
 
@@ -81,14 +85,24 @@ class BaseSolver:
     def supports_interp(self):
         return self._supports_interp
 
-    @property
-    def root_method(self):
-        return self._root_method
-
     @property
     def supports_parallel_solve(self):
         return False
 
+    @property
+    def on_extrapolation(self):
+        return self._on_extrapolation
+
+    @on_extrapolation.setter
+    def on_extrapolation(self, value):
+        if value not in ["warn", "error", "ignore"]:
+            raise ValueError("on_extrapolation must be 'warn', 'raise', or 'ignore'")
+        self._on_extrapolation = value
+
+    @property
+    def root_method(self):
+        return self._root_method
+
     @root_method.setter
     def root_method(self, method):
         if method == "casadi":
@@ -1386,6 +1400,9 @@ class BaseSolver:
         events : dict
             Dictionary of events
         """
+        if self.on_extrapolation == "ignore":
+            return
+
         extrap_events = []
 
         # Add the event dictionary to the solution object
@@ -1423,14 +1440,14 @@ class BaseSolver:
             # no extrapolation events are within the tolerance
             return
 
-        if self._on_extrapolation == "error":
+        if self.on_extrapolation == "error":
             raise pybamm.SolverError(
                 "Solver failed because the following "
                 f"interpolation bounds were exceeded: {extrap_events}. "
                 "You may need to provide additional interpolation points "
                 "outside these bounds."
             )
-        elif self._on_extrapolation == "warn":
+        elif self.on_extrapolation == "warn":
             name = solution.all_models[-1].name
             warnings.warn(
                 f"While solving {name} extrapolation occurred for {extrap_events}",
diff --git a/src/pybamm/solvers/casadi_solver.py b/src/pybamm/solvers/casadi_solver.py
index 745ab12a2..880033341 100644
--- a/src/pybamm/solvers/casadi_solver.py
+++ b/src/pybamm/solvers/casadi_solver.py
@@ -46,6 +46,9 @@ class CasadiSolver(pybamm.BaseSolver):
         the default value is 600 seconds.
     extrap_tol : float, optional
         The tolerance to assert whether extrapolation occurs or not. Default is 0.
+    on_extrapolation : str, optional
+        What to do if the solver is extrapolating. Options are "warn", "error", or "ignore".
+        Default is "error".
     extra_options_setup : dict, optional
         Any options to pass to the CasADi integrator when creating the integrator.
         Please consult `CasADi documentation <https://web.casadi.org/python-api/#integrator>`_ for
@@ -81,19 +84,22 @@ class CasadiSolver(pybamm.BaseSolver):
         max_step_decrease_count=5,
         dt_max=None,
         extrap_tol=None,
+        on_extrapolation=None,
         extra_options_setup=None,
         extra_options_call=None,
         return_solution_if_failed_early=False,
         perturb_algebraic_initial_conditions=None,
         integrators_maxcount=100,
     ):
+        on_extrapolation = on_extrapolation or "error"
         super().__init__(
-            "problem dependent",
-            rtol,
-            atol,
-            root_method,
-            root_tol,
-            extrap_tol,
+            method="problem dependent",
+            rtol=rtol,
+            atol=atol,
+            root_method=root_method,
+            root_tol=root_tol,
+            extrap_tol=extrap_tol,
+            on_extrapolation=on_extrapolation,
         )
         if mode in ["safe", "fast", "fast with events", "safe without grid"]:
             self.mode = mode
@@ -110,8 +116,6 @@ class CasadiSolver(pybamm.BaseSolver):
         self.extra_options_call = extra_options_call or {}
         self.return_solution_if_failed_early = return_solution_if_failed_early
 
-        self._on_extrapolation = "error"
-
         # Decide whether to perturb algebraic initial conditions, True by default for
         # "safe" mode, False by default for other modes
         if perturb_algebraic_initial_conditions is None:
diff --git a/src/pybamm/solvers/idaklu_solver.py b/src/pybamm/solvers/idaklu_solver.py
index 96e2519fb..15f2e1b7a 100644
--- a/src/pybamm/solvers/idaklu_solver.py
+++ b/src/pybamm/solvers/idaklu_solver.py
@@ -51,6 +51,9 @@ class IDAKLUSolver(pybamm.BaseSolver):
         The tolerance for the initial-condition solver (default is 1e-6).
     extrap_tol : float, optional
         The tolerance to assert whether extrapolation occurs or not (default is 0).
+    on_extrapolation : str, optional
+        What to do if the solver is extrapolating. Options are "warn", "error", or "ignore".
+        Default is "warn".
     output_variables : list[str], optional
         List of variables to calculate and return. If none are specified then
         the complete state vector is returned (can be very large) (default is [])
@@ -163,6 +166,7 @@ class IDAKLUSolver(pybamm.BaseSolver):
         root_method="casadi",
         root_tol=1e-6,
         extrap_tol=None,
+        on_extrapolation=None,
         output_variables=None,
         options=None,
     ):
@@ -219,13 +223,14 @@ class IDAKLUSolver(pybamm.BaseSolver):
         self.output_variables = [] if output_variables is None else output_variables
 
         super().__init__(
-            "ida",
-            rtol,
-            atol,
-            root_method,
-            root_tol,
-            extrap_tol,
-            output_variables,
+            method="ida",
+            rtol=rtol,
+            atol=atol,
+            root_method=root_method,
+            root_tol=root_tol,
+            extrap_tol=extrap_tol,
+            output_variables=output_variables,
+            on_extrapolation=on_extrapolation,
         )
         self.name = "IDA KLU solver"
         self._supports_interp = True
diff --git a/src/pybamm/solvers/jax_solver.py b/src/pybamm/solvers/jax_solver.py
index d8f2424ee..beb9581f9 100644
--- a/src/pybamm/solvers/jax_solver.py
+++ b/src/pybamm/solvers/jax_solver.py
@@ -39,6 +39,9 @@ class JaxSolver(pybamm.BaseSolver):
         The absolute tolerance for the solver (default is 1e-6).
     extrap_tol : float, optional
         The tolerance to assert whether extrapolation occurs or not (default is 0).
+    on_extrapolation : str, optional
+        What to do if the solver is extrapolating. Options are "warn", "error", or "ignore".
+        Default is "warn".
     extra_options : dict, optional
         Any options to pass to the solver.
         Please consult `JAX documentation
@@ -53,6 +56,7 @@ class JaxSolver(pybamm.BaseSolver):
         rtol=1e-6,
         atol=1e-6,
         extrap_tol=None,
+        on_extrapolation=None,
         extra_options=None,
     ):
         if not pybamm.has_jax():
@@ -63,7 +67,12 @@ class JaxSolver(pybamm.BaseSolver):
         # note: bdf solver itself calculates consistent initial conditions so can set
         # root_method to none, allow user to override this behavior
         super().__init__(
-            method, rtol, atol, root_method=root_method, extrap_tol=extrap_tol
+            method=method,
+            rtol=rtol,
+            atol=atol,
+            root_method=root_method,
+            extrap_tol=extrap_tol,
+            on_extrapolation=on_extrapolation,
         )
         method_options = ["RK45", "BDF"]
         if method not in method_options:
diff --git a/src/pybamm/solvers/scipy_solver.py b/src/pybamm/solvers/scipy_solver.py
index 5d1ca23d2..46eccfc1b 100644
--- a/src/pybamm/solvers/scipy_solver.py
+++ b/src/pybamm/solvers/scipy_solver.py
@@ -21,6 +21,9 @@ class ScipySolver(pybamm.BaseSolver):
         The absolute tolerance for the solver (default is 1e-6).
     extrap_tol : float, optional
         The tolerance to assert whether extrapolation occurs or not (default is 0).
+    on_extrapolation : str, optional
+        What to do if the solver is extrapolating. Options are "warn", "error", or "ignore".
+        Default is "warn".
     extra_options : dict, optional
         Any options to pass to the solver.
         Please consult `SciPy documentation
@@ -34,6 +37,7 @@ class ScipySolver(pybamm.BaseSolver):
         rtol=1e-6,
         atol=1e-6,
         extrap_tol=None,
+        on_extrapolation=None,
         extra_options=None,
     ):
         super().__init__(
@@ -41,6 +45,7 @@ class ScipySolver(pybamm.BaseSolver):
             rtol=rtol,
             atol=atol,
             extrap_tol=extrap_tol,
+            on_extrapolation=on_extrapolation,
         )
         self._ode_solver = True
         self.extra_options = extra_options or {}


Test patch:
diff --git a/tests/unit/test_solvers/test_base_solver.py b/tests/unit/test_solvers/test_base_solver.py
index 543e27b42..612e55218 100644
--- a/tests/unit/test_solvers/test_base_solver.py
+++ b/tests/unit/test_solvers/test_base_solver.py
@@ -414,3 +414,21 @@ class TestBaseSolver:
             np.testing.assert_allclose(
                 sens_b, exact_diff_b(y, inputs["a"], inputs["b"])
             )
+
+    def test_on_extrapolation_settings(self):
+        # Test setting different on_extrapolation values on BaseSolver
+        base_solver = pybamm.BaseSolver()
+
+        # Test valid values
+        base_solver.on_extrapolation = "warn"
+        assert base_solver.on_extrapolation == "warn"
+        base_solver.on_extrapolation = "error"
+        assert base_solver.on_extrapolation == "error"
+        base_solver.on_extrapolation = "ignore"
+        assert base_solver.on_extrapolation == "ignore"
+
+        # Test invalid value
+        with pytest.raises(
+            ValueError, match="on_extrapolation must be 'warn', 'raise', or 'ignore'"
+        ):
+            base_solver.on_extrapolation = "invalid"
diff --git a/tests/unit/test_solvers/test_idaklu_solver.py b/tests/unit/test_solvers/test_idaklu_solver.py
index 05ed69217..c57735288 100644
--- a/tests/unit/test_solvers/test_idaklu_solver.py
+++ b/tests/unit/test_solvers/test_idaklu_solver.py
@@ -2,7 +2,7 @@ from contextlib import redirect_stdout
 import io
 import pytest
 import numpy as np
-
+import warnings
 import pybamm
 from scipy.sparse import eye
 from tests import get_discretisation_for_testing
@@ -1298,3 +1298,47 @@ class TestIDAKLUSolver:
         np.testing.assert_array_equal(solution.t, t_interp)
         np.testing.assert_allclose(solution.y[0], np.exp(0.1 * solution.t))
         np.testing.assert_allclose(solution.y[-1], 2 * np.exp(0.1 * solution.t))
+
+    def test_interpolant_extrapolate(self):
+        x = np.linspace(0, 2)
+        var = pybamm.Variable("var")
+        rhs = pybamm.FunctionParameter("func", {"var": var})
+
+        model = pybamm.BaseModel()
+        model.rhs[var] = rhs
+        model.initial_conditions[var] = pybamm.Scalar(1)
+
+        # Bug: we need to set the interpolant via parameter values for the extrapolation
+        # to be detected
+        def func(var):
+            return pybamm.Interpolant(x, x, var, interpolator="linear")
+
+        parameter_values = pybamm.ParameterValues({"func": func})
+        parameter_values.process_model(model)
+
+        # Test with on_extrapolation="error"
+        solver = pybamm.IDAKLUSolver(on_extrapolation="error")
+        t_eval = [0, 5]
+
+        with pytest.raises(pybamm.SolverError, match="interpolation bounds"):
+            solver.solve(model, t_eval)
+
+        # Test with on_extrapolation="warn"
+        solver = pybamm.IDAKLUSolver(on_extrapolation="warn")
+        t_eval = [0, 5]
+
+        with warnings.catch_warnings(record=True) as w:
+            warnings.simplefilter("always")
+            solver.solve(model, t_eval)
+            assert len(w) > 0
+            assert "extrapolation occurred" in str(w[0].message)
+
+        # Test with on_extrapolation="ignore"
+        solver = pybamm.IDAKLUSolver(on_extrapolation="ignore")
+        t_eval = [0, 5]
+
+        # Should not raise an error or warning
+        with warnings.catch_warnings(record=True) as w:
+            warnings.simplefilter("always")
+            solver.solve(model, t_eval)
+            assert len(w) == 0