Fix two knock-in defects - This allows computation of EFVs with StrainDesign

straindesign/compute_strain_designs.py

@@ -767,7 +767,7 @@ def compute_strain_designs(model: Model, **kwargs: dict) -> SDSolutions:
 
     # remove ko-costs (and thus knockability) of essential reactions
     [cmp_ko_cost.pop(er) for er in essential_reacs if er in cmp_ko_cost]
-    essential_kis = set(cmp_ki_cost[er] for er in essential_reacs if er in cmp_ki_cost)
+    essential_kis = set(er for er in essential_reacs if er in cmp_ki_cost)
     # Build MILP
     kwargs1 = kwargs
     kwargs1[KOCOST] = cmp_ko_cost

straindesign/strainDesignMILP.py

@@ -275,6 +275,17 @@ def verify_sd(self, sols) -> List:
                             if np.logical_xor(sol[0,z_i],sense==-1) ]
             active_eqs = [i for i in range(self.cont_MILP.z_map_constr_eq.shape[1]) if i not in inactive_eqs]
 
+            # Zeroing a variable whose bounds exclude zero contradicts that bound, so the
+            # region is empty whatever the remaining system does. This has to be tested
+            # here rather than left to the LP: prevent_boundary_knockouts keeps such a
+            # bound (e.g. ATPM >= 3.15) as a row with no z-mapping so that a knockout
+            # contradicts it, but reassign_lb_ub_from_ineq later folds single-variable
+            # rows back into variable bounds, and those vanish together with the column.
+            if any(self.cont_MILP.lb[j] > 0.0 or self.cont_MILP.ub[j] < 0.0 for j in inactive_vars):
+                valid[i] = False
+                continue
+            # Otherwise drop the columns outright. Absence is a stronger statement than an
+            # interval of [0, 0], since it owes nothing to feasibility tolerances.
             lp = MILP_LP(A_ineq=self.cont_MILP.A_ineq[active_ineqs, :][:, active_vars],
                          b_ineq=[self.cont_MILP.b_ineq[i] for i in active_ineqs],
                          A_eq=self.cont_MILP.A_eq[active_eqs, :][:, active_vars],

straindesign/strainDesignProblem.py

@@ -86,9 +86,10 @@ class SDProblem:
             the big-M method by default (with COBRA standard M=1000). M should be chosen 'sufficiently large' 
             to avoid computational artifacts and 'sufficiently small' to avoid numerical issues.
 
-        essential_kis (optional (set)):
-            A set of reactions that are marked as addable and that are essential for at least one of the
-            strain design modules. Providing such "essential knock-ins" may speed up the strain design computation.
+        essential_kis (optional (set of str)):
+            Reaction identifiers that are marked as addable and that are essential for at least one of the
+            strain design modules. Their intervention binaries are fixed to 1, since a solution that omits
+            them cannot satisfy the module they are essential for.
 
     Returns:
         (SDProblem):
@@ -167,7 +168,7 @@ def __init__(self, model: Model, sd_modules: List[SDModule], *args, **kwargs):
         else:
             self.b_ineq = [0.0, float(self.max_cost), np.inf]
         self.z_map_constr_ineq = sparse.csc_matrix((numr, 3))
-        self.lb = [1.0 if r in self.essential_kis else 0.0 for r in model.reactions]
+        self.lb = [1.0 if r.id in self.essential_kis else 0.0 for r in model.reactions]
         self.ub = [1.0 - float(i) for i in self.z_non_targetable]
         self.idx_z = [i for i in range(0, numr)]
         self.c = [0.0] * numr