|
| 1 | +# ___________________________________________________________________________ |
| 2 | +# |
| 3 | +# Pyomo: Python Optimization Modeling Objects |
| 4 | +# Copyright (c) 2008-2025 |
| 5 | +# National Technology and Engineering Solutions of Sandia, LLC |
| 6 | +# Under the terms of Contract DE-NA0003525 with National Technology and |
| 7 | +# Engineering Solutions of Sandia, LLC, the U.S. Government retains certain |
| 8 | +# rights in this software. |
| 9 | +# This software is distributed under the 3-clause BSD License. |
| 10 | +# ___________________________________________________________________________ |
| 11 | + |
| 12 | +import logging |
| 13 | +import re |
| 14 | +import sys |
| 15 | + |
| 16 | +from pyomo.common.collections import ComponentSet, ComponentMap, Bunch |
| 17 | +from pyomo.common.dependencies import attempt_import |
| 18 | +from pyomo.core.base import Suffix, Var, Constraint, SOSConstraint, Objective |
| 19 | +from pyomo.common.errors import ApplicationError |
| 20 | +from pyomo.common.tempfiles import TempfileManager |
| 21 | +from pyomo.common.tee import capture_output |
| 22 | +from pyomo.core.expr.numvalue import is_fixed |
| 23 | +from pyomo.core.expr.numvalue import value |
| 24 | +from pyomo.core.staleflag import StaleFlagManager |
| 25 | +from pyomo.repn import generate_standard_repn |
| 26 | +from pyomo.solvers.plugins.solvers.direct_solver import DirectSolver |
| 27 | +from pyomo.solvers.plugins.solvers.direct_or_persistent_solver import ( |
| 28 | + DirectOrPersistentSolver, |
| 29 | +) |
| 30 | +from pyomo.core.kernel.objective import minimize, maximize |
| 31 | +from pyomo.opt.results.results_ import SolverResults |
| 32 | +from pyomo.opt.results.solution import Solution, SolutionStatus |
| 33 | +from pyomo.opt.results.solver import TerminationCondition, SolverStatus |
| 34 | +from pyomo.opt.base import SolverFactory |
| 35 | +from pyomo.core.base.suffix import Suffix |
| 36 | +import numpy as np |
| 37 | +import time |
| 38 | + |
| 39 | +logger = logging.getLogger('pyomo.solvers') |
| 40 | + |
| 41 | +cuopt, cuopt_available = attempt_import( |
| 42 | + 'cuopt', |
| 43 | + ) |
| 44 | + |
| 45 | +@SolverFactory.register('cuopt_direct', doc='Direct python interface to CUOPT') |
| 46 | +class CUOPTDirect(DirectSolver): |
| 47 | + def __init__(self, **kwds): |
| 48 | + kwds['type'] = 'cuoptdirect' |
| 49 | + super(CUOPTDirect, self).__init__(**kwds) |
| 50 | + self._python_api_exists = True |
| 51 | + |
| 52 | + def _apply_solver(self): |
| 53 | + StaleFlagManager.mark_all_as_stale() |
| 54 | + log_file = None |
| 55 | + if self._log_file: |
| 56 | + log_file = self._log_file |
| 57 | + t0 = time.time() |
| 58 | + self.solution = cuopt.linear_programming.solver.Solve(self._solver_model) |
| 59 | + t1 = time.time() |
| 60 | + self._wallclock_time = t1 - t0 |
| 61 | + return Bunch(rc=None, log=None) |
| 62 | + |
| 63 | + def _add_constraint(self, constraints): |
| 64 | + c_lb, c_ub = [], [] |
| 65 | + matrix_data, matrix_indptr, matrix_indices = [], [0], [] |
| 66 | + for i, con in enumerate(constraints): |
| 67 | + repn = generate_standard_repn(con.body, quadratic=False) |
| 68 | + matrix_data.extend(repn.linear_coefs) |
| 69 | + matrix_indices.extend([self.var_name_dict[str(i)] for i in repn.linear_vars]) |
| 70 | + """for v, c in zip(con.body.linear_vars, con.body.linear_coefs): |
| 71 | + matrix_data.append(value(c)) |
| 72 | + matrix_indices.append(self.var_name_dict[str(v)])""" |
| 73 | + matrix_indptr.append(len(matrix_data)) |
| 74 | + c_lb.append(value(con.lower) if con.lower is not None else -np.inf) |
| 75 | + c_ub.append(value(con.upper) if con.upper is not None else np.inf) |
| 76 | + self._solver_model.set_csr_constraint_matrix(np.array(matrix_data), np.array(matrix_indices), np.array(matrix_indptr)) |
| 77 | + self._solver_model.set_constraint_lower_bounds(np.array(c_lb)) |
| 78 | + self._solver_model.set_constraint_upper_bounds(np.array(c_ub)) |
| 79 | + |
| 80 | + def _add_var(self, variables): |
| 81 | + # Map vriable to index and get var bounds |
| 82 | + var_type_dict = {"Integers": 'I', "Reals": 'C', "Binary": 'I'} # NonNegativeReals ? |
| 83 | + self.var_name_dict = {} |
| 84 | + v_lb, v_ub, v_type = [], [], [] |
| 85 | + |
| 86 | + for i, v in enumerate(variables): |
| 87 | + v_type.append(var_type_dict[str(v.domain)]) |
| 88 | + if v.domain == "Binary": |
| 89 | + v_lb.append(0) |
| 90 | + v_ub.append(1) |
| 91 | + else: |
| 92 | + v_lb.append(v.lb if v.lb is not None else -np.inf) |
| 93 | + v_ub.append(v.ub if v.ub is not None else np.inf) |
| 94 | + self.var_name_dict[str(v)] = i |
| 95 | + self._pyomo_var_to_ndx_map[v] = self._ndx_count |
| 96 | + self._ndx_count += 1 |
| 97 | + |
| 98 | + self._solver_model.set_variable_lower_bounds(np.array(v_lb)) |
| 99 | + self._solver_model.set_variable_upper_bounds(np.array(v_ub)) |
| 100 | + self._solver_model.set_variable_types(np.array(v_type)) |
| 101 | + self._solver_model.set_variable_names(np.array(list(self.var_name_dict.keys()))) |
| 102 | + |
| 103 | + def _set_objective(self, objective): |
| 104 | + repn = generate_standard_repn(objective.expr, quadratic=False) |
| 105 | + obj_coeffs = [0] * len(self.var_name_dict) |
| 106 | + for i, coeff in enumerate(repn.linear_coefs): |
| 107 | + obj_coeffs[self.var_name_dict[str(repn.linear_vars[i])]] = coeff |
| 108 | + self._solver_model.set_objective_coefficients(np.array(obj_coeffs)) |
| 109 | + if objective.sense == maximize: |
| 110 | + self._solver_model.set_maximize(True) |
| 111 | + |
| 112 | + def _set_instance(self, model, kwds={}): |
| 113 | + DirectOrPersistentSolver._set_instance(self, model, kwds) |
| 114 | + self.var_name_dict = None |
| 115 | + self._pyomo_var_to_ndx_map = ComponentMap() |
| 116 | + self._ndx_count = 0 |
| 117 | + |
| 118 | + try: |
| 119 | + self._solver_model = cuopt.linear_programming.DataModel() |
| 120 | + except Exception: |
| 121 | + e = sys.exc_info()[1] |
| 122 | + msg = ( |
| 123 | + "Unable to create CUOPT model. " |
| 124 | + "Have you installed the Python " |
| 125 | + "SDK for CUOPT?\n\n\t" + "Error message: {0}".format(e) |
| 126 | + ) |
| 127 | + self._add_block(model) |
| 128 | + |
| 129 | + def _add_block(self, block): |
| 130 | + self._add_var(block.component_data_objects( |
| 131 | + ctype=Var, descend_into=True, active=True, sort=True) |
| 132 | + ) |
| 133 | + |
| 134 | + for sub_block in block.block_data_objects(descend_into=True, active=True): |
| 135 | + self._add_constraint(sub_block.component_data_objects( |
| 136 | + ctype=Constraint, descend_into=False, active=True, sort=True) |
| 137 | + ) |
| 138 | + obj_counter = 0 |
| 139 | + for obj in sub_block.component_data_objects( |
| 140 | + ctype=Objective, descend_into=False, active=True |
| 141 | + ): |
| 142 | + obj_counter += 1 |
| 143 | + if obj_counter > 1: |
| 144 | + raise ValueError( |
| 145 | + "Solver interface does not support multiple objectives." |
| 146 | + ) |
| 147 | + self._set_objective(obj) |
| 148 | + |
| 149 | + def _postsolve(self): |
| 150 | + extract_duals = False |
| 151 | + extract_slacks = False |
| 152 | + extract_reduced_costs = False |
| 153 | + for suffix in self._suffixes: |
| 154 | + flag = False |
| 155 | + if re.match(suffix, "rc"): |
| 156 | + extract_reduced_costs = True |
| 157 | + flag = True |
| 158 | + if not flag: |
| 159 | + raise RuntimeError( |
| 160 | + "***The cuopt_direct solver plugin cannot extract solution suffix=" |
| 161 | + + suffix |
| 162 | + ) |
| 163 | + |
| 164 | + solution = self.solution |
| 165 | + status = solution.get_termination_status() |
| 166 | + self.results = SolverResults() |
| 167 | + soln = Solution() |
| 168 | + self.results.solver.name = "CUOPT" |
| 169 | + self.results.solver.wallclock_time = self._wallclock_time |
| 170 | + |
| 171 | + prob_type = solution.problem_category |
| 172 | + |
| 173 | + if status in [1]: |
| 174 | + self.results.solver.status = SolverStatus.ok |
| 175 | + self.results.solver.termination_condition = TerminationCondition.optimal |
| 176 | + soln.status = SolutionStatus.optimal |
| 177 | + elif status in [3]: |
| 178 | + self.results.solver.status = SolverStatus.warning |
| 179 | + self.results.solver.termination_condition = TerminationCondition.unbounded |
| 180 | + soln.status = SolutionStatus.unbounded |
| 181 | + elif status in [8]: |
| 182 | + self.results.solver.status = SolverStatus.ok |
| 183 | + self.results.solver.termination_condition = TerminationCondition.feasible |
| 184 | + soln.status = SolutionStatus.feasible |
| 185 | + elif status in [2]: |
| 186 | + self.results.solver.status = SolverStatus.warning |
| 187 | + self.results.solver.termination_condition = TerminationCondition.infeasible |
| 188 | + soln.status = SolutionStatus.infeasible |
| 189 | + elif status in [4]: |
| 190 | + self.results.solver.status = SolverStatus.aborted |
| 191 | + self.results.solver.termination_condition = ( |
| 192 | + TerminationCondition.maxIterations |
| 193 | + ) |
| 194 | + soln.status = SolutionStatus.stoppedByLimit |
| 195 | + elif status in [5]: |
| 196 | + self.results.solver.status = SolverStatus.aborted |
| 197 | + self.results.solver.termination_condition = ( |
| 198 | + TerminationCondition.maxTimeLimit |
| 199 | + ) |
| 200 | + soln.status = SolutionStatus.stoppedByLimit |
| 201 | + elif status in [7]: |
| 202 | + self.results.solver.status = SolverStatus.ok |
| 203 | + self.results.solver.termination_condition = ( |
| 204 | + TerminationCondition.other |
| 205 | + ) |
| 206 | + soln.status = SolutionStatus.other |
| 207 | + else: |
| 208 | + self.results.solver.status = SolverStatus.error |
| 209 | + self.results.solver.termination_condition = TerminationCondition.error |
| 210 | + soln.status = SolutionStatus.error |
| 211 | + |
| 212 | + if self._solver_model.maximize: |
| 213 | + self.results.problem.sense = maximize |
| 214 | + else: |
| 215 | + self.results.problem.sense = minimize |
| 216 | + |
| 217 | + self.results.problem.upper_bound = None |
| 218 | + self.results.problem.lower_bound = None |
| 219 | + try: |
| 220 | + self.results.problem.upper_bound = solution.get_primal_objective() |
| 221 | + self.results.problem.lower_bound = solution.get_primal_objective() |
| 222 | + except Exception as e: |
| 223 | + pass |
| 224 | + |
| 225 | + var_map = self._pyomo_var_to_ndx_map |
| 226 | + primal_solution = solution.get_primal_solution().tolist() |
| 227 | + for i, pyomo_var in enumerate(var_map.keys()): |
| 228 | + pyomo_var.set_value(primal_solution[i], skip_validation=True) |
| 229 | + |
| 230 | + if extract_reduced_costs: |
| 231 | + self._load_rc() |
| 232 | + |
| 233 | + self.results.solution.insert(soln) |
| 234 | + return DirectOrPersistentSolver._postsolve(self) |
| 235 | + |
| 236 | + def warm_start_capable(self): |
| 237 | + return False |
| 238 | + |
| 239 | + def _load_rc(self, vars_to_load=None): |
| 240 | + if not hasattr(self._pyomo_model, 'rc'): |
| 241 | + self._pyomo_model.rc = Suffix(direction=Suffix.IMPORT) |
| 242 | + rc = self._pyomo_model.rc |
| 243 | + var_map = self._pyomo_var_to_ndx_map |
| 244 | + if vars_to_load is None: |
| 245 | + vars_to_load = var_map.keys() |
| 246 | + reduced_costs = self.solution.get_reduced_costs() |
| 247 | + for pyomo_var in vars_to_load: |
| 248 | + rc[pyomo_var] = reduced_costs[var_map[pyomo_var]] |
| 249 | + |
| 250 | + def load_rc(self, vars_to_load): |
| 251 | + """ |
| 252 | + Load the reduced costs into the 'rc' suffix. The 'rc' suffix must live on the parent model. |
| 253 | +
|
| 254 | + Parameters |
| 255 | + ---------- |
| 256 | + vars_to_load: list of Var |
| 257 | + """ |
| 258 | + self._load_rc(vars_to_load) |
0 commit comments