diff --git a/optlang/cplex_interface.py b/optlang/cplex_interface.py index 1af83b21..88b7669a 100644 --- a/optlang/cplex_interface.py +++ b/optlang/cplex_interface.py @@ -711,7 +711,9 @@ def objective(self, value): if self._objective is not None: # Reset previous objective variables = self.objective.variables if len(variables) > 0: - self.problem.objective.set_linear([(variable.name, 0.) for variable in variables]) + name_list = [var.name for var in variables] + index_dict = {n: i for n, i in zip(name_list, self._get_variable_indices(name_list))} + self.problem.objective.set_linear([(index_dict[variable.name], 0.) for variable in variables]) if self.problem.objective.get_num_quadratic_variables() > 0: self.problem.objective.set_quadratic([0. for _ in range(self.problem.variables.get_num())]) super(Model, self.__class__).objective.fset(self, value) @@ -721,7 +723,9 @@ def objective(self, value): # self.problem.objective.set_offset(float(offset)) # Not available prior to 12.6.2 self._objective_offset = offset if linear_coefficients: - self.problem.objective.set_linear([var.name, float(coef)] for var, coef in linear_coefficients.items()) + name_list = [var.name for var in linear_coefficients] + index_dict = {n: i for n, i in zip(name_list, self._get_variable_indices(name_list))} + self.problem.objective.set_linear([index_dict[var.name], float(coef)] for var, coef in linear_coefficients.items()) for key, coef in quadratic_coeffients.items(): if len(key) == 1: @@ -901,3 +905,16 @@ def _get_quadratic_expression(self, quadratic=None): else: pass # Only look at upper triangle return add(terms) + + def _get_variable_indices(self, names): + # Cplex does not keep an index of variable names + # Getting indices thus takes roughly quadratic time + # If many indices are required an alternate and faster method is used, where the full name list must only + # be traversed once + if len(names) < 1000: + return self.problem.variables.get_indices(names) + else: + name_set = set(names) + all_names = self.problem.variables.get_names() + indices = {n: i for i, n in enumerate(all_names) if n in name_set} + return [indices[n] for n in names] diff --git a/optlang/duality.py b/optlang/duality.py index f7a47495..beef5ec1 100644 --- a/optlang/duality.py +++ b/optlang/duality.py @@ -13,8 +13,7 @@ # See the License for the specific language governing permissions and # limitations under the License. -from sympy import Add, Mul - +import optlang # This function is very complex. Should maybe be refactored def convert_linear_problem_to_dual(model, sloppy=False, infinity=None, maintain_standard_form=True, prefix="dual_", dual_model=None): # NOQA @@ -84,6 +83,8 @@ def convert_linear_problem_to_dual(model, sloppy=False, infinity=None, maintain_ if constraint.lb != 0: dual_objective[const_var] = sign * constraint.lb for variable, coef in constraint.expression.as_coefficients_dict().items(): + if variable == 1: + continue coefficients.setdefault(variable.name, {})[const_var] = sign * coef else: if constraint.lb is not None: @@ -105,6 +106,8 @@ def convert_linear_problem_to_dual(model, sloppy=False, infinity=None, maintain_ coefficients_dict = {constraint.expression.args[1]: constraint.expression.args[0]} for variable, coef in coefficients_dict.items(): + if variable == 1: + continue if constraint.lb is not None: coefficients.setdefault(variable.name, {})[lb_var] = -sign * coef if constraint.ub is not None: @@ -131,7 +134,7 @@ def convert_linear_problem_to_dual(model, sloppy=False, infinity=None, maintain_ # Add dual constraints from primal objective primal_objective_dict = model.objective.expression.as_coefficients_dict() for variable in model.variables: - expr = Add(*((coef * dual_var) for dual_var, coef in coefficients[variable.name].items())) + expr = optlang.symbolics.add([(coef * dual_var) for dual_var, coef in coefficients[variable.name].items()]) obj_coef = primal_objective_dict[variable] if maximization: const = model.interface.Constraint(expr, lb=obj_coef, name=prefix + variable.name) @@ -140,7 +143,7 @@ def convert_linear_problem_to_dual(model, sloppy=False, infinity=None, maintain_ dual_model.add(const) # Make dual objective - expr = Add(*((coef * dual_var) for dual_var, coef in dual_objective.items() if coef != 0)) + expr = optlang.symbolics.add([(coef * dual_var) for dual_var, coef in dual_objective.items() if coef != 0]) if maximization: objective = model.interface.Objective(expr, direction="min") else: diff --git a/optlang/glpk_interface.py b/optlang/glpk_interface.py index b7f9952c..6eb2562b 100644 --- a/optlang/glpk_interface.py +++ b/optlang/glpk_interface.py @@ -50,7 +50,7 @@ glp_get_mat_row, glp_get_row_ub, glp_get_row_type, glp_get_row_lb, glp_get_row_name, glp_get_obj_coef, \ glp_get_obj_dir, glp_scale_prob, GLP_SF_AUTO, glp_get_num_int, glp_get_num_bin, glp_mip_col_val, \ glp_mip_obj_val, glp_mip_status, GLP_ETMLIM, glp_adv_basis, glp_read_lp, glp_mip_row_val, \ - get_col_primals, get_col_duals, get_row_primals, get_row_duals + get_col_primals, get_col_duals, get_row_primals, get_row_duals, glp_delete_prob @@ -579,6 +579,10 @@ def __setstate__(self, repr_dict): if repr_dict['glpk_status'] == 'optimal': self.optimize() # since the start is an optimal solution, nothing will happen here + # def __del__(self): # To make sure that the glpk problem is deleted when this is garbage collected + # Gotcha: When objects with a __del__ method are part of a referencing cycle, the entire cycle is never automatically garbage collected + # glp_delete_prob(self.problem) + @property def objective(self): return self._objective diff --git a/optlang/gurobi_interface.py b/optlang/gurobi_interface.py index 94902956..2191be1c 100644 --- a/optlang/gurobi_interface.py +++ b/optlang/gurobi_interface.py @@ -27,6 +27,7 @@ import os import six +from functools import partial from optlang import interface from optlang.util import inheritdocstring, TemporaryFilename from optlang.expression_parsing import parse_optimization_expression @@ -236,7 +237,14 @@ def problem(self, value): @property def primal(self): if self.problem is not None: - return self._internal_constraint.Slack + aux_var = self.problem.problem.getVarByName(self._internal_constraint.getAttr('ConstrName') + '_aux') + if aux_var is not None: + aux_val = aux_var.X + else: + aux_val = 0 + return (self._internal_constraint.RHS - + self._internal_constraint.Slack + + aux_val) # return self._round_primal_to_bounds(primal_from_solver) # Test assertions fail # return primal_from_solver else: @@ -325,9 +333,10 @@ def __init__(self, expression, sloppy=False, *args, **kwargs): @property def value(self): if getattr(self, "problem", None) is not None: + gurobi_problem = self.problem.problem try: - return self.problem.problem.getAttr("ObjVal") + getattr(self.problem, "_objective_offset", 0) - except gurobipy.GurobiError: # TODO: let this check the actual error message + return gurobi_problem.getAttr("ObjVal") + getattr(self.problem, "_objective_offset", 0) + except (gurobipy.GurobiError, AttributeError): # TODO: let this check the actual error message return None else: return None @@ -359,8 +368,9 @@ def _get_expression(self): if self.problem is not None and self._expression_expired and len(self.problem._variables) > 0: grb_obj = self.problem.problem.getObjective() terms = [] + variables = self.problem._variables for i in range(grb_obj.size()): - terms.append(grb_obj.getCoeff(i) * self.problem.variables[grb_obj.getVar(i).getAttr('VarName')]) + terms.append(grb_obj.getCoeff(i) * variables[grb_obj.getVar(i).getAttr('VarName')]) expression = symbolics.add(terms) # TODO implement quadratic objectives self._expression = expression + getattr(self.problem, "_objective_offset", 0) @@ -385,7 +395,8 @@ def lp_method(self): def lp_method(self, value): if value not in _LP_METHODS: raise ValueError("Invalid LP method. Please choose among: " + str(list(_LP_METHODS))) - self.problem.problem.params.Method = _LP_METHODS[value] + if getattr(self, "problem", None) is not None: + self.problem.problem.params.Method = _LP_METHODS[value] @property def presolve(self): @@ -393,12 +404,13 @@ def presolve(self): @presolve.setter def presolve(self, value): - if value is True: - self.problem.problem.params.Presolve = 2 - elif value is False: - self.problem.problem.params.Presolve = 0 - else: - self.problem.problem.params.Presolve = -1 + if getattr(self, "problem", None) is not None: + if value is True: + self.problem.problem.params.Presolve = 2 + elif value is False: + self.problem.problem.params.Presolve = 0 + else: + self.problem.problem.params.Presolve = -1 self._presolve = value @property @@ -427,20 +439,37 @@ def timeout(self, value): self.problem.problem.params.TimeLimit = value self._timeout = value + def __getstate__(self): + return {'presolve': self.presolve, 'verbosity': self.verbosity, 'timeout': self.timeout} + + def __setstate__(self, state): + self.__init__() + for key, val in six.iteritems(state): + setattr(self, key, val) + + def _get_feasibility(self): + return getattr(self.problem.problem.params, "FeasibilityTol") + + def _set_feasibility(self, value): + return setattr(self.problem.problem.params, "FeasibilityTol", value) + + def _get_optimality(self): + return getattr(self.problem.problem.params, "OptimalityTol") + + def _set_optimality(self, value): + return setattr(self.problem.problem.params, "OptimalityTol", value) + + def _get_integrality(self): + return getattr(self.problem.problem.params, "IntFeasTol") + + def _set_integrality(self, value): + return setattr(self.problem.problem.params, "IntFeasTol", value) + def _tolerance_functions(self): return { - "feasibility": ( - lambda: self.problem.problem.params.FeasibilityTol, - lambda x: setattr(self.problem.problem.params, "FeasibilityTol", x) - ), - "optimality": ( - lambda: self.problem.problem.params.OptimalityTol, - lambda x: setattr(self.problem.problem.params, "OptimalityTol", x) - ), - "integrality": ( - lambda: self.problem.problem.params.IntFeasTol, - lambda x: setattr(self.problem.problem.params, "IntFeasTol", x) - ) + "feasibility": (self._get_feasibility, self._set_feasibility), + "optimality": (self._get_optimality, self._set_optimality), + "integrality": (self._get_integrality, self._set_integrality) } diff --git a/optlang/symbolics.py b/optlang/symbolics.py index 33a81e59..db0d9870 100644 --- a/optlang/symbolics.py +++ b/optlang/symbolics.py @@ -72,6 +72,9 @@ def __repr__(self): def __str__(self): return self._name + def __getnewargs__(self): + return (self._name, {}) + def add(*args): if len(args) == 1: args = args[0] diff --git a/optlang/tests/abstract_test_cases.py b/optlang/tests/abstract_test_cases.py index 432d8887..49caca0b 100644 --- a/optlang/tests/abstract_test_cases.py +++ b/optlang/tests/abstract_test_cases.py @@ -731,7 +731,7 @@ def test_objective_expression_includes_constant(self): objective = self.model.objective self.model.objective = self.interface.Objective(objective.expression + 3) self.model.update() - self.assertEqual((self.model.objective.expression - (objective.expression + 3)).expand(), 0) + self.assertEqual((self.model.objective.expression - (objective.expression + 3.)).expand(), 0.) def test_is_integer(self): model = self.model @@ -809,6 +809,19 @@ def test_integer_batch_duals(self): self.assertEqual(model.reduced_costs[x.name], 0) self.assertEqual(model.shadow_prices[c.name], 1) + def test_large_objective(self): + model = self.interface.Model() + model.add([self.interface.Variable(str(i), lb=1) for i in range(1100)]) + model.optimize() + + obj = self.interface.Objective( + optlang.symbolics.add([optlang.symbolics.mul((optlang.symbolics.One, v)) for v in model.variables]), + direction="min" + ) + model.objective = obj + model.optimize() + self.assertAlmostEqual(model.objective.value, len(model.variables)) + @six.add_metaclass(abc.ABCMeta) class AbstractConfigurationTestCase(unittest.TestCase): diff --git a/optlang/tests/test_expression_parsing.py b/optlang/tests/test_expression_parsing.py index 2a0f458b..2b307b82 100644 --- a/optlang/tests/test_expression_parsing.py +++ b/optlang/tests/test_expression_parsing.py @@ -7,6 +7,20 @@ import unittest +def _quad_terms_to_expression(terms): + args = [] + for term, coef in terms.items(): + term = list(term) + args.append(coef * term[0] * term[-1]) + return sum(args) + + +def _compare_term_dicts(test_case, dict1, dict2): + for term, coef1 in dict1.items(): + coef2 = dict2[term] + test_case.assertEqual(coef1 - coef2, 0) + + class ExpressionParsingTestCase(unittest.TestCase): def setUp(self): self.vars = [Variable(name) for name in ["x", "y", "z", "v", "w"]] @@ -22,8 +36,8 @@ def test_parse_linear_expression(self): self.assertEqual(offset_const, 0) self.assertEqual(offset_obj, offset) - self.assertEqual(linear_terms_const, target) - self.assertEqual(linear_terms_obj, target) + _compare_term_dicts(self, linear_terms_const, target) + _compare_term_dicts(self, linear_terms_obj, target) self.assertEqual(quad_terms_const, {}) self.assertEqual(quad_terms_obj, {}) @@ -41,8 +55,10 @@ def test_parse_quadratic_expression(self): self.assertEqual(offset_obj, offset) self.assertEqual(linear_terms_const, {}) self.assertEqual(linear_terms_obj, {}) - self.assertEqual(quad_terms_const, target) - self.assertEqual(quad_terms_obj, target) + _compare_term_dicts(self, quad_terms_const, target) + _compare_term_dicts(self, quad_terms_obj, target) + self.assertEqual((_quad_terms_to_expression(quad_terms_obj) - (expr - offset)).expand(), 0) + self.assertEqual((_quad_terms_to_expression(quad_terms_const) - (expr - offset)).expand(), 0) def test_parse_non_expanded_quadratic_expression(self): x, y, z = self.vars[:3] @@ -57,7 +73,7 @@ def test_parse_non_expanded_quadratic_expression(self): self.assertEqual(offset_const, -4) self.assertEqual(offset_obj, -4 + offset) - self.assertEqual(linear_terms_const, linear_target) - self.assertEqual(linear_terms_obj, linear_target) - self.assertEqual(quad_terms_const, target) - self.assertEqual(quad_terms_obj, target) + _compare_term_dicts(self, linear_terms_const, linear_target) + _compare_term_dicts(self, linear_terms_obj, linear_target) + _compare_term_dicts(self, quad_terms_const, target) + _compare_term_dicts(self, quad_terms_obj, target) diff --git a/optlang/tests/test_scipy_interface.py b/optlang/tests/test_scipy_interface.py index 5b0b022a..ab600f54 100644 --- a/optlang/tests/test_scipy_interface.py +++ b/optlang/tests/test_scipy_interface.py @@ -315,3 +315,6 @@ def test_integer_constraint_dual(self): def test_integer_batch_duals(self): self.skipTest("No duals with scipy") + + def test_large_objective(self): + self.skipTest("Quite slow and not necessary")