Merge branch 'main' into mscontactor_dynamics

idaes/core/util/model_diagnostics.py

@@ -1450,6 +1450,10 @@ def __init__(self, model: _BlockData, **kwargs):
             self._model, scaled=False, equality_constraints_only=True
         )
 
+        # Get list of equality constraint and variable names
+        self._eq_con_list = self.nlp.get_pyomo_equality_constraints()
+        self._var_list = self.nlp.get_pyomo_variables()
+
         if self.jacobian.shape[0] < 2:
             raise ValueError(
                 "Model needs at least 2 equality constraints to perform svd_analysis."
@@ -1562,10 +1566,6 @@ def display_underdetermined_variables_and_constraints(
 
         tol = self.config.size_cutoff_in_singular_vector
 
-        # Get list of equality constraint and variable names
-        eq_con_list = self.nlp.get_pyomo_equality_constraints()
-        var_list = self.nlp.get_pyomo_variables()
-
         if singular_values is None:
             singular_values = range(1, len(self.s) + 1)
 
@@ -1588,11 +1588,11 @@ def display_underdetermined_variables_and_constraints(
             stream.write(f"{TAB}Smallest Singular Value {e}:\n\n")
             stream.write(f"{2 * TAB}Variables:\n\n")
             for v in np.where(abs(self.v[:, e - 1]) > tol)[0]:
-                stream.write(f"{3 * TAB}{var_list[v].name}\n")
+                stream.write(f"{3 * TAB}{self._var_list[v].name}\n")
 
             stream.write(f"\n{2 * TAB}Constraints:\n\n")
             for c in np.where(abs(self.u[:, e - 1]) > tol)[0]:
-                stream.write(f"{3 * TAB}{eq_con_list[c].name}\n")
+                stream.write(f"{3 * TAB}{self._eq_con_list[c].name}\n")
             stream.write("\n")
 
         stream.write("=" * MAX_STR_LENGTH + "\n")
@@ -1620,23 +1620,19 @@ def display_constraints_including_variable(self, variable, stream=None):
                 f"object (got {variable})."
             )
 
-        # Get list of equality constraint and variable names
-        eq_con_list = self.nlp.get_pyomo_equality_constraints()
-        var_list = self.nlp.get_pyomo_variables()
-
         # Get index of variable in Jacobian
         try:
-            var_idx = var_list.index(variable)
-        except (ValueError, PyomoException):
+            var_idx = self.nlp.get_primal_indices([variable])[0]
+        except (KeyError, PyomoException):
             raise AttributeError(f"Could not find {variable.name} in model.")
 
-        nonzeroes = self.jacobian[:, var_idx].nonzero()
+        nonzeros = self.jacobian.getcol(var_idx).nonzero()
 
         # Build a list of all constraints that include var
         cons_w_var = []
-        for i, r in enumerate(nonzeroes[0]):
+        for r in nonzeros[0]:
             cons_w_var.append(
-                f"{eq_con_list[r].name}: {self.jacobian[(r, nonzeroes[1][i])]}"
+                f"{self._eq_con_list[r].name}: {self.jacobian[(r, var_idx)]:.3e}"
             )
 
         # Write the output
@@ -1671,23 +1667,20 @@ def display_variables_in_constraint(self, constraint, stream=None):
                 f"object (got {constraint})."
             )
 
-        # Get list of equality constraint and variable names
-        eq_con_list = self.nlp.get_pyomo_equality_constraints()
-        var_list = self.nlp.get_pyomo_variables()
-
         # Get index of variable in Jacobian
         try:
-            con_idx = eq_con_list.index(constraint)
-        except ValueError:
+            con_idx = self.nlp.get_constraint_indices([constraint])[0]
+        except KeyError:
             raise AttributeError(f"Could not find {constraint.name} in model.")
 
-        nonzeroes = self.jacobian[con_idx, :].nonzero()
+        nonzeros = self.jacobian[con_idx, :].nonzero()
 
         # Build a list of all vars in constraint
         vars_in_cons = []
-        for i, r in enumerate(nonzeroes[0]):
-            c = nonzeroes[1][i]
-            vars_in_cons.append(f"{var_list[c].name}: {self.jacobian[(r, c)]}")
+        for c in nonzeros[1]:
+            vars_in_cons.append(
+                f"{self._var_list[c].name}: {self.jacobian[(con_idx, c)]:.3e}"
+            )
 
         # Write the output
         _write_report_section(
@@ -2281,8 +2274,8 @@ def __init__(self, block_or_jac, solver=None):
             self.jac_eq = get_jacobian(self.block, equality_constraints_only=True)[0]
 
             # Create a list of equality constraint names
-            self.eq_con_list = self.nlp.get_pyomo_equality_constraints()
-            self.var_list = self.nlp.get_pyomo_variables()
+            self._eq_con_list = self.nlp.get_pyomo_equality_constraints()
+            self._var_list = self.nlp.get_pyomo_variables()
 
             self.candidate_eqns = None
 
@@ -2293,7 +2286,7 @@ def __init__(self, block_or_jac, solver=None):
 
             # # TODO: Need to refactor, document, and test support for Jacobian
             # self.jac_eq = block_or_jac
-            # self.eq_con_list = None
+            # self._eq_con_list = None
 
         else:
             raise TypeError("Check the type for 'block_or_jac'")
@@ -2813,11 +2806,11 @@ def underdetermined_variables_and_constraints(self, n_calc=1, tol=0.1, dense=Fal
             )
         print("Column:    Variable")
         for i in np.where(abs(self.v[:, n_calc - 1]) > tol)[0]:
-            print(str(i) + ": " + self.var_list[i].name)
+            print(str(i) + ": " + self._var_list[i].name)
         print("")
         print("Row:    Constraint")
         for i in np.where(abs(self.u[:, n_calc - 1]) > tol)[0]:
-            print(str(i) + ": " + self.eq_con_list[i].name)
+            print(str(i) + ": " + self._eq_con_list[i].name)
 
     def find_candidate_equations(self, verbose=True, tee=False):
         """
@@ -2842,7 +2835,7 @@ def find_candidate_equations(self, verbose=True, tee=False):
         if verbose:
             print("Solving MILP model...")
         ce, ds = self._find_candidate_eqs(
-            self.candidates_milp, self.solver, self.eq_con_list, tee
+            self.candidates_milp, self.solver, self._eq_con_list, tee
         )
 
         if ce is not None:
@@ -2883,7 +2876,7 @@ def find_irreducible_degenerate_sets(self, verbose=True, tee=False):
 
                 # Check if equation 'c' is a major element of an IDS
                 ids_ = self._check_candidate_ids(
-                    self.dh_milp, self.solver, c, self.eq_con_list, tee
+                    self.dh_milp, self.solver, c, self._eq_con_list, tee
                 )
 
                 if ids_ is not None:

idaes/core/util/tests/test_model_diagnostics.py

@@ -467,7 +467,7 @@ def model(self):
         m.b.v2 = Var(units=units.m)
         m.b.v3 = Var(bounds=(0, 5))
         m.b.v4 = Var()
-        m.b.v5 = Var(bounds=(0, 1))
+        m.b.v5 = Var(bounds=(0, 5))
         m.b.v6 = Var()
         m.b.v7 = Var(
             units=units.m, bounds=(0, 1)
@@ -550,7 +550,6 @@ def test_display_variables_at_or_outside_bounds(self, model):
 The following variable(s) have values at or outside their bounds (tol=0.0E+00):
 
     b.v3 (free): value=0.0 bounds=(0, 5)
-    b.v5 (fixed): value=2 bounds=(0, 1)
 
 ====================================================================================
 """
@@ -619,7 +618,6 @@ def test_display_variables_near_bounds(self, model):
 The following variable(s) have values close to their bounds (abs=1.0E-04, rel=1.0E-04):
 
     b.v3: value=0.0 bounds=(0, 5)
-    b.v5: value=2 bounds=(0, 1)
     b.v7: value=1.0000939326524314e-07 bounds=(0, 1)
 
 ====================================================================================
@@ -1009,7 +1007,7 @@ def test_collect_numerical_warnings(self, model):
 
         assert len(warnings) == 2
         assert "WARNING: 1 Constraint with large residuals (>1.0E-05)" in warnings
-        assert "WARNING: 2 Variables at or outside bounds (tol=0.0E+00)" in warnings
+        assert "WARNING: 1 Variable at or outside bounds (tol=0.0E+00)" in warnings
 
         assert len(next_steps) == 2
         assert "display_constraints_with_large_residuals()" in next_steps
@@ -1070,10 +1068,9 @@ def test_collect_numerical_cautions(self, model):
         dt = DiagnosticsToolbox(model=model.b)
 
         cautions = dt._collect_numerical_cautions()
-
         assert len(cautions) == 5
         assert (
-            "Caution: 3 Variables with value close to their bounds (abs=1.0E-04, rel=1.0E-04)"
+            "Caution: 2 Variables with value close to their bounds (abs=1.0E-04, rel=1.0E-04)"
             in cautions
         )
         assert "Caution: 2 Variables with value close to zero (tol=1.0E-08)" in cautions
@@ -1133,7 +1130,6 @@ def test_assert_no_numerical_warnings(self, model):
 
         # Fix numerical issues
         m.b.v3.setlb(-5)
-        m.b.v5.setub(10)
 
         solver = get_solver()
         solver.solve(m)
@@ -1198,12 +1194,12 @@ def test_report_numerical_issues(self, model):
 2 WARNINGS
 
     WARNING: 1 Constraint with large residuals (>1.0E-05)
-    WARNING: 2 Variables at or outside bounds (tol=0.0E+00)
+    WARNING: 1 Variable at or outside bounds (tol=0.0E+00)
 
 ------------------------------------------------------------------------------------
 5 Cautions
 
-    Caution: 3 Variables with value close to their bounds (abs=1.0E-04, rel=1.0E-04)
+    Caution: 2 Variables with value close to their bounds (abs=1.0E-04, rel=1.0E-04)
     Caution: 2 Variables with value close to zero (tol=1.0E-08)
     Caution: 1 Variable with extreme value (<1.0E-04 or >1.0E+04)
     Caution: 1 Variable with None value
@@ -1602,8 +1598,8 @@ def test_display_constraints_including_variable(self):
         expected = """====================================================================================
 The following constraints involve v[1]:
 
-    c1: 1.0
-    c4: 8.0
+    c1: 1.000e+00
+    c4: 8.000e+00
 
 ====================================================================================
 """
@@ -1667,8 +1663,8 @@ def test_display_variables_in_constraint(self):
         expected = """====================================================================================
 The following variables are involved in c1:
 
-    v[1]: 1.0
-    v[2]: 2.0
+    v[1]: 1.000e+00
+    v[2]: 2.000e+00
 
 ====================================================================================
 """

idaes/core/util/tests/test_utility_minimization.py

@@ -127,7 +127,7 @@ def model(self):
             "state_definition": FTPx,
             "state_bounds": {
                 "flow_mol": (0, 100, 1000, pyunits.mol / pyunits.s),
-                "temperature": (273.15, 300, 450, pyunits.K),
+                "temperature": (100, 300, 450, pyunits.K),
                 "pressure": (5e4, 1e5, 1e6, pyunits.Pa),
             },
             "pressure_ref": (1e5, pyunits.Pa),

idaes/models/properties/modular_properties/examples/tests/test_BT_PR.py

@@ -80,7 +80,7 @@ def test_T_sweep(self, m):
         m.fs.obj = Objective(expr=(m.fs.state[1].temperature - 510) ** 2)
         m.fs.state[1].temperature.setub(600)
 
-        for logP in range(8, 13, 1):
+        for logP in [9.5, 10, 10.5, 11, 11.5, 12]:
             m.fs.obj.deactivate()
 
             m.fs.state[1].flow_mol.fix(100)
@@ -115,11 +115,11 @@ def test_P_sweep(self, m):
             assert check_optimal_termination(results)
 
             while m.fs.state[1].pressure.value <= 1e6:
-                m.fs.state[1].pressure.value = m.fs.state[1].pressure.value + 1e5
 
                 results = solver.solve(m)
                 assert check_optimal_termination(results)
-                print(T, m.fs.state[1].pressure.value)
+
+                m.fs.state[1].pressure.value = m.fs.state[1].pressure.value + 1e5
 
     @pytest.mark.component
     def test_T350_P1_x5(self, m):

idaes/models/unit_models/tests/test_equilibrium_reactor.py

@@ -253,13 +253,6 @@ def test_get_performance_contents(self, sapon):
             }
         }
 
-    @pytest.mark.component
-    def test_initialization_error(self, sapon):
-        sapon.fs.unit.outlet.pressure[0].fix(1)
-
-        with pytest.raises(InitializationError):
-            sapon.fs.unit.initialize()
-
 
 class TestInitializers:
     @pytest.fixture

idaes/models/unit_models/tests/test_hx_ntu.py

@@ -505,16 +505,6 @@ def test_conservation(self, model):
             <= 1e-6
         )
 
-    @pytest.mark.component
-    def test_initialization_error(self, model):
-        model.fs.unit.hot_side_outlet.pressure[0].fix(1)
-
-        with pytest.raises(InitializationError):
-            model.fs.unit.initialize()
-
-        # Revert DoF change to avoid contaminating subsequent tests
-        model.fs.unit.hot_side_outlet.pressure[0].unfix()
-
 
 class TestInitializers(object):
     @pytest.fixture(scope="class")

idaes/models_extra/column_models/plate_heat_exchanger.py

@@ -46,6 +46,7 @@
 )
 from pyomo.common.config import ConfigValue, In, Integer
 from pyomo.util.calc_var_value import calculate_variable_from_constraint
+from pyomo.common.deprecation import deprecated
 
 # Import IDAES cores
 from idaes.core import declare_process_block_class
@@ -61,6 +62,13 @@
 _log = idaeslog.getLogger(__name__)
 
 
+@deprecated(
+    "The Plate Heat Exchanger (PHE) model is known to be affected by "
+    "issues causing it to fail to solve on certain platforms starting with Pyomo v6.7.0. "
+    "This might cause the model to be removed in an upcoming IDAES release if these failures are not resolved. "
+    "For more information, see IDAES/idaes-pse#1294",
+    version="2.3.0",
+)
 @declare_process_block_class("PlateHeatExchanger")
 class PlateHeatExchangerData(HeatExchangerNTUData):
     """Plate Heat Exchanger(PHE) Unit Model."""

idaes/models_extra/column_models/tests/test_plate_heat_exchanger.py

@@ -65,6 +65,12 @@ def test_config():
     assert len(m.fs.unit.config) == 9
 
 
+workaround_for_1294 = pytest.mark.xfail(
+    reason="These tests fail with Pyomo 6.7.0. See IDAES/idaes-pse#1294 for details",
+    strict=False,  # the failures only occur for certain platforms, e.g. Windows on GHA
+)
+
+
 # -----------------------------------------------------------------------------
 class TestPHE(object):
     @pytest.fixture(scope="class")
@@ -166,6 +172,7 @@ def test_initialize(self, phe):
             phe, duty=(245000, pyunits.W), optarg={"bound_push": 1e-8, "mu_init": 1e-8}
         )
 
+    @workaround_for_1294
     @pytest.mark.solver
     @pytest.mark.skipif(solver is None, reason="Solver not available")
     @pytest.mark.component
@@ -175,6 +182,7 @@ def test_solve(self, phe):
         # Check for optimal solution
         assert check_optimal_termination(results)
 
+    @workaround_for_1294
     @pytest.mark.solver
     @pytest.mark.skipif(solver is None, reason="Solver not available")
     @pytest.mark.component
@@ -214,6 +222,7 @@ def test_solution(self, phe):
             phe.fs.unit.cold_side_outlet.temperature[0]
         )
 
+    @workaround_for_1294
     @pytest.mark.solver
     @pytest.mark.skipif(solver is None, reason="Solver not available")
     @pytest.mark.component