diff --git a/.github/workflows/test_branches.yml b/.github/workflows/test_branches.yml index e773587ec85..6d803b484b7 100644 --- a/.github/workflows/test_branches.yml +++ b/.github/workflows/test_branches.yml @@ -573,6 +573,13 @@ jobs: $PYTHON_EXE -m pip install --cache-dir cache/pip highspy \ || echo "WARNING: highspy is not available" + - name: Install COPT + if: ${{ ! matrix.slim }} + shell: bash + run: | + $PYTHON_EXE -m pip install --cache-dir cache/pip coptpy \ + || echo "WARNING: coptpy is not available" + - name: Set up coverage tracking run: | if test "${{matrix.TARGET}}" == win; then diff --git a/.github/workflows/test_pr_and_main.yml b/.github/workflows/test_pr_and_main.yml index 2885fd107a8..87cf4bb8890 100644 --- a/.github/workflows/test_pr_and_main.yml +++ b/.github/workflows/test_pr_and_main.yml @@ -603,6 +603,13 @@ jobs: $PYTHON_EXE -m pip install --cache-dir cache/pip highspy \ || echo "WARNING: highspy is not available" + - name: Install COPT + if: ${{ ! matrix.slim }} + shell: bash + run: | + $PYTHON_EXE -m pip install --cache-dir cache/pip coptpy \ + || echo "WARNING: coptpy is not available" + - name: Set up coverage tracking run: | if test "${{matrix.TARGET}}" == win; then diff --git a/pyomo/contrib/appsi/base.py b/pyomo/contrib/appsi/base.py index d279b8b6015..db42cc965ac 100644 --- a/pyomo/contrib/appsi/base.py +++ b/pyomo/contrib/appsi/base.py @@ -1654,7 +1654,7 @@ def license_is_valid(self) -> bool: @property def options(self): - for solver_name in ['gurobi', 'ipopt', 'cplex', 'cbc', 'highs']: + for solver_name in ['gurobi', 'ipopt', 'cplex', 'cbc', 'highs', 'copt']: if hasattr(self, solver_name + '_options'): return getattr(self, solver_name + '_options') raise NotImplementedError('Could not find the correct options') @@ -1662,7 +1662,7 @@ def options(self): @options.setter def options(self, val): found = False - for solver_name in ['gurobi', 'ipopt', 'cplex', 'cbc', 'highs']: + for solver_name in ['gurobi', 'ipopt', 'cplex', 'cbc', 'highs', 'copt']: if hasattr(self, solver_name + '_options'): setattr(self, solver_name + '_options', val) found = True diff --git a/pyomo/contrib/appsi/plugins.py b/pyomo/contrib/appsi/plugins.py index 5333158239e..72e432eb73f 100644 --- a/pyomo/contrib/appsi/plugins.py +++ b/pyomo/contrib/appsi/plugins.py @@ -1,6 +1,6 @@ from pyomo.common.extensions import ExtensionBuilderFactory from .base import SolverFactory -from .solvers import Gurobi, Ipopt, Cbc, Cplex, Highs +from .solvers import Gurobi, Ipopt, Cbc, Cplex, Highs, Copt from .build import AppsiBuilder @@ -21,3 +21,6 @@ def load(): SolverFactory.register( name='appsi_highs', doc='Automated persistent interface to Highs' )(Highs) + SolverFactory.register( + name='appsi_copt', doc='Automated persistent interface to COPT' + )(Copt) diff --git a/pyomo/contrib/appsi/solvers/__init__.py b/pyomo/contrib/appsi/solvers/__init__.py index df58a0cb245..3c303cadeab 100644 --- a/pyomo/contrib/appsi/solvers/__init__.py +++ b/pyomo/contrib/appsi/solvers/__init__.py @@ -3,3 +3,4 @@ from .cbc import Cbc from .cplex import Cplex from .highs import Highs +from .copt import Copt, CoptResults diff --git a/pyomo/contrib/appsi/solvers/copt.py b/pyomo/contrib/appsi/solvers/copt.py new file mode 100644 index 00000000000..03606eef51c --- /dev/null +++ b/pyomo/contrib/appsi/solvers/copt.py @@ -0,0 +1,1223 @@ +import logging +import math +from typing import List, Dict, Optional +from pyomo.common.collections import ComponentMap, OrderedSet +from pyomo.common.log import LogStream +from pyomo.common.dependencies import attempt_import +from pyomo.common.errors import PyomoException +from pyomo.common.tee import capture_output, TeeStream +from pyomo.common.timing import HierarchicalTimer +from pyomo.common.shutdown import python_is_shutting_down +from pyomo.common.config import ConfigValue, NonNegativeInt +from pyomo.core.kernel.objective import minimize, maximize +from pyomo.core.base import SymbolMap, NumericLabeler, TextLabeler +from pyomo.core.base.var import _GeneralVarData +from pyomo.core.base.constraint import _GeneralConstraintData +from pyomo.core.base.sos import _SOSConstraintData +from pyomo.core.base.param import _ParamData +from pyomo.core.expr.numvalue import value, is_constant, native_numeric_types +from pyomo.repn import generate_standard_repn +from pyomo.core.expr.numeric_expr import NPV_MaxExpression, NPV_MinExpression +from pyomo.contrib.appsi.base import ( + PersistentSolver, + Results, + TerminationCondition, + MIPSolverConfig, + PersistentBase, + PersistentSolutionLoader, +) +from pyomo.contrib.appsi.cmodel import cmodel, cmodel_available +from pyomo.core.staleflag import StaleFlagManager +import sys + +logger = logging.getLogger(__name__) + + +def _import_coptpy(): + try: + import coptpy + except ImportError: + Copt._available = Copt.Availability.NotFound + raise + if coptpy.COPT.VERSION_MAJOR < 6: + Copt._available = Copt.Availability.BadVersion + raise ImportError('The APPSI COPT interface requires coptpy>=6.0.0') + return coptpy + + +coptpy, coptpy_available = attempt_import('coptpy', importer=_import_coptpy) + + +class DegreeError(PyomoException): + pass + + +class CoptConfig(MIPSolverConfig): + def __init__( + self, + description=None, + doc=None, + implicit=False, + implicit_domain=None, + visibility=0, + ): + super(CoptConfig, self).__init__( + description=description, + doc=doc, + implicit=implicit, + implicit_domain=implicit_domain, + visibility=visibility, + ) + + self.declare('logfile', ConfigValue(domain=str)) + self.declare('solver_output_logger', ConfigValue()) + self.declare('log_level', ConfigValue(domain=NonNegativeInt)) + + self.logfile = '' + self.solver_output_logger = logger + self.log_level = logging.INFO + + +class CoptSolutionLoader(PersistentSolutionLoader): + def load_vars(self, vars_to_load=None, solution_number=0): + self._assert_solution_still_valid() + self._solver.load_vars( + vars_to_load=vars_to_load, solution_number=solution_number + ) + + def get_primals(self, vars_to_load=None, solution_number=0): + self._assert_solution_still_valid() + return self._solver.get_primals( + vars_to_load=vars_to_load, solution_number=solution_number + ) + + +class CoptResults(Results): + def __init__(self, solver): + super(CoptResults, self).__init__() + self.wallclock_time = None + self.solution_loader = CoptSolutionLoader(solver=solver) + + +class _MutableLowerBound(object): + def __init__(self, expr): + self.var = None + self.expr = expr + + def update(self): + self.var.lb = value(self.expr) + + +class _MutableUpperBound(object): + def __init__(self, expr): + self.var = None + self.expr = expr + + def update(self): + self.var.ub = value(self.expr) + + +class _MutableLinearCoefficient(object): + def __init__(self): + self.expr = None + self.var = None + self.con = None + self.copt_model = None + + def update(self): + self.copt_model.setCoeff(self.con, self.var, value(self.expr)) + + +class _MutableRangeConstant(object): + def __init__(self): + self.lhs_expr = None + self.rhs_expr = None + self.con = None + self.copt_model = None + + def update(self): + lhs_val = value(self.lhs_expr) + rhs_val = value(self.rhs_expr) + self.con.lb = lhs_val + self.con.ub = rhs_val + + +class _MutableConstant(object): + def __init__(self): + self.expr = None + self.con = None + + def update(self): + if isinstance(self.con, coptpy.Constraint): + lb = self.con.lb + ub = self.con.ub + # lb <= con <= ub + if lb > -coptpy.COPT.INFINITY and ub < +coptpy.COPT.INFINITY: + self.con.lb = value(self.expr) + self.con.ub = value(self.expr) + # con >= lb + elif lb > -coptpy.COPT.INFINITY and ub >= +coptpy.COPT.INFINITY: + self.con.lb = value(self.expr) + # con <= ub + elif lb <= -coptpy.COPT.INFINITY and ub < +coptpy.COPT.INFINITY: + self.con.ub = value(self.expr) + else: + raise ValueError( + "Constraint does not has lower/upper bound: {0} \n".format(self.con) + ) + elif isinstance(self.con, coptpy.QConstraint): + self.con.rhs = value(self.expr) + else: + raise TypeError("Unexpected constraint type") + + +class _MutableQuadraticConstraint(object): + def __init__(self, copt_model, copt_con, constant, linear_coefs, quadratic_coefs): + self.con = copt_con + self.copt_model = copt_model + self.constant = constant + self.last_constant_value = value(self.constant.expr) + self.linear_coefs = linear_coefs + self.last_linear_coef_values = [value(i.expr) for i in self.linear_coefs] + self.quadratic_coefs = quadratic_coefs + self.last_quadratic_coef_values = [value(i.expr) for i in self.quadratic_coefs] + + def get_updated_expression(self): + copt_expr = self.copt_model.getQuadRow(self.con) + for ndx, coef in enumerate(self.linear_coefs): + current_coef_value = value(coef.expr) + incremental_coef_value = ( + current_coef_value - self.last_linear_coef_values[ndx] + ) + copt_expr += incremental_coef_value * coef.var + self.last_linear_coef_values[ndx] = current_coef_value + for ndx, coef in enumerate(self.quadratic_coefs): + current_coef_value = value(coef.expr) + incremental_coef_value = ( + current_coef_value - self.last_quadratic_coef_values[ndx] + ) + copt_expr += incremental_coef_value * coef.var1 * coef.var2 + self.last_quadratic_coef_values[ndx] = current_coef_value + return copt_expr + + def get_updated_rhs(self): + return value(self.constant.expr) + + +class _MutableObjective(object): + def __init__(self, copt_model, constant, linear_coefs, quadratic_coefs): + self.copt_model = copt_model + self.constant = constant + self.linear_coefs = linear_coefs + self.quadratic_coefs = quadratic_coefs + self.last_quadratic_coef_values = [value(i.expr) for i in self.quadratic_coefs] + + def get_updated_expression(self): + for ndx, coef in enumerate(self.linear_coefs): + coef.var.obj = value(coef.expr) + self.copt_model.objconst = value(self.constant.expr) + + copt_expr = None + for ndx, coef in enumerate(self.quadratic_coefs): + if value(coef.expr) != self.last_quadratic_coef_values[ndx]: + if copt_expr is None: + copt_expr = self.copt_model.getObjective() + current_coef_value = value(coef.expr) + incremental_coef_value = ( + current_coef_value - self.last_quadratic_coef_values[ndx] + ) + copt_expr += incremental_coef_value * coef.var1 * coef.var2 + self.last_quadratic_coef_values[ndx] = current_coef_value + return copt_expr + + +class _MutableQuadraticCoefficient(object): + def __init__(self): + self.expr = None + self.var1 = None + self.var2 = None + + +class Copt(PersistentBase, PersistentSolver): + """ + Interface to Copt + """ + + _available = None + _num_instances = 0 + _coptenv = None + + def __init__(self, only_child_vars=True): + super(Copt, self).__init__(only_child_vars=only_child_vars) + self._num_instances += 1 + self._config = CoptConfig() + + self._solver_options = dict() + self._solver_model = None + + if coptpy_available and self._coptenv is None: + self._coptenv = coptpy.Envr() + + self._symbol_map = SymbolMap() + self._labeler = None + + self._pyomo_var_to_solver_var_map = dict() + self._pyomo_con_to_solver_con_map = dict() + self._pyomo_sos_to_solver_sos_map = dict() + self._solver_con_to_pyomo_con_map = dict() + + self._mutable_helpers = dict() + self._mutable_bounds = dict() + self._mutable_quadratic_helpers = dict() + self._mutable_objective = None + + self._last_results_object: Optional[CoptResults] = None + + def available(self): + if not coptpy_available: + return self.Availability.NotFound + elif self._available == self.Availability.BadVersion: + return self.Availability.BadVersion + else: + self._available = Copt.Availability.BadLicense + if self._coptenv is not None: + m = self._coptenv.createModel('checklic') + m.setParam("Logging", 0) + try: + # COPT can solve LP up to 10k variables without license + m.addVars(10001) + m.solveLP() + self._available = Copt.Availability.FullLicense + except coptpy.CoptError: + self._available = Copt.Availability.LimitedLicense + return self._available + + def release_license(self): + self._reinit() + + def __del__(self): + if not python_is_shutting_down(): + self._num_instances -= 1 + if self._num_instances == 0: + self.release_license() + + def version(self): + version = ( + coptpy.COPT.VERSION_MAJOR, + coptpy.COPT.VERSION_MINOR, + coptpy.COPT.VERSION_TECHNICAL, + ) + return version + + @property + def config(self) -> CoptConfig: + return self._config + + @config.setter + def config(self, val: CoptConfig): + self._config = val + + @property + def copt_options(self): + """ + A dictionary mapping solver options to values for those options. + These are solver specific. + + Returns + ------- + dict + A dictionary mapping solver options to values for those options. + """ + return self._solver_options + + @copt_options.setter + def copt_options(self, val: Dict): + self._solver_options = val + + @property + def symbol_map(self): + return self._symbol_map + + def _solve(self, timer: HierarchicalTimer): + ostreams = [ + LogStream( + level=self.config.log_level, logger=self.config.solver_output_logger + ) + ] + if self.config.stream_solver: + ostreams.append(sys.stdout) + + with TeeStream(*ostreams) as t: + with capture_output(output=t.STDOUT, capture_fd=False): + config = self.config + options = self.copt_options + + if not config.stream_solver: + self._solver_model.setParam("LogToConsole", 0) + if config.logfile: + self._solver_model.setLogFile(config.logfile) + + if config.time_limit is not None: + self._solver_model.setParam('TimeLimit', config.time_limit) + if config.mip_gap is not None: + self._solver_model.setParam('RelGap', config.mip_gap) + + for key, option in options.items(): + self._solver_model.setParam(key, option) + + timer.start('solve') + self._solver_model.solve() + timer.stop('solve') + + return self._postsolve(timer) + + def solve(self, model, timer: HierarchicalTimer = None) -> Results: + StaleFlagManager.mark_all_as_stale() + + if self._last_results_object is not None: + self._last_results_object.solution_loader.invalidate() + if timer is None: + timer = HierarchicalTimer() + if model is not self._model: + timer.start('set_instance') + self.set_instance(model) + timer.stop('set_instance') + else: + timer.start('update') + self.update(timer=timer) + timer.stop('update') + res = self._solve(timer) + self._last_results_object = res + if self.config.report_timing: + logger.info('\n' + str(timer)) + return res + + def _process_domain_and_bounds( + self, var, var_id, mutable_lbs, mutable_ubs, ndx, coptpy_var + ): + _v, _lb, _ub, _fixed, _domain_interval, _value = self._vars[id(var)] + lb, ub, step = _domain_interval + if lb is None: + lb = -coptpy.COPT.INFINITY + if ub is None: + ub = coptpy.COPT.INFINITY + if step == 0: + vtype = coptpy.COPT.CONTINUOUS + elif step == 1: + if lb == 0 and ub == 1: + vtype = coptpy.COPT.BINARY + else: + vtype = coptpy.COPT.INTEGER + else: + raise ValueError( + f'Unrecognized domain step: {step} (should be either 0 or 1)' + ) + if _fixed: + lb = _value + ub = _value + else: + if _lb is not None: + if not is_constant(_lb): + mutable_bound = _MutableLowerBound(NPV_MaxExpression((_lb, lb))) + if coptpy_var is None: + mutable_lbs[ndx] = mutable_bound + else: + mutable_bound.var = coptpy_var + self._mutable_bounds[var_id, 'lb'] = (var, mutable_bound) + lb = max(value(_lb), lb) + if _ub is not None: + if not is_constant(_ub): + mutable_bound = _MutableUpperBound(NPV_MinExpression((_ub, ub))) + if coptpy_var is None: + mutable_ubs[ndx] = mutable_bound + else: + mutable_bound.var = coptpy_var + self._mutable_bounds[var_id, 'ub'] = (var, mutable_bound) + ub = min(value(_ub), ub) + + return lb, ub, vtype + + def _add_variables(self, variables: List[_GeneralVarData]): + var_names = list() + vtypes = list() + lbs = list() + ubs = list() + mutable_lbs = dict() + mutable_ubs = dict() + for ndx, var in enumerate(variables): + varname = self._symbol_map.getSymbol(var, self._labeler) + lb, ub, vtype = self._process_domain_and_bounds( + var, id(var), mutable_lbs, mutable_ubs, ndx, None + ) + var_names.append(varname) + vtypes.append(vtype) + lbs.append(lb) + ubs.append(ub) + + nvars = len(variables) + + copt_vars = list() + for i in range(nvars): + copt_vars.append( + self._solver_model.addVar( + lb=lbs[i], ub=ubs[i], vtype=vtypes[i], name=var_names[i] + ) + ) + + for ndx, pyomo_var in enumerate(variables): + copt_var = copt_vars[ndx] + self._pyomo_var_to_solver_var_map[id(pyomo_var)] = copt_var + for ndx, mutable_bound in mutable_lbs.items(): + mutable_bound.var = copt_vars[ndx] + for ndx, mutable_bound in mutable_ubs.items(): + mutable_bound.var = copt_vars[ndx] + + def _add_params(self, params: List[_ParamData]): + pass + + def _reinit(self): + saved_config = self.config + saved_options = self.copt_options + saved_update_config = self.update_config + self.__init__(only_child_vars=self._only_child_vars) + self.config = saved_config + self.copt_options = saved_options + self.update_config = saved_update_config + + def set_instance(self, model): + if self._last_results_object is not None: + self._last_results_object.solution_loader.invalidate() + if not self.available(): + c = self.__class__ + raise PyomoException( + f'Solver {c.__module__}.{c.__qualname__} is not available ' + f'({self.available()}).' + ) + self._reinit() + self._model = model + if self.use_extensions and cmodel_available: + self._expr_types = cmodel.PyomoExprTypes() + + if self.config.symbolic_solver_labels: + self._labeler = TextLabeler() + else: + self._labeler = NumericLabeler('x') + + if model.name is not None: + self._solver_model = self._coptenv.createModel(model.name) + else: + self._solver_model = self._coptenv.createModel() + + self.add_block(model) + if self._objective is None: + self.set_objective(None) + + def _get_expr_from_pyomo_expr(self, expr): + mutable_linear_coefficients = list() + mutable_quadratic_coefficients = list() + repn = generate_standard_repn(expr, quadratic=True, compute_values=False) + + degree = repn.polynomial_degree() + if (degree is None) or (degree > 2): + raise DegreeError( + 'CoptAPPSI does not support expressions of degree {0}.'.format(degree) + ) + + if len(repn.quadratic_vars) > 0: + new_expr = coptpy.QuadExpr(0.0) + else: + new_expr = coptpy.LinExpr(0.0) + + if len(repn.linear_vars) > 0: + linear_coef_vals = list() + for ndx, coef in enumerate(repn.linear_coefs): + if not is_constant(coef): + mutable_linear_coefficient = _MutableLinearCoefficient() + mutable_linear_coefficient.expr = coef + mutable_linear_coefficient.var = self._pyomo_var_to_solver_var_map[ + id(repn.linear_vars[ndx]) + ] + mutable_linear_coefficients.append(mutable_linear_coefficient) + linear_coef_vals.append(value(coef)) + new_expr += coptpy.LinExpr( + [self._pyomo_var_to_solver_var_map[id(i)] for i in repn.linear_vars], + linear_coef_vals, + ) + + for ndx, v in enumerate(repn.quadratic_vars): + x, y = v + copt_x = self._pyomo_var_to_solver_var_map[id(x)] + copt_y = self._pyomo_var_to_solver_var_map[id(y)] + coef = repn.quadratic_coefs[ndx] + if not is_constant(coef): + mutable_quadratic_coefficient = _MutableQuadraticCoefficient() + mutable_quadratic_coefficient.expr = coef + mutable_quadratic_coefficient.var1 = copt_x + mutable_quadratic_coefficient.var2 = copt_y + mutable_quadratic_coefficients.append(mutable_quadratic_coefficient) + coef_val = value(coef) + new_expr += coef_val * copt_x * copt_y + + return ( + new_expr, + repn.constant, + mutable_linear_coefficients, + mutable_quadratic_coefficients, + ) + + def _add_constraints(self, cons: List[_GeneralConstraintData]): + for con in cons: + conname = self._symbol_map.getSymbol(con, self._labeler) + ( + copt_expr, + repn_constant, + mutable_linear_coefficients, + mutable_quadratic_coefficients, + ) = self._get_expr_from_pyomo_expr(con.body) + + if ( + copt_expr.__class__ in {coptpy.LinExpr, coptpy.Var} + or copt_expr.__class__ in native_numeric_types + ): + if con.equality: + rhs_expr = con.lower - repn_constant + rhs_val = value(rhs_expr) + coptpy_con = self._solver_model.addConstr( + copt_expr == rhs_val, name=conname + ) + + if not is_constant(rhs_expr): + mutable_constant = _MutableConstant() + mutable_constant.expr = rhs_expr + mutable_constant.con = coptpy_con + self._mutable_helpers[con] = [mutable_constant] + elif con.has_lb() and con.has_ub(): + lhs_expr = con.lower - repn_constant + rhs_expr = con.upper - repn_constant + lhs_val = value(lhs_expr) + rhs_val = value(rhs_expr) + coptpy_con = self._solver_model.addBoundConstr( + copt_expr, lhs_val, rhs_val, name=conname + ) + + if not is_constant(lhs_expr) or not is_constant(rhs_expr): + mutable_range_constant = _MutableRangeConstant() + mutable_range_constant.lhs_expr = lhs_expr + mutable_range_constant.rhs_expr = rhs_expr + mutable_range_constant.con = coptpy_con + mutable_range_constant.copt_model = self._solver_model + self._mutable_helpers[con] = [mutable_range_constant] + elif con.has_lb(): + rhs_expr = con.lower - repn_constant + rhs_val = value(rhs_expr) + coptpy_con = self._solver_model.addConstr( + copt_expr >= rhs_val, name=conname + ) + if not is_constant(rhs_expr): + mutable_constant = _MutableConstant() + mutable_constant.expr = rhs_expr + mutable_constant.con = coptpy_con + self._mutable_helpers[con] = [mutable_constant] + elif con.has_ub(): + rhs_expr = con.upper - repn_constant + rhs_val = value(rhs_expr) + coptpy_con = self._solver_model.addConstr( + copt_expr <= rhs_val, name=conname + ) + if not is_constant(rhs_expr): + mutable_constant = _MutableConstant() + mutable_constant.expr = rhs_expr + mutable_constant.con = coptpy_con + self._mutable_helpers[con] = [mutable_constant] + else: + raise ValueError( + "Constraint does not have a lower or an upper bound: {0} \n".format( + con + ) + ) + + for tmp in mutable_linear_coefficients: + tmp.con = coptpy_con + tmp.copt_model = self._solver_model + if len(mutable_linear_coefficients) > 0: + if con not in self._mutable_helpers: + self._mutable_helpers[con] = mutable_linear_coefficients + else: + self._mutable_helpers[con].extend(mutable_linear_coefficients) + elif copt_expr.__class__ is coptpy.QuadExpr: + if con.equality: + rhs_expr = con.lower - repn_constant + rhs_val = value(rhs_expr) + coptpy_con = self._solver_model.addQConstr( + copt_expr == rhs_val, name=conname + ) + elif con.has_lb() and con.has_ub(): + raise NotImplementedError( + 'Quadratic range constraints are not supported' + ) + elif con.has_lb(): + rhs_expr = con.lower - repn_constant + rhs_val = value(rhs_expr) + coptpy_con = self._solver_model.addQConstr( + copt_expr >= rhs_val, name=conname + ) + elif con.has_ub(): + rhs_expr = con.upper - repn_constant + rhs_val = value(rhs_expr) + coptpy_con = self._solver_model.addQConstr( + copt_expr <= rhs_val, name=conname + ) + else: + raise ValueError( + "Constraint does not have a lower or an upper bound: {0} \n".format( + con + ) + ) + if ( + len(mutable_linear_coefficients) > 0 + or len(mutable_quadratic_coefficients) > 0 + or not is_constant(repn_constant) + ): + mutable_constant = _MutableConstant() + mutable_constant.expr = rhs_expr + mutable_quadratic_constraint = _MutableQuadraticConstraint( + self._solver_model, + coptpy_con, + mutable_constant, + mutable_linear_coefficients, + mutable_quadratic_coefficients, + ) + self._mutable_quadratic_helpers[con] = mutable_quadratic_constraint + else: + raise ValueError( + 'Unrecognized COPT expression type: ' + str(copt_expr.__class__) + ) + + self._pyomo_con_to_solver_con_map[con] = coptpy_con + self._solver_con_to_pyomo_con_map[id(coptpy_con)] = con + + def _add_sos_constraints(self, cons: List[_SOSConstraintData]): + for con in cons: + self._symbol_map.getSymbol(con, self._labeler) + + level = con.level + if level == 1: + sos_type = coptpy.COPT.SOS_TYPE1 + elif level == 2: + sos_type = coptpy.COPT.SOS_TYPE2 + else: + raise ValueError( + "Solver does not support SOS level {0} constraints".format(level) + ) + + copt_vars = [] + weights = [] + + for v, w in con.get_items(): + v_id = id(v) + copt_vars.append(self._pyomo_var_to_solver_var_map[v_id]) + weights.append(w) + + coptpy_con = self._solver_model.addSOS(sos_type, copt_vars, weights) + self._pyomo_sos_to_solver_sos_map[con] = coptpy_con + + def _remove_constraints(self, cons: List[_GeneralConstraintData]): + for con in cons: + solver_con = self._pyomo_con_to_solver_con_map[con] + self._solver_model.remove(solver_con) + self._symbol_map.removeSymbol(con) + del self._pyomo_con_to_solver_con_map[con] + del self._solver_con_to_pyomo_con_map[id(solver_con)] + self._mutable_helpers.pop(con, None) + self._mutable_quadratic_helpers.pop(con, None) + + def _remove_sos_constraints(self, cons: List[_SOSConstraintData]): + for con in cons: + solver_sos_con = self._pyomo_sos_to_solver_sos_map[con] + self._solver_model.remove(solver_sos_con) + self._symbol_map.removeSymbol(con) + del self._pyomo_sos_to_solver_sos_map[con] + + def _remove_variables(self, variables: List[_GeneralVarData]): + for var in variables: + v_id = id(var) + solver_var = self._pyomo_var_to_solver_var_map[v_id] + self._solver_model.remove(solver_var) + self._symbol_map.removeSymbol(var) + del self._pyomo_var_to_solver_var_map[v_id] + self._mutable_bounds.pop(v_id, None) + + def _remove_params(self, params: List[_ParamData]): + pass + + def _update_variables(self, variables: List[_GeneralVarData]): + for var in variables: + var_id = id(var) + if var_id not in self._pyomo_var_to_solver_var_map: + raise ValueError( + 'The Var provided to update_var needs to be added first: {0}'.format( + var + ) + ) + self._mutable_bounds.pop((var_id, 'lb'), None) + self._mutable_bounds.pop((var_id, 'ub'), None) + coptpy_var = self._pyomo_var_to_solver_var_map[var_id] + lb, ub, vtype = self._process_domain_and_bounds( + var, var_id, None, None, None, coptpy_var + ) + coptpy_var.lb = lb + coptpy_var.ub = ub + coptpy_var.vtype = vtype + + def update_params(self): + for con, helpers in self._mutable_helpers.items(): + for helper in helpers: + helper.update() + for k, (v, helper) in self._mutable_bounds.items(): + helper.update() + + for con, helper in self._mutable_quadratic_helpers.items(): + copt_con = helper.con + new_copt_expr = helper.get_updated_expression() + new_rhs = helper.get_updated_rhs() + new_sense = copt_con.sense + pyomo_con = self._solver_con_to_pyomo_con_map[id(copt_con)] + name = self._symbol_map.getSymbol(pyomo_con, self._labeler) + self._solver_model.remove(copt_con) + new_con = self._solver_model.addQConstr( + new_copt_expr, new_sense, new_rhs, name=name + ) + self._pyomo_con_to_solver_con_map[id(pyomo_con)] = new_con + del self._solver_con_to_pyomo_con_map[id(copt_con)] + self._solver_con_to_pyomo_con_map[id(new_con)] = pyomo_con + helper.con = new_con + + helper = self._mutable_objective + pyomo_obj = self._objective + new_copt_expr = helper.get_updated_expression() + if new_copt_expr is not None: + if pyomo_obj.sense == minimize: + sense = coptpy.COPT.MINIMIZE + else: + sense = coptpy.COPT.MAXIMIZE + self._solver_model.setObjective(new_copt_expr, sense=sense) + + def _set_objective(self, obj): + if obj is None: + sense = coptpy.COPT.MINIMIZE + copt_expr = 0 + repn_constant = 0 + mutable_linear_coefficients = list() + mutable_quadratic_coefficients = list() + else: + if obj.sense == minimize: + sense = coptpy.COPT.MINIMIZE + elif obj.sense == maximize: + sense = coptpy.COPT.MAXIMIZE + else: + raise ValueError( + 'Objective sense is not recognized: {0}'.format(obj.sense) + ) + + ( + copt_expr, + repn_constant, + mutable_linear_coefficients, + mutable_quadratic_coefficients, + ) = self._get_expr_from_pyomo_expr(obj.expr) + + mutable_constant = _MutableConstant() + mutable_constant.expr = repn_constant + mutable_objective = _MutableObjective( + self._solver_model, + mutable_constant, + mutable_linear_coefficients, + mutable_quadratic_coefficients, + ) + self._mutable_objective = mutable_objective + + self._solver_model.setObjective(copt_expr + value(repn_constant), sense=sense) + + def _postsolve(self, timer: HierarchicalTimer): + config = self.config + + status = self._solver_model.status + + results = CoptResults(self) + results.wallclock_time = self._solver_model.SolvingTime + + if status == coptpy.COPT.UNSTARTED: + results.termination_condition = TerminationCondition.unknown + elif status == coptpy.COPT.OPTIMAL: + results.termination_condition = TerminationCondition.optimal + elif status == coptpy.COPT.INFEASIBLE: + results.termination_condition = TerminationCondition.infeasible + elif status == coptpy.COPT.UNBOUNDED: + results.termination_condition = TerminationCondition.unbounded + elif status == coptpy.COPT.INF_OR_UNB: + results.termination_condition = TerminationCondition.infeasibleOrUnbounded + elif status == coptpy.COPT.NUMERICAL: + results.termination_condition = TerminationCondition.error + elif status == coptpy.COPT.NODELIMIT: + results.termination_condition = TerminationCondition.maxIterations + elif status == coptpy.COPT.IMPRECISE: + results.termination_condition = TerminationCondition.optimal + elif status == coptpy.COPT.TIMEOUT: + results.termination_condition = TerminationCondition.maxTimeLimit + elif status == coptpy.COPT.INTERRUPTED: + results.termination_condition = TerminationCondition.interrupted + else: + results.termination_condition = TerminationCondition.unknown + + results.best_feasible_objective = None + results.best_objective_bound = None + if self._objective is not None: + if self._solver_model.ismip and self._solver_model.hasmipsol: + results.best_feasible_objective = self._solver_model.objval + results.best_objective_bound = self._solver_model.bestbnd + elif ( + self._solver_model.ismip == 0 + and results.termination_condition == TerminationCondition.optimal + ): + results.best_feasible_objective = self._solver_model.lpobjval + results.best_objective_bound = self._solver_model.lpobjval + if results.best_feasible_objective is not None and not math.isfinite( + results.best_feasible_objective + ): + results.best_feasible_objective = None + + timer.start('load solution') + if config.load_solution: + if self._solver_model.haslpsol or self._solver_model.hasmipsol: + if results.termination_condition != TerminationCondition.optimal: + logger.warning( + 'Loading a feasible but suboptimal solution. ' + 'Please set load_solution=False and check ' + 'results.termination_condition and ' + 'results.found_feasible_solution() before loading a solution.' + ) + self.load_vars() + else: + raise RuntimeError( + 'A feasible solution was not found, so no solution can be loaded.' + 'Please set opt.config.load_solution=False and check ' + 'results.termination_condition and ' + 'results.best_feasible_objective before loading a solution.' + ) + timer.stop('load solution') + + return results + + def _load_suboptimal_mip_solution(self, vars_to_load, solution_number): + if self.get_model_attr("IsMIP") == 0: + raise ValueError( + 'Cannot obtain suboptimal solutions for a continuous model' + ) + var_map = self._pyomo_var_to_solver_var_map + ref_vars = self._referenced_variables + copt_vars_to_load = [var_map[pyomo_var] for pyomo_var in vars_to_load] + vals = self._solver_model.getPoolSolution(solution_number, copt_vars_to_load) + res = ComponentMap() + for var_id, val in zip(vars_to_load, vals): + using_cons, using_sos, using_obj = ref_vars[var_id] + if using_cons or using_sos or (using_obj is not None): + res[self._vars[var_id][0]] = val + return res + + def load_vars(self, vars_to_load=None, solution_number=0): + for v, val in self.get_primals( + vars_to_load=vars_to_load, solution_number=solution_number + ).items(): + v.set_value(val, skip_validation=True) + StaleFlagManager.mark_all_as_stale(delayed=True) + + def get_primals(self, vars_to_load=None, solution_number=0): + if self._solver_model.haslpsol == 0 and self._solver_model.hasmipsol == 0: + raise RuntimeError( + 'Solver does not currently have a valid solution. Please ' + 'check the termination condition.' + ) + + var_map = self._pyomo_var_to_solver_var_map + ref_vars = self._referenced_variables + if vars_to_load is None: + vars_to_load = self._pyomo_var_to_solver_var_map.keys() + else: + vars_to_load = [id(v) for v in vars_to_load] + + if solution_number != 0 and self._solver_model.poolsols > 0: + return self._load_suboptimal_mip_solution( + vars_to_load=vars_to_load, solution_number=solution_number + ) + else: + copt_vars_to_load = [var_map[pyomo_var_id] for pyomo_var_id in vars_to_load] + vals = self._solver_model.getInfo("Value", copt_vars_to_load) + + res = ComponentMap() + for var_id, val in zip(vars_to_load, vals): + using_cons, using_sos, using_obj = ref_vars[var_id] + if using_cons or using_sos or (using_obj is not None): + res[self._vars[var_id][0]] = val + return res + + def get_reduced_costs(self, vars_to_load=None): + if self._solver_model.status != coptpy.COPT.OPTIMAL: + raise RuntimeError( + 'Solver does not currently have valid reduced costs. Please ' + 'check the termination condition.' + ) + + var_map = self._pyomo_var_to_solver_var_map + ref_vars = self._referenced_variables + res = ComponentMap() + if vars_to_load is None: + vars_to_load = self._pyomo_var_to_solver_var_map.keys() + else: + vars_to_load = [id(v) for v in vars_to_load] + + copt_vars_to_load = [var_map[pyomo_var_id] for pyomo_var_id in vars_to_load] + vals = self._solver_model.getInfo("RedCost", copt_vars_to_load) + + for var_id, val in zip(vars_to_load, vals): + using_cons, using_sos, using_obj = ref_vars[var_id] + if using_cons or using_sos or (using_obj is not None): + res[self._vars[var_id][0]] = val + + return res + + def get_duals(self, cons_to_load=None): + if self._solver_model.Status != coptpy.COPT.OPTIMAL: + raise RuntimeError( + 'Solver does not currently have valid duals. Please ' + 'check the termination condition.' + ) + # TODO: Cannot get duals for quadratic constraints so far with COPT + if self._solver_model.qconstrs > 0: + raise RuntimeError( + "Dual solution for quadratic constraints is not available in COPT" + ) + + con_map = self._pyomo_con_to_solver_con_map + reverse_con_map = self._solver_con_to_pyomo_con_map + dual = dict() + + if cons_to_load is None: + cons_to_load = con_map.keys() + + copt_cons_to_load = [con_map[pyomo_con_id] for pyomo_con_id in cons_to_load] + vals = list() + for copt_con in copt_cons_to_load: + vals.append(copt_con.dual) + + for copt_con, val in zip(copt_cons_to_load, vals): + pyomo_con = reverse_con_map[id(copt_con)] + dual[pyomo_con] = val + + return dual + + def get_slacks(self, cons_to_load=None): + # NOTE: Slacks of COPT are activities of constraints + if self._solver_model.haslpsol == 0 and self._solver_model.hasmipsol == 0: + raise RuntimeError( + 'Solver does not currently have valid slacks. Please ' + 'check the termination condition.' + ) + + con_map = self._pyomo_con_to_solver_con_map + reverse_con_map = self._solver_con_to_pyomo_con_map + slack = dict() + + if cons_to_load is None: + cons_to_load = con_map.keys() + + copt_cons_to_load = [con_map[pyomo_con_id] for pyomo_con_id in cons_to_load] + vals = list() + for copt_con in copt_cons_to_load: + val = 0.0 + if isinstance(copt_con, coptpy.Constraint): + if self._solver_model.ismip: + copt_row = self._solver_model.getRow(copt_con) + slack = copt_row.getValue() + else: + slack = copt_con.slack + lb = copt_con.lb + ub = copt_con.ub + # lb <= copt_con <= ub + if lb > -coptpy.COPT.INFINITY and ub < +coptpy.COPT.INFINITY: + if lb == ub: + val = ub - slack + else: + val = lb - slack if slack - lb > ub - slack else ub - slack + # copt_con >= lb + elif lb > -coptpy.COPT.INFINITY and ub >= +coptpy.COPT.INFINITY: + val = lb - slack + # copt_con <= ub + elif lb <= -coptpy.COPT.INFINITY and ub < +coptpy.COPT.INFINITY: + val = ub - slack + # free copt_con + else: + val = slack + elif isinstance(copt_con, coptpy.QConstraint): + if self._solver_model.ismip: + copt_row = self._solver_model.getQuadRow(copt_con) + slack = copt_row.getValue() + else: + slack = copt_con.slack + rhs = copt_con.rhs + val = rhs - slack + else: + raise TypeError("Unexpected constraint type for computing slack") + vals.append(val) + + for copt_con, val in zip(copt_cons_to_load, vals): + pyomo_con = reverse_con_map[id(copt_con)] + slack[pyomo_con] = val + + return slack + + def update(self, timer: HierarchicalTimer = None): + super(Copt, self).update(timer=timer) + + def get_model_attr(self, attr): + """ + Get the value of an attribute on the COPT model. + + Parameters + ---------- + attr: str + The attribute to get. See COPT documentation for descriptions of + the attributes. + """ + return getattr(self._solver_model, attr) + + def write(self, filename): + """ + Write the model to a file. + + Parameters + ---------- + filename: str + Name of the file to which the model should be written. + """ + self._solver_model.write(filename) + + def set_linear_constraint_attr(self, con, attr, val): + """ + Set the value of information on a COPT linear constraint. + + Parameters + ---------- + con: pyomo.core.base.constraint._GeneralConstraintData + The pyomo constraint for which the corresponding COPT constraint attribute + should be modified. + attr: str + The information to be modified. See the COPT documentation. + val: any + See COPT documentation for acceptable values. + """ + setattr(self._pyomo_con_to_solver_con_map[con], attr, val) + + def set_var_attr(self, var, attr, val): + """ + Set the value of information on a COPT variable. + + Parameters + ---------- + var: pyomo.core.base.var._GeneralVarData + The pyomo var for which the corresponding COPT variable information + should be modified. + attr: str + The information to be modified. See the COPT documentation. + val: any + See COPT documentation for acceptable values. + """ + setattr(self._pyomo_var_to_solver_var_map[id(var)], attr, val) + + def get_var_attr(self, var, attr): + """ + Get the value of information on a COPT variable. + + Parameters + ---------- + var: pyomo.core.base.var._GeneralVarData + The pyomo var for which the corresponding COPT variable information + should be retrieved. + attr: str + The information to get. See the COPT documentation. + """ + return getattr(self._pyomo_var_to_solver_var_map[id(var)], attr) + + def get_linear_constraint_attr(self, con, attr): + """ + Get the value of information on a COPT linear constraint. + + Parameters + ---------- + con: pyomo.core.base.constraint._GeneralConstraintData + The pyomo constraint for which the corresponding COPT constraint information + should be retrieved. + attr: str + The information to get. See the COPT documentation. + """ + return getattr(self._pyomo_con_to_solver_con_map[con], attr) + + def get_sos_attr(self, con, attr): + """ + Get the value of information on a COPT sos constraint. + + Parameters + ---------- + con: pyomo.core.base.sos._SOSConstraintData + The pyomo SOS constraint for which the corresponding COPT SOS constraint + information should be retrieved. + attr: str + The information to get. See the COPT documentation. + """ + return getattr(self._pyomo_sos_to_solver_sos_map[con], attr) + + def get_quadratic_constraint_attr(self, con, attr): + """ + Get the value of information on a COPT quadratic constraint. + + Parameters + ---------- + con: pyomo.core.base.constraint._GeneralConstraintData + The pyomo constraint for which the corresponding COPT constraint information + should be retrieved. + attr: str + The information to get. See the COPT documentation. + """ + return getattr(self._pyomo_con_to_solver_con_map[con], attr) + + def set_copt_param(self, param, val): + """ + Set a COPT parameter. + + Parameters + ---------- + param: str + The COPT parameter to set. Options include any COPT parameter. + Please see the COPT documentation for options. + val: any + The value to set the parameter to. See COPT documentation for possible values. + """ + self._solver_model.setParam(param, val) + + def get_copt_param_info(self, param): + """ + Get information about a COPT parameter. + + Parameters + ---------- + param: str + The COPT parameter to get info for. See COPT documentation for possible options. + + Returns + ------- + A 5-tuple containing the parameter name, current value, default value, + minimum value and maximum value. + """ + return self._solver_model.getParamInfo(param) + + def reset(self): + self._solver_model.reset() diff --git a/pyomo/contrib/appsi/solvers/tests/test_persistent_solvers.py b/pyomo/contrib/appsi/solvers/tests/test_persistent_solvers.py index 33f6877aaf8..fb9a6e0692f 100644 --- a/pyomo/contrib/appsi/solvers/tests/test_persistent_solvers.py +++ b/pyomo/contrib/appsi/solvers/tests/test_persistent_solvers.py @@ -6,7 +6,7 @@ parameterized = parameterized.parameterized from pyomo.contrib.appsi.base import TerminationCondition, Results, PersistentSolver from pyomo.contrib.appsi.cmodel import cmodel_available -from pyomo.contrib.appsi.solvers import Gurobi, Ipopt, Cplex, Cbc, Highs +from pyomo.contrib.appsi.solvers import Gurobi, Ipopt, Cplex, Cbc, Highs, Copt from typing import Type from pyomo.core.expr.numeric_expr import LinearExpression import os @@ -25,11 +25,18 @@ ('cplex', Cplex), ('cbc', Cbc), ('highs', Highs), + ('copt', Copt), +] +mip_solvers = [ + ('gurobi', Gurobi), + ('cplex', Cplex), + ('cbc', Cbc), + ('highs', Highs), + ('copt', Copt), ] -mip_solvers = [('gurobi', Gurobi), ('cplex', Cplex), ('cbc', Cbc), ('highs', Highs)] nlp_solvers = [('ipopt', Ipopt)] -qcp_solvers = [('gurobi', Gurobi), ('ipopt', Ipopt), ('cplex', Cplex)] -miqcqp_solvers = [('gurobi', Gurobi), ('cplex', Cplex)] +qcp_solvers = [('gurobi', Gurobi), ('ipopt', Ipopt), ('cplex', Cplex), ('copt', Copt)] +miqcqp_solvers = [('gurobi', Gurobi), ('cplex', Cplex), ('copt', Copt)] only_child_vars_options = [True, False] @@ -457,7 +464,7 @@ def test_results_infeasible( opt.config.load_solution = False res = opt.solve(m) self.assertNotEqual(res.termination_condition, TerminationCondition.optimal) - if opt_class is Ipopt: + if opt_class in (Copt, Ipopt): acceptable_termination_conditions = { TerminationCondition.infeasible, TerminationCondition.unbounded, diff --git a/pyomo/contrib/doe/examples/fim_doe_tutorial.ipynb b/pyomo/contrib/doe/examples/fim_doe_tutorial.ipynb index 12d5a610db4..36ec42fbe49 100644 --- a/pyomo/contrib/doe/examples/fim_doe_tutorial.ipynb +++ b/pyomo/contrib/doe/examples/fim_doe_tutorial.ipynb @@ -87,6 +87,7 @@ "if \"google.colab\" in sys.modules:\n", " !wget \"https://raw.githubusercontent.com/IDAES/idaes-pse/main/scripts/colab_helper.py\"\n", " import colab_helper\n", + "\n", " colab_helper.install_idaes()\n", " colab_helper.install_ipopt()\n", "\n", diff --git a/pyomo/contrib/iis/iis.py b/pyomo/contrib/iis/iis.py index bd192d04eb3..39b9692f137 100644 --- a/pyomo/contrib/iis/iis.py +++ b/pyomo/contrib/iis/iis.py @@ -1,7 +1,7 @@ """ This module contains functions for computing an irreducible infeasible set for a Pyomo MILP or LP using a specified commercial solver, one of CPLEX, -Gurobi, or Xpress. +Gurobi, Xpress or COPT. """ import abc @@ -24,8 +24,8 @@ def write_iis(pyomo_model, iis_file_name, solver=None, logger=logger): iis_file_name:str A file name to write the IIS to, e.g., infeasible_model.ilp solver:str - Specify the solver to use, one of "cplex", "gurobi", or "xpress". - If None, the tool will use the first solver available. + Specify the solver to use, one of "cplex", "gurobi", "xpress" or + "copt". If None, the tool will use the first solver available. logger:logging.Logger A logger for messages. Uses pyomo.contrib.iis logger by default. @@ -127,10 +127,20 @@ def write(self, file_name): return _remove_suffix(file_name, ".lp") + ".lp" +class CoptIIS(_IISBase): + def compute(self): + self._solver._solver_model.computeIIS() + + def write(self, file_name): + self._solver._solver_model.writeIIS(file_name) + return file_name + + _solver_map = { "cplex_persistent": CplexConflict, "gurobi_persistent": GurobiIIS, "xpress_persistent": XpressIIS, + "copt_persistent": CoptIIS, } diff --git a/pyomo/contrib/iis/tests/test_iis.py b/pyomo/contrib/iis/tests/test_iis.py index b1b675d5081..ecc09ee1f67 100644 --- a/pyomo/contrib/iis/tests/test_iis.py +++ b/pyomo/contrib/iis/tests/test_iis.py @@ -43,11 +43,19 @@ def test_write_iis_gurobi(self): def test_write_iis_xpress(self): _test_iis("xpress") + @unittest.skipUnless( + pyo.SolverFactory("copt_persistent").available(exception_flag=False), + "COPT not available", + ) + def test_write_iis_copt(self): + _test_iis("copt") + @unittest.skipUnless( ( pyo.SolverFactory("cplex_persistent").available(exception_flag=False) or pyo.SolverFactory("gurobi_persistent").available(exception_flag=False) or pyo.SolverFactory("xpress_persistent").available(exception_flag=False) + or pyo.SolverFactory("copt_persistent").available(exception_flag=False) ), "Persistent solver not available", ) @@ -75,11 +83,19 @@ def test_exception_gurobi_not_available(self): def test_exception_xpress_not_available(self): self._assert_raises_unavailable_solver("xpress") + @unittest.skipIf( + pyo.SolverFactory("copt_persistent").available(exception_flag=False), + "COPT available", + ) + def test_exception_copt_not_available(self): + self._assert_raises_unavailable_solver("copt") + @unittest.skipIf( ( pyo.SolverFactory("cplex_persistent").available(exception_flag=False) or pyo.SolverFactory("gurobi_persistent").available(exception_flag=False) or pyo.SolverFactory("xpress_persistent").available(exception_flag=False) + or pyo.SolverFactory("copt_persistent").available(exception_flag=False) ), "Persistent solver available", ) diff --git a/pyomo/contrib/mindtpy/config_options.py b/pyomo/contrib/mindtpy/config_options.py index ed0c86baae9..c1c0876b203 100644 --- a/pyomo/contrib/mindtpy/config_options.py +++ b/pyomo/contrib/mindtpy/config_options.py @@ -536,8 +536,10 @@ def _add_subsolver_configs(CONFIG): 'gams', 'gurobi_persistent', 'cplex_persistent', + 'copt_persistent', 'appsi_cplex', 'appsi_gurobi', + 'appsi_copt', # 'appsi_highs', TODO: feasibility pump now fails with appsi_highs #2951 ] ), @@ -618,8 +620,10 @@ def _add_subsolver_configs(CONFIG): 'gams', 'gurobi_persistent', 'cplex_persistent', + 'copt_persistent', 'appsi_cplex', 'appsi_gurobi', + 'appsi_copt', # 'appsi_highs', ] ), diff --git a/pyomo/solvers/plugins/solvers/COPT.py b/pyomo/solvers/plugins/solvers/COPT.py new file mode 100644 index 00000000000..303b5ea8290 --- /dev/null +++ b/pyomo/solvers/plugins/solvers/COPT.py @@ -0,0 +1,55 @@ +# ___________________________________________________________________________ +# +# Pyomo: Python Optimization Modeling Objects +# Copyright (c) 2008-2022 +# National Technology and Engineering Solutions of Sandia, LLC +# Under the terms of Contract DE-NA0003525 with National Technology and +# Engineering Solutions of Sandia, LLC, the U.S. Government retains certain +# rights in this software. +# This software is distributed under the 3-clause BSD License. +# ___________________________________________________________________________ + +import logging + +from pyomo.opt.base import OptSolver +from pyomo.opt.base.solvers import SolverFactory + +logger = logging.getLogger('pyomo.solvers') + + +@SolverFactory.register('copt', doc='The COPT solver') +class COPT(OptSolver): + """ + The COPT LP/MIP solver + """ + + def __new__(cls, *args, **kwds): + mode = kwds.pop('solver_io', 'python') + if mode is None: + mode = 'python' + if mode == 'lp' or mode == 'mps': + # COPT provides command line tool 'copt_cmd' but not integrated into Pyomo + logger.error('LP/MPS mode not supported by Pyomo interface of COPT') + return + if mode in ['python', 'direct']: + opt = SolverFactory('copt_direct', **kwds) + if opt is None: + logger.error('Python API for COPT is not installed') + return + return opt + if mode == 'persistent': + opt = SolverFactory('copt_persistent', **kwds) + if opt is None: + logger.error('Python API for COPT is not installed') + return + return opt + + if mode == 'os': + opt = SolverFactory('_ossolver', **kwds) + elif mode == 'nl': + opt = SolverFactory('asl', **kwds) + else: + logger.error('Unknown IO type: %s' % mode) + return + opt.set_options('solver=coptampl') + return opt diff --git a/pyomo/solvers/plugins/solvers/__init__.py b/pyomo/solvers/plugins/solvers/__init__.py index c5fbfa97e42..9792c62518c 100644 --- a/pyomo/solvers/plugins/solvers/__init__.py +++ b/pyomo/solvers/plugins/solvers/__init__.py @@ -30,3 +30,6 @@ import pyomo.solvers.plugins.solvers.mosek_persistent import pyomo.solvers.plugins.solvers.xpress_direct import pyomo.solvers.plugins.solvers.xpress_persistent +import pyomo.solvers.plugins.solvers.COPT +import pyomo.solvers.plugins.solvers.copt_direct +import pyomo.solvers.plugins.solvers.copt_persistent diff --git a/pyomo/solvers/plugins/solvers/copt_direct.py b/pyomo/solvers/plugins/solvers/copt_direct.py new file mode 100644 index 00000000000..09f7ee546a4 --- /dev/null +++ b/pyomo/solvers/plugins/solvers/copt_direct.py @@ -0,0 +1,826 @@ +# ___________________________________________________________________________ +# +# Pyomo: Python Optimization Modeling Objects +# Copyright (c) 2008-2022 +# National Technology and Engineering Solutions of Sandia, LLC +# Under the terms of Contract DE-NA0003525 with National Technology and +# Engineering Solutions of Sandia, LLC, the U.S. Government retains certain +# rights in this software. +# This software is distributed under the 3-clause BSD License. +# ___________________________________________________________________________ + +import logging +import re +import sys + +from pyomo.common.collections import ComponentSet, ComponentMap, Bunch +from pyomo.common.dependencies import attempt_import +from pyomo.common.errors import ApplicationError +from pyomo.common.tempfiles import TempfileManager +from pyomo.core.expr.numvalue import value, is_fixed +from pyomo.core.staleflag import StaleFlagManager +from pyomo.repn import generate_standard_repn +from pyomo.solvers.plugins.solvers.direct_solver import DirectSolver +from pyomo.solvers.plugins.solvers.direct_or_persistent_solver import ( + DirectOrPersistentSolver, +) +from pyomo.core.kernel.objective import minimize, maximize +from pyomo.opt.results.results_ import SolverResults +from pyomo.opt.results.solution import Solution, SolutionStatus +from pyomo.opt.results.solver import TerminationCondition, SolverStatus +from pyomo.opt.base import SolverFactory +from pyomo.core.base.suffix import Suffix + +logger = logging.getLogger('pyomo.solvers') + + +class DegreeError(ValueError): + pass + + +def _parse_coptpy_version(coptpy, avail): + if not avail: + return + coptpy_major = coptpy.COPT.VERSION_MAJOR + coptpy_minor = coptpy.COPT.VERSION_MINOR + coptpy_tech = coptpy.COPT.VERSION_TECHNICAL + CoptDirect._version = (coptpy_major, coptpy_minor, coptpy_tech) + CoptDirect._name = "COPT %s.%s.%s" % CoptDirect._version + while len(CoptDirect._version) < 4: + CoptDirect._version += (0,) + CoptDirect._version = CoptDirect._version[:4] + CoptDirect._version_major = CoptDirect._version[0] + + +coptpy, coptpy_available = attempt_import( + 'coptpy', catch_exceptions=(Exception,), callback=_parse_coptpy_version +) + + +@SolverFactory.register('copt_direct', doc='Direct python interface to COPT') +class CoptDirect(DirectSolver): + _name = None + _version = 0 + _version_major = 0 + _coptenv = None + + def __init__(self, **kwds): + if 'type' not in kwds: + kwds['type'] = 'copt_direct' + + super(CoptDirect, self).__init__(**kwds) + + self._python_api_exists = True + + self._pyomo_var_to_solver_var_map = ComponentMap() + self._solver_var_to_pyomo_var_map = ComponentMap() + self._pyomo_con_to_solver_con_map = dict() + self._solver_con_to_pyomo_con_map = ComponentMap() + + self._max_obj_degree = 2 + self._max_constraint_degree = 2 + + self._capabilities.linear = True + self._capabilities.quadratic_objective = True + self._capabilities.quadratic_constraint = True + self._capabilities.integer = True + self._capabilities.sos1 = True + self._capabilities.sos2 = True + + if coptpy_available and self._coptenv is None: + self._coptenv = coptpy.Envr() + self._coptmodel_name = "coptprob" + self._solver_model = None + + def available(self, exception_flag=True): + if not coptpy_available: + if exception_flag: + raise ApplicationError( + "No Python bindings available for %d solver plugin" % (type(self),) + ) + return False + else: + return True + + def _apply_solver(self): + StaleFlagManager.mark_all_as_stale() + + if not self._tee: + self._solver_model.setParam('Logging', 0) + + if self._keepfiles: + self._solver_model.setLogFile(self._log_file) + print("Solver log file: " + self._log_file) + + for key, option in self.options.items(): + if key.lower() == "writeprob": + self._solver_model.write(option) + else: + self._solver_model.setParam(key, option) + + self._solver_model.solve() + + return Bunch(rc=None, log=None) + + def _get_expr_from_pyomo_repn(self, repn, max_degree=2): + referenced_vars = ComponentSet() + + degree = repn.polynomial_degree() + if (degree is None) or (degree > max_degree): + raise DegreeError( + 'CoptDirect does not support expressions of degree {0}.'.format(degree) + ) + + if len(repn.quadratic_vars) > 0: + new_expr = coptpy.QuadExpr(0.0) + else: + new_expr = coptpy.LinExpr(0.0) + + if len(repn.linear_vars) > 0: + referenced_vars.update(repn.linear_vars) + new_expr += coptpy.LinExpr( + [self._pyomo_var_to_solver_var_map[i] for i in repn.linear_vars], + repn.linear_coefs, + ) + + for i, v in enumerate(repn.quadratic_vars): + x, y = v + new_expr.addTerm( + repn.quadratic_coefs[i], + self._pyomo_var_to_solver_var_map[x], + self._pyomo_var_to_solver_var_map[y], + ) + referenced_vars.add(x) + referenced_vars.add(y) + + new_expr += repn.constant + + return new_expr, referenced_vars + + def _get_expr_from_pyomo_expr(self, expr, max_degree=2): + if max_degree == 2: + repn = generate_standard_repn(expr, quadratic=True) + else: + repn = generate_standard_repn(expr, quadratic=False) + + try: + copt_expr, referenced_vars = self._get_expr_from_pyomo_repn( + repn, max_degree + ) + except DegreeError as e: + msg = e.args[0] + msg += '\nexpr: {0}'.format(expr) + raise DegreeError(msg) + + return copt_expr, referenced_vars + + def _copt_lb_ub_from_var(self, var): + if var.is_fixed(): + val = var.value + return val, val + if var.has_lb(): + lb = value(var.lb) + else: + lb = -coptpy.COPT.INFINITY + if var.has_ub(): + ub = value(var.ub) + else: + ub = +coptpy.COPT.INFINITY + return lb, ub + + def _copt_vtype_from_var(self, var): + if var.is_binary(): + vtype = coptpy.COPT.BINARY + elif var.is_integer(): + vtype = coptpy.COPT.INTEGER + elif var.is_continuous(): + vtype = coptpy.COPT.CONTINUOUS + else: + raise ValueError( + 'Variable domain type is not recognized for {0}'.format(var.domain) + ) + return vtype + + def _add_var(self, var): + varname = self._symbol_map.getSymbol(var, self._labeler) + vtype = self._copt_vtype_from_var(var) + lb, ub = self._copt_lb_ub_from_var(var) + + coptpy_var = self._solver_model.addVar(lb=lb, ub=ub, vtype=vtype, name=varname) + + self._pyomo_var_to_solver_var_map[var] = coptpy_var + self._solver_var_to_pyomo_var_map[coptpy_var] = var + self._referenced_variables[var] = 0 + + def _add_constraint(self, con): + if not con.active: + return None + if is_fixed(con.body): + if self._skip_trivial_constraints: + return None + + conname = self._symbol_map.getSymbol(con, self._labeler) + + if con._linear_canonical_form: + copt_expr, referenced_vars = self._get_expr_from_pyomo_repn( + con.canonical_form(), self._max_constraint_degree + ) + else: + copt_expr, referenced_vars = self._get_expr_from_pyomo_expr( + con.body, self._max_constraint_degree + ) + + if con.has_lb(): + if not is_fixed(con.lower): + raise ValueError( + "Lower bound of constraint {0} is not constant.".format(con) + ) + + if con.has_ub(): + if not is_fixed(con.upper): + raise ValueError( + "Upper bound of constraint {0} is not constant.".format(con) + ) + + if con.equality: + coptpy_con = self._solver_model.addQConstr( + copt_expr == value(con.lower), name=conname + ) + elif con.has_lb() and con.has_ub(): + coptpy_con = self._solver_model.addBoundConstr( + copt_expr, value(con.lower), value(con.upper), name=conname + ) + elif con.has_lb(): + coptpy_con = self._solver_model.addQConstr( + copt_expr >= value(con.lower), name=conname + ) + elif con.has_ub(): + coptpy_con = self._solver_model.addQConstr( + copt_expr <= value(con.upper), name=conname + ) + else: + raise ValueError( + "Constraint does not has lower/upper bound: {0} \n".format(con) + ) + + for var in referenced_vars: + self._referenced_variables[var] += 1 + self._vars_referenced_by_con[con] = referenced_vars + + self._pyomo_con_to_solver_con_map[con] = coptpy_con + self._solver_con_to_pyomo_con_map[coptpy_con] = con + + def _add_sos_constraint(self, con): + if not con.active: + return None + + self._symbol_map.getSymbol(con, self._labeler) + + level = con.level + if level == 1: + sos_type = coptpy.COPT.SOS_TYPE1 + elif level == 2: + sos_type = coptpy.COPT.SOS_TYPE2 + else: + raise ValueError( + "Solver does not support SOS level {0} constraints".format(level) + ) + + copt_vars = [] + weights = [] + + self._vars_referenced_by_con[con] = ComponentSet() + + if hasattr(con, 'get_items'): + sos_items = list(con.get_items()) + else: + sos_items = list(con.items()) + + for v, w in sos_items: + self._vars_referenced_by_con[con].add(v) + copt_vars.append(self._pyomo_var_to_solver_var_map[v]) + self._referenced_variables[v] += 1 + weights.append(w) + + coptpy_con = self._solver_model.addSOS(sos_type, copt_vars, weights) + self._pyomo_con_to_solver_con_map[con] = coptpy_con + self._solver_con_to_pyomo_con_map[coptpy_con] = con + + def _set_objective(self, obj): + if self._objective is not None: + for var in self._vars_referenced_by_obj: + self._referenced_variables[var] -= 1 + self._vars_referenced_by_obj = ComponentSet() + self._objective = None + + if obj.active is False: + raise ValueError('Cannot add inactive objective to solver.') + + if obj.sense == minimize: + sense = coptpy.COPT.MINIMIZE + elif obj.sense == maximize: + sense = coptpy.COPT.MAXIMIZE + else: + raise ValueError('Objective sense is not recognized: {0}'.format(obj.sense)) + + copt_expr, referenced_vars = self._get_expr_from_pyomo_expr( + obj.expr, self._max_obj_degree + ) + + for var in referenced_vars: + self._referenced_variables[var] += 1 + + self._solver_model.setObjective(copt_expr, sense=sense) + self._objective = obj + self._vars_referenced_by_obj = referenced_vars + + def _add_block(self, block): + DirectOrPersistentSolver._add_block(self, block) + + def _set_instance(self, model, kwds={}): + DirectOrPersistentSolver._set_instance(self, model, kwds) + + self._pyomo_con_to_solver_con_map = dict() + self._solver_con_to_pyomo_con_map = ComponentMap() + self._pyomo_var_to_solver_var_map = ComponentMap() + self._solver_var_to_pyomo_var_map = ComponentMap() + + if self._solver_model is not None: + self._solver_model.clear() + self._solver_model = None + try: + if model.name is not None: + self._solver_model = self._coptenv.createModel(model.name) + self._coptmodel_name = model.name + else: + self._solver_model = self._coptenv.createModel() + except Exception: + e = sys.exc_info()[1] + msg = ( + "Unable to create COPT model. Have you installed the Python bindings for COPT?\n\n\t" + + "Error message: {0}".format(e) + ) + raise Exception(msg) + + self._add_block(model) + + for var, n_ref in self._referenced_variables.items(): + if n_ref != 0: + if var.fixed: + if not self._output_fixed_variable_bounds: + raise ValueError( + "Encountered a fixed variable (%s) inside an active objective or constraint " + "expression on model %s, which is usually indicative of a preprocessing error." + "Use the IO-option 'output_fixed_variable_bounds=True' to suppress this error" + "and fix the variable by overwriting its bounds in the COPT instance." + % (var.name, self._pyomo_model.name) + ) + + def _postsolve(self): + extract_duals = False + extract_slacks = False + extract_reduced_costs = False + + for suffix in self._suffixes: + flag = False + if re.match(suffix, "dual"): + extract_duals = True + flag = True + if re.match(suffix, "slack"): + extract_slacks = True + flag = True + if re.match(suffix, "rc"): + extract_reduced_costs = True + flag = True + if not flag: + raise RuntimeError( + "***The copt_direct solver plugin cannot extract solution suffix=" + + suffix + ) + + if self._solver_model.ismip: + if extract_reduced_costs: + logger.warning("Cannot get reduced costs for MIP.") + extract_reduced_costs = False + if extract_duals: + logger.warning("Cannot get duals for MIP.") + extract_duals = False + else: + if self._solver_model.lpstatus != coptpy.COPT.OPTIMAL: + extract_reduced_costs = False + extract_duals = False + extract_slacks = False + + self.results = SolverResults() + soln = Solution() + + self.results.solver.name = self._name + self.results.solver.wallclock_time = self._solver_model.SolvingTime + + status = self._solver_model.status + if status == coptpy.COPT.UNSTARTED: + self.results.solver.status = SolverStatus.unknown + self.results.solver.termination_message = "Model was not solved yet." + self.results.solver.termination_condition = TerminationCondition.unknown + soln.status = SolutionStatus.unknown + elif status == coptpy.COPT.OPTIMAL: + self.results.solver.status = SolverStatus.ok + self.results.solver.termination_message = ( + "Model was solved to optimality within tolerances." + ) + self.results.solver.termination_condition = TerminationCondition.optimal + soln.status = SolutionStatus.optimal + elif status == coptpy.COPT.INFEASIBLE: + self.results.solver.status = SolverStatus.warning + self.results.solver.termination_message = ( + "Model was proven to be infeasible." + ) + self.results.solver.termination_condition = TerminationCondition.infeasible + soln.status = SolutionStatus.infeasible + elif status == coptpy.COPT.UNBOUNDED: + self.results.solver.status = SolverStatus.warning + self.results.solver.termination_message = ( + "Model was proven to be unbounded." + ) + self.results.solver.termination_condition = TerminationCondition.unbounded + soln.status = SolutionStatus.unbounded + elif status == coptpy.COPT.INF_OR_UNB: + self.results.solver.status = SolverStatus.warning + self.results.solver.termination_message = ( + "Model was proven to be infeasible or unbounded." + ) + self.results.solver.termination_condition = ( + TerminationCondition.infeasibleOrUnbounded + ) + soln.status = SolutionStatus.unsure + elif status == coptpy.COPT.NUMERICAL: + self.results.solver.status = SolverStatus.error + self.results.solver.termination_message = ( + "Optimization was terminated due to numerical difficulties." + ) + self.results.solver.termination_condition = TerminationCondition.error + soln.status = SolutionStatus.error + elif status == coptpy.COPT.NODELIMIT: + self.results.solver.status = SolverStatus.aborted + self.results.solver.termination_message = ( + "Optimization terminated because the node limit was reached." + ) + self.results.solver.termination_condition = ( + TerminationCondition.maxEvaluations + ) + soln.status = SolutionStatus.stoppedByLimit + elif status == coptpy.COPT.IMPRECISE: + self.results.solver.status = SolverStatus.warning + self.results.solver.termination_message = "Unable to satisfy optimality tolerances and returns a sub-optimal solution." + self.results.solver.termination_condition = TerminationCondition.other + soln.status = SolutionStatus.feasible + elif status == coptpy.COPT.TIMEOUT: + self.results.solver.status = SolverStatus.aborted + self.results.solver.termination_message = ( + "Optimization terminated because the time limit was reached." + ) + self.results.solver.termination_condition = ( + TerminationCondition.maxTimeLimit + ) + soln.status = SolutionStatus.stoppedByLimit + elif status == coptpy.COPT.UNFINISHED: + self.results.solver.status = SolverStatus.error + self.results.solver.termination_message = ( + "Optimization was terminated unexpectedly." + ) + self.results.solver.termination_condition = TerminationCondition.error + soln.status = SolutionStatus.error + elif status == coptpy.COPT.INTERRUPTED: + self.results.solver.status = SolverStatus.aborted + self.results.solver.termination_message = ( + "Optimization was terminated by the user." + ) + self.results.solver.termination_condition = ( + TerminationCondition.userInterrupt + ) + soln.status = SolutionStatus.stoppedByLimit + else: + self.results.solver.status = SolverStatus.error + self.results.solver.termination_message = ( + "Unknown COPT status " + "(" + str(status) + ")" + ) + self.results.solver.termination_condition = TerminationCondition.error + self.status = SolutionStatus.error + + self.results.problem.name = self._coptmodel_name + + if self._solver_model.objsense == coptpy.COPT.MINIMIZE: + self.results.problem.sense = minimize + elif self._solver_model.objsense == coptpy.COPT.MAXIMIZE: + self.results.problem.sense = maximize + else: + raise RuntimeError( + 'Unrecognized COPT objective sense: {0}'.format( + self._solver_model.objsense + ) + ) + + self.results.problem.upper_bound = None + self.results.problem.lower_bound = None + if self._solver_model.ismip == 0: + self.results.problem.upper_bound = self._solver_model.lpobjval + self.results.problem.lower_bound = self._solver_model.lpobjval + elif self._solver_model.objsense == coptpy.COPT.MINIMIZE: + self.results.problem.upper_bound = self._solver_model.objval + self.results.problem.lower_bound = self._solver_model.bestbnd + elif self._solver_model.objsense == coptpy.COPT.MAXIMIZE: + self.results.problem.upper_bound = self._solver_model.bestbnd + self.results.problem.lower_bound = self._solver_model.objval + else: + raise RuntimeError( + 'Unrecognized COPT objective sense: {0}'.format( + self._solver_model.objsense + ) + ) + + try: + soln.gap = ( + self.results.problem.upper_bound - self.results.problem.lower_bound + ) + except TypeError: + soln.gap = None + + self.results.problem.number_of_constraints = ( + self._solver_model.rows + + self._solver_model.qconstrs + + self._solver_model.soss + ) + self.results.problem.number_of_nonzeros = self._solver_model.elems + self.results.problem.number_of_variables = self._solver_model.cols + self.results.problem.number_of_binary_variables = self._solver_model.bins + self.results.problem.number_of_integer_variables = self._solver_model.ints + self.results.problem.number_of_continuous_variables = ( + self._solver_model.cols - self._solver_model.ints - self._solver_model.bins + ) + self.results.problem.number_of_objectives = 1 + self.results.problem.number_of_solutions = ( + self._solver_model.haslpsol or self._solver_model.hasmipsol + ) + + if self._save_results: + if self._solver_model.haslpsol or self._solver_model.hasmipsol: + soln_variables = soln.variable + soln_constraints = soln.constraint + + var_map = self._pyomo_var_to_solver_var_map + vars_to_load = var_map.keys() + copt_vars = [var_map[pyomo_var] for pyomo_var in vars_to_load] + + var_vals = self._solver_model.getInfo('Value', copt_vars) + names = [] + for copt_var in copt_vars: + names.append(copt_var.name) + for copt_var, val, name in zip(copt_vars, var_vals, names): + pyomo_var = self._solver_var_to_pyomo_var_map[copt_var] + if self._referenced_variables[pyomo_var] > 0: + soln_variables[name] = {'Value': val} + + if extract_reduced_costs: + vals = self._solver_model.getInfo('RedCost', copt_vars) + for copt_var, val, name in zip(copt_vars, vals, names): + pyomo_var = self._solver_var_to_pyomo_var_map[copt_var] + if self._referenced_variables[pyomo_var] > 0: + soln_variables[name]['Rc'] = val + + if extract_duals or extract_slacks: + copt_cons = self._solver_model.getConstrs() + con_names = [] + for copt_con in copt_cons: + con_names.append(copt_con.name) + for name in con_names: + soln_constraints[name] = {} + + if self._solver_model.qconstrs > 0: + copt_q_cons = self._solver_model.getQConstrs() + q_con_names = [] + for copt_q_con in copt_q_cons: + q_con_names.append(copt_q_con.name) + for name in q_con_names: + soln_constraints[name] = {} + + if extract_duals: + vals = self._solver_model.getInfo('Dual', copt_cons) + for val, name in zip(vals, con_names): + soln_constraints[name]['Dual'] = val + # TODO: Get duals for quadratic constraints + if self._solver_model.qconstrs > 0: + raise RuntimeError( + "Dual solution for quadratic constraints is not available in COPT" + ) + + if extract_slacks: + # NOTE: Slacks in COPT are activities of constraints + if self._solver_model.ismip: + l_vals = list() + for con in copt_cons: + row = self._solver_model.getRow(con) + l_vals.append(row.getValue()) + else: + l_vals = self._solver_model.getInfo('Slack', copt_cons) + l_lb = self._solver_model.getInfo("LB", copt_cons) + l_ub = self._solver_model.getInfo("UB", copt_cons) + + for val, lb, ub, name in zip(l_vals, l_lb, l_ub, con_names): + # lb <= l_con <= ub + if lb > -coptpy.COPT.INFINITY and ub < +coptpy.COPT.INFINITY: + if lb == ub: + soln_constraints[name]['Slack'] = ub - val + else: + soln_constraints[name]['Slack'] = ( + lb - val if val - lb > ub - val else ub - val + ) + # l_con >= lb + elif lb > -coptpy.COPT.INFINITY and ub >= +coptpy.COPT.INFINITY: + soln_constraints[name]['Slack'] = lb - val + # l_con <= ub + elif lb <= -coptpy.COPT.INFINITY and ub < +coptpy.COPT.INFINITY: + soln_constraints[name]['Slack'] = ub - val + # free l_con + else: + soln_constraints[name]['Slack'] = val + + if self._solver_model.qconstrs > 0: + if self._solver_model.ismip: + q_vals = list() + for q_con in copt_q_cons: + q_row = self._solver_model.getQuadRow(q_con) + q_vals.append(q_row.getValue()) + else: + q_vals = self._solver_model.getInfo('Slack', copt_q_cons) + q_rhs = list() + for q_con in copt_q_cons: + q_rhs.append(q_con.rhs) + + for val, rhs, name in zip(q_vals, q_rhs, q_con_names): + soln_constraints[name]['Slack'] = rhs - val + elif self._load_solutions: + if self._solver_model.haslpsol or self._solver_model.hasmipsol: + self.load_vars() + if extract_reduced_costs: + self._load_rc() + if extract_duals: + self._load_duals() + if extract_slacks: + self._load_slacks() + self.results.solution.insert(soln) + + TempfileManager.pop(remove=not self._keepfiles) + + return DirectOrPersistentSolver._postsolve(self) + + def warm_start_capable(self): + return True + + def _warm_start(self): + for pyomo_var, coptpy_var in self._pyomo_var_to_solver_var_map.items(): + if pyomo_var.value is not None: + self._solver_model.setMipStart(coptpy_var, value(pyomo_var)) + self._solver_model.loadMipStart() + + def _load_vars(self, vars_to_load=None): + var_map = self._pyomo_var_to_solver_var_map + ref_vars = self._referenced_variables + + if vars_to_load is None: + vars_to_load = var_map.keys() + + copt_vars_to_load = [var_map[pyomo_var] for pyomo_var in vars_to_load] + vals = self._solver_model.getInfo('Value', copt_vars_to_load) + + for var, val in zip(vars_to_load, vals): + if ref_vars[var] > 0: + var.set_value(val, skip_validation=True) + + def _load_rc(self, vars_to_load=None): + if not hasattr(self._pyomo_model, 'rc'): + self._pyomo_model.rc = Suffix(direction=Suffix.IMPORT) + + rc = self._pyomo_model.rc + var_map = self._pyomo_var_to_solver_var_map + ref_vars = self._referenced_variables + + if vars_to_load is None: + vars_to_load = var_map.keys() + + copt_vars_to_load = [var_map[pyomo_var] for pyomo_var in vars_to_load] + vals = self._solver_model.getInfo('RedCost', copt_vars_to_load) + + for var, val in zip(vars_to_load, vals): + if ref_vars[var] > 0: + rc[var] = val + + def _load_duals(self, cons_to_load=None): + # TODO: Dual solution for quadratic constraints are not available + if self._solver_model.qconstrs > 0: + raise RuntimeError( + "Dual solution for quadratic constraints is not available in COPT" + ) + if not hasattr(self._pyomo_model, 'dual'): + self._pyomo_model.dual = Suffix(direction=Suffix.IMPORT) + + dual = self._pyomo_model.dual + con_map = self._pyomo_con_to_solver_con_map + reverse_con_map = self._solver_con_to_pyomo_con_map + + if cons_to_load is None: + cons_to_load = con_map.keys() + + copt_cons_to_load = [con_map[pyomo_con] for pyomo_con in cons_to_load] + vals = list() + for copt_con in copt_cons_to_load: + vals.append(copt_con.dual) + + for copt_con, val in zip(copt_cons_to_load, vals): + pyomo_con = reverse_con_map[copt_con] + dual[pyomo_con] = val + + def _load_slacks(self, cons_to_load=None): + # NOTE: Slacks in COPT are activities of constraints + if not hasattr(self._pyomo_model, 'slack'): + self._pyomo_model.slack = Suffix(direction=Suffix.IMPORT) + + slack = self._pyomo_model.slack + con_map = self._pyomo_con_to_solver_con_map + reverse_con_map = self._solver_con_to_pyomo_con_map + + if cons_to_load is None: + cons_to_load = con_map.keys() + + copt_cons_to_load = [con_map[pyomo_con] for pyomo_con in cons_to_load] + vals = list() + for copt_con in copt_cons_to_load: + val = 0.0 + if isinstance(copt_con, coptpy.Constraint): + if self._solver_model.ismip: + copt_row = self._solver_model.getRow(copt_con) + slack = copt_row.getValue() + else: + slack = copt_con.slack + lb = copt_con.lb + ub = copt_con.ub + # lb <= copt_con <= ub + if lb > -coptpy.COPT.INFINITY and ub < +coptpy.COPT.INFINITY: + if lb == ub: + val = ub - slack + else: + val = lb - slack if slack - lb > ub - slack else ub - slack + # copt_con >= lb + elif lb > -coptpy.COPT.INFINITY and ub >= +coptpy.COPT.INFINITY: + val = lb - slack + # copt_con <= ub + elif lb <= -coptpy.COPT.INFINITY and ub < +coptpy.COPT.INFINITY: + val = ub - slack + # free copt_con + else: + val = slack + elif isinstance(copt_con, coptpy.QConstraint): + if self._solver_model.ismip: + copt_row = self._solver_model.getQuadRow(copt_con) + slack = copt_row.getValue() + else: + slack = copt_con.slack + rhs = copt_con.rhs + val = rhs - slack + else: + raise TypeError("Unexpected constraint type for computing slack") + vals.append(val) + + for copt_con, val in zip(copt_cons_to_load, vals): + pyomo_con = reverse_con_map[copt_con] + slack[pyomo_con] = val + + def load_duals(self, cons_to_load=None): + """ + Load the duals into the 'dual' suffix. The 'dual' suffix must live on the parent model. + + Parameters + ---------- + cons_to_load: list of Constraint + """ + self._load_duals(cons_to_load) + + def load_rc(self, vars_to_load): + """ + Load the reduced costs into the 'rc' suffix. The 'rc' suffix must live on the parent model. + + Parameters + ---------- + vars_to_load: list of Var + """ + self._load_rc(vars_to_load) + + def load_slacks(self, cons_to_load=None): + """ + Load the values of the slack variables into the 'slack' suffix. The 'slack' suffix must live + on the parent model. + + Parameters + ---------- + cons_to_load: list of Constraint + """ + self._load_slacks(cons_to_load) diff --git a/pyomo/solvers/plugins/solvers/copt_persistent.py b/pyomo/solvers/plugins/solvers/copt_persistent.py new file mode 100644 index 00000000000..1f991563932 --- /dev/null +++ b/pyomo/solvers/plugins/solvers/copt_persistent.py @@ -0,0 +1,252 @@ +# ___________________________________________________________________________ +# +# Pyomo: Python Optimization Modeling Objects +# Copyright (c) 2008-2022 +# National Technology and Engineering Solutions of Sandia, LLC +# Under the terms of Contract DE-NA0003525 with National Technology and +# Engineering Solutions of Sandia, LLC, the U.S. Government retains certain +# rights in this software. +# This software is distributed under the 3-clause BSD License. +# ___________________________________________________________________________ + +from pyomo.solvers.plugins.solvers.copt_direct import CoptDirect, coptpy +from pyomo.solvers.plugins.solvers.persistent_solver import PersistentSolver +from pyomo.opt.base import SolverFactory + + +@SolverFactory.register('copt_persistent', doc='Persistent python interface to COPT') +class CoptPersistent(PersistentSolver, CoptDirect): + """ + A class that provides a persistent interface to COPT. Persistent interface is similar to direct + interface excepts that it 'remember' their model, thus it allows incremental changes to the model. + Note that users are responsible for notifying the persistent solver interface when changes were + made to the corresponding pyomo model. + + Keyword Arguments + ----------------- + model: ConcreteModel + Passing a model to the constructor is equivalent to calling the set_instance method. + type: str + String indicating the class type of the solver instance. + name: str + String representing either the class type of the solver instance or an assigned name. + doc: str + Documentation for the solver + options: dict + Dictionary of solver options + """ + + def __init__(self, **kwds): + kwds['type'] = 'copt_persistent' + CoptDirect.__init__(self, **kwds) + + self._pyomo_model = kwds.pop('model', None) + if self._pyomo_model is not None: + self.set_instance(self._pyomo_model, **kwds) + + def _remove_constraint(self, solver_con): + self._solver_model.remove(solver_con) + + def _remove_sos_constraint(self, solver_sos_con): + self._solver_model.remove(solver_sos_con) + + def _remove_var(self, solver_var): + self._solver_model.remove(solver_var) + + def _warm_start(self): + CoptDirect._warm_start(self) + + def update_var(self, var): + """ + Update a single variable in the solver's model. + This will update bounds, fix/unfix the variable as needed, and update the variable type. + + Parameters + ---------- + var: Var (scalar Var or single _VarData) + + """ + if var not in self._pyomo_var_to_solver_var_map: + raise ValueError( + 'The Var provided to update_var needs to be added first: {0}'.format( + var + ) + ) + + coptpy_var = self._pyomo_var_to_solver_var_map[var] + vtype = self._copt_vtype_from_var(var) + lb, ub = self._copt_lb_ub_from_var(var) + + coptpy_var.lb = lb + coptpy_var.ub = ub + coptpy_var.vtype = vtype + + def write(self, filename): + """ + Write the model to a file. + + Parameters + ---------- + filename: str + Name of the file to which the model should be written. + """ + self._solver_model.write(filename) + + def set_linear_constraint_attr(self, con, attr, val): + """ + Set the value of information on a COPT linear constraint. + + Parameters + ---------- + con: pyomo.core.base.constraint._GeneralConstraintData + The pyomo constraint for which the corresponding COPT constraint information should be modified. + attr: str + The information to be modified. See COPT documentation for details. + val: any + See COPT documentation for acceptable values. + """ + setattr(self._pyomo_con_to_solver_con_map[con], attr, val) + + def set_var_attr(self, var, attr, val): + """ + Set the value of information on a COPT variable. + + Parameters + ---------- + con: pyomo.core.base.var._GeneralVarData + The pyomo var for which the corresponding COPT variable information should be modified. + attr: str + The information to be modified. See COPT documentation for details. + val: any + See COPT documentation for acceptable values. + """ + setattr(self._pyomo_var_to_solver_var_map[var], attr, val) + + def get_model_attr(self, attr): + """ + Get the value of an attribute on the COPT model. + + Parameters + ---------- + attr: str + The attribute to get. See COPT documentation for descriptions of the attributes. + """ + return getattr(self._solver_model, attr) + + def get_var_attr(self, var, attr): + """ + Get the value of information on a COPT variable. + + Parameters + ---------- + var: pyomo.core.base.var._GeneralVarData + The pyomo var for which the corresponding COPT variable attribute should be retrieved. + attr: str + The information to get. See COPT documentation for details. + """ + return getattr(self._pyomo_var_to_solver_var_map[var], attr) + + def get_linear_constraint_attr(self, con, attr): + """ + Get the value of information on a COPT linear constraint. + + Parameters + ---------- + con: pyomo.core.base.constraint._GeneralConstraintData + The pyomo constraint for which the corresponding COPT constraint information should be retrieved. + attr: str + The attribute to get. See COPT documentation for details. + """ + return getattr(self._pyomo_con_to_solver_con_map[con], attr) + + def get_sos_attr(self, con, attr): + """ + Get the value of information on a COPT sos constraint. + + Parameters + ---------- + con: pyomo.core.base.sos._SOSConstraintData + The pyomo SOS constraint for which the corresponding COPT SOS constraint information + should be retrieved. + attr: str + The information to get. See COPT documentation for details. + """ + return getattr(self._pyomo_con_to_solver_con_map[con], attr) + + def get_quadratic_constraint_attr(self, con, attr): + """ + Get the value of information on a COPT quadratic constraint. + + Parameters + ---------- + con: pyomo.core.base.constraint._GeneralConstraintData + The pyomo constraint for which the corresponding COPT quadratic constraint information + should be retrieved. + attr: str + The information to get. See COPT documentation for details. + """ + return getattr(self._pyomo_con_to_solver_con_map[con], attr) + + def set_copt_param(self, param, val): + """ + Set a COPT parameter. + + Parameters + ---------- + param: str + The COPT parameter to set. Options include any COPT parameter. + Please see the COPT documentation for options. + val: any + The value to set the parameter to. See COPT documentation for possible values. + """ + self._solver_model.setParam(param, val) + + def get_copt_param_info(self, param): + """ + Get information about a COPT parameter. + + Parameters + ---------- + param: str + The COPT parameter to get information for. See COPT documentation for possible options. + + Returns + ------- + A 5-tuple containing the parameter name, current value, default value, minimum value and maximum value. + """ + return self._solver_model.getParamInfo(param) + + def _add_column(self, var, obj_coef, constraints, coefficients): + """ + Add a column to the solver's model + + This will add the Pyomo variable var to the solver's model, and put the coefficients on the + associated constraints in the solver model. If the obj_coef is not zero, it will add + obj_coef*var to the objective of the solver's model. + + Parameters + ---------- + var: Var (scalar Var or single _VarData) + obj_coef: float + constraints: list of solver constraints + coefficients: list of coefficients to put on var in the associated constraint + """ + varname = self._symbol_map.getSymbol(var, self._labeler) + vtype = self._copt_vtype_from_var(var) + lb, ub = self._copt_lb_ub_from_var(var) + + coptpy_var = self._solver_model.addVar( + obj=obj_coef, + lb=lb, + ub=ub, + vtype=vtype, + name=varname, + column=coptpy.Column(constrs=constraints, coeffs=coefficients), + ) + + self._pyomo_var_to_solver_var_map[var] = coptpy_var + self._solver_var_to_pyomo_var_map[coptpy_var] = var + self._referenced_variables[var] = len(coefficients) + + def reset(self): + self._solver_model.reset() diff --git a/pyomo/solvers/tests/checks/test_copt_persistent.py b/pyomo/solvers/tests/checks/test_copt_persistent.py new file mode 100644 index 00000000000..fe47a1b6163 --- /dev/null +++ b/pyomo/solvers/tests/checks/test_copt_persistent.py @@ -0,0 +1,320 @@ +# ___________________________________________________________________________ +# +# Pyomo: Python Optimization Modeling Objects +# Copyright (c) 2008-2022 +# National Technology and Engineering Solutions of Sandia, LLC +# Under the terms of Contract DE-NA0003525 with National Technology and +# Engineering Solutions of Sandia, LLC, the U.S. Government retains certain +# rights in this software. +# This software is distributed under the 3-clause BSD License. +# ___________________________________________________________________________ + +import pyomo.common.unittest as unittest +import pyomo.environ as pyo + +coptpy_available = False +try: + import coptpy + + coptpy_available = True +except: + pass + + +class TestCoptPersistent(unittest.TestCase): + @unittest.skipIf(not coptpy_available, "coptpy is not available") + def test_basics(self): + m = pyo.ConcreteModel() + m.x = pyo.Var(bounds=(-10, 10)) + m.y = pyo.Var() + m.obj = pyo.Objective(expr=m.x**2 + m.y**2) + m.c1 = pyo.Constraint(expr=m.y >= 2 * m.x + 1) + + opt = pyo.SolverFactory('copt_persistent') + opt.set_instance(m) + + self.assertEqual(opt.get_model_attr('Cols'), 2) + self.assertEqual(opt.get_model_attr('Rows'), 1) + self.assertEqual(opt.get_model_attr('QConstrs'), 0) + self.assertEqual(opt.get_var_attr(m.x, 'LB'), -10) + self.assertEqual(opt.get_var_attr(m.x, 'UB'), 10) + + res = opt.solve() + self.assertAlmostEqual(m.x.value, -0.4) + self.assertAlmostEqual(m.y.value, 0.2) + opt.load_duals() + self.assertAlmostEqual(m.dual[m.c1], -0.4) + del m.dual + + m.c2 = pyo.Constraint(expr=m.y >= -m.x + 1) + opt.add_constraint(m.c2) + self.assertEqual(opt.get_model_attr('Cols'), 2) + self.assertEqual(opt.get_model_attr('Rows'), 2) + self.assertEqual(opt.get_model_attr('QConstrs'), 0) + + res = opt.solve(save_results=False, load_solutions=False) + self.assertAlmostEqual(m.x.value, -0.4) + self.assertAlmostEqual(m.y.value, 0.2) + opt.load_vars() + self.assertAlmostEqual(m.x.value, 0) + self.assertAlmostEqual(m.y.value, 1) + + opt.remove_constraint(m.c2) + m.del_component(m.c2) + self.assertEqual(opt.get_model_attr('Cols'), 2) + self.assertEqual(opt.get_model_attr('Rows'), 1) + self.assertEqual(opt.get_model_attr('QConstrs'), 0) + + self.assertEqual(opt.get_copt_param_info('FeasTol')[1], 1e-6) + res = opt.solve(options={'FeasTol': '1e-7'}) + self.assertEqual(opt.get_copt_param_info('FeasTol')[1], 1e-7) + self.assertAlmostEqual(m.x.value, -0.4) + self.assertAlmostEqual(m.y.value, 0.2) + + m.x.setlb(-5) + m.x.setub(5) + opt.update_var(m.x) + self.assertEqual(opt.get_var_attr(m.x, 'LB'), -5) + self.assertEqual(opt.get_var_attr(m.x, 'UB'), 5) + + m.x.fix(0) + opt.update_var(m.x) + self.assertEqual(opt.get_var_attr(m.x, 'LB'), 0) + self.assertEqual(opt.get_var_attr(m.x, 'UB'), 0) + + m.x.unfix() + opt.update_var(m.x) + self.assertEqual(opt.get_var_attr(m.x, 'LB'), -5) + self.assertEqual(opt.get_var_attr(m.x, 'UB'), 5) + + m.c2 = pyo.Constraint(expr=m.y >= m.x**2) + opt.add_constraint(m.c2) + self.assertEqual(opt.get_model_attr('Cols'), 2) + self.assertEqual(opt.get_model_attr('Rows'), 1) + self.assertEqual(opt.get_model_attr('QConstrs'), 1) + + opt.remove_constraint(m.c2) + m.del_component(m.c2) + self.assertEqual(opt.get_model_attr('Cols'), 2) + self.assertEqual(opt.get_model_attr('Rows'), 1) + self.assertEqual(opt.get_model_attr('QConstrs'), 0) + + m.z = pyo.Var() + opt.add_var(m.z) + self.assertEqual(opt.get_model_attr('Cols'), 3) + opt.remove_var(m.z) + del m.z + self.assertEqual(opt.get_model_attr('Cols'), 2) + + @unittest.skipIf(not coptpy_available, "coptpy is not available") + def test_update1(self): + m = pyo.ConcreteModel() + m.x = pyo.Var() + m.y = pyo.Var() + m.z = pyo.Var() + m.obj = pyo.Objective(expr=m.z) + m.c1 = pyo.Constraint(expr=m.z >= m.x**2 + m.y**2) + + opt = pyo.SolverFactory('copt_persistent') + opt.set_instance(m) + self.assertEqual(opt._solver_model.getAttr('QConstrs'), 1) + + opt.remove_constraint(m.c1) + self.assertEqual(opt._solver_model.getAttr('QConstrs'), 0) + + opt.add_constraint(m.c1) + self.assertEqual(opt._solver_model.getAttr('QConstrs'), 1) + + @unittest.skipIf(not coptpy_available, "coptpy is not available") + def test_update2(self): + m = pyo.ConcreteModel() + m.x = pyo.Var() + m.y = pyo.Var() + m.z = pyo.Var() + m.obj = pyo.Objective(expr=m.z) + m.c2 = pyo.Constraint(expr=m.x + m.y == 1) + + opt = pyo.SolverFactory('copt_persistent') + opt.set_instance(m) + self.assertEqual(opt._solver_model.getAttr('Rows'), 1) + + opt.remove_constraint(m.c2) + self.assertEqual(opt._solver_model.getAttr('Rows'), 0) + + opt.add_constraint(m.c2) + self.assertEqual(opt._solver_model.getAttr('Rows'), 1) + + @unittest.skipIf(not coptpy_available, "coptpy is not available") + def test_update3(self): + m = pyo.ConcreteModel() + m.x = pyo.Var() + m.y = pyo.Var() + m.z = pyo.Var() + m.obj = pyo.Objective(expr=m.z) + m.c1 = pyo.Constraint(expr=m.z >= m.x**2 + m.y**2) + + opt = pyo.SolverFactory('copt_persistent') + opt.set_instance(m) + self.assertEqual(opt._solver_model.getAttr('QConstrs'), 1) + m.c2 = pyo.Constraint(expr=m.y >= m.x**2) + opt.add_constraint(m.c2) + self.assertEqual(opt._solver_model.getAttr('QConstrs'), 2) + opt.remove_constraint(m.c2) + self.assertEqual(opt._solver_model.getAttr('QConstrs'), 1) + + @unittest.skipIf(not coptpy_available, "coptpy is not available") + def test_update4(self): + m = pyo.ConcreteModel() + m.x = pyo.Var() + m.y = pyo.Var() + m.z = pyo.Var() + m.obj = pyo.Objective(expr=m.z) + m.c1 = pyo.Constraint(expr=m.z >= m.x + m.y) + + opt = pyo.SolverFactory('copt_persistent') + opt.set_instance(m) + self.assertEqual(opt._solver_model.getAttr('Rows'), 1) + m.c2 = pyo.Constraint(expr=m.y >= m.x) + opt.add_constraint(m.c2) + self.assertEqual(opt._solver_model.getAttr('Rows'), 2) + opt.remove_constraint(m.c2) + self.assertEqual(opt._solver_model.getAttr('Rows'), 1) + + @unittest.skipIf(not coptpy_available, "coptpy is not available") + def test_update5(self): + m = pyo.ConcreteModel() + m.a = pyo.Set(initialize=[1, 2, 3], ordered=True) + m.x = pyo.Var(m.a, within=pyo.Binary) + m.y = pyo.Var(within=pyo.Binary) + m.obj = pyo.Objective(expr=m.y) + m.c1 = pyo.SOSConstraint(var=m.x, sos=1) + + opt = pyo.SolverFactory('copt_persistent') + opt.set_instance(m) + self.assertEqual(opt._solver_model.getAttr('Soss'), 1) + + opt.remove_sos_constraint(m.c1) + self.assertEqual(opt._solver_model.getAttr('Soss'), 0) + + opt.add_sos_constraint(m.c1) + self.assertEqual(opt._solver_model.getAttr('Soss'), 1) + + @unittest.skipIf(not coptpy_available, "coptpy is not available") + def test_update6(self): + m = pyo.ConcreteModel() + m.a = pyo.Set(initialize=[1, 2, 3], ordered=True) + m.x = pyo.Var(m.a, within=pyo.Binary) + m.y = pyo.Var(within=pyo.Binary) + m.obj = pyo.Objective(expr=m.y) + m.c1 = pyo.SOSConstraint(var=m.x, sos=1) + + opt = pyo.SolverFactory('copt_persistent') + opt.set_instance(m) + self.assertEqual(opt._solver_model.getAttr('Soss'), 1) + m.c2 = pyo.SOSConstraint(var=m.x, sos=2) + opt.add_sos_constraint(m.c2) + self.assertEqual(opt._solver_model.getAttr('Soss'), 2) + opt.remove_sos_constraint(m.c2) + self.assertEqual(opt._solver_model.getAttr('Soss'), 1) + + @unittest.skipIf(not coptpy_available, "coptpy is not available") + def test_update7(self): + m = pyo.ConcreteModel() + m.x = pyo.Var() + m.y = pyo.Var() + + opt = pyo.SolverFactory('copt_persistent') + opt.set_instance(m) + self.assertEqual(opt._solver_model.getAttr('Cols'), 2) + + opt.remove_var(m.x) + self.assertEqual(opt._solver_model.getAttr('Cols'), 1) + + opt.add_var(m.x) + self.assertEqual(opt._solver_model.getAttr('Cols'), 2) + + opt.remove_var(m.x) + opt.add_var(m.x) + opt.remove_var(m.x) + self.assertEqual(opt._solver_model.getAttr('Cols'), 1) + + @unittest.skipIf(not coptpy_available, "coptpy is not available") + def test_quadratic_constraint_attr(self): + m = pyo.ConcreteModel() + m.x = pyo.Var() + m.y = pyo.Var() + m.c = pyo.Constraint(expr=m.y >= m.x**2) + + opt = pyo.SolverFactory('copt_persistent') + opt.set_instance(m) + self.assertEqual(opt.get_quadratic_constraint_attr(m.c, 'rhs'), 0) + + @unittest.skipIf(not coptpy_available, "coptpy is not available") + def test_add_column(self): + m = pyo.ConcreteModel() + m.x = pyo.Var(within=pyo.NonNegativeReals) + m.c = pyo.Constraint(expr=(0, m.x, 1)) + m.obj = pyo.Objective(expr=-m.x) + + opt = pyo.SolverFactory('copt_persistent') + opt.set_instance(m) + opt.solve() + self.assertAlmostEqual(m.x.value, 1) + + m.y = pyo.Var(within=pyo.NonNegativeReals) + + opt.add_column(m, m.y, -3, [m.c], [2]) + opt.solve() + + self.assertAlmostEqual(m.x.value, 0) + self.assertAlmostEqual(m.y.value, 0.5) + + @unittest.skipIf(not coptpy_available, "coptpy is not available") + def test_add_column_exceptions(self): + m = pyo.ConcreteModel() + m.x = pyo.Var() + m.c = pyo.Constraint(expr=(0, m.x, 1)) + m.ci = pyo.Constraint([1, 2], rule=lambda m, i: (0, m.x, i + 1)) + m.cd = pyo.Constraint(expr=(0, -m.x, 1)) + m.cd.deactivate() + m.obj = pyo.Objective(expr=-m.x) + + opt = pyo.SolverFactory('copt_persistent') + + # set_instance not called + self.assertRaises(RuntimeError, opt.add_column, m, m.x, 0, [m.c], [1]) + + opt.set_instance(m) + + m2 = pyo.ConcreteModel() + m2.y = pyo.Var() + m2.c = pyo.Constraint(expr=(0, m.x, 1)) + + # different model than attached to opt + self.assertRaises(RuntimeError, opt.add_column, m2, m2.y, 0, [], []) + # pyomo var attached to different model + self.assertRaises(RuntimeError, opt.add_column, m, m2.y, 0, [], []) + + z = pyo.Var() + # pyomo var floating + self.assertRaises(RuntimeError, opt.add_column, m, z, -2, [m.c, z], [1]) + + m.y = pyo.Var() + # len(coefficients) == len(constraints) + self.assertRaises(RuntimeError, opt.add_column, m, m.y, -2, [m.c], [1, 2]) + self.assertRaises(RuntimeError, opt.add_column, m, m.y, -2, [m.c, z], [1]) + + # add indexed constraint + self.assertRaises(AttributeError, opt.add_column, m, m.y, -2, [m.ci], [1]) + # add something not a _ConstraintData + self.assertRaises(AttributeError, opt.add_column, m, m.y, -2, [m.x], [1]) + + # constraint not on solver model + self.assertRaises(KeyError, opt.add_column, m, m.y, -2, [m2.c], [1]) + + # inactive constraint + self.assertRaises(KeyError, opt.add_column, m, m.y, -2, [m.cd], [1]) + + opt.add_var(m.y) + # var already in solver model + self.assertRaises(RuntimeError, opt.add_column, m, m.y, -2, [m.c], [1]) diff --git a/pyomo/solvers/tests/solvers.py b/pyomo/solvers/tests/solvers.py index 6bbfe08c7c7..edfe2709976 100644 --- a/pyomo/solvers/tests/solvers.py +++ b/pyomo/solvers/tests/solvers.py @@ -426,6 +426,35 @@ def test_solver_cases(*args): import_suffixes=['dual'], ) + # + # COPT + # + _copt_capabilities = set( + [ + 'linear', + 'integer', + 'quadratic_objective', + 'quadratic_constraint', + 'sos1', + 'sos2', + ] + ) + + # TODO: Enable 'dual' suffix tests when duals for quad rows available too + _test_solver_cases['copt_direct', 'python'] = initialize( + name='copt_direct', + io='python', + capabilities=_copt_capabilities, + import_suffixes=['slack', 'rc'], + ) + + _test_solver_cases['copt_persistent', 'python'] = initialize( + name='copt_persistent', + io='python', + capabilities=_copt_capabilities, + import_suffixes=['slack', 'rc'], + ) + logging.disable(logging.NOTSET) #