diff --git a/ecoli/experiments/metabolism_redux_sim.py b/ecoli/experiments/metabolism_redux_sim.py index f15c3e557..1bd706879 100644 --- a/ecoli/experiments/metabolism_redux_sim.py +++ b/ecoli/experiments/metabolism_redux_sim.py @@ -210,6 +210,57 @@ def test_ecoli_with_metabolism_redux_div( ) +@pytest.mark.slow +def test_ecoli_metabolism_redux_solver_fallback(monkeypatch): + """Run the *entire* MetabolismRedux process inside a short EcoliSim and + verify the cvxpy solver fallback recovers a transient primary-solver + (GLOP) failure, so the simulation runs to completion instead of losing the + seed. + + This is the in-sim counterpart to the standalone-LP unit tests in + ``ecoli/processes/test_metabolism_redux_solver_fallback.py`` and resolves + the ``# TODO (Cyrus) Add test for entire process`` in + ``ecoli/processes/metabolism_redux.py``. It also speaks to the review + question "do the other solvers fail on the same simulations GLOP fails + on?": here GLOP is made to raise on a *real* in-sim FBA problem and a + fallback solver reaches an optimal solution on that identical problem, so + the tick -- and the whole run -- survive. + """ + import cvxpy as cp + from cvxpy.error import SolverError + + real_solve = cp.Problem.solve + injected = {"glop_failures": 0} + + def flaky_glop_solve(self, solver=None, **kwargs): + # Inject a single transient GLOP failure on the first GLOP solve of + # the run; every other solve -- including the fallback retry -- is real. + if solver == cp.GLOP and injected["glop_failures"] == 0: + injected["glop_failures"] += 1 + raise SolverError("simulated transient GLOP failure (test)") + return real_solve(self, solver=solver, **kwargs) + + monkeypatch.setattr(cp.Problem, "solve", flaky_glop_solve) + + sim = EcoliSim.from_file(CONFIG_DIR_PATH + "metabolism_redux.json") + sim.max_duration = 4 + sim.divide = False + sim.progress_bar = False + sim.log_updates = False + sim.emitter = "timeseries" + sim.build_ecoli() + + # Would raise if the injected GLOP failure were not caught and retried by + # solve_with_fallback; a completed run means a fallback solver handled the + # identical problem GLOP choked on. + sim.run() + + # Guard against a vacuous pass: the fallback path was actually exercised. + assert injected["glop_failures"] == 1, ( + "expected exactly one injected GLOP failure to be triggered in-sim" + ) + + @pytest.mark.slow def test_ecoli_with_metabolism_classic( filename="metabolism_redux_classic", diff --git a/ecoli/processes/metabolism_redux.py b/ecoli/processes/metabolism_redux.py index 822c3e5fa..6045f51f4 100644 --- a/ecoli/processes/metabolism_redux.py +++ b/ecoli/processes/metabolism_redux.py @@ -2,6 +2,7 @@ MetabolismRedux """ +import logging import numpy as np import numpy.typing as npt import time @@ -19,6 +20,7 @@ from ecoli.processes.registries import topology_registry import cvxpy as cp +from cvxpy.error import SolverError from typing import Iterable, Mapping from dataclasses import dataclass @@ -28,6 +30,106 @@ VOLUME_UNITS = units.L MASS_UNITS = units.g TIME_UNITS = units.s + +log = logging.getLogger(__name__) + +# Fallback solvers to try, in order, if the primary solver (GLOP by default) +# raises cvxpy.error.SolverError or fails to reach an optimal solution. +# GLOP failures are typically transient/numerical (e.g. an ill-conditioned +# basis at a particular tick/seed) rather than problem infeasibility, so +# retrying with a different LP-capable solver that is installed in this +# environment (see `cvxpy.installed_solvers()`) usually recovers a solution +# instead of losing the whole simulation seed. +DEFAULT_FALLBACK_SOLVERS: tuple = (cp.PDLP, cp.CLARABEL, cp.SCS) + + +def solve_with_fallback( + problem: "cp.Problem", + primary_solver=cp.GLOP, + fallback_solvers: Iterable = DEFAULT_FALLBACK_SOLVERS, + **solve_kwargs, +) -> "cp.Problem": + """Solve a cvxpy Problem in place, falling back to alternate solvers if + the primary solver errors out or fails to find an optimal solution. + + This does NOT change the problem/objective/constraints in any way -- it + only changes which underlying solver cvxpy dispatches to. `primary_solver` + is tried first (fast path); on `cvxpy.error.SolverError` or a non-optimal + status, each solver in `fallback_solvers` is tried in turn, logging a + warning each time a fallback is attempted. If every solver fails, the + original exception from the primary solver is re-raised (or a ValueError + for a non-optimal status), with the list of attempted solvers noted. + + Returns the same `problem` object (now populated with the solution from + whichever solver succeeded). + """ + solvers_to_try = [primary_solver, *fallback_solvers] + first_error: Optional[BaseException] = None + + for i, solver in enumerate(solvers_to_try): + is_primary = i == 0 + try: + problem.solve(solver=solver, **solve_kwargs) + except SolverError as e: + if first_error is None: + first_error = e + if not is_primary: + log.warning( + "Fallback solver %s also failed while solving network " + "flow metabolism model: %s", + solver, + e, + ) + continue + log.warning( + "Primary solver %s raised SolverError (%s); retrying with " + "fallback solvers %s.", + solver, + e, + list(fallback_solvers), + ) + continue + + if problem.status == "optimal": + if not is_primary: + log.warning( + "Network flow metabolism model solved successfully with " + "fallback solver %s after primary solver %s failed.", + solver, + primary_solver, + ) + return problem + + # Solved without raising, but did not reach an optimal status. + status_error = ValueError( + f"Solver {solver} returned non-optimal status " + f"'{problem.status}' for network flow metabolism model." + ) + if first_error is None: + first_error = status_error + if is_primary: + log.warning( + "Primary solver %s returned non-optimal status '%s'; " + "retrying with fallback solvers %s.", + solver, + problem.status, + list(fallback_solvers), + ) + else: + log.warning( + "Fallback solver %s also returned non-optimal status '%s'.", + solver, + problem.status, + ) + + # Every solver either raised or failed to reach an optimal solution. + raise ValueError( + "Network flow model of metabolism did not converge to an optimal " + f"solution with any solver (tried {solvers_to_try}). " + f"Original error from primary solver {primary_solver}: {first_error}" + ) from first_error + + CONC_UNITS = COUNTS_UNITS / VOLUME_UNITS CONVERSION_UNITS = MASS_UNITS * TIME_UNITS / VOLUME_UNITS GDCW_BASIS = units.mmol / units.g / units.h @@ -1058,21 +1160,7 @@ def solve( p = cp.Problem(cp.Minimize(loss), constr) - try: - p.solve(solver=solver, verbose=False) - except cp.error.SolverError: - p.solve( - solver=solver, - verbose=True, - ) - raise ValueError( - "Network flow model of metabolism did not converge to a solution." - ) - if p.status != "optimal": - raise ValueError( - "Network flow model of metabolism did not " - "converge to an optimal solution." - ) + solve_with_fallback(p, primary_solver=solver, verbose=False) velocities = np.array(v.value) dm_dt = np.array(dm.value) @@ -1127,7 +1215,12 @@ def test_network_flow_model(): ) -# TODO (Cyrus) Add test for entire process +# The entire process is exercised in-sim by +# test_ecoli_metabolism_redux_solver_fallback in +# ecoli/experiments/metabolism_redux_sim.py, which runs MetabolismRedux inside +# a short EcoliSim and asserts the solver fallback recovers a transient GLOP +# failure. The standalone solver-fallback LP is covered by +# ecoli/processes/test_metabolism_redux_solver_fallback.py. if __name__ == "__main__": test_network_flow_model() diff --git a/ecoli/processes/test_metabolism_redux_solver_fallback.py b/ecoli/processes/test_metabolism_redux_solver_fallback.py new file mode 100644 index 000000000..c8c6bc70a --- /dev/null +++ b/ecoli/processes/test_metabolism_redux_solver_fallback.py @@ -0,0 +1,97 @@ +""" +Unit tests for the cvxpy solver fallback in +`ecoli.processes.metabolism_redux.solve_with_fallback`. + +These use a trivial standalone cvxpy LP (not the full MetabolismRedux +network flow model) so they run fast and don't require a whole-cell +simulation. They verify that when the primary solver (standing in for +GLOP) raises `cvxpy.error.SolverError`, the solve transparently falls +back to another solver and still returns an optimal solution, without +altering the problem's objective or constraints. See +`ecoli/processes/metabolism_redux.py` for the production usage inside +`NetworkFlowModel.solve`. +""" + +import cvxpy as cp +import pytest +from cvxpy.error import SolverError + +from ecoli.processes.metabolism_redux import solve_with_fallback + + +def test_solve_with_fallback_recovers_from_primary_solver_error(monkeypatch): + """A SolverError from the primary solver should trigger a retry with a + fallback solver, and the fallback should still find the correct + optimum (i.e. the problem itself is untouched).""" + x = cp.Variable() + problem = cp.Problem(cp.Minimize(cp.abs(x - 3)), [x >= 0, x <= 10]) + + real_solve = cp.Problem.solve + calls = [] + + def fake_solve(self, solver=None, **kwargs): + calls.append(solver) + if solver == "FAKE_PRIMARY": + raise SolverError("simulated GLOP failure") + return real_solve(self, solver=solver, **kwargs) + + monkeypatch.setattr(cp.Problem, "solve", fake_solve) + + solve_with_fallback( + problem, + primary_solver="FAKE_PRIMARY", + fallback_solvers=(cp.CLARABEL,), + ) + + assert problem.status == "optimal" + assert x.value is not None + assert abs(x.value - 3) < 1e-4 + # Primary was tried first and failed; fallback was then used. + assert calls == ["FAKE_PRIMARY", cp.CLARABEL] + + +def test_solve_with_fallback_tries_multiple_fallbacks_in_order(monkeypatch): + """If the first fallback also errors, the next fallback in the list + should be attempted before giving up.""" + x = cp.Variable() + problem = cp.Problem(cp.Minimize(cp.abs(x - 3)), [x >= 0, x <= 10]) + + real_solve = cp.Problem.solve + calls = [] + + def fake_solve(self, solver=None, **kwargs): + calls.append(solver) + if solver in ("FAKE_PRIMARY", "FAKE_FALLBACK_1"): + raise SolverError(f"simulated failure for {solver}") + return real_solve(self, solver=solver, **kwargs) + + monkeypatch.setattr(cp.Problem, "solve", fake_solve) + + solve_with_fallback( + problem, + primary_solver="FAKE_PRIMARY", + fallback_solvers=("FAKE_FALLBACK_1", cp.CLARABEL), + ) + + assert problem.status == "optimal" + assert calls == ["FAKE_PRIMARY", "FAKE_FALLBACK_1", cp.CLARABEL] + + +def test_solve_with_fallback_raises_with_context_when_all_solvers_fail(monkeypatch): + """If every solver (primary + all fallbacks) fails, the original error + should be re-raised (wrapped) with context about what was tried, + instead of silently swallowing the failure.""" + x = cp.Variable() + problem = cp.Problem(cp.Minimize(cp.abs(x)), [x >= 0]) + + def always_fail(self, solver=None, **kwargs): + raise SolverError(f"simulated failure for {solver}") + + monkeypatch.setattr(cp.Problem, "solve", always_fail) + + with pytest.raises(ValueError, match="did not converge"): + solve_with_fallback( + problem, + primary_solver="FAKE_PRIMARY", + fallback_solvers=("FAKE_FALLBACK",), + )