Feat/956 compute iis highs

doc/release_notes.rst

@@ -271,6 +271,7 @@ Version 0.7.0
 
 * OETC: ``Model.solve()`` forwards solver options to the handler; ``OetcSettings.from_env()`` reads ``OETC_*``.
 * SCIP supports quadratic problems on Windows.
+* ``Model.compute_infeasibilities`` now also supports HiGHS, via ``Highs.getIis``.
 
 **Performance**
 

linopy/model.py

@@ -2348,9 +2348,9 @@ def compute_infeasibilities(self) -> list[int]:
         """
         Compute a set of infeasible constraints.
 
-        This function requires that the model was solved with `gurobi` or `xpress`
-        and the termination condition was infeasible. The solver must have detected
-        the infeasibility during the solve process.
+        This function requires that the model was solved with `gurobi`, `xpress`
+        or `highs` and the termination condition was infeasible. The solver must
+        have detected the infeasibility during the solve process.
 
         Returns
         -------
@@ -2383,6 +2383,16 @@ def compute_infeasibilities(self) -> list[int]:
             except ImportError:
                 pass
 
+        # Check for HiGHS
+        if "highs" in available_solvers:
+            try:
+                import highspy
+
+                if solver_model is not None and isinstance(solver_model, highspy.Highs):
+                    return self._compute_infeasibilities_highs(solver_model)
+            except ImportError:
+                pass
+
         # If we get here, either the solver doesn't support IIS or no solver model is available
         if solver_model is None:
             # Check if this is a supported solver without a stored model
@@ -2398,12 +2408,12 @@ def compute_infeasibilities(self) -> list[int]:
                 # This is an unsupported solver
                 raise NotImplementedError(
                     f"Computing infeasibilities is not supported for '{solver_name}' solver. "
-                    "Only Gurobi and Xpress solvers support IIS computation."
+                    "Only Gurobi, Xpress and HiGHS solvers support IIS computation."
                 )
         else:
             # We have a solver model but it's not a supported type
             raise NotImplementedError(
-                "Computing infeasibilities is only supported for Gurobi and Xpress solvers. "
+                "Computing infeasibilities is only supported for Gurobi, Xpress and HiGHS solvers. "
                 f"Current solver model type: {type(solver_model).__name__}"
             )
 
@@ -2529,12 +2539,51 @@ def _extract_iis_constraints(self, solver_model: Any, iis_num: int) -> list[Any]
 
         return miisrow
 
+    def _compute_infeasibilities_highs(self, solver_model: Any) -> list[int]:
+        """Compute infeasibilities for the HiGHS solver."""
+        if not hasattr(solver_model, "getIis"):
+            raise NotImplementedError(
+                "Computing infeasibilities requires a `highspy` version that "
+                "supports `Highs.getIis` (HiGHS IIS computation). "
+                "Please upgrade the `highspy` package."
+            )
+
+        import highspy
+
+        solver = self.solver
+        assert solver is not None
+        if "iis_strategy" not in solver.solver_options:
+            solver_model.setOptionValue(
+                "iis_strategy",
+                int(highspy.IisStrategy.kIisStrategyFromLp)
+                | int(highspy.IisStrategy.kIisStrategyIrreducible),
+            )
+        status, iis = solver_model.getIis()
+        if status == highspy.HighsStatus.kError or not iis.valid_:
+            raise RuntimeError(
+                "HiGHS failed to compute an irreducible infeasible subsystem (IIS)."
+            )
+
+        row_index = np.asarray(iis.row_index_, dtype=np.intp)
+        if not len(row_index):
+            return []
+
+        if solver.io_api == "direct":
+            clabels = self.constraints.label_index.clabels
+        else:
+            from linopy.solvers import _names_to_labels
+
+            clabels = _names_to_labels(solver_model.getLp().row_names_)
+
+        labels = {int(clabels[pos]) for pos in row_index if clabels[pos] >= 0}
+        return sorted(labels)
+
     def format_infeasibilities(self, display_max_terms: int | None = None) -> str:
         """
         Return a string representation of infeasible constraints.
 
-        This function requires that the model was solved using `gurobi` or `xpress`
-        and the termination condition was infeasible.
+        This function requires that the model was solved using `gurobi`, `xpress`
+        or `highs` and the termination condition was infeasible.
 
         Parameters
         ----------
@@ -2573,8 +2622,8 @@ def compute_set_of_infeasible_constraints(self) -> Dataset:
         """
         Compute a set of infeasible constraints.
 
-        This function requires that the model was solved with `gurobi` or `xpress` and the
-        termination condition was infeasible.
+        This function requires that the model was solved with `gurobi`, `xpress` or `highs`
+        and the termination condition was infeasible.
 
         Returns
         -------

linopy/solvers.py

@@ -1575,6 +1575,7 @@ class Highs(Solver[None]):
             SolverFeature.SOLUTION_FILE_NOT_NEEDED,
             SolverFeature.SEMI_CONTINUOUS_VARIABLES,
             SolverFeature.MIP_DUAL_BOUND_REPORT,
+            SolverFeature.IIS_COMPUTATION,
         }
     )
     supports_persistent_update: ClassVar[bool] = True

test/test_infeasibility.py

@@ -75,7 +75,7 @@ def multi_dimensional_infeasible_model(self) -> Model:
 
         return m
 
-    @pytest.mark.parametrize("solver", ["gurobi", "xpress"])
+    @pytest.mark.parametrize("solver", ["gurobi", "xpress", "highs"])
     def test_simple_infeasibility_detection(
         self, simple_infeasible_model: Model, solver: str
     ) -> None:
@@ -98,7 +98,7 @@ def test_simple_infeasibility_detection(
         assert isinstance(formatted, str)
         assert formatted
 
-    @pytest.mark.parametrize("solver", ["gurobi", "xpress"])
+    @pytest.mark.parametrize("solver", ["gurobi", "xpress", "highs"])
     def test_complex_infeasibility_detection(
         self, complex_infeasible_model: Model, solver: str
     ) -> None:
@@ -121,7 +121,7 @@ def test_complex_infeasibility_detection(
         # We expect at least 2 constraints to be involved
         assert len(labels) >= 2
 
-    @pytest.mark.parametrize("solver", ["gurobi", "xpress"])
+    @pytest.mark.parametrize("solver", ["gurobi", "xpress", "highs"])
     def test_multi_dimensional_infeasibility(
         self, multi_dimensional_infeasible_model: Model, solver: str
     ) -> None:
@@ -152,7 +152,7 @@ def test_unsolved_model_error(self) -> None:
         ):
             m.compute_infeasibilities()
 
-    @pytest.mark.parametrize("solver", ["gurobi", "xpress"])
+    @pytest.mark.parametrize("solver", ["gurobi", "xpress", "highs"])
     def test_no_solver_model_error(self, solver: str) -> None:
         """Test error when solver model is not available after solving."""
         if solver not in available_solvers:
@@ -173,7 +173,7 @@ def test_no_solver_model_error(self, solver: str) -> None:
         with pytest.raises(ValueError, match="No solver model available"):
             m.compute_infeasibilities()
 
-    @pytest.mark.parametrize("solver", ["gurobi", "xpress"])
+    @pytest.mark.parametrize("solver", ["gurobi", "xpress", "highs"])
     def test_feasible_model_iis(self, solver: str) -> None:
         """Test IIS computation on a feasible model."""
         if solver not in available_solvers:
@@ -220,7 +220,7 @@ def test_unsupported_solver_error(self) -> None:
             with pytest.raises(NotImplementedError):
                 m.compute_infeasibilities()
 
-    @pytest.mark.parametrize("solver", ["gurobi", "xpress"])
+    @pytest.mark.parametrize("solver", ["gurobi", "xpress", "highs"])
     def test_deprecated_method(
         self, simple_infeasible_model: Model, solver: str
     ) -> None:
@@ -246,7 +246,7 @@ def test_deprecated_method(
         # Check that it contains constraint labels
         assert len(subset) > 0
 
-    @pytest.mark.parametrize("solver", ["gurobi", "xpress"])
+    @pytest.mark.parametrize("solver", ["gurobi", "xpress", "highs"])
     def test_masked_constraint_infeasibility(
         self, solver: str, capsys: pytest.CaptureFixture[str]
     ) -> None:

test/test_optimization.py

@@ -698,7 +698,7 @@ def test_infeasible_model(
     assert status == "warning"
     assert "infeasible" in condition
 
-    if solver in ["gurobi", "xpress"]:
+    if solver in ["gurobi", "xpress", "highs"]:
         # ignore deprecated warning
         with pytest.warns(DeprecationWarning):
             model.compute_set_of_infeasible_constraints()