complete transition to hybrid solver

docs/source/installing.ipynb

@@ -10,7 +10,7 @@
     "\n",
     "By default MICOM will also install a linear and quadratic programming solver that will work for large problems. For this it leverages\n",
     "a custom hybrid solver that combines [HIGHS](https://highs.dev/) and [OSQP](https://osqp.org/). Those will be installed along with MICOM\n",
-    "automatically on Windows and Linux, but require one additional step on MacOS (see below).\n",
+    "automatically.\n",
     "\n",
     "If you have several supported solvers installed you may later specify one with the `solver` argument in `micom.Community` or\n",
     "`micom.workflows.build`.\n",
@@ -28,10 +28,10 @@
     "Now install the CPLEX python package into your activated environment:\n",
     "\n",
     "```bash\n",
-    "pip install ibm/cplex/python/3.8/x86-64_linux\n",
+    "pip install ibm/cplex/python/3.10/x86-64_linux\n",
     "```\n",
     "\n",
-    "Substitute `3.8` with your Python version. Substitute `x86-64_linux` with the folder corresponding to your system (there will only be one subfolder in that directory).\n",
+    "Substitute `3.10` with your Python version. Substitute `x86-64_linux` with the folder corresponding to your system (there will only be one subfolder in that directory).\n",
     "\n",
     "**Gurobi**\n",
     "\n",
@@ -49,6 +49,11 @@
     "grbgetkey YOUR-LICENSE-KEY\n",
     "```"
    ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": []
   }
  ],
  "metadata": {

micom/community.py

@@ -139,7 +139,7 @@ def __init__(
         if not solver:
             solver = [
                 s
-                for s in ["cplex", "gurobi", "osqp", "glpk"]
+                for s in ["cplex", "gurobi", "hybrid", "glpk"]
                 if s in cobra.util.solver.solvers
             ][0]
         logger.info("using the %s solver." % solver)
@@ -499,13 +499,6 @@ def optimize(
         solution by default but only growth rates which is much quicker. You can
         request a full solution by setting the `fluxes` and `pFBA` arguments.
 
-        Note
-        ----
-        for OSQP this will use an additional quadratic regularization in order
-        to convert all LP problems to QPs and stabilize the solution. This may not
-        perform well if the community growth rate is larger 1000 (which it should never
-        be as this would be a doubling time of a few seconds).
-
         Parameters
         ----------
         slim : boolean, optional
@@ -883,7 +876,7 @@ def to_pickle(self, filename):
 
         """
         with open(filename, mode="wb") as out:
-            pickle.dump(self, out, protocol=4)
+            pickle.dump(self, out, protocol=pickle.HIGHEST_PROTOCOL)
 
     @cobra.Model.solver.setter
     def solver(self, s):

micom/problems.py

@@ -91,9 +91,7 @@ def cooperative_tradeoff(community, min_growth, fraction, fluxes, pfba, atol, rt
             com.variables.community_objective.lb = fr * min_growth
             com.variables.community_objective.ub = min_growth
             sol = solve(community, fluxes=fluxes, pfba=pfba, atol=atol, rtol=rtol)
-            # OSQP is better with QPs then LPs
-            # so it won't get better with the crossover
-            if not pfba and sol.status != OPTIMAL and solver != "osqp":
+            if not pfba and sol.status != OPTIMAL:
                 sol = crossover(
                     com, sol, fluxes=fluxes
                 )

micom/solution.py

@@ -114,10 +114,6 @@ def __init__(self, community, slim=False, reactions=None, metabolites=None):
         )
         self.members.index.name = "compartments"
         self.growth_rate = sum(community.abundances * gcs)
-        # Save estimated accuracy for osqp
-        if interface_to_str(community.problem) == "osqp":
-            self.primal_residual = community.solver.problem.info.pri_res
-            self.dual_residual = community.solver.problem.info.dua_res
 
     def _repr_html_(self):
         if self.status in good:
@@ -178,8 +174,6 @@ def add_pfba_objective(community, atol=1e-6, rtol=1e-6):
     community.objective = Zero
     community.objective_direction = "min"
     community.objective.set_linear_coefficients(dict.fromkeys(variables, 1.0))
-    if interface_to_str(community.solver.interface) == "osqp":
-        community.objective += 1e-6 * community.variables.community_objective ** 2
     if community.modification is None:
         community.modification = "pFBA"
     else:
@@ -189,15 +183,6 @@ def add_pfba_objective(community, atol=1e-6, rtol=1e-6):
 def solve(community, fluxes=True, pfba=True, raise_error=False, atol=1e-6, rtol=1e-6):
     """Get all fluxes stratified by taxa."""
     solver_name = interface_to_str(community.solver.interface)
-    term = None
-    if solver_name == "osqp" and community.objective.is_Linear:
-        # This improves OSQP by soooo much
-        term = (
-            1e-6
-            * community.solver.problem.direction
-            * community.variables.community_objective ** 2
-        )
-        community.objective += term
     community.solver.optimize()
     status = community.solver.status
     if status in good:
@@ -216,10 +201,6 @@ def solve(community, fluxes=True, pfba=True, raise_error=False, atol=1e-6, rtol=
             sol = CommunitySolution(community)
         else:
             sol = CommunitySolution(community, slim=True)
-        if term:
-            correction = 1e-6 * community.variables.community_objective.primal ** 2
-            sol.objective_value -= community.solver.problem.direction * correction
-            community.objective -= term
         return sol
     logger.warning("solver encountered an error %s" % status)
     return None
@@ -236,7 +217,7 @@ def reset_solver(community):
         community.solver.problem.reset()
     elif interface == "glpk":
         glp_adv_basis(community.solver.problem, 0)
-    elif interface == "osqp":
+    elif interface == "hybrid":
         community.solver.problem.reset()
 
 
@@ -268,7 +249,7 @@ def crossover(community, sol, fluxes=False):
         com.objective = com.scale * com.variables.community_objective
         for sp in com.taxa:
             const = com.constraints["objective_" + sp]
-            const.ub = gcs[sp]
+            const.ub = max(const.lb, gcs[sp])
         logger.info("finding closest feasible solution")
         s = com.optimize()
         if s is None:

micom/util.py

@@ -112,8 +112,9 @@ def serialize_models(files, dir="."):
         model = load_model(f)
         logger.info("serializing {}".format(f))
         pickle.dump(
-            model, open(path.join(dir, fname + ".pickle"), "wb"), protocol=2
-        )  # required for Python 2 compat
+            model, open(path.join(dir, fname + ".pickle"), "wb"),
+            protocol=pickle.HIGHEST_PROTOCOL
+        )
 
 
 def chr_or_input(m):

setup.cfg

@@ -12,7 +12,7 @@ url = https://github.com/micom-dev/micom
 author = Christian Diener
 author_email = [email protected]
 license = Apache License 2.0
-classifiers = 
+classifiers =
 	Development Status :: 5 - Production/Stable
 	Intended Audience :: Science/Research
 	Topic :: Scientific/Engineering :: Bio-Informatics
@@ -23,7 +23,7 @@ classifiers =
 	Programming Language :: Python :: 3.9
 	Programming Language :: Python :: 3.10
 	Programming Language :: Python :: 3.11
-keywords = 
+keywords =
 	microbiome
 	microbiota
 	modeling
@@ -35,26 +35,27 @@ keywords =
 packages = find:
 zip_safe = True
 python_requires = >=3.7
-install_requires = 
-	cobra>=0.26.2
+install_requires =
+	cobra>=0.29.0
 	jinja2>=2.10.0
-	scikit-learn>=0.24.1
+	scikit-learn>=1.4.2
 	scipy>=1.8.1
 	symengine>=0.6.1
 	osqp>=0.6.2
+	optlang>=1.8.1
 	highspy>=1.7.1.dev1
-tests_require = 
+tests_require =
 	coverage
 	pytest
 	pytest-cov
 
 [options.packages.find]
-exclude = 
+exclude =
 	docs
 	tests
 
 [options.package_data]
-micom = 
+micom =
 	data/*.csv
 	data/*.gz
 	data/templates/*.*
@@ -69,7 +70,7 @@ search = __version__ = "{current_version}"
 replace = __version__ = "{new_version}"
 
 [tool:pytest]
-filterwarnings = 
+filterwarnings =
 	ignore::DeprecationWarning
 	ignore::FutureWarning
 
@@ -81,23 +82,23 @@ source = micom
 branch = True
 
 [coverage:report]
-exclude_lines = 
+exclude_lines =
 	pragma: no cover
-	
+
 	def __repr__
 	if self\.debug
-	
+
 	raise AssertionError
 	raise NotImplementedError
-	
+
 	if 0:
 	if __name__ == .__main__.:
 ignore_errors = True
 
 [flake8]
 exclude = tests/*, docs/*, versioneer.py, _version.py
 max-line-length = 88
-ignore = 
+ignore =
 	E203
 	W503
 	D202