From a036c557d82861b4ceb677cf2582a91eb8f06067 Mon Sep 17 00:00:00 2001 From: Andrew Lee Date: Fri, 6 Oct 2023 13:01:38 -0400 Subject: [PATCH] Better unit testing of degeenracy hunter methods --- idaes/core/util/model_diagnostics.py | 50 +++++++++------- .../core/util/tests/test_model_diagnostics.py | 60 +++++++++++++++++++ 2 files changed, 89 insertions(+), 21 deletions(-) diff --git a/idaes/core/util/model_diagnostics.py b/idaes/core/util/model_diagnostics.py index 2cdc44f2c5..f5594bfe65 100644 --- a/idaes/core/util/model_diagnostics.py +++ b/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 diff --git a/idaes/core/util/tests/test_model_diagnostics.py b/idaes/core/util/tests/test_model_diagnostics.py index a742014ac2..af6d6131c1 100644 --- a/idaes/core/util/tests/test_model_diagnostics.py +++ b/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")