From 08f1040c897013fc11655eb73137da7cc9aadce1 Mon Sep 17 00:00:00 2001 From: Marc Berliner <34451391+MarcBerliner@users.noreply.github.com> Date: Thu, 20 Aug 2026 09:40:47 -0400 Subject: [PATCH 1/3] feat: solve the composite electrode SOH by bracketed rootfind Each stoichiometry was a variable in one nine-variable algebraic system handed to a Newton solver, which is a basin-of-attraction problem: worn and lithium-poor states sat outside it. Now each is an expression, found by inverting the open-circuit potential it sits on with a `pybamm.Brent`. The OCPs diverge outside [0, 1], so every inversion brackets a root and converges from any target, including targets outside the physical range, whose answer is non-physical but exact. Nothing is left for a solver to iterate, so the stoichiometries are read back through one CasADi function rather than per-variable post-processing. The capacities reach that function as inputs, so only the parameters its graph substitutes can invalidate it: a caller ageing a cell reuses it, and a caller swapping an open-circuit potential does not. Over 8000 targets -- a thousand states of charge and a thousand voltages, three hysteresis directions, four wear states: main branch failures 147/8000 0/8000 build 12.2 ms 18.5 ms first solve 58.9 ms 14.8 ms Main's failures are silent: it returns a state violating lithium conservation rather than raising. Repeat calls through a reused simulation, reading the stoichiometries back: main branch new target 18.8 ms 7.2 ms aged capacities 15.3 ms 15.6 ms no reuse 75.3 ms 91.2 ms Reading the nine variables back through ProcessedVariable is still slower than main, which the evaluator bypasses only for get_initial_stoichiometries_composite. Co-Authored-By: Claude Opus 5 --- CHANGELOG.md | 1 + .../lithium_ion/electrode_soh_composite.py | 637 ++++++++++-------- .../test_electrode_soh_composite_sweeps.py | 258 +++++++ .../test_lithium_ion/test_electrode_soh.py | 197 ++++++ 4 files changed, 803 insertions(+), 290 deletions(-) create mode 100644 packages/pybamm/tests/integration/test_models/test_full_battery_models/test_lithium_ion/test_electrode_soh_composite_sweeps.py diff --git a/CHANGELOG.md b/CHANGELOG.md index a17810776c..b8d3cd75e1 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -2,6 +2,7 @@ ## Features +- `ElectrodeSOHComposite` now solves each stoichiometry by bracketed rootfind, converging from any target. ([#5730](https://github.com/pybamm-team/PyBaMM/pull/5730)) - 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/models/full_battery_models/lithium_ion/electrode_soh_composite.py b/packages/pybamm/src/pybamm/models/full_battery_models/lithium_ion/electrode_soh_composite.py index 2eb9d7fa9b..2ed1d9cd1b 100644 --- a/packages/pybamm/src/pybamm/models/full_battery_models/lithium_ion/electrode_soh_composite.py +++ b/packages/pybamm/src/pybamm/models/full_battery_models/lithium_ion/electrode_soh_composite.py @@ -6,6 +6,9 @@ import warnings from typing import Any +import casadi +import numpy as np + import pybamm from .electrode_soh import _ElectrodeSOH, get_esoh_default_solver @@ -15,6 +18,14 @@ get_lithiation_delithiation, ) +# `LithiumIonParameters.U` clips its stoichiometry and adds an asymptote, so every OCP +# diverges outside [0, 1] and both brackets straddle a root. The inner one is wider so +# an inversion can always reach the potential the outer solve asks it for. +_STOICH_LO, _STOICH_HI = -10.0, 10.0 +_BRACKET_LO, _BRACKET_HI = -5.0, 5.0 + +_ABSTOL, _MAX_ITER = 1e-14, 200 + def _get_primary_only_options( options: dict | pybamm.BatteryModelOptions | None, @@ -74,6 +85,19 @@ def _get_primary_only_options( return options_dict +def _initialization_method(initial_value: float | str) -> str: + """``"voltage"`` for a string ending in ``"V"``, ``"SOC"`` for a float.""" + if isinstance(initial_value, str) and initial_value.endswith("V"): + return "voltage" + if isinstance(initial_value, float): + return "SOC" + raise ValueError( + "Invalid initial value. Expected a float between 0 and 1 " + "(for SOC) or a string ending in 'V' (for voltage), got " + f"{initial_value!r} of type {type(initial_value).__name__}" + ) + + def _get_stoich_variables(options): """Create stoichiometry variables for composite electrodes.""" variables = { @@ -300,234 +324,195 @@ def __init__( pybamm.citations.register("Mohtat2019") super().__init__(name) param = pybamm.LithiumIonParameters(options) + neg_composite = check_if_composite(options, "negative") + pos_composite = check_if_composite(options, "positive") Q_Li = pybamm.InputParameter("Q_Li") - is_negative_composite = check_if_composite(options, "negative") - is_positive_composite = check_if_composite(options, "positive") - variables = _get_stoich_variables(options) - x_100_1 = variables["x_100_1"] - y_100_1 = variables["y_100_1"] - x_0_1 = variables["x_0_1"] - y_0_1 = variables["y_0_1"] - V_max = param.ocp_soc_100 - V_min = param.ocp_soc_0 - if is_negative_composite: - x_100_2 = variables["x_100_2"] - x_0_2 = variables["x_0_2"] - self.algebraic[x_100_2] = param.n.sec.U( - x_100_2, - param.T_ref, - get_lithiation_delithiation( - get_equilibrium_direction("100", "negative", options, "secondary"), - "negative", - options, - phase="secondary", - ), - ) - param.n.prim.U( - x_100_1, - param.T_ref, - get_lithiation_delithiation( - get_equilibrium_direction("100", "negative", options, "primary"), - "negative", - options, - phase="primary", - ), - ) - self.algebraic[x_0_2] = param.n.sec.U( - x_0_2, - param.T_ref, - get_lithiation_delithiation( - get_equilibrium_direction("0", "negative", options, "secondary"), - "negative", - options, - phase="secondary", - ), - ) - param.n.prim.U( - x_0_1, - param.T_ref, - get_lithiation_delithiation( - get_equilibrium_direction("0", "negative", options, "primary"), - "negative", - options, - phase="primary", - ), - ) - if is_positive_composite: - y_100_2 = variables["y_100_2"] - y_0_2 = variables["y_0_2"] - self.algebraic[y_100_2] = param.p.sec.U( - y_100_2, - param.T_ref, - get_lithiation_delithiation( - get_equilibrium_direction("100", "positive", options, "secondary"), - "positive", - options, - phase="secondary", - ), - ) - param.p.prim.U( - y_100_1, - param.T_ref, - get_lithiation_delithiation( - get_equilibrium_direction("100", "positive", options, "primary"), - "positive", - options, - phase="primary", - ), - ) - self.algebraic[y_0_2] = param.p.prim.U( - y_0_1, - param.T_ref, - get_lithiation_delithiation( - get_equilibrium_direction("0", "positive", options, "primary"), - "positive", - options, - phase="primary", - ), - ) - param.p.sec.U( - y_0_2, - param.T_ref, - get_lithiation_delithiation( - get_equilibrium_direction("0", "positive", options, "secondary"), - "positive", - options, - phase="secondary", - ), - ) - self.algebraic[x_100_1] = ( - param.p.prim.U( - y_100_1, - param.T_ref, - get_lithiation_delithiation( - get_equilibrium_direction("100", "positive", options, "primary"), - "positive", - options, - phase="primary", - ), + Q_n = [pybamm.InputParameter("Q_n_1")] + Q_p = [pybamm.InputParameter("Q_p_1")] + if neg_composite: + Q_n.append(pybamm.InputParameter("Q_n_2")) + if pos_composite: + Q_p.append(pybamm.InputParameter("Q_p_2")) + + T_ref = param.T_ref + T_init = param.T_init if initialization_method == "voltage" else T_ref + + def ocps(electrode, soc): + """(open-circuit potential, branch, name) per phase, at an SOC level.""" + side = param.n if electrode == "negative" else param.p + phases = ["primary"] + if check_if_composite(options, electrode): + phases.append("secondary") + out = [] + for phase in phases: + U = (side.prim if phase == "primary" else side.sec).U + if soc == "init": + branch = get_lithiation_delithiation( + direction, electrode, options, phase=phase + ) + else: + branch = get_lithiation_delithiation( + get_equilibrium_direction(soc, electrode, options, phase), + electrode, + options, + phase=phase, + ) + out.append((U, branch, f"{electrode} {phase} {soc}")) + return out + + def invert(U, branch, target, T, name): + """Stoichiometry at a given potential, by Brent.""" + sto = pybamm._BrentUnknown(name) + return pybamm._Brent( + U(sto, T, branch) - target, + sto, + (_STOICH_LO, _STOICH_HI), + abstol=_ABSTOL, + max_iter=_MAX_ITER, ) - - param.n.prim.U( - x_100_1, - param.T_ref, - get_lithiation_delithiation( - get_equilibrium_direction("100", "negative", options, "primary"), - "negative", - options, - phase="primary", - ), + + def bracket(pairs, offsets, T): + """Potentials for which every inversion in `pairs` brackets.""" + los = [ + U(pybamm.Scalar(_BRACKET_HI), T, b) - off + for (U, b, _), off in zip(pairs, offsets, strict=True) + ] + his = [ + U(pybamm.Scalar(_BRACKET_LO), T, b) - off + for (U, b, _), off in zip(pairs, offsets, strict=True) + ] + lo, hi = los[0], his[0] + for value in los[1:]: + lo = pybamm.maximum(lo, value) + for value in his[1:]: + hi = pybamm.minimum(hi, value) + return lo, hi + + def limits(soc, V_target, closure): + """A limit state: one unknown, the shared negative potential.""" + neg, pos = ocps("negative", soc), ocps("positive", soc) + + # built twice: once against the unknown potential, once against the solved + def states(U_n, suffix): + return ( + [invert(U, b, U_n, T_ref, f"{n} {suffix}") for U, b, n in neg], + [ + invert(U, b, U_n + V_target, T_ref, f"{n} {suffix}") + for U, b, n in pos + ], + ) + + U_n = pybamm._BrentUnknown(f"negative potential {soc}") + lo, hi = bracket(neg + pos, [0] * len(neg) + [V_target] * len(pos), T_ref) + solved = pybamm._Brent( + closure(*states(U_n, "guess")), + U_n, + (lo, hi), + abstol=_ABSTOL, + max_iter=_MAX_ITER, ) - - V_max + return states(solved, "solved") + + def total_lithium(x, y): + return sum(q * s for q, s in zip(Q_n + Q_p, x + y, strict=True)) + + x_100, y_100 = limits( + "100", param.ocp_soc_100, lambda x, y: Q_Li - total_lithium(x, y) ) - self.algebraic[x_0_1] = ( - param.p.prim.U( - y_0_1, - param.T_ref, - get_lithiation_delithiation( - get_equilibrium_direction("0", "positive", options, "primary"), - "positive", - options, - phase="primary", - ), - ) - - param.n.prim.U( - x_0_1, - param.T_ref, - get_lithiation_delithiation( - get_equilibrium_direction("0", "negative", options, "primary"), - "negative", - options, - phase="primary", - ), - ) - - V_min + x_0, y_0 = limits( + "0", + param.ocp_soc_0, + lambda x, y: ( + -sum(q * (a - b) for q, a, b in zip(Q_p, y_100, y, strict=True)) + - sum(q * (a - b) for q, a, b in zip(Q_n, x_100, x, strict=True)) + ), ) - self.algebraic[y_0_1] = _get_electrode_capacity_equation( - options, "positive" - ) - _get_electrode_capacity_equation(options, "negative") - self.algebraic[y_100_1] = Q_Li - _get_cyclable_lithium_equation(options) - x_init_1 = variables["x_init_1"] - y_init_1 = variables["y_init_1"] - if initialization_method == "voltage": - V_init = pybamm.InputParameter("V_init") - self.algebraic[x_init_1] = ( - param.p.prim.U( - y_init_1, - param.T_init, - get_lithiation_delithiation( - direction, "positive", options, phase="primary" - ), - ) - - param.n.prim.U( - x_init_1, - param.T_init, - get_lithiation_delithiation( - direction, "negative", options, phase="primary" - ), + neg_init, pos_init = ocps("negative", "init"), ocps("positive", "init") + + def init_states(x_init_1, suffix): + U_n = neg_init[0][0](x_init_1, T_init, neg_init[0][1]) + x = [x_init_1] + [ + invert(U, b, U_n, T_init, f"{n} {suffix}") for U, b, n in neg_init[1:] + ] + if initialization_method == "voltage": + U_p = U_n + pybamm.InputParameter("V_init") + return x, [ + invert(U, b, U_p, T_init, f"{n} {suffix}") for U, b, n in pos_init + ] + # No cell-voltage relation: lithium conservation fixes the positive side. + remaining = Q_Li - sum(q * s for q, s in zip(Q_n, x, strict=True)) + if not pos_composite: + return x, [remaining / Q_p[0]] + U_p = pybamm._BrentUnknown(f"positive potential {suffix}") + lo, hi = bracket(pos_init, [0] * len(pos_init), T_init) + solved = pybamm._Brent( + sum( + q * invert(U, b, U_p, T_init, f"{n} {suffix} guess") + for q, (U, b, n) in zip(Q_p, pos_init, strict=True) ) - - V_init - ) - self.algebraic[y_init_1] = ( - _get_cyclable_lithium_equation(options, "init") - Q_Li + - remaining, + U_p, + (lo, hi), + abstol=_ABSTOL, + max_iter=_MAX_ITER, ) + return x, [ + invert(U, b, solved, T_init, f"{n} {suffix} solved") + for U, b, n in pos_init + ] + + unknown = pybamm._BrentUnknown("x_init_1") + x_guess, y_guess = init_states(unknown, "guess") + if initialization_method == "voltage": + closure = total_lithium(x_guess, y_guess) - Q_Li elif initialization_method == "SOC": - soc_init = pybamm.InputParameter("SOC_init") - negative_soc = x_init_1 * pybamm.InputParameter("Q_n_1") - if is_negative_composite: - x_init_2 = variables["x_init_2"] - negative_soc += x_init_2 * pybamm.InputParameter("Q_n_2") - - negative_0_soc = x_0_1 * pybamm.InputParameter("Q_n_1") - if is_negative_composite: - negative_0_soc += x_0_2 * pybamm.InputParameter("Q_n_2") - - negative_100_soc = x_100_1 * pybamm.InputParameter("Q_n_1") - if is_negative_composite: - negative_100_soc += x_100_2 * pybamm.InputParameter("Q_n_2") - self.algebraic[x_init_1] = ( - (negative_soc - negative_0_soc) / (negative_100_soc - negative_0_soc) - ) - soc_init - self.algebraic[y_init_1] = ( - _get_cyclable_lithium_equation(options, "init") - Q_Li - ) + + def charge(x): + return sum(q * s for q, s in zip(Q_n, x, strict=True)) + + closure = (charge(x_guess) - charge(x_0)) / ( + charge(x_100) - charge(x_0) + ) - pybamm.InputParameter("SOC_init") else: - raise ValueError("Invalid initialization method") - T = param.T_init if initialization_method == "voltage" else param.T_ref - if is_positive_composite: - y_init_2 = variables["y_init_2"] - self.algebraic[y_init_2] = param.p.prim.U( - y_init_1, - T, - get_lithiation_delithiation( - direction, "positive", options, phase="primary" - ), - ) - param.p.sec.U( - y_init_2, - T, - get_lithiation_delithiation( - direction, "positive", options, phase="secondary" - ), - ) - if is_negative_composite: - x_init_2 = variables["x_init_2"] - self.algebraic[x_init_2] = param.n.prim.U( - x_init_1, - T, - get_lithiation_delithiation( - direction, "negative", options, phase="primary" - ), - ) - param.n.sec.U( - x_init_2, - T, - get_lithiation_delithiation( - direction, "negative", options, phase="secondary" - ), + raise pybamm.OptionError( + f"Invalid initialization method '{initialization_method}', " + "expected 'voltage' or 'SOC'" ) - self.variables.update(variables) - if initialization_method == "SOC": - soc_init = pybamm.InputParameter("SOC_init") - else: - soc_init = (V_init - V_min) / (V_max - V_min) - self.initial_conditions.update(_get_initial_conditions(options, soc_init)) + # Solve for x_init_1 rather than the potential: the residual is sensitive to + # it, and a non-physical target still has an exact answer. + x_init_1 = pybamm._Brent( + closure, + unknown, + (_BRACKET_LO, _BRACKET_HI), + abstol=_ABSTOL, + max_iter=_MAX_ITER, + ) + x_init, y_init = init_states(x_init_1, "solved") + + # the stoichiometries are expressions, so the solver needs a state of its own; + # the placeholder is kept out of `variables` + placeholder = pybamm.Variable("ESOH placeholder") + self.algebraic = {placeholder: placeholder} + self.initial_conditions = {placeholder: pybamm.Scalar(0)} + + for index, value in enumerate(x_100): + self.variables[f"x_100_{index + 1}"] = value + for index, value in enumerate(y_100): + self.variables[f"y_100_{index + 1}"] = value + for index, value in enumerate(x_0): + self.variables[f"x_0_{index + 1}"] = value + for index, value in enumerate(y_0): + self.variables[f"y_0_{index + 1}"] = value + for index, value in enumerate(x_init): + self.variables[f"x_init_{index + 1}"] = value + for index, value in enumerate(y_init): + self.variables[f"y_init_{index + 1}"] = value + + # set by `_esoh_evaluator`, keyed on the parameter values + self._evaluator: tuple | None = None @property def default_solver(self): @@ -803,10 +788,11 @@ def solve_full( inputs : dict, optional Additional inputs initial_conditions : dict, optional - Dictionary of initial conditions for variables (e.g., from split solve) + Accepted and ignored. Each stoichiometry is found by a bracketed + rootfind, which needs a bracket rather than a starting guess. esoh_sim : :class:`pybamm.Simulation`, optional - A pre-built simulation wrapping an :class:`ElectrodeSOHComposite` model - to reuse across calls. If not provided, a new one is created. + A pre-built simulation wrapping an :class:`ElectrodeSOHComposite` model. + Passing one back reuses its compiled evaluator across calls. Returns ------- @@ -831,28 +817,20 @@ def solve_full( Q_Li = parameter_values.evaluate(param.Q_Li_particles_init, inputs=inputs) - if isinstance(initial_value, str) and initial_value.endswith("V"): + initialization_method = _initialization_method(initial_value) + if initialization_method == "voltage": V_init = float(initial_value[:-1]) - initialization_method = "voltage" - elif isinstance(initial_value, float): - initialization_method = "SOC" - if initial_value > 1: - warnings.warn( - message=f"Initial SoC {initial_value} is greater than 1", - category=UserWarning, - stacklevel=2, - ) - elif initial_value < 0: - warnings.warn( - message=f"Initial SoC {initial_value} is less than 0", - category=UserWarning, - stacklevel=2, - ) - else: - raise ValueError( - "Invalid initial value. Expected a float between 0 and 1 " - "(for SOC) or a string ending in 'V' (for voltage), got " - f"{initial_value!r} of type {type(initial_value).__name__}" + elif initial_value > 1: + warnings.warn( + message=f"Initial SoC {initial_value} is greater than 1", + category=UserWarning, + stacklevel=2, + ) + elif initial_value < 0: + warnings.warn( + message=f"Initial SoC {initial_value} is less than 0", + category=UserWarning, + stacklevel=2, ) all_inputs = {**inputs, **Qs, "Q_Li": Q_Li} @@ -865,20 +843,114 @@ def solve_full( model = ElectrodeSOHComposite( options, direction, initialization_method=initialization_method ) - esoh_sim = pybamm.Simulation( - model, - parameter_values=parameter_values, - solver=get_esoh_default_solver(tol), + else: + model = esoh_sim.model + + names, input_names, function = _esoh_evaluator(model, parameter_values) + try: + values = np.asarray( + function(*[all_inputs[name] for name in input_names]) + ).reshape(-1) + except RuntimeError as error: + # the native rootfinder reports failure as RuntimeError; callers, and the + # fallback in get_initial_stoichiometries_composite, expect a SolverError + raise pybamm.SolverError( + f"Composite electrode SOH solve failed: {error}" + ) from error + if not np.all(np.isfinite(values)): + raise pybamm.SolverError( + "Composite electrode SOH solve returned a non-finite stoichiometry" ) + return dict(zip(names, values, strict=True)) - if initial_conditions is not None: - esoh_sim.build() - esoh_sim.built_model.set_initial_conditions_from( - initial_conditions, inputs=all_inputs - ) - sol = esoh_sim.solve([0], inputs=all_inputs) - return {var: sol[var].entries[0] for var in sol.all_models[0].variables} +def _unique_nodes(roots): + """Every symbol reachable from ``roots``, visited once. + + Unlike ``pre_order``, which re-yields a shared node once per path to it. + """ + seen: set = set() + nodes = [] + stack = list(roots) + while stack: + symbol = stack.pop() + if id(symbol) in seen: + continue + seen.add(id(symbol)) + nodes.append(symbol) + stack.extend(symbol.children) + return nodes + + +def _parameter_fingerprint(parameter_values, names): + """A comparable summary of the parameters in ``names``. + + Numbers compare by value and everything else by identity, so replacing an OCP + function counts as a change but re-reading the same one does not. + """ + return tuple( + (name, value) + if isinstance(value := parameter_values[name], (int, float)) + else (name, id(value)) + for name in sorted(names) + if name in parameter_values + ) + + +def _esoh_evaluator(model, parameter_values): + """Map capacities and target straight to stoichiometries. + + The model defines each stoichiometry as an expression, not a state, so one CasADi + function replaces a solve. Cached on ``model`` per parameter set. + + Returns + ------- + tuple + ``(names, input_names, function)``, where ``function`` maps the inputs in + ``input_names`` order to the stoichiometries in ``names`` order. + """ + cached = model._evaluator + if cached is not None: + baked, fingerprint, names, input_names, function = cached + if _parameter_fingerprint(parameter_values, baked) == fingerprint: + return names, input_names, function + + names = sorted(model.variables) + nodes = _unique_nodes(model.variables[name] for name in names) + input_names = sorted( + {s.name for s in nodes if isinstance(s, pybamm.InputParameter)} + ) + # Only the parameters the graph substitutes can invalidate it. The capacities + # arrive as InputParameter, so a caller ageing a cell does not rebuild anything. + baked = frozenset( + s.name + for s in nodes + if isinstance(s, pybamm.Parameter | pybamm.FunctionParameter) + ) + symbols = {name: casadi.MX.sym(name) for name in input_names} + time = casadi.MX.sym("t") + state = casadi.MX.sym("y", 1) + # one conversion cache: the variables are views of a single shared graph + converted: dict = {} + expressions = [ + parameter_values.process_symbol(model.variables[name]).to_casadi( + time, state, inputs=symbols, casadi_symbols=converted + ) + for name in names + ] + function = casadi.Function( + "electrode_soh_composite", + list(symbols.values()), + [casadi.vertcat(*expressions)], + ) + model._evaluator = ( + baked, + _parameter_fingerprint(parameter_values, baked), + names, + input_names, + function, + ) + return names, input_names, function def get_initial_stoichiometries_composite( @@ -953,6 +1025,10 @@ def get_initial_stoichiometries_composite( "Only `cyclable lithium capacity` is supported for composite electrodes" ) + # A value that is neither an SOC nor a voltage is the caller's mistake, not a + # solve that failed, so it must not reach the fallback. + _initialization_method(initial_value) + try: return ElectrodeSOHComposite.solve_full( initial_value, @@ -964,45 +1040,26 @@ def get_initial_stoichiometries_composite( inputs=inputs, esoh_sim=esoh_sim, ) - except (pybamm.SolverError, ValueError) as first_error: - if try_split_solve: - try: - split_results = ElectrodeSOHComposite.solve_split( - initial_value, - parameter_values, - direction=direction, - param=param, - options=options, - tol=tol, - inputs=inputs, - ) - - try: - return ElectrodeSOHComposite.solve_full( - initial_value, - parameter_values, - direction=direction, - param=param, - options=options, - tol=tol, - inputs=inputs, - initial_conditions=split_results, - ) - except (pybamm.SolverError, ValueError) as retry_error: - raise ValueError( - f"Failed to solve composite electrode SOH. " - f"Initial full solve error: {first_error}. " - f"Retry with split solve initial conditions also failed: " - f"{retry_error}" - ) from retry_error - - except (pybamm.SolverError, ValueError) as split_error: - raise ValueError( - f"Failed to solve composite electrode SOH. " - f"Full solve error: {first_error}. " - f"Split solve error: {split_error}" - ) from split_error - else: - raise ValueError( + except (pybamm.SolverError, ValueError, RuntimeError) as first_error: + if not try_split_solve: + raise pybamm.SolverError( f"Failed to solve composite electrode SOH: {first_error}" ) from first_error + # The split solve reaches the answer a different way, one electrode at a + # time, so its result is the fallback rather than a guess to re-solve from. + try: + return ElectrodeSOHComposite.solve_split( + initial_value, + parameter_values, + direction=direction, + param=param, + options=options, + tol=tol, + inputs=inputs, + ) + except (pybamm.SolverError, ValueError, RuntimeError) as split_error: + raise pybamm.SolverError( + f"Failed to solve composite electrode SOH. " + f"Full solve error: {first_error}. " + f"Split solve error: {split_error}" + ) from split_error diff --git a/packages/pybamm/tests/integration/test_models/test_full_battery_models/test_lithium_ion/test_electrode_soh_composite_sweeps.py b/packages/pybamm/tests/integration/test_models/test_full_battery_models/test_lithium_ion/test_electrode_soh_composite_sweeps.py new file mode 100644 index 0000000000..24d0b25a29 --- /dev/null +++ b/packages/pybamm/tests/integration/test_models/test_full_battery_models/test_lithium_ion/test_electrode_soh_composite_sweeps.py @@ -0,0 +1,258 @@ +# +# Sweep the composite electrode SOH solver densely over its targets +# +from __future__ import annotations + +import numpy as np +import pytest + +import pybamm +from pybamm.models.full_battery_models.lithium_ion import electrode_soh_composite as esc +from pybamm.models.full_battery_models.lithium_ion.util import ( + get_lithiation_delithiation, +) + +# Dense enough that a bracket that fails only on a sliver of the range is still hit. +SWEEP = 1000 + +OPTIONS = { + "particle phases": ("2", "1"), + "open-circuit potential": (("single", "current sigmoid"), "single"), +} + +# (Q_n_1, Q_n_2, Q_p_1, Q_Li) multipliers on the nominal capacities +WEAR = { + "nominal": (1.0, 1.0, 1.0, 1.0), + "very low lithium": (1.0, 1.0, 1.0, 0.60), + "lost secondary": (1.0, 0.50, 1.0, 0.90), + "worn all": (0.70, 1.0, 0.90, 0.75), +} + + +class Case: + """A composite SOH solver bound to one direction and initialisation method.""" + + def __init__(self, direction, method): + self.options = pybamm.BatteryModelOptions(OPTIONS) + self.parameter_values = pybamm.ParameterValues("Chen2020_composite") + self.param = pybamm.LithiumIonParameters(self.options) + self.method = method + model = pybamm.lithium_ion.ElectrodeSOHComposite( + self.options, direction, initialization_method=method + ) + self.names, self.input_names, self.function = esc._esoh_evaluator( + model, self.parameter_values + ) + + # mirror the temperature and branch choice in ElectrodeSOHComposite.__init__ + reference = self.parameter_values.evaluate(self.param.T_ref) + temperature = ( + self.parameter_values.evaluate(self.param.T_init) + if method == "voltage" + else reference + ) + + def potential(side, phase, electrode): + stoichiometry = pybamm.InputParameter("sto") + branch = get_lithiation_delithiation( + direction, electrode, self.options, phase=phase + ) + processed = self.parameter_values.process_symbol( + (side.prim if phase == "primary" else side.sec).U( + stoichiometry, temperature, branch + ) + ) + return lambda value: float( + np.asarray(processed.evaluate(inputs={"sto": value})).reshape(-1)[0] + ) + + self.U_n = potential(self.param.n, "primary", "negative") + self.U_n2 = potential(self.param.n, "secondary", "negative") + self.U_p = potential(self.param.p, "primary", "positive") + + def capacities(self, scales): + evaluate = self.parameter_values.evaluate + return { + "Q_n_1": evaluate(self.param.n.prim.Q_init) * scales[0], + "Q_n_2": evaluate(self.param.n.sec.Q_init) * scales[1], + "Q_p_1": evaluate(self.param.p.prim.Q_init) * scales[2], + "Q_Li": evaluate(self.param.Q_Li_particles_init) * scales[3], + } + + def solve(self, capacities, target): + key = "V_init" if self.method == "voltage" else "SOC_init" + inputs = {**capacities, key: float(target)} + values = np.asarray( + self.function(*[inputs[name] for name in self.input_names]) + ).reshape(-1) + return dict(zip(self.names, values, strict=True)) + + def lithium(self, capacities, state): + return ( + capacities["Q_n_1"] * state["x_init_1"] + + capacities["Q_n_2"] * state["x_init_2"] + + capacities["Q_p_1"] * state["y_init_1"] + ) + + def state_of_charge(self, capacities, state): + def charge(tag): + return ( + capacities["Q_n_1"] * state[f"x_{tag}_1"] + + capacities["Q_n_2"] * state[f"x_{tag}_2"] + ) + + return (charge("init") - charge("0")) / (charge("100") - charge("0")) + + +class TestCompositeElectrodeSOHSweeps: + """Dense sweeps checking physics rather than stored answers.""" + + @pytest.mark.parametrize("direction", ["discharge", "charge", None]) + @pytest.mark.parametrize("wear", list(WEAR)) + def test_a_thousand_states_of_charge(self, direction, wear): + case = Case(direction, "SOC") + capacities = case.capacities(WEAR[wear]) + for target in np.linspace(0.0, 1.0, SWEEP): + state = case.solve(capacities, target) + assert all(np.isfinite(v) for v in state.values()), target + assert case.lithium(capacities, state) == pytest.approx( + capacities["Q_Li"], rel=1e-9 + ), target + assert case.state_of_charge(capacities, state) == pytest.approx( + target, abs=1e-7 + ) + assert case.U_n(state["x_init_1"]) == pytest.approx( + case.U_n2(state["x_init_2"]), abs=1e-6 + ), target + + @pytest.mark.parametrize("direction", ["discharge", "charge", None]) + @pytest.mark.parametrize("wear", list(WEAR)) + def test_a_thousand_voltages(self, direction, wear): + case = Case(direction, "voltage") + capacities = case.capacities(WEAR[wear]) + for target in np.linspace(2.5, 4.2, SWEEP): + state = case.solve(capacities, target) + assert all(np.isfinite(v) for v in state.values()), target + assert case.lithium(capacities, state) == pytest.approx( + capacities["Q_Li"], rel=1e-9 + ), target + voltage = case.U_p(state["y_init_1"]) - case.U_n(state["x_init_1"]) + assert voltage == pytest.approx(target, abs=1e-7) + + def test_repeating_a_solve_returns_identical_bits(self): + # the rootfinder caches its last solve, which must not change an answer + case = Case("discharge", "SOC") + capacities = case.capacities(WEAR["nominal"]) + for target in np.linspace(0.0, 1.0, SWEEP): + first = case.solve(capacities, target) + assert case.solve(capacities, target) == first + + def test_the_answer_moves_smoothly_with_the_target(self): + # a bracket that flipped to another root would show up as a jump + case = Case("discharge", "SOC") + capacities = case.capacities(WEAR["nominal"]) + targets = np.linspace(0.0, 1.0, SWEEP) + x_init = np.array([case.solve(capacities, t)["x_init_1"] for t in targets]) + steps = np.diff(x_init) + assert np.all(steps > 0), "x_init_1 must increase with state of charge" + assert np.max(steps) < 50 * np.median(steps), "discontinuity in x_init_1" + + def test_a_state_of_charge_round_trips_through_its_voltage(self): + by_soc = Case("discharge", "SOC") + by_voltage = Case("discharge", "voltage") + soc_capacities = by_soc.capacities(WEAR["nominal"]) + voltage_capacities = by_voltage.capacities(WEAR["nominal"]) + for target in np.linspace(0.05, 0.95, SWEEP): + state = by_soc.solve(soc_capacities, target) + voltage = by_voltage.U_p(state["y_init_1"]) - by_voltage.U_n( + state["x_init_1"] + ) + back = by_voltage.solve(voltage_capacities, voltage) + assert back["x_init_1"] == pytest.approx(state["x_init_1"], abs=1e-6) + + @pytest.mark.parametrize( + ("method", "targets"), + [ + ("SOC", np.linspace(-5.0, 5.0, SWEEP)), + ("voltage", np.linspace(0.5, 6.0, SWEEP)), + ], + ) + def test_a_target_outside_the_window_still_solves_exactly(self, method, targets): + """A non-physical target has an exact answer, and must be given it.""" + case = Case("discharge", method) + capacities = case.capacities(WEAR["nominal"]) + for target in targets: + state = case.solve(capacities, target) + assert all(np.isfinite(v) for v in state.values()), target + assert case.lithium(capacities, state) == pytest.approx( + capacities["Q_Li"], rel=1e-9 + ), target + if method == "voltage": + voltage = case.U_p(state["y_init_1"]) - case.U_n(state["x_init_1"]) + assert voltage == pytest.approx(target, abs=1e-6), target + else: + assert case.state_of_charge(capacities, state) == pytest.approx( + target, abs=1e-7 + ), target + + def test_a_non_physical_target_returns_a_non_physical_stoichiometry(self): + case = Case("discharge", "SOC") + capacities = case.capacities(WEAR["nominal"]) + assert case.solve(capacities, -3.0)["x_init_1"] < 0 + assert case.solve(capacities, 3.0)["x_init_1"] > 1 + + +class TestCompositeElectrodeSOHReuse: + """Ageing a cell must not rebuild the evaluator. + + The capacities reach the compiled function as inputs, so only the parameters its + graph substitutes can invalidate it. + """ + + # both fall as a cell ages, and both feed the capacities the solve is given + FADE = ( + "Initial concentration in positive electrode [mol.m-3]", + "Primary: Negative electrode active material volume fraction", + ) + + @staticmethod + def _reused_simulation(parameter_values): + model = pybamm.lithium_ion.ElectrodeSOHComposite( + pybamm.BatteryModelOptions(OPTIONS), + "discharge", + initialization_method="SOC", + ) + return pybamm.Simulation(model, parameter_values=parameter_values) + + def _call(self, parameter_values, simulation): + return pybamm.lithium_ion.get_initial_stoichiometries_composite( + 0.5, + parameter_values, + direction="discharge", + options=OPTIONS, + esoh_sim=simulation, + ) + + def test_fading_capacities_reuse_the_compiled_evaluator(self): + parameter_values = pybamm.ParameterValues("Chen2020_composite") + simulation = self._reused_simulation(parameter_values) + self._call(parameter_values, simulation) + built = simulation.model._evaluator[-1] + + nominal = {key: parameter_values[key] for key in self.FADE} + for step in range(1, 11): + for key in self.FADE: + parameter_values[key] = nominal[key] * (1 - 0.002 * step) + state = self._call(parameter_values, simulation) + assert np.isfinite(state["x_init_1"]) + assert simulation.model._evaluator[-1] is built + + def test_a_new_open_circuit_potential_does_rebuild(self): + parameter_values = pybamm.ParameterValues("Chen2020_composite") + simulation = self._reused_simulation(parameter_values) + self._call(parameter_values, simulation) + built = simulation.model._evaluator[-1] + + parameter_values["Primary: Negative electrode OCP [V]"] = lambda sto: 1.5 - sto + self._call(parameter_values, simulation) + assert simulation.model._evaluator[-1] is not built diff --git a/packages/pybamm/tests/unit/test_models/test_full_battery_models/test_lithium_ion/test_electrode_soh.py b/packages/pybamm/tests/unit/test_models/test_full_battery_models/test_lithium_ion/test_electrode_soh.py index 460a11179a..d268eca254 100644 --- a/packages/pybamm/tests/unit/test_models/test_full_battery_models/test_lithium_ion/test_electrode_soh.py +++ b/packages/pybamm/tests/unit/test_models/test_full_battery_models/test_lithium_ion/test_electrode_soh.py @@ -4,9 +4,13 @@ import contextlib +import numpy as np import pytest import pybamm +from pybamm.models.full_battery_models.lithium_ion.util import ( + get_lithiation_delithiation, +) # Fixture for TestElectrodeSOHMSMR, TestCalculateTheoreticalEnergy and TestGetInitialOCPMSMR class. @@ -1163,3 +1167,196 @@ def test_min_max_ocp(self, options): ) assert Up_100 - Un_100 == pytest.approx(4.2) assert Up_0 - Un_0 == pytest.approx(2.8) + + +class TestElectrodeSOHCompositeHardCases: + """Composite-electrode states checked against the physics, not stored answers. + + The requested state must be the one that comes back, and lithium conserved. + """ + + OPTIONS = { + "particle phases": ("2", "1"), + "open-circuit potential": (("single", "current sigmoid"), "single"), + } + + # (Q_n_1, Q_n_2, Q_p_1, Q_Li) multipliers on the nominal capacities + NOMINAL = (1.0, 1.0, 1.0, 1.0) + LOW_LI = (1.0, 1.0, 1.0, 0.85) + VERY_LOW_LI = (1.0, 1.0, 1.0, 0.6) + WORN_NEGATIVE = (0.85, 0.85, 1.0, 0.85) + LOST_SECONDARY = (1.0, 0.5, 1.0, 0.9) + WORN_ALL = (0.7, 1.0, 0.9, 0.75) + + CASES = [ + ("voltage", NOMINAL, 3.35), + ("voltage", WORN_ALL, 2.6416666666666666), + ("SOC", NOMINAL, 0.5833333333333334), + ("SOC", LOW_LI, 0.0), + ("SOC", LOW_LI, 0.08333333333333333), + ("SOC", VERY_LOW_LI, 0.0), + ("SOC", VERY_LOW_LI, 0.08333333333333333), + ("SOC", WORN_NEGATIVE, 0.08333333333333333), + ("SOC", LOST_SECONDARY, 0.0), + ("SOC", LOST_SECONDARY, 0.08333333333333333), + ] + + @staticmethod + def _setup(scales): + options = pybamm.BatteryModelOptions(TestElectrodeSOHCompositeHardCases.OPTIONS) + parameter_values = pybamm.ParameterValues("Chen2020_composite") + param = pybamm.LithiumIonParameters(options) + scale_n1, scale_n2, scale_p1, scale_li = scales + capacities = { + "Q_Li": parameter_values.evaluate(param.Q_Li_particles_init) * scale_li, + "Q_n_1": parameter_values.evaluate(param.n.prim.Q_init) * scale_n1, + "Q_n_2": parameter_values.evaluate(param.n.sec.Q_init) * scale_n2, + "Q_p_1": parameter_values.evaluate(param.p.prim.Q_init) * scale_p1, + } + return options, parameter_values, param, capacities + + @staticmethod + def _ocp(parameter_values, potential, stoichiometry, temperature, branch): + sto = pybamm.InputParameter("sto") + processed = parameter_values.process_symbol(potential(sto, temperature, branch)) + return float( + np.asarray(processed.evaluate(inputs={"sto": stoichiometry})).reshape(-1)[0] + ) + + @pytest.mark.parametrize(("initialization_method", "scales", "target"), CASES) + def test_reaches_the_requested_state(self, initialization_method, scales, target): + options, parameter_values, param, capacities = self._setup(scales) + model = pybamm.lithium_ion.ElectrodeSOHComposite( + options, initialization_method=initialization_method + ) + key = "V_init" if initialization_method == "voltage" else "SOC_init" + inputs = {**capacities, key: target} + sim = pybamm.Simulation(model, parameter_values=parameter_values) + solution = sim.solve([0], inputs=inputs) + state = {name: float(solution[name](0)) for name in model.variables} + + def branch(electrode, phase): + return get_lithiation_delithiation(None, electrode, options, phase=phase) + + if initialization_method == "voltage": + voltage = self._ocp( + parameter_values, + param.p.prim.U, + state["y_init_1"], + param.T_init, + branch("positive", "primary"), + ) - self._ocp( + parameter_values, + param.n.prim.U, + state["x_init_1"], + param.T_init, + branch("negative", "primary"), + ) + assert voltage == pytest.approx(target, abs=1e-08) + else: + + def charge(tag): + return ( + capacities["Q_n_1"] * state[f"x_{tag}_1"] + + capacities["Q_n_2"] * state[f"x_{tag}_2"] + ) + + soc = (charge("init") - charge("0")) / (charge("100") - charge("0")) + assert soc == pytest.approx(target, abs=1e-08) + + lithium = ( + capacities["Q_n_1"] * state["x_init_1"] + + capacities["Q_n_2"] * state["x_init_2"] + + capacities["Q_p_1"] * state["y_init_1"] + ) + assert lithium == pytest.approx(capacities["Q_Li"], rel=1e-10) + + +class TestElectrodeSOHCompositeFallback: + """The split solve is the fallback when the bracketed solve cannot answer. + + Nothing in the parameter sweep fails, so these force the first attempt to. + """ + + OPTIONS = { + "particle phases": ("2", "1"), + "open-circuit potential": (("single", "current sigmoid"), "single"), + } + + def _call(self, **kwargs): + return pybamm.lithium_ion.get_initial_stoichiometries_composite( + 0.5, + pybamm.ParameterValues("Chen2020_composite"), + direction="discharge", + options=self.OPTIONS, + **kwargs, + ) + + @pytest.mark.parametrize( + "error", + [ + pybamm.SolverError("full solve failed"), + RuntimeError("rootfinder process failed"), + ValueError("bad bracket"), + ], + ) + def test_a_failed_full_solve_falls_back_to_the_split_solve( + self, monkeypatch, error + ): + def fail(*args, **kwargs): + raise error + + expected = pybamm.lithium_ion.ElectrodeSOHComposite.solve_split( + 0.5, + pybamm.ParameterValues("Chen2020_composite"), + direction="discharge", + options=self.OPTIONS, + ) + monkeypatch.setattr( + pybamm.lithium_ion.ElectrodeSOHComposite, "solve_full", fail + ) + from_split = self._call(try_split_solve=True) + assert from_split.keys() == expected.keys() + for name, value in expected.items(): + assert from_split[name] == pytest.approx(value, abs=1e-10), name + + def test_the_fallback_can_be_turned_off(self, monkeypatch): + def fail(*args, **kwargs): + raise RuntimeError("rootfinder process failed") + + monkeypatch.setattr( + pybamm.lithium_ion.ElectrodeSOHComposite, "solve_full", fail + ) + with pytest.raises(pybamm.SolverError, match="Failed to solve composite"): + self._call(try_split_solve=False) + + def test_both_failing_reports_both_errors(self, monkeypatch): + def fail_full(*args, **kwargs): + raise pybamm.SolverError("full is broken") + + def fail_split(*args, **kwargs): + raise pybamm.SolverError("split is broken") + + monkeypatch.setattr( + pybamm.lithium_ion.ElectrodeSOHComposite, "solve_full", fail_full + ) + monkeypatch.setattr( + pybamm.lithium_ion.ElectrodeSOHComposite, "solve_split", fail_split + ) + with pytest.raises( + pybamm.SolverError, match=r"full is broken.*split is broken" + ): + self._call(try_split_solve=True) + + def test_a_non_finite_stoichiometry_is_a_failure_not_a_result(self, monkeypatch): + # a solve that returns NaN must reach the fallback, not be handed back + def nan_evaluator(model, parameter_values): + names = sorted(model.variables) + return names, [], lambda *a: np.full(len(names), np.nan) + + monkeypatch.setattr( + pybamm.models.full_battery_models.lithium_ion.electrode_soh_composite, + "_esoh_evaluator", + nan_evaluator, + ) + assert np.isfinite(self._call(try_split_solve=True)["x_init_1"]) From c2e203c464061df3f8b45279d2f163b5f8860d9c Mon Sep 17 00:00:00 2001 From: Marc Berliner <34451391+MarcBerliner@users.noreply.github.com> Date: Mon, 24 Aug 2026 12:25:25 -0400 Subject: [PATCH 2/3] fix: take the split solve on Windows, where the plugin cannot be reached The full solve evaluates the stoichiometries as one expression holding rootfinds nested about three deep -- one composite solve makes on the order of a million of them -- which only the native "brent" plugin gets through quickly. Windows cannot load that plugin, and neither Python-side route is usable there: walking the tree could not finish a single solve in fourteen minutes, and driving the compiled residual from a CasADi callback still pays 830k Python round trips. So Windows takes the path it took before this branch. solve_split is algebraic throughout -- _ElectrodeSOH through a Simulation, no rootfind node anywhere in it -- and get_initial_stoichiometries_composite already falls back to it, so raising from solve_full is all it takes. 47 of the 63 tests in the file pass that way in 35s; the 16 that exercise the full solve itself skip, since what they test is genuinely absent. pybamm.is_windows() joins is_macos_intel() in util, so the check has one home. Co-Authored-By: Claude Opus 5 --- packages/pybamm/src/pybamm/__init__.py | 1 + .../lithium_ion/electrode_soh_composite.py | 9 +++++++ packages/pybamm/src/pybamm/util.py | 7 ++++++ .../test_lithium_ion/test_electrode_soh.py | 25 +++++++++++++------ 4 files changed, 34 insertions(+), 8 deletions(-) diff --git a/packages/pybamm/src/pybamm/__init__.py b/packages/pybamm/src/pybamm/__init__.py index 740fca3e26..405783a81f 100644 --- a/packages/pybamm/src/pybamm/__init__.py +++ b/packages/pybamm/src/pybamm/__init__.py @@ -11,6 +11,7 @@ from .util import ( get_parameters_filepath, has_jax, + is_windows, raise_jax_not_found, import_optional_dependency, ) diff --git a/packages/pybamm/src/pybamm/models/full_battery_models/lithium_ion/electrode_soh_composite.py b/packages/pybamm/src/pybamm/models/full_battery_models/lithium_ion/electrode_soh_composite.py index 2ed1d9cd1b..51c5adce3c 100644 --- a/packages/pybamm/src/pybamm/models/full_battery_models/lithium_ion/electrode_soh_composite.py +++ b/packages/pybamm/src/pybamm/models/full_battery_models/lithium_ion/electrode_soh_composite.py @@ -833,6 +833,15 @@ def solve_full( stacklevel=2, ) + # The stoichiometries are one expression holding rootfinds nested several + # deep, which only the native plugin evaluates quickly enough; on Windows it + # cannot be reached, so the caller's split solve answers instead. + if pybamm.is_windows(): + raise pybamm.SolverError( + "the full composite electrode SOH solve needs the 'brent' CasADi " + "rootfinder plugin, which is not available on Windows" + ) + all_inputs = {**inputs, **Qs, "Q_Li": Q_Li} if initialization_method == "voltage": all_inputs["V_init"] = V_init diff --git a/packages/pybamm/src/pybamm/util.py b/packages/pybamm/src/pybamm/util.py index d7ecc7043b..46de9976f1 100644 --- a/packages/pybamm/src/pybamm/util.py +++ b/packages/pybamm/src/pybamm/util.py @@ -368,6 +368,13 @@ def is_macos_intel(): return sys.platform == "darwin" and platform.machine() == "x86_64" +def is_windows(): + """Check if running on Windows.""" + import sys + + return sys.platform == "win32" + + def raise_jax_not_found(): """Raise an appropriate error when JAX is not available.""" if is_macos_intel(): diff --git a/packages/pybamm/tests/unit/test_models/test_full_battery_models/test_lithium_ion/test_electrode_soh.py b/packages/pybamm/tests/unit/test_models/test_full_battery_models/test_lithium_ion/test_electrode_soh.py index d268eca254..d9e9a1f627 100644 --- a/packages/pybamm/tests/unit/test_models/test_full_battery_models/test_lithium_ion/test_electrode_soh.py +++ b/packages/pybamm/tests/unit/test_models/test_full_battery_models/test_lithium_ion/test_electrode_soh.py @@ -289,6 +289,12 @@ def ocv(x, y, negative_branch, positive_branch): assert ocv(x_0, y_0, None, None) == pytest.approx(V_min + 0.1, abs=1e-6) +needs_full_solve = pytest.mark.skipif( + pybamm.is_windows(), + reason="the full composite solve needs the brent plugin, which Windows cannot load", +) + + class TestElectrodeSOHComposite: @staticmethod def _check_phases_equal(results, xy, soc): @@ -377,6 +383,7 @@ def _get_params_and_options(composite_electrode): "positive", # positive-only composite ], ) + @needs_full_solve def test_half_cell_with_same_ocp_curves(self, composite_electrode, initial_value): pvals, options = self._get_params_and_options(composite_electrode) # Use composite ESOH helper to compute initial stoichiometries at a voltage @@ -435,6 +442,7 @@ def test_chen2020_composite_defaults(self, initial_value): - param.n.prim.U(results["x_init_1"], param.T_ref) ) == pytest.approx(V_target, abs=1e-05) + @needs_full_solve def test_chen2020_composite_default_solve(self): pvals = pybamm.ParameterValues("Chen2020_composite") options = {"particle phases": ("2", "1")} @@ -1060,15 +1068,16 @@ def test_error_warning(self): "5 A", parameter_values_composite, options=options_composite ) - with pytest.warns(UserWarning, match=r"is greater than 1"): - pybamm.lithium_ion.ElectrodeSOHComposite.solve_full( - 1.001, parameter_values_composite, options=options_composite - ) + if not pybamm.is_windows(): + with pytest.warns(UserWarning, match=r"is greater than 1"): + pybamm.lithium_ion.ElectrodeSOHComposite.solve_full( + 1.001, parameter_values_composite, options=options_composite + ) - with pytest.warns(UserWarning, match=r"is less than 0"): - pybamm.lithium_ion.ElectrodeSOHComposite.solve_full( - -0.001, parameter_values_composite, options=options_composite - ) + with pytest.warns(UserWarning, match=r"is less than 0"): + pybamm.lithium_ion.ElectrodeSOHComposite.solve_full( + -0.001, parameter_values_composite, options=options_composite + ) with pytest.warns(UserWarning, match=r"is greater than 1"): pybamm.lithium_ion.ElectrodeSOHComposite.solve_split( 1.001, parameter_values_composite, options=options_composite From 3a8b1da600553a6fc50195101f5756a68fa275a6 Mon Sep 17 00:00:00 2001 From: Marc Berliner <34451391+MarcBerliner@users.noreply.github.com> Date: Mon, 24 Aug 2026 14:47:02 -0400 Subject: [PATCH 3/3] fix: skip the model-through-solver test on Windows too test_reaches_the_requested_state solves ElectrodeSOHComposite as a model and reads its variables back, so the rootfinds travel into IDAKLU as a serialised function. On Windows they are a CasADi Callback by then, which the solver's own CasADi cannot deserialise -- "not found 'CallbackInternal'" -- so this is the same gap as the other full-solve tests, one layer further in. Co-Authored-By: Claude Opus 5 --- .../test_lithium_ion/test_electrode_soh.py | 2 ++ 1 file changed, 2 insertions(+) diff --git a/packages/pybamm/tests/unit/test_models/test_full_battery_models/test_lithium_ion/test_electrode_soh.py b/packages/pybamm/tests/unit/test_models/test_full_battery_models/test_lithium_ion/test_electrode_soh.py index d9e9a1f627..08b0a0f2e7 100644 --- a/packages/pybamm/tests/unit/test_models/test_full_battery_models/test_lithium_ion/test_electrode_soh.py +++ b/packages/pybamm/tests/unit/test_models/test_full_battery_models/test_lithium_ion/test_electrode_soh.py @@ -1232,6 +1232,8 @@ def _ocp(parameter_values, potential, stoichiometry, temperature, branch): np.asarray(processed.evaluate(inputs={"sto": stoichiometry})).reshape(-1)[0] ) + # solves the model itself, whose variables hold the rootfinds + @needs_full_solve @pytest.mark.parametrize(("initialization_method", "scales", "target"), CASES) def test_reaches_the_requested_state(self, initialization_method, scales, target): options, parameter_values, param, capacities = self._setup(scales)