fix: better handle infeasible optimizations

src/memote/suite/tests/test_biomass.py

@@ -151,8 +151,8 @@ def test_biomass_default_production(model, reaction_id):
 
     """
     ann = test_biomass_default_production.annotation
-    ann["data"][reaction_id] = helpers.run_fba(model, reaction_id)
-    outcome = ann["data"][reaction_id] > 1e-07
+    ann["data"][reaction_id] = helpers.get_biomass_flux(model, reaction_id)
+    outcome = ann["data"][reaction_id] > model.tolerance
     ann["metric"][reaction_id] = 1.0 - float(outcome)
     ann["message"][reaction_id] = wrapper.fill(
         """Using the biomass reaction {} this is the growth rate (1/h) that
@@ -188,8 +188,8 @@ def test_biomass_open_production(model, reaction_id):
     """
     ann = test_biomass_open_production.annotation
     helpers.open_boundaries(model)
-    ann["data"][reaction_id] = helpers.run_fba(model, reaction_id)
-    outcome = ann["data"][reaction_id] > 1e-07
+    ann["data"][reaction_id] = helpers.get_biomass_flux(model, reaction_id)
+    outcome = ann["data"][reaction_id] > model.tolerance
     ann["metric"][reaction_id] = 1.0 - float(outcome)
     ann["message"][reaction_id] = wrapper.fill(
         """Using the biomass reaction {} this is the growth rate that can be
@@ -408,7 +408,7 @@ def test_fast_growth_default(model, reaction_id):
 
     """
     ann = test_fast_growth_default.annotation
-    outcome = helpers.run_fba(model, reaction_id) > 2.81
+    outcome = helpers.get_biomass_flux(model, reaction_id) > 2.81
     ann["data"][reaction_id] = outcome
     ann["metric"][reaction_id] = 1.0 - float(outcome)
     if ann["data"][reaction_id]:

src/memote/support/biomass.py

@@ -184,8 +184,7 @@ def find_blocked_biomass_precursors(reaction, model):
             dm_rxn = model.add_boundary(
                 precursor, type="safe-demand", reaction_id="safe_demand", lb=0, ub=ub
             )
-            flux = helpers.run_fba(model, dm_rxn.id, direction="max")
-            if np.isnan(flux) or abs(flux) < 1e-08:
+            if helpers.get_biomass_flux(model, dm_rxn.id) <= 0.0:
                 blocked_precursors.append(precursor)
     return blocked_precursors
 

src/memote/support/helpers.py

@@ -24,6 +24,7 @@
 from collections import defaultdict
 from operator import attrgetter, itemgetter
 
+import cobra
 import numpy as np
 import pandas as pd
 from cobra.exceptions import Infeasible
@@ -630,10 +631,7 @@ def run_fba(model, rxn_id, direction="max", single_value=True):
     model.objective = model.reactions.get_by_id(rxn_id)
     model.objective_direction = direction
     if single_value:
-        try:
-            return model.slim_optimize()
-        except Infeasible:
-            return np.nan
+        return model.slim_optimize()
     else:
         try:
             solution = model.optimize()
@@ -642,6 +640,17 @@ def run_fba(model, rxn_id, direction="max", single_value=True):
             return np.nan
 
 
+def get_biomass_flux(model: cobra.Model, rxn_id: str) -> float:
+    """"""
+    model.objective = model.reactions.get_by_id(rxn_id)
+    model.objective_direction = "max"
+    flux = model.slim_optimize(error_value=0.0)
+    if abs(flux) < model.tolerance:
+        return 0.0
+    else:
+        return flux
+
+
 def close_boundaries_sensibly(model):
     """
     Return a cobra model with all boundaries closed and changed constraints.