Fix GLOA convex affine cut slope

pyomo/contrib/gdpopt/gloa.py

@@ -195,31 +195,49 @@ def _add_cuts_to_discrete_problem(
                 aff_utils_blocks[parent_block] = aff_utils
                 aff_utils.GDPopt_aff_cons = Constraint(NonNegativeIntegers)
             aff_cuts = aff_utils.GDPopt_aff_cons
-            cut_body = sum(
+            concave_cut_body = sum(
                 ccSlope[var] * (var - var.value)
                 for var in vars_in_constr
                 if not var.fixed
             )
-            if not is_potentially_variable(cut_body):
+            convex_cut_body = sum(
+                cvSlope[var] * (var - var.value)
+                for var in vars_in_constr
+                if not var.fixed
+            )
+            concave_cut_is_variable = is_potentially_variable(concave_cut_body)
+            convex_cut_is_variable = is_potentially_variable(convex_cut_body)
+
+            if not concave_cut_is_variable:
                 if (
-                    cut_body + ccStart >= lb_int - config.constraint_tolerance
-                    and cut_body + cvStart <= ub_int + config.constraint_tolerance
+                    value(concave_cut_body + ccStart)
+                    < lb_int - config.constraint_tolerance
                 ):
-                    # We won't add them, but nothing is wrong--they hold
-                    config.logger.debug("Affine cut is trivially True.")
-                else:
-                    # something went wrong.
                     raise DeveloperError("One of the affine cuts is trivially False.")
+            if not convex_cut_is_variable:
+                if (
+                    value(convex_cut_body + cvStart)
+                    > ub_int + config.constraint_tolerance
+                ):
+                    raise DeveloperError("One of the affine cuts is trivially False.")
+
+            if not concave_cut_is_variable and not convex_cut_is_variable:
+                # We won't add them, but nothing is wrong--they hold
+                config.logger.debug("Affine cut is trivially True.")
             else:
-                concave_cut = cut_body + ccStart >= lb_int
-                convex_cut = cut_body + cvStart <= ub_int
                 idx = len(aff_cuts)
-                aff_cuts[idx] = concave_cut
-                aff_cuts[idx + 1] = convex_cut
-                _add_bigm_constraint_to_transformed_model(m, aff_cuts[idx], aff_cuts)
-                _add_bigm_constraint_to_transformed_model(
-                    m, aff_cuts[idx + 1], aff_cuts
-                )
-                counter += 2
+                if concave_cut_is_variable:
+                    aff_cuts[idx] = concave_cut_body + ccStart >= lb_int
+                    _add_bigm_constraint_to_transformed_model(
+                        m, aff_cuts[idx], aff_cuts
+                    )
+                    idx += 1
+                    counter += 1
+                if convex_cut_is_variable:
+                    aff_cuts[idx] = convex_cut_body + cvStart <= ub_int
+                    _add_bigm_constraint_to_transformed_model(
+                        m, aff_cuts[idx], aff_cuts
+                    )
+                    counter += 1
 
                 config.logger.debug("Added %s affine cuts" % counter)

pyomo/contrib/gdpopt/tests/test_gdpopt.py

@@ -312,6 +312,41 @@ def test_gloa_cut_generation_ignores_deactivated_constraints(self):
         self.assertEqual(c2.lower, 0)
         self.assertIsNone(c2.upper)
 
+    @unittest.skipIf(not mcpp_available(), "MC++ is not available")
+    def test_gloa_affine_cut_uses_convex_slope_for_upper_cut(self):
+        m = ConcreteModel()
+        m.P = Var(bounds=(0, 10), initialize=0.05)
+        m.Q = Var(bounds=(1, 10), initialize=1.000002190520329)
+        m.F = Var(bounds=(0, 10), initialize=0.050000453446344)
+        m.QP = Var(bounds=(0, 10), initialize=1.0)
+        m.c = Constraint(expr=m.P * m.Q - m.F * m.QP == 0)
+
+        m.GDPopt_utils = Block()
+        util_block = m.GDPopt_utils
+        util_block.algebraic_variable_list = [m.P, m.Q, m.F, m.QP]
+        util_block.global_constraint_list = [m.c]
+        util_block.constraints_by_disjunct = {}
+
+        config = Bunch(
+            logger=logging.getLogger('pyomo.contrib.gdpopt.tests'),
+            integer_tolerance=1e-6,
+            constraint_tolerance=1e-6,
+        )
+
+        SolverFactory('gdpopt.gloa')._add_cuts_to_discrete_problem(
+            util_block, util_block, None, config, Bunch()
+        )
+
+        feasible_point = [(m.P, 0.05), (m.Q, 1.1), (m.F, 0.055), (m.QP, 1.0)]
+        for var, val in feasible_point:
+            var.set_value(val)
+
+        aff_cuts = m.GDPopt_aff.GDPopt_aff_cons
+        self.assertEqual(len(aff_cuts), 2)
+        convex_cut = aff_cuts[1]
+        self.assertLessEqual(value(convex_cut.body), value(convex_cut.upper) + 1e-12)
+        self.assertAlmostEqual(value(convex_cut.body), -0.5)
+
     def test_complain_when_no_algorithm_specified(self):
         m = self.get_GDP_on_block()
         with self.assertRaisesRegex(

pyomo/contrib/mindtpy/tests/unit_test.py

@@ -7,14 +7,25 @@
 # software.  This software is distributed under the 3-clause BSD License.
 # ____________________________________________________________________________________
 
-import pyomo.common.unittest as unittest
-from pyomo.contrib.mindtpy.util import set_var_valid_value
+import logging
 
-from pyomo.environ import Var, Integers, ConcreteModel, Integers
+import pyomo.common.unittest as unittest
+from pyomo.common.collections import Bunch
+from pyomo.contrib.mcpp.pyomo_mcpp import mcpp_available
 from pyomo.contrib.mindtpy.algorithm_base_class import _MindtPyAlgorithm
 from pyomo.contrib.mindtpy.config_options import _get_MindtPy_OA_config
+from pyomo.contrib.mindtpy.cut_generation import add_affine_cuts
 from pyomo.contrib.mindtpy.tests.MINLP5_simple import SimpleMINLP5
-from pyomo.contrib.mindtpy.util import add_var_bound
+from pyomo.contrib.mindtpy.util import add_var_bound, set_var_valid_value
+from pyomo.environ import (
+    Block,
+    ConcreteModel,
+    Constraint,
+    ConstraintList,
+    Integers,
+    Var,
+    value,
+)
 
 
 class UnitTestMindtPy(unittest.TestCase):
@@ -94,6 +105,33 @@ def test_add_var_bound(self):
             solver_object.working_model.y.upper, solver_object.config.integer_var_bound
         )
 
+    @unittest.skipIf(not mcpp_available(), "MC++ is not available")
+    def test_goa_affine_cut_uses_convex_slope_for_upper_cut(self):
+        m = ConcreteModel()
+        m.P = Var(bounds=(0, 10), initialize=0.05)
+        m.Q = Var(bounds=(1, 10), initialize=1.000002190520329)
+        m.F = Var(bounds=(0, 10), initialize=0.050000453446344)
+        m.QP = Var(bounds=(0, 10), initialize=1.0)
+        m.c = Constraint(expr=m.P * m.Q - m.F * m.QP == 0)
+
+        m.MindtPy_utils = Block()
+        m.MindtPy_utils.nonlinear_constraint_list = [m.c]
+        m.MindtPy_utils.cuts = Block()
+        m.MindtPy_utils.cuts.aff_cuts = ConstraintList()
+
+        config = Bunch(logger=logging.getLogger('pyomo.contrib.mindtpy.tests'))
+        add_affine_cuts(m, config, Bunch())
+
+        feasible_point = [(m.P, 0.05), (m.Q, 1.1), (m.F, 0.055), (m.QP, 1.0)]
+        for var, val in feasible_point:
+            var.set_value(val)
+
+        aff_cuts = list(m.MindtPy_utils.cuts.aff_cuts.values())
+        self.assertEqual(len(aff_cuts), 2)
+        convex_cut = aff_cuts[1]
+        self.assertLessEqual(value(convex_cut.body), value(convex_cut.upper) + 1e-12)
+        self.assertAlmostEqual(value(convex_cut.body), -0.5)
+
 
 if __name__ == '__main__':
     unittest.main()