Better unit testing of degeenracy hunter methods

idaes/core/util/model_diagnostics.py

@@ -1779,6 +1779,19 @@ def eq_abs_upper(b, c):
 
         self.candidates_milp = m_dh
 
+    def _identify_candidates(self):
+        eq_con_list = self.nlp.get_pyomo_equality_constraints()
+
+        for i in self.candidates_milp.C:
+            # Check if constraint is included
+            if self.candidates_milp.abs_nu[i]() > self.config.m_small * MMULT:
+                # If it is, save the value of nu
+                if eq_con_list is None:
+                    name = i
+                else:
+                    name = eq_con_list[i]
+                self.degenerate_set[name] = self.candidates_milp.nu[i]()
+
     def _solve_candidates_milp(self, tee: bool = False):
         """Solve MILP to generate set of candidate equations
 
@@ -1789,23 +1802,13 @@ def _solve_candidates_milp(self, tee: bool = False):
         """
         _log.info("Solving Candidates MILP model.")
 
-        eq_con_list = self.nlp.get_pyomo_equality_constraints()
-
         results = self._get_solver().solve(self.candidates_milp, tee=tee)
 
         self.degenerate_set = {}
 
         if check_optimal_termination(results):
             # We found a degenerate set
-            for i in self.candidates_milp.C:
-                # Check if constraint is included
-                if self.candidates_milp.abs_nu[i]() > self.config.m_small * MMULT:
-                    # If it is, save the value of nu
-                    if eq_con_list is None:
-                        name = i
-                    else:
-                        name = eq_con_list[i]
-                    self.degenerate_set[name] = self.candidates_milp.nu[i]()
+            self._identify_candidates()
         else:
             _log.debug(
                 "Solver did not return an optimal termination condition for "
@@ -1867,6 +1870,20 @@ def eq_upper(m_dh, c):
 
         self.ids_milp = m_dh
 
+    def _get_ids(self):
+        # Create empty dictionary
+        ids_ = {}
+
+        eq_con_list = self.nlp.get_pyomo_equality_constraints()
+
+        for i in self.ids_milp.C:
+            # Check if constraint is included
+            if self.ids_milp.y[i]() > YTOL:
+                # If it is, save the value of nu
+                ids_[eq_con_list[i]] = self.ids_milp.nu[i]()
+
+        return ids_
+
     def _solve_ids_milp(self, cons: Constraint, tee: bool = False):
         """Solve MILP to check if equation 'cons' is a significant component
         in an irreducible degenerate set
@@ -1891,24 +1908,15 @@ def _solve_ids_milp(self, cons: Constraint, tee: bool = False):
 
         self.ids_milp.nu[cons_idx].unfix()
 
-        # Create empty dictionary
-        ids_ = {}
-
         if check_optimal_termination(results):
             # We found an irreducible degenerate set
-            for i in self.ids_milp.C:
-                # Check if constraint is included
-                if self.ids_milp.y[i]() > YTOL:
-                    # If it is, save the value of nu
-                    ids_[eq_con_list[i]] = self.ids_milp.nu[i]()
+            return self._get_ids()
         else:
             raise ValueError(
                 f"Solver did not return an optimal termination condition for "
                 f"IDS MILP with constraint {cons.name}."
             )
 
-        return ids_
-
     def find_irreducible_degenerate_sets(self, tee=False):
         """
         Compute irreducible degenerate sets

idaes/core/util/tests/test_model_diagnostics.py

@@ -1756,6 +1756,30 @@ def test_prepare_candidates_milp(self, model):
 
         assert isinstance(dh.candidates_milp, ConcreteModel)
 
+    @pytest.mark.unit
+    def test_identify_candidates(self, model):
+        dh = DegeneracyHunter2(model)
+        dh._prepare_candidates_milp()
+
+        dh.candidates_milp.nu[0].set_value(-1e-05)
+        dh.candidates_milp.nu[1].set_value(1e-05)
+
+        dh.candidates_milp.y_pos[0].set_value(0)
+        dh.candidates_milp.y_pos[1].set_value(1)
+
+        dh.candidates_milp.y_neg[0].set_value(-0)
+        dh.candidates_milp.y_neg[1].set_value(-0)
+
+        dh.candidates_milp.abs_nu[0].set_value(1e-05)
+        dh.candidates_milp.abs_nu[1].set_value(1e-05)
+
+        dh._identify_candidates()
+
+        assert dh.degenerate_set == {
+            model.con2: -1e-05,
+            model.con5: 1e-05,
+        }
+
     @pytest.mark.solver
     @pytest.mark.component
     @pytest.mark.skipif(not solver_available, reason="SCIP is not available")
@@ -1764,6 +1788,18 @@ def test_solve_candidates_milp(self, model):
         dh._prepare_candidates_milp()
         dh._solve_candidates_milp()
 
+        assert value(dh.candidates_milp.nu[0]) == pytest.approx(-1e-05, rel=1e-5)
+        assert value(dh.candidates_milp.nu[1]) == pytest.approx(1e-05, rel=1e-5)
+
+        assert value(dh.candidates_milp.y_pos[0]) == pytest.approx(0, abs=1e-5)
+        assert value(dh.candidates_milp.y_pos[1]) == pytest.approx(1, rel=1e-5)
+
+        assert value(dh.candidates_milp.y_neg[0]) == pytest.approx(-0, abs=1e-5)
+        assert value(dh.candidates_milp.y_neg[1]) == pytest.approx(-0, abs=1e-5)
+
+        assert value(dh.candidates_milp.abs_nu[0]) == pytest.approx(1e-05, rel=1e-5)
+        assert value(dh.candidates_milp.abs_nu[1]) == pytest.approx(1e-05, rel=1e-5)
+
         assert dh.degenerate_set == {
             model.con2: -1e-05,
             model.con5: 1e-05,
@@ -1776,6 +1812,24 @@ def test_prepare_ids_milp(self, model):
 
         assert isinstance(dh.ids_milp, ConcreteModel)
 
+    @pytest.mark.unit
+    def test_solve_ids_milp(self, model):
+        dh = DegeneracyHunter2(model)
+        dh._prepare_ids_milp()
+
+        dh.ids_milp.nu[0].set_value(1)
+        dh.ids_milp.nu[1].set_value(-1)
+
+        dh.ids_milp.y[0].set_value(1)
+        dh.ids_milp.y[1].set_value(1)
+
+        ids_ = dh._get_ids()
+
+        assert ids_ == {
+            model.con2: 1,
+            model.con5: -1,
+        }
+
     @pytest.mark.solver
     @pytest.mark.component
     @pytest.mark.skipif(not solver_available, reason="SCIP is not available")
@@ -1789,6 +1843,12 @@ def test_solve_ids_milp(self, model):
             model.con5: -1,
         }
 
+        assert value(dh.ids_milp.nu[0]) == pytest.approx(1, rel=1e-5)
+        assert value(dh.ids_milp.nu[1]) == pytest.approx(-1, rel=1e-5)
+
+        assert value(dh.ids_milp.y[0]) == pytest.approx(1, rel=1e-5)
+        assert value(dh.ids_milp.y[1]) == pytest.approx(1, rel=1e-5)
+
     @pytest.mark.solver
     @pytest.mark.component
     @pytest.mark.skipif(not solver_available, reason="SCIP is not available")