From 6890840c6dd46cb60373dbc445ce5a23f1f6c0df Mon Sep 17 00:00:00 2001 From: "Alinson S. Xavier" Date: Thu, 21 Jan 2021 09:15:14 -0600 Subject: [PATCH] InternalSolver: Better specify and test infeasibility --- miplearn/instance.py | 12 ++++- miplearn/solvers/gurobi.py | 28 ++++++----- miplearn/solvers/internal.py | 43 +++++++---------- miplearn/solvers/pyomo/base.py | 24 +++++++--- miplearn/solvers/tests/__init__.py | 46 +++++++++++++++---- .../solvers/tests/test_internal_solver.py | 28 +++++++++-- miplearn/solvers/tests/test_lazy_cb.py | 4 +- .../solvers/tests/test_learning_solver.py | 12 ++--- miplearn/types.py | 2 +- 9 files changed, 134 insertions(+), 65 deletions(-) diff --git a/miplearn/instance.py b/miplearn/instance.py index e1b657d..bf8379f 100644 --- a/miplearn/instance.py +++ b/miplearn/instance.py @@ -10,6 +10,7 @@ from typing import Any, List import numpy as np from miplearn.types import TrainingSample +import pyomo.environ as pe class Instance(ABC): @@ -30,7 +31,7 @@ class Instance(ABC): @abstractmethod def to_model(self) -> Any: """ - Returns a concrete Pyomo model corresponding to this instance. + Returns the optimization model corresponding to this instance. """ pass @@ -163,3 +164,12 @@ class Instance(ABC): data = json.dumps(self.__dict__, indent=2).encode("utf-8") with gzip.GzipFile(filename, "w") as f: f.write(data) + + +class PyomoInstance(Instance, ABC): + @abstractmethod + def to_model(self) -> pe.ConcreteModel: + """ + Returns the concrete Pyomo model corresponding to this instance. + """ + pass diff --git a/miplearn/solvers/gurobi.py b/miplearn/solvers/gurobi.py index ac17a04..f60f3d0 100644 --- a/miplearn/solvers/gurobi.py +++ b/miplearn/solvers/gurobi.py @@ -128,8 +128,11 @@ class GurobiSolver(InternalSolver): for (idx, var) in vardict.items(): var.vtype = self.GRB.BINARY log = streams[0].getvalue() + opt_value = None + if not self.is_infeasible(): + opt_value = self.model.objVal return { - "Optimal value": self.model.objVal, + "Optimal value": opt_value, "Log": log, } @@ -173,14 +176,15 @@ class GurobiSolver(InternalSolver): if not should_repeat: break log = streams[0].getvalue() - if self.model.modelSense == 1: - sense = "min" - lb = self.model.objBound - ub = self.model.objVal - else: - sense = "max" - lb = self.model.objVal - ub = self.model.objBound + ub, lb = None, None + sense = "min" if self.model.modelSense == 1 else "max" + if self.model.solCount > 0: + if self.model.modelSense == 1: + lb = self.model.objBound + ub = self.model.objVal + else: + lb = self.model.objVal + ub = self.model.objBound ws_value = self._extract_warm_start_value(log) stats: MIPSolveStats = { "Lower bound": lb, @@ -194,8 +198,10 @@ class GurobiSolver(InternalSolver): } return stats - def get_solution(self) -> Dict: + def get_solution(self) -> Optional[Dict]: self._raise_if_callback() + if self.model.solCount == 0: + return None solution: Dict = {} for (varname, vardict) in self._all_vars.items(): solution[varname] = {} @@ -228,7 +234,7 @@ class GurobiSolver(InternalSolver): var = self._all_vars[var_name][index] return self._get_value(var) - def is_infeasible(self): + def is_infeasible(self) -> bool: return self.model.status in [self.GRB.INFEASIBLE, self.GRB.INF_OR_UNBD] def get_dual(self, cid): diff --git a/miplearn/solvers/internal.py b/miplearn/solvers/internal.py index 785821d..b6d0103 100644 --- a/miplearn/solvers/internal.py +++ b/miplearn/solvers/internal.py @@ -4,7 +4,7 @@ import logging from abc import ABC, abstractmethod -from typing import Any, Dict, List +from typing import Any, Dict, List, Optional from miplearn.instance import Instance from miplearn.types import ( @@ -39,15 +39,13 @@ class InternalSolver(ABC): Solves the LP relaxation of the currently loaded instance. After this method finishes, the solution can be retrieved by calling `get_solution`. + This method should not permanently modify the problem. That is, subsequent + calls to `solve` should solve the original MIP, not the LP relaxation. + Parameters ---------- - tee: bool + tee If true, prints the solver log to the screen. - - Returns - ------- - dict - A dictionary of solver statistics. """ pass @@ -64,34 +62,27 @@ class InternalSolver(ABC): Parameters ---------- - iteration_cb: () -> Bool + iteration_cb: By default, InternalSolver makes a single call to the native `solve` method and returns the result. If an iteration callback is provided instead, InternalSolver enters a loop, where `solve` and `iteration_cb` - are called alternatively. To stop the loop, `iteration_cb` should - return False. Any other result causes the solver to loop again. - lazy_cb: (internal_solver, model) -> None + are called alternatively. To stop the loop, `iteration_cb` should return + False. Any other result causes the solver to loop again. + lazy_cb: This function is called whenever the solver finds a new candidate - solution and can be used to add lazy constraints to the model. Only - the following operations within the callback are allowed: - - Querying the value of a variable, through `get_value(var, idx)` - - Querying if a constraint is satisfied, through `is_constraint_satisfied(cobj)` - - Adding a new constraint to the problem, through `add_constraint` + solution and can be used to add lazy constraints to the model. Only the + following operations within the callback are allowed: + - Querying the value of a variable + - Querying if a constraint is satisfied + - Adding a new constraint to the problem Additional operations may be allowed by specific subclasses. - tee: Bool + tee If true, prints the solver log to the screen. - - Returns - ------- - dict - A dictionary of solver statistics containing the following keys: - "Lower bound", "Upper bound", "Wallclock time", "Nodes", "Sense", - "Log" and "Warm start value". """ pass @abstractmethod - def get_solution(self) -> Dict: + def get_solution(self) -> Optional[Dict]: """ Returns current solution found by the solver. @@ -201,7 +192,7 @@ class InternalSolver(ABC): pass @abstractmethod - def set_constraint_rhs(self, cid: str, rhs: str) -> None: + def set_constraint_rhs(self, cid: str, rhs: float) -> None: pass @abstractmethod diff --git a/miplearn/solvers/pyomo/base.py b/miplearn/solvers/pyomo/base.py index de5bfe0..b19409a 100644 --- a/miplearn/solvers/pyomo/base.py +++ b/miplearn/solvers/pyomo/base.py @@ -11,6 +11,7 @@ from typing import Any, List, Dict, Optional import pyomo from pyomo import environ as pe from pyomo.core import Var, Constraint +from pyomo.opt import TerminationCondition from miplearn.instance import Instance from miplearn.solvers import RedirectOutput @@ -44,6 +45,7 @@ class BasePyomoSolver(InternalSolver): self._obj_sense = None self._varname_to_var = {} self._cname_to_constr = {} + self._termination_condition = None for (key, value) in params.items(): self._pyomo_solver.options[key] = value @@ -65,8 +67,11 @@ class BasePyomoSolver(InternalSolver): for var in self._bin_vars: var.domain = pyomo.core.base.set_types.Binary self._pyomo_solver.update_var(var) + opt_value = None + if not self.is_infeasible(): + opt_value = results["Problem"][0]["Lower bound"] return { - "Optimal value": results["Problem"][0]["Lower bound"], + "Optimal value": opt_value, "Log": streams[0].getvalue(), } @@ -100,9 +105,14 @@ class BasePyomoSolver(InternalSolver): log = streams[0].getvalue() node_count = self._extract_node_count(log) ws_value = self._extract_warm_start_value(log) + self._termination_condition = results["Solver"][0]["Termination condition"] + lb, ub = None, None + if not self.is_infeasible(): + lb = results["Problem"][0]["Lower bound"] + ub = results["Problem"][0]["Upper bound"] stats: MIPSolveStats = { - "Lower bound": results["Problem"][0]["Lower bound"], - "Upper bound": results["Problem"][0]["Upper bound"], + "Lower bound": lb, + "Upper bound": ub, "Wallclock time": total_wallclock_time, "Sense": self._obj_sense, "Log": log, @@ -112,7 +122,9 @@ class BasePyomoSolver(InternalSolver): } return stats - def get_solution(self) -> Dict: + def get_solution(self) -> Optional[Dict]: + if self.is_infeasible(): + return None solution: Dict = {} for var in self.model.component_objects(Var): solution[str(var)] = {} @@ -276,8 +288,8 @@ class BasePyomoSolver(InternalSolver): def set_constraint_rhs(self, cid, rhs): raise Exception("Not implemented") - def is_infeasible(self): - raise Exception("Not implemented") + def is_infeasible(self) -> bool: + return self._termination_condition == TerminationCondition.infeasible def get_dual(self, cid): raise Exception("Not implemented") diff --git a/miplearn/solvers/tests/__init__.py b/miplearn/solvers/tests/__init__.py index 5d65d2a..c101c0c 100644 --- a/miplearn/solvers/tests/__init__.py +++ b/miplearn/solvers/tests/__init__.py @@ -3,8 +3,11 @@ # Released under the modified BSD license. See COPYING.md for more details. from inspect import isclass -from typing import List, Callable +from typing import List, Callable, Any +from pyomo import environ as pe + +from miplearn.instance import Instance, PyomoInstance from miplearn.problems.knapsack import KnapsackInstance, GurobiKnapsackInstance from miplearn.solvers.gurobi import GurobiSolver from miplearn.solvers.internal import InternalSolver @@ -13,28 +16,55 @@ from miplearn.solvers.pyomo.gurobi import GurobiPyomoSolver from miplearn.solvers.pyomo.xpress import XpressPyomoSolver -def _get_instance(solver): - def _is_subclass_or_instance(obj, parent_class): - return isinstance(obj, parent_class) or ( - isclass(obj) and issubclass(obj, parent_class) - ) +class InfeasiblePyomoInstance(PyomoInstance): + def to_model(self) -> pe.ConcreteModel: + model = pe.ConcreteModel() + model.x = pe.Var(domain=pe.Binary) + model.OBJ = pe.Objective(expr=model.x, sense=pe.maximize) + model.eq = pe.Constraint(expr=model.x >= 2) + return model + + +class InfeasibleGurobiInstance(Instance): + def to_model(self) -> Any: + import gurobipy as gp + from gurobipy import GRB + + model = gp.Model() + x = model.addVars(1, vtype=GRB.BINARY, name="x") + model.addConstr(x[0] >= 2) + model.setObjective(x[0]) + return model + + +def _is_subclass_or_instance(obj, parent_class): + return isinstance(obj, parent_class) or ( + isclass(obj) and issubclass(obj, parent_class) + ) + +def _get_knapsack_instance(solver): if _is_subclass_or_instance(solver, BasePyomoSolver): return KnapsackInstance( weights=[23.0, 26.0, 20.0, 18.0], prices=[505.0, 352.0, 458.0, 220.0], capacity=67.0, ) - if _is_subclass_or_instance(solver, GurobiSolver): return GurobiKnapsackInstance( weights=[23.0, 26.0, 20.0, 18.0], prices=[505.0, 352.0, 458.0, 220.0], capacity=67.0, ) - assert False +def _get_infeasible_instance(solver): + if _is_subclass_or_instance(solver, BasePyomoSolver): + return InfeasiblePyomoInstance() + if _is_subclass_or_instance(solver, GurobiSolver): + return InfeasibleGurobiInstance() + + def _get_internal_solvers() -> List[Callable[[], InternalSolver]]: return [GurobiPyomoSolver, GurobiSolver, XpressPyomoSolver] diff --git a/miplearn/solvers/tests/test_internal_solver.py b/miplearn/solvers/tests/test_internal_solver.py index 1b70fd8..ee64f64 100644 --- a/miplearn/solvers/tests/test_internal_solver.py +++ b/miplearn/solvers/tests/test_internal_solver.py @@ -11,7 +11,11 @@ import pyomo.environ as pe from miplearn.solvers import RedirectOutput from miplearn.solvers.gurobi import GurobiSolver from miplearn.solvers.pyomo.base import BasePyomoSolver -from miplearn.solvers.tests import _get_instance, _get_internal_solvers +from miplearn.solvers.tests import ( + _get_knapsack_instance, + _get_internal_solvers, + _get_infeasible_instance, +) logger = logging.getLogger(__name__) @@ -30,7 +34,7 @@ def test_redirect_output(): def test_internal_solver_warm_starts(): for solver_class in _get_internal_solvers(): logger.info("Solver: %s" % solver_class) - instance = _get_instance(solver_class) + instance = _get_knapsack_instance(solver_class) model = instance.to_model() solver = solver_class() solver.set_instance(instance, model) @@ -82,7 +86,7 @@ def test_internal_solver(): for solver_class in _get_internal_solvers(): logger.info("Solver: %s" % solver_class) - instance = _get_instance(solver_class) + instance = _get_knapsack_instance(solver_class) model = instance.to_model() solver = solver_class() solver.set_instance(instance, model) @@ -158,10 +162,26 @@ def test_internal_solver(): assert round(stats["Lower bound"]) == 1179.0 +def test_infeasible_instance(): + for solver_class in _get_internal_solvers(): + instance = _get_infeasible_instance(solver_class) + solver = solver_class() + solver.set_instance(instance) + stats = solver.solve() + + assert solver.get_solution() is None + assert stats["Upper bound"] is None + assert stats["Lower bound"] is None + + stats = solver.solve_lp() + assert solver.get_solution() is None + assert stats["Optimal value"] is None + + def test_iteration_cb(): for solver_class in _get_internal_solvers(): logger.info("Solver: %s" % solver_class) - instance = _get_instance(solver_class) + instance = _get_knapsack_instance(solver_class) solver = solver_class() solver.set_instance(instance) count = 0 diff --git a/miplearn/solvers/tests/test_lazy_cb.py b/miplearn/solvers/tests/test_lazy_cb.py index 0d7b3eb..9fb2ad0 100644 --- a/miplearn/solvers/tests/test_lazy_cb.py +++ b/miplearn/solvers/tests/test_lazy_cb.py @@ -5,14 +5,14 @@ import logging from miplearn.solvers.gurobi import GurobiSolver -from miplearn.solvers.tests import _get_instance +from miplearn.solvers.tests import _get_knapsack_instance logger = logging.getLogger(__name__) def test_lazy_cb(): solver = GurobiSolver() - instance = _get_instance(solver) + instance = _get_knapsack_instance(solver) model = instance.to_model() def lazy_cb(cb_solver, cb_model): diff --git a/miplearn/solvers/tests/test_learning_solver.py b/miplearn/solvers/tests/test_learning_solver.py index aa8d199..02b112e 100644 --- a/miplearn/solvers/tests/test_learning_solver.py +++ b/miplearn/solvers/tests/test_learning_solver.py @@ -9,7 +9,7 @@ import os from miplearn.solvers.gurobi import GurobiSolver from miplearn.solvers.learning import LearningSolver -from miplearn.solvers.tests import _get_instance, _get_internal_solvers +from miplearn.solvers.tests import _get_knapsack_instance, _get_internal_solvers logger = logging.getLogger(__name__) @@ -18,7 +18,7 @@ def test_learning_solver(): for mode in ["exact", "heuristic"]: for internal_solver in _get_internal_solvers(): logger.info("Solver: %s" % internal_solver) - instance = _get_instance(internal_solver) + instance = _get_knapsack_instance(internal_solver) solver = LearningSolver( solver=internal_solver, mode=mode, @@ -50,7 +50,7 @@ def test_learning_solver(): def test_solve_without_lp(): for internal_solver in _get_internal_solvers(): logger.info("Solver: %s" % internal_solver) - instance = _get_instance(internal_solver) + instance = _get_knapsack_instance(internal_solver) solver = LearningSolver( solver=internal_solver, solve_lp_first=False, @@ -62,7 +62,7 @@ def test_solve_without_lp(): def test_parallel_solve(): for internal_solver in _get_internal_solvers(): - instances = [_get_instance(internal_solver) for _ in range(10)] + instances = [_get_knapsack_instance(internal_solver) for _ in range(10)] solver = LearningSolver(solver=internal_solver) results = solver.parallel_solve(instances, n_jobs=3) assert len(results) == 10 @@ -76,7 +76,7 @@ def test_solve_fit_from_disk(): # Create instances and pickle them filenames = [] for k in range(3): - instance = _get_instance(internal_solver) + instance = _get_knapsack_instance(internal_solver) with tempfile.NamedTemporaryFile(suffix=".pkl", delete=False) as file: filenames += [file.name] pickle.dump(instance, file) @@ -114,7 +114,7 @@ def test_solve_fit_from_disk(): def test_simulate_perfect(): internal_solver = GurobiSolver - instance = _get_instance(internal_solver) + instance = _get_knapsack_instance(internal_solver) with tempfile.NamedTemporaryFile(suffix=".pkl", delete=False) as tmp: pickle.dump(instance, tmp) tmp.flush() diff --git a/miplearn/types.py b/miplearn/types.py index 6cd701f..4455af8 100644 --- a/miplearn/types.py +++ b/miplearn/types.py @@ -22,7 +22,7 @@ TrainingSample = TypedDict( LPSolveStats = TypedDict( "LPSolveStats", { - "Optimal value": float, + "Optimal value": Optional[float], "Log": str, }, )