diff --git a/pyomo/contrib/gdpopt/enumerate.py b/pyomo/contrib/gdpopt/enumerate.py index aae1b7bc0f9..d372ceb0067 100644 --- a/pyomo/contrib/gdpopt/enumerate.py +++ b/pyomo/contrib/gdpopt/enumerate.py @@ -7,6 +7,7 @@ # software. This software is distributed under the 3-clause BSD License. # ____________________________________________________________________________________ +import math from itertools import product from pyomo.common.collections import ComponentSet @@ -73,7 +74,9 @@ def solve(self, model, **kwds): def _discrete_solution_iterator( self, disjunctions, non_indicator_boolean_vars, discrete_var_list, config ): - discrete_var_values = [range(v.lb, v.ub + 1) for v in discrete_var_list] + discrete_var_values = [ + range(math.ceil(v.lb), math.floor(v.ub) + 1) for v in discrete_var_list + ] # we will calculate all the possible indicator_var realizations, and # then multiply those out by all the boolean var realizations and all # the integer var realizations. @@ -110,8 +113,6 @@ def _log_current_state(self, logger, subproblem_type, primal_improved=False): ) def _solve_gdp(self, original_model, config): - logger = config.logger - util_block = self.original_util_block # From preprocessing to make sure this *is* a GDP, we already have # lists of: @@ -124,19 +125,34 @@ def _solve_gdp(self, original_model, config): subproblem, subproblem_util_block = get_subproblem(original_model, util_block) - discrete_solns = list( - self._discrete_solution_iterator( - subproblem_util_block.disjunction_list, - subproblem_util_block.non_indicator_boolean_variable_list, - subproblem_util_block.discrete_variable_list, - config, - ) + disjunctions = subproblem_util_block.disjunction_list + non_indicator_boolean_vars = ( + subproblem_util_block.non_indicator_boolean_variable_list + ) + discrete_vars = subproblem_util_block.discrete_variable_list + + for v in discrete_vars: + if v.lb is None or v.ub is None: + raise ValueError( + f"GDPopt enumeration requires finite bounds on integer variable {v.name}." + ) + + self.num_discrete_solns = math.prod( + len(disjunction.disjuncts) for disjunction in disjunctions ) - self.num_discrete_solns = len(discrete_solns) - for soln in discrete_solns: - # We will interrupt based on time limit or iteration limit: + if config.force_subproblem_nlp: + self.num_discrete_solns *= 2 ** len(non_indicator_boolean_vars) * math.prod( + max(0, math.floor(v.ub) - math.ceil(v.lb) + 1) for v in discrete_vars + ) + + if self.reached_time_limit(config) or self.reached_iteration_limit(config): + return + for soln in self._discrete_solution_iterator( + disjunctions, non_indicator_boolean_vars, discrete_vars, config + ): if self.reached_time_limit(config) or self.reached_iteration_limit(config): - break + return + self.iteration += 1 with time_code(self.timing, 'nlp'): @@ -159,26 +175,23 @@ def _solve_gdp(self, original_model, config): # the whole problem is unbounded, we can stop self._update_primal_bound_to_unbounded(config) self._log_current_state(config.logger, 'subproblem', True) - break + return else: # Just log where we are self._log_current_state(config.logger, 'subproblem') - if self.iteration == self.num_discrete_solns: - # We can terminate optimally or declare infeasibility: We have - # enumerated all solutions, so our incumbent is optimal (or - # locally optimal, depending on how we solved the subproblems) - # if it exists, and if not then there is no solution. - if self.incumbent_boolean_soln is None: - self._update_dual_bound_to_infeasible() - self._load_infeasible_termination_status(config) - else: # the incumbent is optimal - self._update_bounds(dual=self.primal_bound(), force_update=True) - self._log_current_state(config.logger, '') - config.logger.info( - 'GDPopt exiting--all discrete solutions have been ' - 'enumerated.' - ) - self.pyomo_results.solver.termination_condition = tc.optimal - break + # We can terminate optimally or declare infeasibility: We have + # enumerated all solutions, so our incumbent is optimal (or + # locally optimal, depending on how we solved the subproblems) + # if it exists, and if not then there is no solution. + if self.incumbent_boolean_soln is None: + self._update_dual_bound_to_infeasible() + self._load_infeasible_termination_status(config) + else: # the incumbent is optimal + self._update_bounds(dual=self.primal_bound(), force_update=True) + self._log_current_state(config.logger, '') + config.logger.info( + 'GDPopt exiting--all discrete solutions have been enumerated.' + ) + self.pyomo_results.solver.termination_condition = tc.optimal diff --git a/pyomo/contrib/gdpopt/tests/test_enumerate.py b/pyomo/contrib/gdpopt/tests/test_enumerate.py index f0860fe577f..66044640d7f 100644 --- a/pyomo/contrib/gdpopt/tests/test_enumerate.py +++ b/pyomo/contrib/gdpopt/tests/test_enumerate.py @@ -25,6 +25,71 @@ import pyomo.gdp.tests.models as models +class _ExpiredEnumerationSolver(GDP_Enumeration_Solver): + """Enumeration solver that enters its GDP solve after the time limit.""" + + def _solve_gdp(self, original_model, config): + """Expire the active timer before running the real GDP solve.""" + self.timing.main_timer_start_time -= config.time_limit + return super()._solve_gdp(original_model, config) + + def _discrete_solution_iterator(self, *args): + """Fail instead of requesting a discrete solution.""" + raise AssertionError('discrete solutions were enumerated') + + +class TestGDPoptEnumerateUnit(unittest.TestCase): + def test_large_space_not_enumerated_after_time_limit(self): + m = ConcreteModel() + m.x = Var(bounds=(-1, 1)) + m.i = Var(domain=Integers, bounds=(0.5, 10**20)) + m.obj = Objective(expr=m.x + m.i) + + def disjunction_rule(m, _): + return [[m.x <= 0], [m.x >= 0]] + + m.disjunctions = Disjunction(range(32), rule=disjunction_rule) + solver = _ExpiredEnumerationSolver() + results = solver.solve(m, force_subproblem_nlp=True, time_limit=1) + + self.assertEqual(solver.num_discrete_solns, 2**32 * 10**20) + self.assertEqual( + results.solver.termination_condition, TerminationCondition.maxTimeLimit + ) + + def test_noninteger_integer_bounds(self): + m = ConcreteModel() + m.i = Var(domain=Integers, bounds=(0.5, 3.5)) + solver = GDP_Enumeration_Solver(force_subproblem_nlp=True) + + solutions = solver._discrete_solution_iterator([], [], [m.i], solver.config) + + self.assertEqual([solution[2] for solution in solutions], [(1,), (2,), (3,)]) + + def test_unbounded_integer_variable(self): + m = ConcreteModel() + m.i = Var(domain=Integers) + m.disjunction = Disjunction(expr=[[m.i >= 0], [m.i <= 0]]) + m.obj = Objective(expr=m.i) + + with self.assertRaisesRegex(ValueError, 'finite bounds'): + GDP_Enumeration_Solver().solve(m, force_subproblem_nlp=True) + + def test_completion(self): + m = ConcreteModel() + m.disjunction = Disjunction( + expr=[[Constraint.Infeasible], [Constraint.Infeasible]] + ) + m.obj = Objective(expr=0) + + results = GDP_Enumeration_Solver().solve(m, iterlim=2) + + self.assertEqual(results.solver.iterations, 2) + self.assertEqual( + results.solver.termination_condition, TerminationCondition.infeasible + ) + + @unittest.skipUnless(SolverFactory('gurobi').available(), 'Gurobi not available') @unittest.skipUnless(SolverFactory('gurobi').license_is_valid(), 'Gurobi not licensed') class TestGDPoptEnumerate(unittest.TestCase):