diff --git a/miplearn/solvers/gurobi.py b/miplearn/solvers/gurobi.py index d8c0829..f95c052 100644 --- a/miplearn/solvers/gurobi.py +++ b/miplearn/solvers/gurobi.py @@ -26,6 +26,7 @@ from miplearn.types import ( UserCutCallback, Solution, VariableName, + Category, ) logger = logging.getLogger(__name__) @@ -67,6 +68,7 @@ class GurobiSolver(InternalSolver): self.lazy_cb_frequency = lazy_cb_frequency self._bin_vars: List["gurobipy.Var"] = [] self._varname_to_var: Dict[str, "gurobipy.Var"] = {} + self._original_vtype: Dict["gurobipy.Var", str] = {} self._dirty = True self._has_lp_solution = False self._has_mip_solution = False @@ -101,6 +103,7 @@ class GurobiSolver(InternalSolver): def _update_vars(self) -> None: assert self.model is not None self._varname_to_var.clear() + self._original_vtype = {} self._bin_vars.clear() for var in self.model.getVars(): assert var.varName not in self._varname_to_var, ( @@ -112,6 +115,7 @@ class GurobiSolver(InternalSolver): "Only binary and continuous variables are currently supported. " "Variable {var.varName} has type {var.vtype}." ) + self._original_vtype[var] = var.vtype if var.vtype == "B": self._bin_vars.append(var) @@ -433,7 +437,7 @@ class GurobiSolver(InternalSolver): var.lower_bound = gp_var.lb var.upper_bound = gp_var.ub var.obj_coeff = gp_var.obj - var.type = gp_var.vtype + var.type = self._original_vtype[gp_var] if self._has_lp_solution: var.reduced_cost = gp_var.rc @@ -587,8 +591,9 @@ class GurobiTestInstanceKnapsack(PyomoTestInstanceKnapsack): model = gp.Model("Knapsack") n = len(self.weights) x = model.addVars(n, vtype=GRB.BINARY, name="x") + z = model.addVar(vtype=GRB.CONTINUOUS, name="z", ub=self.capacity) model.addConstr( - gp.quicksum(x[i] * self.weights[i] for i in range(n)) <= self.capacity, + gp.quicksum(x[i] * self.weights[i] for i in range(n)) == z, "eq_capacity", ) model.setObjective( diff --git a/miplearn/solvers/learning.py b/miplearn/solvers/learning.py index 5420b41..6a70696 100644 --- a/miplearn/solvers/learning.py +++ b/miplearn/solvers/learning.py @@ -253,7 +253,6 @@ class LearningSolver: # ------------------------------------------------------- logger.info("Extracting features (after-mip)...") features = FeaturesExtractor(self.internal_solver).extract(instance) - features.lp_solve = lp_stats features.mip_solve = mip_stats instance.features_after_mip.append(features) diff --git a/miplearn/solvers/pyomo/base.py b/miplearn/solvers/pyomo/base.py index dca2062..d805c90 100644 --- a/miplearn/solvers/pyomo/base.py +++ b/miplearn/solvers/pyomo/base.py @@ -12,7 +12,7 @@ import numpy as np import pyomo from overrides import overrides from pyomo import environ as pe -from pyomo.core import Var, Suffix +from pyomo.core import Var, Suffix, Objective from pyomo.core.base import _GeneralVarData from pyomo.core.base.constraint import ConstraintList from pyomo.core.expr.numeric_expr import SumExpression, MonomialTermExpression @@ -64,6 +64,7 @@ class BasePyomoSolver(InternalSolver): self._termination_condition: str = "" self._has_lp_solution = False self._has_mip_solution = False + self._obj: Dict[str, float] = {} for (key, value) in params.items(): self._pyomo_solver.options[key] = value @@ -223,6 +224,9 @@ class BasePyomoSolver(InternalSolver): self._all_vars += [var[idx]] if var[idx].domain == pyomo.core.base.set_types.Binary: self._bin_vars += [var[idx]] + for obj in self.model.component_objects(Objective): + self._obj = self._parse_pyomo_expr(obj.expr) + break def _update_constrs(self) -> None: assert self.model is not None @@ -371,11 +375,42 @@ class BasePyomoSolver(InternalSolver): for var in self.model.component_objects(pyomo.core.Var): for idx in var: varname = f"{var}[{idx}]" + if idx is None: + varname = str(var) variables[varname] = self._parse_pyomo_variable(var[idx]) return variables def _parse_pyomo_variable(self, var: pyomo.core.Var) -> Variable: - return Variable() + # Variable type + vtype: Optional[str] = None + if var.domain == pyomo.core.Binary: + vtype = "B" + elif var.domain in [ + pyomo.core.Reals, + pyomo.core.NonNegativeReals, + pyomo.core.NonPositiveReals, + pyomo.core.NegativeReals, + pyomo.core.PositiveReals, + ]: + vtype = "C" + if vtype is None: + raise Exception(f"unknown variable domain: {var.domain}") + + # Bounds + lb, ub = var.bounds + + # Objective coefficient + obj_coeff = 0.0 + if var.name in self._obj: + obj_coeff = self._obj[var.name] + + return Variable( + value=var.value, + type=vtype, + lower_bound=float(lb), + upper_bound=float(ub), + obj_coeff=obj_coeff, + ) @overrides def get_constraints(self) -> Dict[str, Constraint]: @@ -408,10 +443,10 @@ class BasePyomoSolver(InternalSolver): assert ( (not has_lb) or (not has_ub) or pyomo_constr.upper() == pyomo_constr.lower() ), "range constraints not supported" - if has_lb: + if not has_ub: constr.sense = ">" constr.rhs = pyomo_constr.lower() - elif has_ub: + elif not has_lb: constr.sense = "<" constr.rhs = pyomo_constr.upper() else: @@ -419,22 +454,7 @@ class BasePyomoSolver(InternalSolver): constr.rhs = pyomo_constr.upper() # Extract LHS - lhs = {} - if isinstance(pyomo_constr.body, SumExpression): - for term in pyomo_constr.body._args_: - if isinstance(term, MonomialTermExpression): - lhs[term._args_[1].name] = term._args_[0] - elif isinstance(term, _GeneralVarData): - lhs[term.name] = 1.0 - else: - raise Exception(f"Unknown term type: {term.__class__.__name__}") - elif isinstance(pyomo_constr.body, _GeneralVarData): - lhs[pyomo_constr.body.name] = 1.0 - else: - raise Exception( - f"Unknown expression type: {pyomo_constr.body.__class__.__name__}" - ) - constr.lhs = lhs + constr.lhs = self._parse_pyomo_expr(pyomo_constr.body) # Extract solution attributes if self._has_lp_solution: @@ -446,6 +466,22 @@ class BasePyomoSolver(InternalSolver): # Build constraint return constr + def _parse_pyomo_expr(self, expr): + lhs = {} + if isinstance(expr, SumExpression): + for term in expr._args_: + if isinstance(term, MonomialTermExpression): + lhs[term._args_[1].name] = float(term._args_[0]) + elif isinstance(term, _GeneralVarData): + lhs[term.name] = 1.0 + else: + raise Exception(f"Unknown term type: {term.__class__.__name__}") + elif isinstance(expr, _GeneralVarData): + lhs[expr.name] = 1.0 + else: + raise Exception(f"Unknown expression type: {expr.__class__.__name__}") + return lhs + @overrides def are_callbacks_supported(self) -> bool: return False @@ -453,23 +489,20 @@ class BasePyomoSolver(InternalSolver): @overrides def get_constraint_attrs(self) -> List[str]: return [ - "category", "dual_value", "lazy", "lhs", "rhs", "sense", "slack", - "user_features", ] @overrides def get_variable_attrs(self) -> List[str]: return [ # "basis_status", - # "category", - # "lower_bound", - # "obj_coeff", + "lower_bound", + "obj_coeff", # "reduced_cost", # "sa_lb_down", # "sa_lb_up", @@ -477,10 +510,9 @@ class BasePyomoSolver(InternalSolver): # "sa_obj_up", # "sa_ub_down", # "sa_ub_up", - # "type", - # "upper_bound", - # "user_features", - # "value", + "type", + "upper_bound", + "value", ] @@ -529,12 +561,13 @@ class PyomoTestInstanceKnapsack(Instance): model = pe.ConcreteModel() items = range(len(self.weights)) model.x = pe.Var(items, domain=pe.Binary) + model.z = pe.Var(domain=pe.Reals, bounds=(0, self.capacity)) model.OBJ = pe.Objective( expr=sum(model.x[v] * self.prices[v] for v in items), sense=pe.maximize, ) model.eq_capacity = pe.Constraint( - expr=sum(model.x[v] * self.weights[v] for v in items) <= self.capacity + expr=sum(model.x[v] * self.weights[v] for v in items) == model.z ) return model @@ -552,3 +585,9 @@ class PyomoTestInstanceKnapsack(Instance): self.weights[item], self.prices[item], ] + + @overrides + def get_variable_category(self, var_name: VariableName) -> Optional[Category]: + if var_name.startswith("x"): + return "default" + return None diff --git a/miplearn/solvers/tests/__init__.py b/miplearn/solvers/tests/__init__.py index 9901dc0..3e806b6 100644 --- a/miplearn/solvers/tests/__init__.py +++ b/miplearn/solvers/tests/__init__.py @@ -117,6 +117,12 @@ def run_basic_usage_tests(solver: InternalSolver) -> None: type="B", upper_bound=1.0, ), + "z": Variable( + lower_bound=0.0, + obj_coeff=0.0, + type="C", + upper_bound=67.0, + ), }, ), ) @@ -127,9 +133,9 @@ def run_basic_usage_tests(solver: InternalSolver) -> None: { "eq_capacity": Constraint( lazy=False, - lhs={"x[0]": 23.0, "x[1]": 26.0, "x[2]": 20.0, "x[3]": 18.0}, - rhs=67.0, - sense="<", + lhs={"x[0]": 23.0, "x[1]": 26.0, "x[2]": 20.0, "x[3]": 18.0, "z": -1.0}, + rhs=0.0, + sense="=", ) }, ) @@ -161,7 +167,7 @@ def run_basic_usage_tests(solver: InternalSolver) -> None: sa_obj_up=inf, sa_ub_down=0.913043, sa_ub_up=2.043478, - type="C", + type="B", upper_bound=1.0, value=1.0, ), @@ -176,7 +182,7 @@ def run_basic_usage_tests(solver: InternalSolver) -> None: sa_obj_up=570.869565, sa_ub_down=0.923077, sa_ub_up=inf, - type="C", + type="B", upper_bound=1.0, value=0.923077, ), @@ -191,7 +197,7 @@ def run_basic_usage_tests(solver: InternalSolver) -> None: sa_obj_up=inf, sa_ub_down=0.9, sa_ub_up=2.2, - type="C", + type="B", upper_bound=1.0, value=1.0, ), @@ -206,10 +212,25 @@ def run_basic_usage_tests(solver: InternalSolver) -> None: sa_obj_up=243.692308, sa_ub_down=0.0, sa_ub_up=inf, - type="C", + type="B", upper_bound=1.0, value=0.0, ), + "z": Variable( + basis_status="U", + lower_bound=0.0, + obj_coeff=0.0, + reduced_cost=13.538462, + sa_lb_down=-inf, + sa_lb_up=67.0, + sa_obj_down=-13.538462, + sa_obj_up=inf, + sa_ub_down=43.0, + sa_ub_up=69.0, + type="C", + upper_bound=67.0, + value=67.0, + ), }, ), ) @@ -221,15 +242,21 @@ def run_basic_usage_tests(solver: InternalSolver) -> None: solver, { "eq_capacity": Constraint( + basis_status="N", + dual_value=13.538462, lazy=False, - lhs={"x[0]": 23.0, "x[1]": 26.0, "x[2]": 20.0, "x[3]": 18.0}, - rhs=67.0, - sense="<", + lhs={ + "x[0]": 23.0, + "x[1]": 26.0, + "x[2]": 20.0, + "x[3]": 18.0, + "z": -1.0, + }, + rhs=0.0, + sa_rhs_down=-24.0, + sa_rhs_up=1.9999999999999987, + sense="=", slack=0.0, - dual_value=13.538462, - sa_rhs_down=43.0, - sa_rhs_up=69.0, - basis_status="N", ) }, ), @@ -238,9 +265,6 @@ def run_basic_usage_tests(solver: InternalSolver) -> None: # Solve MIP mip_stats = solver.solve( tee=True, - iteration_cb=None, - lazy_cb=None, - user_cut_cb=None, ) assert not solver.is_infeasible() assert mip_stats.mip_log is not None @@ -289,6 +313,13 @@ def run_basic_usage_tests(solver: InternalSolver) -> None: upper_bound=1.0, value=1.0, ), + "z": Variable( + lower_bound=0.0, + obj_coeff=0.0, + type="C", + upper_bound=67.0, + value=61.0, + ), }, ), ) @@ -298,11 +329,12 @@ def run_basic_usage_tests(solver: InternalSolver) -> None: _round_constraints(solver.get_constraints()), { "eq_capacity": Constraint( - lhs={"x[0]": 23.0, "x[1]": 26.0, "x[2]": 20.0, "x[3]": 18.0}, - rhs=67.0, - sense="<", - slack=6.0, - ), + lazy=False, + lhs={"x[0]": 23.0, "x[1]": 26.0, "x[2]": 20.0, "x[3]": 18.0, "z": -1.0}, + rhs=0.0, + sense="=", + slack=0.0, + ) }, ) @@ -317,11 +349,13 @@ def run_basic_usage_tests(solver: InternalSolver) -> None: _round_constraints(solver.get_constraints()), { "eq_capacity": Constraint( - lhs={"x[0]": 23.0, "x[1]": 26.0, "x[2]": 20.0, "x[3]": 18.0}, - rhs=67.0, - sense="<", + lazy=False, + lhs={"x[0]": 23.0, "x[1]": 26.0, "x[2]": 20.0, "x[3]": 18.0, "z": -1.0}, + rhs=0.0, + sense="=", ), "cut": Constraint( + lazy=False, lhs={"x[0]": 1.0}, rhs=0.0, sense="<", diff --git a/tests/solvers/test_learning_solver.py b/tests/solvers/test_learning_solver.py index ed0c981..90f1fe2 100644 --- a/tests/solvers/test_learning_solver.py +++ b/tests/solvers/test_learning_solver.py @@ -92,7 +92,7 @@ def test_parallel_solve( for instance in instances: data = instance.training_data[0] assert data.solution is not None - assert len(data.solution.keys()) == 4 + assert len(data.solution.keys()) == 5 def test_solve_fit_from_disk( diff --git a/tests/test_features.py b/tests/test_features.py index c9d5830..bb0d4f8 100644 --- a/tests/test_features.py +++ b/tests/test_features.py @@ -41,7 +41,7 @@ def test_knapsack() -> None: sa_obj_up=inf, sa_ub_down=0.913043, sa_ub_up=2.043478, - type="C", + type="B", upper_bound=1.0, user_features=[23.0, 505.0], value=1.0, @@ -59,7 +59,7 @@ def test_knapsack() -> None: sa_obj_up=570.869565, sa_ub_down=0.923077, sa_ub_up=inf, - type="C", + type="B", upper_bound=1.0, user_features=[26.0, 352.0], value=0.923077, @@ -86,7 +86,7 @@ def test_knapsack() -> None: sa_obj_up=inf, sa_ub_down=0.9, sa_ub_up=2.2, - type="C", + type="B", upper_bound=1.0, user_features=[20.0, 458.0], value=1.0, @@ -104,12 +104,30 @@ def test_knapsack() -> None: sa_obj_up=243.692308, sa_ub_down=0.0, sa_ub_up=inf, - type="C", + type="B", upper_bound=1.0, user_features=[18.0, 220.0], value=0.0, alvarez_2017=[1.0, 0.143322, 0.0, 0.0, 1.0, -1.0, 46.051702, 3.16515], ), + "z": Variable( + basis_status="U", + category=None, + lower_bound=0.0, + obj_coeff=0.0, + reduced_cost=13.538462, + sa_lb_down=-inf, + sa_lb_up=67.0, + sa_obj_down=-13.538462, + sa_obj_up=inf, + sa_ub_down=43.0, + sa_ub_up=69.0, + type="C", + upper_bound=67.0, + user_features=None, + value=67.0, + alvarez_2017=[0.0, 0.0, 0.0, 0.0, 1.0, -1.0, 0.0, 0.0], + ), }, ) assert_equals( @@ -120,11 +138,11 @@ def test_knapsack() -> None: category="eq_capacity", dual_value=13.538462, lazy=False, - lhs={"x[0]": 23.0, "x[1]": 26.0, "x[2]": 20.0, "x[3]": 18.0}, - rhs=67.0, - sa_rhs_down=43.0, - sa_rhs_up=69.0, - sense="<", + lhs={"x[0]": 23.0, "x[1]": 26.0, "x[2]": 20.0, "x[3]": 18.0, "z": -1.0}, + rhs=0.0, + sa_rhs_down=-24.0, + sa_rhs_up=1.9999999999999987, + sense="=", slack=0.0, user_features=[0.0], )