implemented Solution.infeasible

src/csnlp/core/solutions.py

@@ -6,6 +6,7 @@
 from itertools import product as _product
 from typing import TYPE_CHECKING
 from typing import Any as _Any
+from typing import Optional
 from typing import Protocol as _Protocol
 from typing import TypeVar as _TypeVar
 from typing import Union
@@ -136,6 +137,79 @@ def success(self) -> bool:
         return self.stats["success"]
 
     @_cached_property
+    def infeasible(self) -> Optional[bool]:
+        r"""Gets whether the solver status indicates infeasibility. If ``False``, it
+        does not imply feasibility as the solver or its CasADi interface may have not
+        detect it.
+
+        For different solvers, the infeasibility status is stored in different ways.
+        Here is a list of what I gathered so far. The solvers are grouped based on the
+        type of problem they solve. An (F) next to the solver's name indicates that the
+        the solver will crash the program if ``"error_on_fail": True`` and the solver
+        fails. The ``return_status`` and ``unified_return_status`` can both be found in
+        the solver's stats, or in this solution object.
+
+        * NLPs
+
+          - **fatrop**: unclear; better to return ``None`` for now
+          - **ipopt**: ``return_status == "Infeasible_Problem_Detected"``
+          - **qrsqp** (F): ``return_status == "Search_Direction_Becomes_Too_Small"``
+            (dubious)
+          - **sqpmethod** (F): ``return_status == "Search_Direction_Becomes_Too_Small"``
+            (dubious)
+
+        * QPs
+
+          - **ipqp** (F): no clear way to detect infeasibility; return ``None`` for now
+          - **osqp** (F): ``unified_return_status == "SOLVER_RET_INFEASIBLE"`` or
+            ``"infeasible" in return_status``
+          - **proxqp** (F): ``return_status == "PROXQP_PRIMAL_INFEASIBLE"`` or
+            ``return_status == "PROXQP_DUAL_INFEASIBLE"``
+          - **qpoases** (F): ``"infeasib" in return_status``
+          - **qrqp** (F): ``return_status == "Failed to calculate search direction"``
+
+        * LPs
+
+          - **clp** (F): ``return_status == "primal infeasible"`` or
+            ``return_status == "dual infeasible"``
+
+        * Mixed-Iteger Problems (MIPs)
+
+          - **bonmin** (F): ``return_status == "INFEASIBLE"``
+          - **cbc** (F): ``"not feasible" in return_status``
+          - **gurobi** (F): ``return_status == "INFEASIBLE"``
+        """
+        solver_plugin = self.solver_plugin
+        # NLPs
+        if solver_plugin == "ipopt":
+            return self.status == "Infeasible_Problem_Detected"
+        if solver_plugin in ("qrsqp", "sqpmethod"):
+            return self.status == "Search_Direction_Becomes_Too_Small"
+        # QPs
+        if solver_plugin == "osqp":
+            return (
+                self.unified_return_status == "SOLVER_RET_INFEASIBLE"
+                or "infeasible" in self.status
+            )
+        if solver_plugin == "proxqp":
+            return (
+                self.status == "PROXQP_PRIMAL_INFEASIBLE"
+                or self.status == "PROXQP_DUAL_INFEASIBLE"
+            )
+        if solver_plugin == "qpoases":
+            return "infeasib" in self.status
+        if solver_plugin == "qrqp":
+            return self.status == "Failed to calculate search direction"
+        # LPs
+        if solver_plugin == "clp":
+            return "infeasible" in self.status
+        # MIPs
+        if solver_plugin in ("bonmin", "gurobi"):
+            return self.status == "INFEASIBLE"
+        if solver_plugin == "cbc":
+            return "not feasible" in self.status
+        return None
+
     @property
     def barrier_parameter(self) -> float:
         """Gets the IPOPT barrier parameter at the optimal solution"""

tests/test_core.py

@@ -9,6 +9,7 @@
 from casadi.tools import entry, struct_MX, struct_SX
 from parameterized import parameterized
 
+from csnlp import Nlp
 from csnlp.core.cache import invalidate_cache
 from csnlp.core.data import array2cs, cs2array, find_index_in_vector
 from csnlp.core.debug import NlpDebug, NlpDebugEntry
@@ -346,6 +347,94 @@ def test_cmp_key__returns_correct_solution(self):
                 min_index, _ = min(sols, key=lambda sol: DummySolution.cmp_key(sol[1]))
                 self.assertEqual(min_index, case["expected"])
 
+    @parameterized.expand(
+        product(
+            ("MX", "SX"),
+            (True, False),
+            [
+                # (
+                #     "nlp",
+                #     "sqpmethod",
+                #     {
+                #         "error_on_fail": False,
+                #         "print_time": False,
+                #         "print_status": False,
+                #         "print_header": False,
+                #         "print_iteration": False,
+                #         "qpsol_options": {
+                #             "error_on_fail": False, "printLevel": "none"
+                #         },
+                #     },
+                # ),
+                (
+                    "nlp",
+                    "ipopt",
+                    {
+                        "print_time": False,
+                        "ipopt": {"print_level": 0, "sb": "yes"},
+                    },
+                ),
+                # (
+                #     "conic",
+                #     "osqp",
+                #     {
+                #         "error_on_fail": False,
+                #         "print_time": False,
+                #         "osqp": {"verbose": False},
+                #     },
+                # ),
+                ("conic", "qpoases", {"error_on_fail": False, "printLevel": "none"}),
+                # ("conic", "proxqp", {"error_on_fail": False}),
+                (
+                    "conic",
+                    "qrqp",
+                    {
+                        "error_on_fail": False,
+                        "print_time": False,
+                        "print_header": False,
+                        "print_info": False,
+                        "print_iter": False,
+                    },
+                ),
+                ("conic", "clp", {"error_on_fail": False}),
+                (
+                    "nlp",
+                    "bonmin",
+                    {
+                        "print_time": False,
+                        "bonmin": {
+                            "fp_log_level": 0,
+                            "lp_log_level": 0,
+                            "milp_log_level": 0,
+                            "nlp_log_level": 0,
+                            "oa_cuts_log_level": 0,
+                            "oa_log_level": 0,
+                        },
+                    },
+                ),
+                # ("conic", "cbc", {"error_on_fail": False}),
+            ],
+        )
+    )
+    def test_infeasible(self, sym_type, is_feas, solver_data):
+        solver_type, solver, solver_options = solver_data
+
+        prob = Nlp(sym_type=sym_type)
+        discrete = solver in ("bonmin", "cbc", "gurobi", "knitro")
+        x, _, _ = prob.variable("x", discrete=discrete)
+        lb = -abs(np.random.randn())
+        ub = abs(np.random.randn())
+        if not is_feas:
+            lb, ub = ub, lb
+        prob.constraint("lb", x, ">=", lb)
+        prob.constraint("ub", x, "<=", ub)
+        prob.minimize(x)
+        prob.init_solver(solver_options, solver, solver_type)
+        sol = prob.solve()
+
+        self.assertEqual(sol.success, is_feas)
+        self.assertEqual(not sol.infeasible, is_feas)
+
 
 class TestData(unittest.TestCase):
     @parameterized.expand(product([cs.MX, cs.SX], [(1, 1), (3, 1), (1, 3), (3, 3)]))