always-cautious/ortools-clone
1
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

fix model_builder

Browse Source
This commit is contained in:
Laurent Perron
2022-09-12 13:21:53 +02:00
parent 81b31ff24b
commit 9c8c2ba602
16 changed files with 99 additions and 2468 deletions
Download Patch File Download Diff File Expand all files Collapse all files

34
ortools/linear_solver/BUILD.bazel
View File

@@ -79,8 +79,7 @@ cc_library(
"linear_solver_callback.cc",
"linear_solver.cc",
"lpi_glop.cpp",
"model_validator.cc",
"pdlp_proto_solver.cc",
"pdlp_interface.cc",
"sat_interface.cc",
"scip_callback.cc",
"scip_interface.cc",
@@ -103,17 +102,11 @@ cc_library(
hdrs = [
"glop_interface.cc",
"glop_utils.h",
"gurobi_proto_solver.h",
"linear_expr.h",
"linear_solver.h",
"linear_solver_callback.h",
"model_validator.h",
"pdlp_proto_solver.h",
"sat_proto_solver.h",
"sat_solver_utils.h",
"scip_callback.h",
"scip_helper_macros.h",
"scip_proto_solver.h",
],
copts = [
"-DUSE_PDLP",
@@ -123,6 +116,7 @@ cc_library(
":linear_solver_cc_proto",
":scip_with_glop",
":model_exporter",
":model_validator",
"@com_google_absl//absl/status",
"@com_google_absl//absl/status:statusor",
"@com_google_absl//absl/strings",
@@ -159,6 +153,30 @@ cc_library(
}),
)
cc_library(
name = "model_validator",
srcs = ["model_validator.cc"],
hdrs = ["model_validator.h"],
visibility = ["//visibility:public"],
deps = [
":linear_solver_cc_proto",
"//ortools/base",
"//ortools/base:accurate_sum",
"//ortools/base:map_util",
"//ortools/port:file",
"//ortools/port:proto_utils",
"//ortools/util:fp_utils",
"//ortools/util:lazy_mutable_copy",
"@com_google_absl//absl/container:flat_hash_map",
"@com_google_absl//absl/container:flat_hash_set",
"@com_google_absl//absl/status",
"@com_google_absl//absl/status:statusor",
"@com_google_absl//absl/strings",
"@com_google_absl//absl/strings:str_format",
"@com_google_absl//absl/types:optional",
],
)
copy_file(
name = "lpi_glop",
src = "@scip//:src/lpi/lpi_glop.cpp",

591
ortools/linear_solver/gurobi_proto_solver.cc
View File

@@ -1,591 +0,0 @@
// Copyright 2010-2022 Google LLC
// Licensed under the Apache License, Version 2.0 (the "License");
// you may not use this file except in compliance with the License.
// You may obtain a copy of the License at
//
// http://www.apache.org/licenses/LICENSE-2.0
//
// Unless required by applicable law or agreed to in writing, software
// distributed under the License is distributed on an "AS IS" BASIS,
// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
// See the License for the specific language governing permissions and
// limitations under the License.
#include "ortools/linear_solver/gurobi_proto_solver.h"
#include <algorithm>
#include <cmath>
#include <limits>
#include <memory>
#include <numeric>
#include <string>
#include <vector>
#include "absl/status/status.h"
#include "absl/status/statusor.h"
#include "absl/strings/str_cat.h"
#include "absl/strings/str_format.h"
#include "absl/strings/str_join.h"
#include "absl/strings/str_split.h"
#include "absl/types/optional.h"
#include "ortools/base/cleanup.h"
#include "ortools/base/status_macros.h"
#include "ortools/base/timer.h"
#include "ortools/gurobi/environment.h"
#include "ortools/linear_solver/linear_solver.h"
#include "ortools/linear_solver/linear_solver.pb.h"
#include "ortools/linear_solver/model_validator.h"
#include "ortools/util/lazy_mutable_copy.h"
namespace operations_research {
namespace {
constexpr int GRB_OK = 0;
inline absl::Status GurobiCodeToUtilStatus(int error_code,
const char* source_file,
int source_line,
const char* statement,
GRBenv* const env) {
if (error_code == GRB_OK) return absl::OkStatus();
return absl::InvalidArgumentError(absl::StrFormat(
"Gurobi error code %d (file '%s', line %d) on '%s': %s", error_code,
source_file, source_line, statement, GRBgeterrormsg(env)));
}
int AddIndicatorConstraint(const MPGeneralConstraintProto& gen_cst,
GRBmodel* gurobi_model,
std::vector<int>* tmp_variables,
std::vector<double>* tmp_coefficients) {
CHECK(gurobi_model != nullptr);
CHECK(tmp_variables != nullptr);
CHECK(tmp_coefficients != nullptr);
const auto& ind_cst = gen_cst.indicator_constraint();
MPConstraintProto cst = ind_cst.constraint();
if (cst.lower_bound() > -std::numeric_limits<double>::infinity()) {
int status = GRBaddgenconstrIndicator(
gurobi_model, gen_cst.name().c_str(), ind_cst.var_index(),
ind_cst.var_value(), cst.var_index_size(),
cst.mutable_var_index()->mutable_data(),
cst.mutable_coefficient()->mutable_data(),
cst.upper_bound() == cst.lower_bound() ? GRB_EQUAL : GRB_GREATER_EQUAL,
cst.lower_bound());
if (status != GRB_OK) return status;
}
if (cst.upper_bound() < std::numeric_limits<double>::infinity() &&
cst.lower_bound() != cst.upper_bound()) {
return GRBaddgenconstrIndicator(gurobi_model, gen_cst.name().c_str(),
ind_cst.var_index(), ind_cst.var_value(),
cst.var_index_size(),
cst.mutable_var_index()->mutable_data(),
cst.mutable_coefficient()->mutable_data(),
GRB_LESS_EQUAL, cst.upper_bound());
}
return GRB_OK;
}
int AddSosConstraint(const MPSosConstraint& sos_cst, GRBmodel* gurobi_model,
std::vector<int>* tmp_variables,
std::vector<double>* tmp_weights) {
CHECK(gurobi_model != nullptr);
CHECK(tmp_variables != nullptr);
CHECK(tmp_weights != nullptr);
tmp_variables->resize(sos_cst.var_index_size(), 0);
for (int v = 0; v < sos_cst.var_index_size(); ++v) {
(*tmp_variables)[v] = sos_cst.var_index(v);
}
tmp_weights->resize(sos_cst.var_index_size(), 0);
if (sos_cst.weight_size() == sos_cst.var_index_size()) {
for (int w = 0; w < sos_cst.weight_size(); ++w) {
(*tmp_weights)[w] = sos_cst.weight(w);
}
} else {
DCHECK_EQ(sos_cst.weight_size(), 0);
// Gurobi requires variable weights in their SOS constraints.
std::iota(tmp_weights->begin(), tmp_weights->end(), 1);
}
std::vector<int> types = {sos_cst.type() == MPSosConstraint::SOS1_DEFAULT
? GRB_SOS_TYPE1
: GRB_SOS_TYPE2};
std::vector<int> begins = {0};
return GRBaddsos(gurobi_model, /*numsos=*/1,
/*nummembers=*/sos_cst.var_index_size(),
/*types=*/types.data(),
/*beg=*/begins.data(), /*ind=*/tmp_variables->data(),
/*weight*/ tmp_weights->data());
}
int AddQuadraticConstraint(const MPGeneralConstraintProto& gen_cst,
GRBmodel* gurobi_model) {
CHECK(gurobi_model != nullptr);
constexpr double kInfinity = std::numeric_limits<double>::infinity();
CHECK(gen_cst.has_quadratic_constraint());
const MPQuadraticConstraint& quad_cst = gen_cst.quadratic_constraint();
auto addqconstr = [](GRBmodel* gurobi_model, MPQuadraticConstraint quad_cst,
char sense, double rhs, const std::string& name) {
return GRBaddqconstr(
gurobi_model,
/*numlnz=*/quad_cst.var_index_size(),
/*lind=*/quad_cst.mutable_var_index()->mutable_data(),
/*lval=*/quad_cst.mutable_coefficient()->mutable_data(),
/*numqnz=*/quad_cst.qvar1_index_size(),
/*qrow=*/quad_cst.mutable_qvar1_index()->mutable_data(),
/*qcol=*/quad_cst.mutable_qvar2_index()->mutable_data(),
/*qval=*/quad_cst.mutable_qcoefficient()->mutable_data(),
/*sense=*/sense,
/*rhs=*/rhs,
/*QCname=*/name.c_str());
};
if (quad_cst.has_lower_bound() && quad_cst.lower_bound() > -kInfinity) {
const int grb_status =
addqconstr(gurobi_model, gen_cst.quadratic_constraint(),
GRB_GREATER_EQUAL, quad_cst.lower_bound(),
gen_cst.has_name() ? gen_cst.name() + "_lb" : "");
if (grb_status != GRB_OK) return grb_status;
}
if (quad_cst.has_upper_bound() && quad_cst.upper_bound() < kInfinity) {
const int grb_status =
addqconstr(gurobi_model, gen_cst.quadratic_constraint(), GRB_LESS_EQUAL,
quad_cst.upper_bound(),
gen_cst.has_name() ? gen_cst.name() + "_ub" : "");
if (grb_status != GRB_OK) return grb_status;
}
return GRB_OK;
}
int AddAndConstraint(const MPGeneralConstraintProto& gen_cst,
GRBmodel* gurobi_model, std::vector<int>* tmp_variables) {
CHECK(gurobi_model != nullptr);
CHECK(tmp_variables != nullptr);
auto and_cst = gen_cst.and_constraint();
return GRBaddgenconstrAnd(
gurobi_model,
/*name=*/gen_cst.name().c_str(),
/*resvar=*/and_cst.resultant_var_index(),
/*nvars=*/and_cst.var_index_size(),
/*vars=*/and_cst.mutable_var_index()->mutable_data());
}
int AddOrConstraint(const MPGeneralConstraintProto& gen_cst,
GRBmodel* gurobi_model, std::vector<int>* tmp_variables) {
CHECK(gurobi_model != nullptr);
CHECK(tmp_variables != nullptr);
auto or_cst = gen_cst.or_constraint();
return GRBaddgenconstrOr(gurobi_model,
/*name=*/gen_cst.name().c_str(),
/*resvar=*/or_cst.resultant_var_index(),
/*nvars=*/or_cst.var_index_size(),
/*vars=*/or_cst.mutable_var_index()->mutable_data());
}
int AddMinConstraint(const MPGeneralConstraintProto& gen_cst,
GRBmodel* gurobi_model, std::vector<int>* tmp_variables) {
CHECK(gurobi_model != nullptr);
CHECK(tmp_variables != nullptr);
auto min_cst = gen_cst.min_constraint();
return GRBaddgenconstrMin(
gurobi_model,
/*name=*/gen_cst.name().c_str(),
/*resvar=*/min_cst.resultant_var_index(),
/*nvars=*/min_cst.var_index_size(),
/*vars=*/min_cst.mutable_var_index()->mutable_data(),
/*constant=*/min_cst.has_constant()
? min_cst.constant()
: std::numeric_limits<double>::infinity());
}
int AddMaxConstraint(const MPGeneralConstraintProto& gen_cst,
GRBmodel* gurobi_model, std::vector<int>* tmp_variables) {
CHECK(gurobi_model != nullptr);
CHECK(tmp_variables != nullptr);
auto max_cst = gen_cst.max_constraint();
return GRBaddgenconstrMax(
gurobi_model,
/*name=*/gen_cst.name().c_str(),
/*resvar=*/max_cst.resultant_var_index(),
/*nvars=*/max_cst.var_index_size(),
/*vars=*/max_cst.mutable_var_index()->mutable_data(),
/*constant=*/max_cst.has_constant()
? max_cst.constant()
: -std::numeric_limits<double>::infinity());
}
} // namespace
absl::Status SetSolverSpecificParameters(const std::string& parameters,
GRBenv* gurobi) {
if (parameters.empty()) return absl::OkStatus();
std::vector<std::string> error_messages;
for (absl::string_view line : absl::StrSplit(parameters, '\n')) {
// Comment tokens end at the next new-line, or the end of the string.
// The first character must be '#'
if (line[0] == '#') continue;
for (absl::string_view token :
absl::StrSplit(line, ',', absl::SkipWhitespace())) {
if (token.empty()) continue;
std::vector<std::string> key_value =
absl::StrSplit(token, absl::ByAnyChar(" ="), absl::SkipWhitespace());
// If one parameter fails, we keep processing the list of parameters.
if (key_value.size() != 2) {
const std::string current_message =
absl::StrCat("Cannot parse parameter '", token,
"'. Expected format is 'ParameterName value' or "
"'ParameterName=value'");
error_messages.push_back(current_message);
continue;
}
const int gurobi_code =
GRBsetparam(gurobi, key_value[0].c_str(), key_value[1].c_str());
if (gurobi_code != GRB_OK) {
const std::string current_message = absl::StrCat(
"Error setting parameter '", key_value[0], "' to value '",
key_value[1], "': ", GRBgeterrormsg(gurobi));
error_messages.push_back(current_message);
continue;
}
VLOG(2) << absl::StrCat("Set parameter '", key_value[0], "' to value '",
key_value[1]);
}
}
if (error_messages.empty()) return absl::OkStatus();
return absl::InvalidArgumentError(absl::StrJoin(error_messages, "\n"));
}
absl::StatusOr<MPSolutionResponse> GurobiSolveProto(
const MPModelRequest& request, GRBenv* gurobi_env) {
MPSolutionResponse response;
const absl::optional<LazyMutableCopy<MPModelProto>> optional_model =
ExtractValidMPModelOrPopulateResponseStatus(request, &response);
if (!optional_model) return response;
const MPModelProto& model = optional_model->get();
// We set `gurobi_env` to point to a new environment if no existing one is
// provided. We must make sure that we free this environment when we exit this
// function.
bool gurobi_env_was_created = false;
auto gurobi_env_deleter = absl::MakeCleanup([&]() {
if (gurobi_env_was_created && gurobi_env != nullptr) {
GRBfreeenv(gurobi_env);
}
});
if (gurobi_env == nullptr) {
ASSIGN_OR_RETURN(gurobi_env, GetGurobiEnv());
gurobi_env_was_created = true;
}
GRBmodel* gurobi_model = nullptr;
auto gurobi_model_deleter = absl::MakeCleanup([&]() {
const int error_code = GRBfreemodel(gurobi_model);
LOG_IF(DFATAL, error_code != GRB_OK)
<< "GRBfreemodel failed with error " << error_code << ": "
<< GRBgeterrormsg(gurobi_env);
});
// `gurobi_env` references ther GRBenv argument.
#define RETURN_IF_GUROBI_ERROR(x) \
RETURN_IF_ERROR( \
GurobiCodeToUtilStatus(x, __FILE__, __LINE__, #x, gurobi_env));
RETURN_IF_GUROBI_ERROR(GRBnewmodel(gurobi_env, &gurobi_model,
model.name().c_str(),
/*numvars=*/0,
/*obj=*/nullptr,
/*lb=*/nullptr,
/*ub=*/nullptr,
/*vtype=*/nullptr,
/*varnames=*/nullptr));
GRBenv* const model_env = GRBgetenv(gurobi_model);
if (request.has_solver_specific_parameters()) {
const auto parameters_status = SetSolverSpecificParameters(
request.solver_specific_parameters(), model_env);
if (!parameters_status.ok()) {
response.set_status(MPSOLVER_MODEL_INVALID_SOLVER_PARAMETERS);
response.set_status_str(
std::string(parameters_status.message())); // NOLINT
return response;
}
}
if (request.solver_time_limit_seconds() > 0) {
RETURN_IF_GUROBI_ERROR(GRBsetdblparam(model_env, GRB_DBL_PAR_TIMELIMIT,
request.solver_time_limit_seconds()));
}
RETURN_IF_GUROBI_ERROR(
GRBsetintparam(model_env, GRB_INT_PAR_OUTPUTFLAG,
request.enable_internal_solver_output()));
const int variable_size = model.variable_size();
bool has_integer_variables = false;
{
std::vector<double> obj_coeffs(variable_size, 0);
std::vector<double> lb(variable_size);
std::vector<double> ub(variable_size);
std::vector<char> ctype(variable_size);
std::vector<const char*> varnames(variable_size);
for (int v = 0; v < variable_size; ++v) {
const MPVariableProto& variable = model.variable(v);
obj_coeffs[v] = variable.objective_coefficient();
lb[v] = variable.lower_bound();
ub[v] = variable.upper_bound();
ctype[v] = variable.is_integer() && SolverTypeIsMip(request.solver_type())
? GRB_INTEGER
: GRB_CONTINUOUS;
if (variable.is_integer()) has_integer_variables = true;
if (!variable.name().empty()) varnames[v] = variable.name().c_str();
}
RETURN_IF_GUROBI_ERROR(
GRBaddvars(gurobi_model, variable_size, 0, nullptr, nullptr, nullptr,
/*obj=*/obj_coeffs.data(),
/*lb=*/lb.data(), /*ub=*/ub.data(), /*vtype=*/ctype.data(),
/*varnames=*/const_cast<char**>(varnames.data())));
// Set solution hints if any.
for (int i = 0; i < model.solution_hint().var_index_size(); ++i) {
RETURN_IF_GUROBI_ERROR(GRBsetdblattrelement(
gurobi_model, GRB_DBL_ATTR_START, model.solution_hint().var_index(i),
model.solution_hint().var_value(i)));
}
}
{
std::vector<int> ct_variables;
std::vector<double> ct_coefficients;
for (int c = 0; c < model.constraint_size(); ++c) {
const MPConstraintProto& constraint = model.constraint(c);
const int size = constraint.var_index_size();
ct_variables.resize(size, 0);
ct_coefficients.resize(size, 0);
for (int i = 0; i < size; ++i) {
ct_variables[i] = constraint.var_index(i);
ct_coefficients[i] = constraint.coefficient(i);
}
// Using GRBaddrangeconstr for constraints that don't require it adds
// a slack which is not always removed by presolve.
if (constraint.lower_bound() == constraint.upper_bound()) {
RETURN_IF_GUROBI_ERROR(GRBaddconstr(
gurobi_model, /*numnz=*/size, /*cind=*/ct_variables.data(),
/*cval=*/ct_coefficients.data(),
/*sense=*/GRB_EQUAL, /*rhs=*/constraint.lower_bound(),
/*constrname=*/constraint.name().c_str()));
} else if (constraint.lower_bound() ==
-std::numeric_limits<double>::infinity()) {
RETURN_IF_GUROBI_ERROR(GRBaddconstr(
gurobi_model, /*numnz=*/size, /*cind=*/ct_variables.data(),
/*cval=*/ct_coefficients.data(),
/*sense=*/GRB_LESS_EQUAL, /*rhs=*/constraint.upper_bound(),
/*constrname=*/constraint.name().c_str()));
} else if (constraint.upper_bound() ==
std::numeric_limits<double>::infinity()) {
RETURN_IF_GUROBI_ERROR(GRBaddconstr(
gurobi_model, /*numnz=*/size, /*cind=*/ct_variables.data(),
/*cval=*/ct_coefficients.data(),
/*sense=*/GRB_GREATER_EQUAL, /*rhs=*/constraint.lower_bound(),
/*constrname=*/constraint.name().c_str()));
} else {
RETURN_IF_GUROBI_ERROR(GRBaddrangeconstr(
gurobi_model, /*numnz=*/size, /*cind=*/ct_variables.data(),
/*cval=*/ct_coefficients.data(),
/*lower=*/constraint.lower_bound(),
/*upper=*/constraint.upper_bound(),
/*constrname=*/constraint.name().c_str()));
}
}
for (const auto& gen_cst : model.general_constraint()) {
switch (gen_cst.general_constraint_case()) {
case MPGeneralConstraintProto::kIndicatorConstraint: {
RETURN_IF_GUROBI_ERROR(AddIndicatorConstraint(
gen_cst, gurobi_model, &ct_variables, &ct_coefficients));
break;
}
case MPGeneralConstraintProto::kSosConstraint: {
RETURN_IF_GUROBI_ERROR(AddSosConstraint(gen_cst.sos_constraint(),
gurobi_model, &ct_variables,
&ct_coefficients));
break;
}
case MPGeneralConstraintProto::kQuadraticConstraint: {
RETURN_IF_GUROBI_ERROR(AddQuadraticConstraint(gen_cst, gurobi_model));
break;
}
case MPGeneralConstraintProto::kAbsConstraint: {
RETURN_IF_GUROBI_ERROR(GRBaddgenconstrAbs(
gurobi_model,
/*name=*/gen_cst.name().c_str(),
/*resvar=*/gen_cst.abs_constraint().resultant_var_index(),
/*argvar=*/gen_cst.abs_constraint().var_index()));
break;
}
case MPGeneralConstraintProto::kAndConstraint: {
RETURN_IF_GUROBI_ERROR(
AddAndConstraint(gen_cst, gurobi_model, &ct_variables));
break;
}
case MPGeneralConstraintProto::kOrConstraint: {
RETURN_IF_GUROBI_ERROR(
AddOrConstraint(gen_cst, gurobi_model, &ct_variables));
break;
}
case MPGeneralConstraintProto::kMinConstraint: {
RETURN_IF_GUROBI_ERROR(
AddMinConstraint(gen_cst, gurobi_model, &ct_variables));
break;
}
case MPGeneralConstraintProto::kMaxConstraint: {
RETURN_IF_GUROBI_ERROR(
AddMaxConstraint(gen_cst, gurobi_model, &ct_variables));
break;
}
default:
return absl::UnimplementedError(
absl::StrFormat("General constraints of type %i not supported.",
gen_cst.general_constraint_case()));
}
}
}
RETURN_IF_GUROBI_ERROR(GRBsetintattr(gurobi_model, GRB_INT_ATTR_MODELSENSE,
model.maximize() ? -1 : 1));
RETURN_IF_GUROBI_ERROR(GRBsetdblattr(gurobi_model, GRB_DBL_ATTR_OBJCON,
model.objective_offset()));
if (model.has_quadratic_objective()) {
MPQuadraticObjective qobj = model.quadratic_objective();
if (qobj.coefficient_size() > 0) {
RETURN_IF_GUROBI_ERROR(
GRBaddqpterms(gurobi_model, /*numqnz=*/qobj.coefficient_size(),
/*qrow=*/qobj.mutable_qvar1_index()->mutable_data(),
/*qcol=*/qobj.mutable_qvar2_index()->mutable_data(),
/*qval=*/qobj.mutable_coefficient()->mutable_data()));
}
}
RETURN_IF_GUROBI_ERROR(GRBupdatemodel(gurobi_model));
const absl::Time time_before = absl::Now();
UserTimer user_timer;
user_timer.Start();
RETURN_IF_GUROBI_ERROR(GRBoptimize(gurobi_model));
const absl::Duration solving_duration = absl::Now() - time_before;
user_timer.Stop();
VLOG(1) << "Finished solving in GurobiSolveProto(), walltime = "
<< solving_duration << ", usertime = " << user_timer.GetDuration();
response.mutable_solve_info()->set_solve_wall_time_seconds(
absl::ToDoubleSeconds(solving_duration));
response.mutable_solve_info()->set_solve_user_time_seconds(
absl::ToDoubleSeconds(user_timer.GetDuration()));
int optimization_status = 0;
RETURN_IF_GUROBI_ERROR(
GRBgetintattr(gurobi_model, GRB_INT_ATTR_STATUS, &optimization_status));
int solution_count = 0;
RETURN_IF_GUROBI_ERROR(
GRBgetintattr(gurobi_model, GRB_INT_ATTR_SOLCOUNT, &solution_count));
switch (optimization_status) {
case GRB_OPTIMAL:
response.set_status(MPSOLVER_OPTIMAL);
break;
case GRB_INF_OR_UNBD:
DLOG(INFO) << "Gurobi solve returned GRB_INF_OR_UNBD, which we treat as "
"INFEASIBLE even though it may mean UNBOUNDED.";
response.set_status_str(
"The model may actually be unbounded: Gurobi returned "
"GRB_INF_OR_UNBD");
ABSL_FALLTHROUGH_INTENDED;
case GRB_INFEASIBLE:
response.set_status(MPSOLVER_INFEASIBLE);
break;
case GRB_UNBOUNDED:
response.set_status(MPSOLVER_UNBOUNDED);
break;
default: {
if (solution_count > 0) {
response.set_status(MPSOLVER_FEASIBLE);
} else {
response.set_status(MPSOLVER_NOT_SOLVED);
response.set_status_str(
absl::StrFormat("Gurobi status code %d", optimization_status));
}
break;
}
}
if (solution_count > 0 && (response.status() == MPSOLVER_FEASIBLE ||
response.status() == MPSOLVER_OPTIMAL)) {
double objective_value = 0;
RETURN_IF_GUROBI_ERROR(
GRBgetdblattr(gurobi_model, GRB_DBL_ATTR_OBJVAL, &objective_value));
response.set_objective_value(objective_value);
double best_objective_bound = 0;
const int error = GRBgetdblattr(gurobi_model, GRB_DBL_ATTR_OBJBOUND,
&best_objective_bound);
if (response.status() == MPSOLVER_OPTIMAL &&
error == GRB_ERROR_DATA_NOT_AVAILABLE) {
// If the presolve deletes all variables, there's no best bound.
response.set_best_objective_bound(objective_value);
} else {
RETURN_IF_GUROBI_ERROR(error);
response.set_best_objective_bound(best_objective_bound);
}
response.mutable_variable_value()->Resize(variable_size, 0);
RETURN_IF_GUROBI_ERROR(
GRBgetdblattrarray(gurobi_model, GRB_DBL_ATTR_X, 0, variable_size,
response.mutable_variable_value()->mutable_data()));
// NOTE, GurobiSolveProto() is exposed to external clients via MPSolver API,
// which assumes the solution values of integer variables are rounded to
// integer values.
auto round_values_of_integer_variables_fn =
[&](google::protobuf::RepeatedField<double>* values) {
for (int v = 0; v < variable_size; ++v) {
if (model.variable(v).is_integer()) {
(*values)[v] = std::round((*values)[v]);
}
}
};
round_values_of_integer_variables_fn(response.mutable_variable_value());
if (!has_integer_variables && model.general_constraint_size() == 0) {
response.mutable_dual_value()->Resize(model.constraint_size(), 0);
RETURN_IF_GUROBI_ERROR(GRBgetdblattrarray(
gurobi_model, GRB_DBL_ATTR_PI, 0, model.constraint_size(),
response.mutable_dual_value()->mutable_data()));
}
const int additional_solutions = std::min(
solution_count, std::min(request.populate_additional_solutions_up_to(),
std::numeric_limits<int32_t>::max() - 1) +
1);
for (int i = 1; i < additional_solutions; ++i) {
RETURN_IF_GUROBI_ERROR(
GRBsetintparam(model_env, GRB_INT_PAR_SOLUTIONNUMBER, i));
MPSolution* solution = response.add_additional_solutions();
solution->mutable_variable_value()->Resize(variable_size, 0);
double objective_value = 0;
RETURN_IF_GUROBI_ERROR(GRBgetdblattr(
gurobi_model, GRB_DBL_ATTR_POOLOBJVAL, &objective_value));
solution->set_objective_value(objective_value);
RETURN_IF_GUROBI_ERROR(GRBgetdblattrarray(
gurobi_model, GRB_DBL_ATTR_XN, 0, variable_size,
solution->mutable_variable_value()->mutable_data()));
round_values_of_integer_variables_fn(solution->mutable_variable_value());
}
}
#undef RETURN_IF_GUROBI_ERROR
return response;
}
} // namespace operations_research

52
ortools/linear_solver/gurobi_proto_solver.h
View File

@@ -1,52 +0,0 @@
// Copyright 2010-2022 Google LLC
// Licensed under the Apache License, Version 2.0 (the "License");
// you may not use this file except in compliance with the License.
// You may obtain a copy of the License at
//
// http://www.apache.org/licenses/LICENSE-2.0
//
// Unless required by applicable law or agreed to in writing, software
// distributed under the License is distributed on an "AS IS" BASIS,
// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
// See the License for the specific language governing permissions and
// limitations under the License.
#ifndef OR_TOOLS_LINEAR_SOLVER_GUROBI_PROTO_SOLVER_H_
#define OR_TOOLS_LINEAR_SOLVER_GUROBI_PROTO_SOLVER_H_
#include <string>
#include "absl/status/status.h"
#include "absl/status/statusor.h"
#include "ortools/gurobi/environment.h"
#include "ortools/linear_solver/linear_solver.pb.h"
namespace operations_research {
// Solves the input request.
//
// By default this function creates a new primary Gurobi environment, but an
// existing one can be passed as parameter. This can be useful with single-use
// Gurobi licenses since it is not possible to create a second environment if
// one already exists with those licenses.
//
// Please note though that the provided environment should not be actively used
// by another thread at the same time.
absl::StatusOr<MPSolutionResponse> GurobiSolveProto(
const MPModelRequest& request, GRBenv* gurobi_env = nullptr);
// Set parameters specified in the string. The format of the string is a series
// of tokens separated by either '\n' or by ',' characters.
// Any token whose first character is a '#' or has zero length is skiped.
// Comment tokens (i.e. those starting with #) can contain ',' characters.
// Any other token has the form:
// parameter_name(separator)value
// where (separator) is either '=' or ' '.
// A valid string can look-like:
// "#\n# Gurobi-specific parameters, still part of the
// comment\n\nThreads=1\nPresolve 2,SolutionLimit=100" This function will
// process each and every token, even if an intermediate token is unrecognized.
absl::Status SetSolverSpecificParameters(const std::string& parameters,
GRBenv* gurobi);
} // namespace operations_research
#endif // OR_TOOLS_LINEAR_SOLVER_GUROBI_PROTO_SOLVER_H_

135
ortools/linear_solver/pdlp_proto_solver.cc
View File

@@ -1,135 +0,0 @@
// Copyright 2010-2022 Google LLC
// Licensed under the Apache License, Version 2.0 (the "License");
// you may not use this file except in compliance with the License.
// You may obtain a copy of the License at
//
// http://www.apache.org/licenses/LICENSE-2.0
//
// Unless required by applicable law or agreed to in writing, software
// distributed under the License is distributed on an "AS IS" BASIS,
// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
// See the License for the specific language governing permissions and
// limitations under the License.
#include "ortools/linear_solver/pdlp_proto_solver.h"
#include <atomic>
#include <optional>
#include <string>
#include <utility>
#include "absl/status/statusor.h"
#include "ortools/base/logging.h"
#include "ortools/base/status_macros.h"
#include "ortools/linear_solver/linear_solver.pb.h"
#include "ortools/linear_solver/model_validator.h"
#include "ortools/pdlp/iteration_stats.h"
#include "ortools/pdlp/primal_dual_hybrid_gradient.h"
#include "ortools/pdlp/quadratic_program.h"
#include "ortools/pdlp/solve_log.pb.h"
#include "ortools/pdlp/solvers.pb.h"
#include "ortools/port/proto_utils.h"
#include "ortools/util/lazy_mutable_copy.h"
namespace operations_research {
absl::StatusOr<MPSolutionResponse> PdlpSolveProto(
const MPModelRequest& request, const bool relax_integer_variables,
const std::atomic<bool>* interrupt_solve) {
pdlp::PrimalDualHybridGradientParams params;
if (request.enable_internal_solver_output()) {
params.set_verbosity_level(3);
} else {
params.set_verbosity_level(0);
}
MPSolutionResponse error_response;
if (!ProtobufTextFormatMergeFromString(request.solver_specific_parameters(),
&params)) {
error_response.set_status(
MPSolverResponseStatus::MPSOLVER_MODEL_INVALID_SOLVER_PARAMETERS);
return error_response;
}
if (interrupt_solve != nullptr && interrupt_solve->load() == true) {
error_response.set_status(MPSolverResponseStatus::MPSOLVER_NOT_SOLVED);
return error_response;
}
if (request.has_solver_time_limit_seconds()) {
params.mutable_termination_criteria()->set_time_sec_limit(
request.solver_time_limit_seconds());
}
const absl::optional<LazyMutableCopy<MPModelProto>> optional_model =
ExtractValidMPModelOrPopulateResponseStatus(request, &error_response);
if (!optional_model) {
LOG_IF(WARNING, request.enable_internal_solver_output())
<< "Failed to extract a valid model from protocol buffer. Status: "
<< ProtoEnumToString<MPSolverResponseStatus>(error_response.status())
<< " (" << error_response.status()
<< "): " << error_response.status_str();
return error_response;
}
ASSIGN_OR_RETURN(
pdlp::QuadraticProgram qp,
pdlp::QpFromMpModelProto(optional_model->get(), relax_integer_variables));
const double objective_scaling_factor = qp.objective_scaling_factor;
pdlp::SolverResult pdhg_result =
pdlp::PrimalDualHybridGradient(std::move(qp), params, interrupt_solve);
// PDLP's statuses don't map very cleanly to MPSolver statuses. Do the best
// we can for now.
MPSolutionResponse response;
switch (pdhg_result.solve_log.termination_reason()) {
case pdlp::TERMINATION_REASON_OPTIMAL:
response.set_status(MPSOLVER_OPTIMAL);
break;
case pdlp::TERMINATION_REASON_NUMERICAL_ERROR:
response.set_status(MPSOLVER_ABNORMAL);
break;
case pdlp::TERMINATION_REASON_PRIMAL_INFEASIBLE:
response.set_status(MPSOLVER_INFEASIBLE);
break;
case pdlp::TERMINATION_REASON_INTERRUPTED_BY_USER:
response.set_status(MPSOLVER_CANCELLED_BY_USER);
break;
default:
response.set_status(MPSOLVER_NOT_SOLVED);
}
if (pdhg_result.solve_log.has_termination_string()) {
response.set_status_str(pdhg_result.solve_log.termination_string());
}
const std::optional<pdlp::ConvergenceInformation> convergence_information =
pdlp::GetConvergenceInformation(pdhg_result.solve_log.solution_stats(),
pdhg_result.solve_log.solution_type());
if (convergence_information.has_value()) {
response.set_objective_value(convergence_information->primal_objective());
}
// variable_value and dual_value are supposed to be set iff 'status' is
// OPTIMAL or FEASIBLE. However, we set them in all cases.
for (const double v : pdhg_result.primal_solution) {
response.add_variable_value(v);
}
// QpFromMpModelProto converts maximization problems to minimization problems
// for PDLP by negating the objective and setting objective_scaling_factor to
// -1. This maintains the same set of primal solutions. Dual solutions need to
// be negated if objective_scaling_factor is -1.
for (const double v : pdhg_result.dual_solution) {
response.add_dual_value(objective_scaling_factor * v);
}
for (const double v : pdhg_result.reduced_costs) {
response.add_reduced_cost(objective_scaling_factor * v);
}
response.set_solver_specific_info(pdhg_result.solve_log.SerializeAsString());
return response;
}
} // namespace operations_research

45
ortools/linear_solver/pdlp_proto_solver.h
View File

@@ -1,45 +0,0 @@
// Copyright 2010-2022 Google LLC
// Licensed under the Apache License, Version 2.0 (the "License");
// you may not use this file except in compliance with the License.
// You may obtain a copy of the License at
//
// http://www.apache.org/licenses/LICENSE-2.0
//
// Unless required by applicable law or agreed to in writing, software
// distributed under the License is distributed on an "AS IS" BASIS,
// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
// See the License for the specific language governing permissions and
// limitations under the License.
#ifndef OR_TOOLS_LINEAR_SOLVER_PDLP_PROTO_SOLVER_H_
#define OR_TOOLS_LINEAR_SOLVER_PDLP_PROTO_SOLVER_H_
#include <atomic>
#include "absl/status/statusor.h"
#include "ortools/linear_solver/linear_solver.pb.h"
namespace operations_research {
// Uses pdlp::PrimalDualHybridGradient to solve the problem specified by the
// MPModelRequest. Users of this interface should be aware of the size
// limitations of MPModelProto (see, e.g., large_linear_program.proto).
//
// The optional interrupt_solve can be used to interrupt the solve early. The
// solver will periodically check its value and stop if it holds true.
//
// If relax_integer_variables is true, integrality constraints are relaxed
// before solving. If false, integrality constraints result in an error. The
// solver_specific_info field in the MPSolutionResponse contains a serialized
// SolveLog.
//
// Returns an error if the conversion from MPModelProto to
// pdlp::QuadraticProgram fails. The lack of an error does not imply success.
// Check the SolveLog's termination_reason for more refined status details.
absl::StatusOr<MPSolutionResponse> PdlpSolveProto(
const MPModelRequest& request, bool relax_integer_variables = false,
const std::atomic<bool>* interrupt_solve = nullptr);
} // namespace operations_research
#endif // OR_TOOLS_LINEAR_SOLVER_PDLP_PROTO_SOLVER_H_

24
ortools/linear_solver/proto_solver/BUILD.bazel
View File

@@ -36,14 +36,7 @@ cc_library(
"-DUSE_SCIP",
],
deps = [
"//ortools/linear_solver:linear_solver_cc_proto",
"//ortools/linear_solver:scip_with_glop",
"//ortools/linear_solver:model_exporter",
"@com_google_absl//absl/status",
"@com_google_absl//absl/status:statusor",
"@com_google_absl//absl/strings",
"@com_google_absl//absl/synchronization",
"@com_google_absl//absl/types:optional",
"//ortools/base",
"//ortools/base:accurate_sum",
"//ortools/base:dynamic_library",
"//ortools/base:hash",
@@ -51,13 +44,16 @@ cc_library(
"//ortools/base:status_macros",
"//ortools/base:stl_util",
"//ortools/base:timer",
"//ortools/base",
"//ortools/bop:bop_parameters_cc_proto",
"//ortools/bop:integral_solver",
"//ortools/glop:lp_solver",
"//ortools/glop:parameters_cc_proto",
"//ortools/gscip:legacy_scip_params",
"//ortools/gurobi:environment",
"//ortools/linear_solver:linear_solver_cc_proto",
"//ortools/linear_solver:model_exporter",
"//ortools/linear_solver:model_validator",
"//ortools/linear_solver:scip_with_glop",
"//ortools/pdlp:primal_dual_hybrid_gradient",
"//ortools/pdlp:solve_log_cc_proto",
"//ortools/pdlp:solvers_cc_proto",
@@ -68,8 +64,10 @@ cc_library(
"//ortools/sat:lp_utils",
"//ortools/util:fp_utils",
"//ortools/util:lazy_mutable_copy",
] + select({
":with_glpk": ["@glpk//:glpk"],
"//conditions:default": [],
}),
"@com_google_absl//absl/status",
"@com_google_absl//absl/status:statusor",
"@com_google_absl//absl/strings",
"@com_google_absl//absl/synchronization",
"@com_google_absl//absl/types:optional",
],
)

5
ortools/linear_solver/proto_solver/gurobi_proto_solver.cc
View File

@@ -32,7 +32,6 @@
#include "ortools/base/status_macros.h"
#include "ortools/base/timer.h"
#include "ortools/gurobi/environment.h"
#include "ortools/linear_solver/linear_solver.h"
#include "ortools/linear_solver/linear_solver.pb.h"
#include "ortools/linear_solver/model_validator.h"
#include "ortools/util/lazy_mutable_copy.h"
@@ -339,7 +338,9 @@ absl::StatusOr<MPSolutionResponse> GurobiSolveProto(
obj_coeffs[v] = variable.objective_coefficient();
lb[v] = variable.lower_bound();
ub[v] = variable.upper_bound();
ctype[v] = variable.is_integer() && SolverTypeIsMip(request.solver_type())
ctype[v] = variable.is_integer() &&
request.solver_type() ==
MPModelRequest::GUROBI_MIXED_INTEGER_PROGRAMMING
? GRB_INTEGER
: GRB_CONTINUOUS;
if (variable.is_integer()) has_integer_variables = true;

56
ortools/linear_solver/python/BUILD.bazel Normal file
View File

@@ -0,0 +1,56 @@
# Copyright 2010-2022 Google LLC
# Licensed under the Apache License, Version 2.0 (the "License");
# you may not use this file except in compliance with the License.
# You may obtain a copy of the License at
#
# http://www.apache.org/licenses/LICENSE-2.0
#
# Unless required by applicable law or agreed to in writing, software
# distributed under the License is distributed on an "AS IS" BASIS,
# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
# See the License for the specific language governing permissions and
# limitations under the License.
# Python wrapper for model_builder.
load("@ortools_deps//:requirements.bzl", "requirement")
load("@pybind11_bazel//:build_defs.bzl", "pybind_extension")
load("@rules_python//python:defs.bzl", "py_library")
pybind_extension(
name = "pywrap_model_builder_helper",
srcs = ["pywrap_model_builder_helper.cc"],
visibility = ["//visibility:public"],
deps = [
"//ortools/linear_solver:linear_solver_cc_proto",
"//ortools/linear_solver:model_exporter",
"//ortools/linear_solver/wrappers:model_builder_helper",
"@com_google_absl//absl/strings",
"@eigen//:eigen3",
],
)
py_library(
name = "model_builder_helper",
srcs = ["model_builder_helper.py"],
data = [
":pywrap_model_builder_helper.so",
],
visibility = ["//visibility:public"],
deps = [
requirement("numpy"),
"//ortools/linear_solver:linear_solver_py_pb2",
],
)
py_library(
name = "model_builder",
srcs = ["model_builder.py"],
data = [
":pywrap_model_builder_helper.so",
],
visibility = ["//visibility:public"],
deps = [
":model_builder_helper",
"//ortools/linear_solver:linear_solver_py_pb2",
],
)

2
ortools/linear_solver/samples/code_samples.bzl
View File

@@ -13,6 +13,8 @@
"""Helper macro to compile and test code samples."""
load("@ortools_deps//:requirements.bzl", "requirement")
def code_sample_cc(name):
native.cc_binary(
name = name + "_cc",

409
ortools/linear_solver/sat_proto_solver.cc
View File

@@ -1,409 +0,0 @@
// Copyright 2010-2022 Google LLC
// Licensed under the Apache License, Version 2.0 (the "License");
// you may not use this file except in compliance with the License.
// You may obtain a copy of the License at
//
// http://www.apache.org/licenses/LICENSE-2.0
//
// Unless required by applicable law or agreed to in writing, software
// distributed under the License is distributed on an "AS IS" BASIS,
// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
// See the License for the specific language governing permissions and
// limitations under the License.
#include "ortools/linear_solver/sat_proto_solver.h"
#include <cmath>
#include <cstdint>
#include <functional>
#include <memory>
#include <string>
#include <type_traits>
#include <utility>
#include <vector>
#include "absl/status/statusor.h"
#include "ortools/linear_solver/linear_solver.pb.h"
#include "ortools/linear_solver/model_validator.h"
#include "ortools/linear_solver/sat_solver_utils.h"
#include "ortools/port/proto_utils.h"
#include "ortools/sat/cp_model.pb.h"
#include "ortools/sat/cp_model_solver.h"
#include "ortools/sat/lp_utils.h"
#include "ortools/sat/parameters_validation.h"
#include "ortools/sat/sat_parameters.pb.h"
#include "ortools/util/logging.h"
#include "ortools/util/time_limit.h"
namespace operations_research {
namespace {
#if defined(PROTOBUF_INTERNAL_IMPL)
using google::protobuf::Message;
#else
using google::protobuf::Message;
#endif
// Proto-lite disables some features of protos (see
// go/abp-libraries/proto2-lite) and messages inherit from MessageLite directly
// instead of inheriting from Message (which is itself a specialization of
// MessageLite).
constexpr bool kProtoLiteSatParameters =
!std::is_base_of<Message, sat::SatParameters>::value;
MPSolverResponseStatus ToMPSolverResponseStatus(sat::CpSolverStatus status,
bool has_objective) {
switch (status) {
case sat::CpSolverStatus::UNKNOWN:
return MPSOLVER_NOT_SOLVED;
case sat::CpSolverStatus::MODEL_INVALID:
return MPSOLVER_MODEL_INVALID;
case sat::CpSolverStatus::FEASIBLE:
return MPSOLVER_FEASIBLE;
case sat::CpSolverStatus::INFEASIBLE:
return MPSOLVER_INFEASIBLE;
case sat::CpSolverStatus::OPTIMAL:
return MPSOLVER_OPTIMAL;
default: {
}
}
return MPSOLVER_ABNORMAL;
}
sat::CpSolverStatus FromMPSolverResponseStatus(MPSolverResponseStatus status) {
switch (status) {
case MPSolverResponseStatus::MPSOLVER_OPTIMAL:
return sat::OPTIMAL;
case MPSolverResponseStatus::MPSOLVER_INFEASIBLE:
return sat::INFEASIBLE;
case MPSolverResponseStatus::MPSOLVER_MODEL_INVALID:
return sat::MODEL_INVALID;
default: {
}
}
return sat::UNKNOWN;
}
MPSolutionResponse InfeasibleResponse(SolverLogger& logger,
std::string message) {
SOLVER_LOG(&logger, "Infeasible model detected in sat_solve_proto.\n",
message);
// This is needed for our benchmark scripts.
if (logger.LoggingIsEnabled()) {
sat::CpSolverResponse cp_response;
cp_response.set_status(sat::CpSolverStatus::INFEASIBLE);
SOLVER_LOG(&logger, CpSolverResponseStats(cp_response));
}
MPSolutionResponse response;
response.set_status(MPSolverResponseStatus::MPSOLVER_INFEASIBLE);
response.set_status_str(message);
return response;
}
MPSolutionResponse ModelInvalidResponse(SolverLogger& logger,
std::string message) {
SOLVER_LOG(&logger, "Invalid model/parameters in sat_solve_proto.\n",
message);
// This is needed for our benchmark scripts.
if (logger.LoggingIsEnabled()) {
sat::CpSolverResponse cp_response;
cp_response.set_status(sat::CpSolverStatus::MODEL_INVALID);
SOLVER_LOG(&logger, CpSolverResponseStats(cp_response));
}
MPSolutionResponse response;
response.set_status(MPSolverResponseStatus::MPSOLVER_MODEL_INVALID);
response.set_status_str(message);
return response;
}
} // namespace
absl::StatusOr<MPSolutionResponse> SatSolveProto(
MPModelRequest request, std::atomic<bool>* interrupt_solve,
std::function<void(const std::string&)> logging_callback,
std::function<void(const MPSolution&)> solution_callback) {
sat::SatParameters params;
params.set_log_search_progress(request.enable_internal_solver_output());
// Set it now so that it can be overwritten by the solver specific parameters.
if (request.has_solver_specific_parameters()) {
// See EncodeSatParametersAsString() documentation.
if (kProtoLiteSatParameters) {
if (!params.MergeFromString(request.solver_specific_parameters())) {
return absl::InvalidArgumentError(
"solver_specific_parameters is not a valid binary stream of the "
"SatParameters proto");
}
} else {
if (!ProtobufTextFormatMergeFromString(
request.solver_specific_parameters(), &params)) {
return absl::InvalidArgumentError(
"solver_specific_parameters is not a valid textual representation "
"of the SatParameters proto");
}
}
}
if (request.has_solver_time_limit_seconds()) {
params.set_max_time_in_seconds(request.solver_time_limit_seconds());
}
// TODO(user): We do not support all the parameters here. In particular the
// logs before the solver is called will not be appended to the response. Fix
// that, and remove code duplication for the logger config. One way should be
// to not touch/configure anything if the logger is already created while
// calling SolveCpModel() and call a common config function from here or from
// inside Solve()?
SolverLogger logger;
if (logging_callback != nullptr) {
logger.AddInfoLoggingCallback(logging_callback);
}
logger.EnableLogging(params.log_search_progress());
logger.SetLogToStdOut(params.log_to_stdout());
// Model validation and delta handling.
MPSolutionResponse response;
if (!ExtractValidMPModelInPlaceOrPopulateResponseStatus(&request,
&response)) {
// Note that the ExtractValidMPModelInPlaceOrPopulateResponseStatus() can
// also close trivial model (empty or trivially infeasible). So this is not
// always the MODEL_INVALID status.
//
// The logging is only needed for our benchmark script, so we use UNKNOWN
// here, but we could log the proper status instead.
if (logger.LoggingIsEnabled()) {
sat::CpSolverResponse cp_response;
cp_response.set_status(FromMPSolverResponseStatus(response.status()));
SOLVER_LOG(&logger, CpSolverResponseStats(cp_response));
}
return response;
}
// We start by some extra validation since our code do not accept any kind
// of input.
MPModelProto* const mp_model = request.mutable_model();
if (!sat::MPModelProtoValidationBeforeConversion(params, *mp_model,
&logger)) {
return ModelInvalidResponse(logger, "Extra CP-SAT validation failed.");
}
{
const std::string error = sat::ValidateParameters(params);
if (!error.empty()) {
return ModelInvalidResponse(
logger, absl::StrCat("Invalid CP-SAT parameters: ", error));
}
}
// This is good to do before any presolve.
if (!sat::MakeBoundsOfIntegerVariablesInteger(params, mp_model, &logger)) {
return InfeasibleResponse(logger,
"An integer variable has an empty domain");
}
// Coefficients really close to zero can cause issues.
// We remove them right away according to our parameters.
RemoveNearZeroTerms(params, mp_model, &logger);
// Note(user): the LP presolvers API is a bit weird and keep a reference to
// the given GlopParameters, so we need to make sure it outlive them.
const glop::GlopParameters glop_params;
std::vector<std::unique_ptr<glop::Preprocessor>> for_postsolve;
if (!params.enumerate_all_solutions()) {
const glop::ProblemStatus status =
ApplyMipPresolveSteps(glop_params, mp_model, &for_postsolve, &logger);
switch (status) {
case glop::ProblemStatus::INIT:
// Continue with the solve.
break;
case glop::ProblemStatus::PRIMAL_INFEASIBLE:
return InfeasibleResponse(
logger, "Problem proven infeasible during MIP presolve");
case glop::ProblemStatus::INVALID_PROBLEM:
return ModelInvalidResponse(
logger, "Problem detected invalid during MIP presolve");
default:
// TODO(user): We put the INFEASIBLE_OR_UNBOUNBED case here since there
// is no return status that exactly matches it.
if (params.log_search_progress()) {
// This is needed for our benchmark scripts.
sat::CpSolverResponse cp_response;
cp_response.set_status(sat::CpSolverStatus::UNKNOWN);
LOG(INFO) << CpSolverResponseStats(cp_response);
}
response.set_status(MPSolverResponseStatus::MPSOLVER_UNKNOWN_STATUS);
if (status == glop::ProblemStatus::INFEASIBLE_OR_UNBOUNDED) {
response.set_status_str(
"Problem proven infeasible or unbounded during MIP presolve");
}
return response;
}
}
// We need to do that before the automatic detection of integers.
RemoveNearZeroTerms(params, mp_model, &logger);
SOLVER_LOG(&logger, "");
SOLVER_LOG(&logger, "Scaling to pure integer problem.");
const int num_variables = mp_model->variable_size();
std::vector<double> var_scaling(num_variables, 1.0);
if (params.mip_automatically_scale_variables()) {
var_scaling = sat::DetectImpliedIntegers(mp_model, &logger);
if (!sat::MakeBoundsOfIntegerVariablesInteger(params, mp_model, &logger)) {
return InfeasibleResponse(
logger, "A detected integer variable has an empty domain");
}
}
if (params.mip_var_scaling() != 1.0) {
const std::vector<double> other_scaling = sat::ScaleContinuousVariables(
params.mip_var_scaling(), params.mip_max_bound(), mp_model);
for (int i = 0; i < var_scaling.size(); ++i) {
var_scaling[i] *= other_scaling[i];
}
}
// Abort if one only want to solve pure-IP model and we don't have one.
if (params.only_solve_ip()) {
bool all_integer = true;
for (const MPVariableProto& var : mp_model->variable()) {
if (!var.is_integer()) {
all_integer = false;
break;
}
}
if (!all_integer) {
return ModelInvalidResponse(
logger,
"The model contains non-integer variables but the parameter "
"'only_solve_ip' was set. Change this parameter if you "
"still want to solve a more constrained version of the original MIP "
"where non-integer variables can only take a finite set of values.");
}
}
sat::CpModelProto cp_model;
if (!ConvertMPModelProtoToCpModelProto(params, *mp_model, &cp_model,
&logger)) {
return ModelInvalidResponse(logger,
"Failed to convert model into CP-SAT model");
}
DCHECK_EQ(cp_model.variables().size(), var_scaling.size());
DCHECK_EQ(cp_model.variables().size(), mp_model->variable().size());
// Copy and scale the hint if there is one.
if (request.model().has_solution_hint()) {
auto* cp_model_hint = cp_model.mutable_solution_hint();
const int size = request.model().solution_hint().var_index().size();
for (int i = 0; i < size; ++i) {
const int var = request.model().solution_hint().var_index(i);
if (var >= var_scaling.size()) continue;
// To handle weird hint input values, we cap any large value to +/-
// mip_max_bound() which is also the min/max value of any variable once
// scaled.
double value =
request.model().solution_hint().var_value(i) * var_scaling[var];
if (std::abs(value) > params.mip_max_bound()) {
value = value > 0 ? params.mip_max_bound() : -params.mip_max_bound();
}
cp_model_hint->add_vars(var);
cp_model_hint->add_values(static_cast<int64_t>(std::round(value)));
}
}
// We no longer need the request. Reclaim its memory.
const int old_num_variables = mp_model->variable().size();
const int old_num_constraints = mp_model->constraint().size();
request.Clear();
// Configure model.
sat::Model sat_model;
sat_model.Register<SolverLogger>(&logger);
sat_model.Add(NewSatParameters(params));
if (interrupt_solve != nullptr) {
sat_model.GetOrCreate<TimeLimit>()->RegisterExternalBooleanAsLimit(
interrupt_solve);
}
auto post_solve = [&](const sat::CpSolverResponse& cp_response) {
MPSolution mp_solution;
mp_solution.set_objective_value(cp_response.objective_value());
// Postsolve the bound shift and scaling.
glop::ProblemSolution glop_solution((glop::RowIndex(old_num_constraints)),
(glop::ColIndex(old_num_variables)));
for (int v = 0; v < glop_solution.primal_values.size(); ++v) {
glop_solution.primal_values[glop::ColIndex(v)] =
static_cast<double>(cp_response.solution(v)) / var_scaling[v];
}
for (int i = for_postsolve.size(); --i >= 0;) {
for_postsolve[i]->RecoverSolution(&glop_solution);
}
for (int v = 0; v < glop_solution.primal_values.size(); ++v) {
mp_solution.add_variable_value(
glop_solution.primal_values[glop::ColIndex(v)]);
}
return mp_solution;
};
if (solution_callback != nullptr) {
sat_model.Add(sat::NewFeasibleSolutionObserver(
[&](const sat::CpSolverResponse& cp_response) {
solution_callback(post_solve(cp_response));
}));
}
// Solve.
const sat::CpSolverResponse cp_response =
sat::SolveCpModel(cp_model, &sat_model);
// Convert the response.
//
// TODO(user): Implement the row and column status.
response.mutable_solve_info()->set_solve_wall_time_seconds(
cp_response.wall_time());
response.mutable_solve_info()->set_solve_user_time_seconds(
cp_response.user_time());
response.set_status(
ToMPSolverResponseStatus(cp_response.status(), cp_model.has_objective()));
if (response.status() == MPSOLVER_FEASIBLE ||
response.status() == MPSOLVER_OPTIMAL) {
response.set_objective_value(cp_response.objective_value());
response.set_best_objective_bound(cp_response.best_objective_bound());
MPSolution post_solved_solution = post_solve(cp_response);
*response.mutable_variable_value() =
std::move(*post_solved_solution.mutable_variable_value());
}
// Copy and postsolve any additional solutions.
//
// TODO(user): Remove the postsolve hack of copying to a response.
for (int i = 0; i < cp_response.additional_solutions().size(); ++i) {
sat::CpSolverResponse temp;
*temp.mutable_solution() = cp_response.additional_solutions(i).values();
MPSolution post_solved_solution = post_solve(temp);
*(response.add_additional_solutions()->mutable_variable_value()) =
std::move(*post_solved_solution.mutable_variable_value());
}
return response;
}
std::string EncodeSatParametersAsString(const sat::SatParameters& parameters) {
if (kProtoLiteSatParameters) {
// Here we use SerializeToString() instead of SerializeAsString() since the
// later ignores errors and returns an empty string instead (which can be a
// valid value when no fields are set).
std::string bytes;
CHECK(parameters.SerializeToString(&bytes));
return bytes;
}
return ProtobufShortDebugString(parameters);
}
} // namespace operations_research

75
ortools/linear_solver/sat_proto_solver.h
View File

@@ -1,75 +0,0 @@
// Copyright 2010-2022 Google LLC
// Licensed under the Apache License, Version 2.0 (the "License");
// you may not use this file except in compliance with the License.
// You may obtain a copy of the License at
//
// http://www.apache.org/licenses/LICENSE-2.0
//
// Unless required by applicable law or agreed to in writing, software
// distributed under the License is distributed on an "AS IS" BASIS,
// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
// See the License for the specific language governing permissions and
// limitations under the License.
#ifndef OR_TOOLS_LINEAR_SOLVER_SAT_PROTO_SOLVER_H_
#define OR_TOOLS_LINEAR_SOLVER_SAT_PROTO_SOLVER_H_
#include <functional>
#include <string>
#include "absl/status/statusor.h"
#include "ortools/linear_solver/linear_solver.pb.h"
#include "ortools/sat/sat_parameters.pb.h"
#include "ortools/util/logging.h"
namespace operations_research {
// Solve the input MIP model with the SAT solver.
//
// If possible, std::move the request into this function call to avoid a copy.
//
// If you need to change the solver parameters, please use the
// EncodeSatParametersAsString() function below to set the request's
// solver_specific_parameters field.
//
// The optional interrupt_solve can be used to interrupt the solve early. It
// must only be set to true, never reset to false. It is also used internally by
// the solver that will set it to true for its own internal logic. As a
// consequence the caller should ignore the stored value and should not use the
// same atomic for different concurrent calls.
//
// The optional logging_callback will be called when the SAT parameter
// log_search_progress is set to true. Passing a callback will disable the
// default logging to INFO. Note though that by default the SAT parameter
// log_to_stdout is true so even with a callback, the logs will appear on stdout
// too unless log_to_stdout is set to false. The enable_internal_solver_output
// in the request will act as the SAT parameter log_search_progress.
//
// The optional solution_callback will be called on each intermediate solution
// found by the solver. The solver may call solution_callback from multiple
// threads, but it will ensure that at most one thread executes
// solution_callback at a time.
absl::StatusOr<MPSolutionResponse> SatSolveProto(
MPModelRequest request, std::atomic<bool>* interrupt_solve = nullptr,
std::function<void(const std::string&)> logging_callback = nullptr,
std::function<void(const MPSolution&)> solution_callback = nullptr);
// Returns a string that should be used in MPModelRequest's
// solver_specific_parameters field to encode the SAT parameters.
//
// The returned string's content depends on the version of the proto library
// that is linked in the binary.
//
// By default it will contain the textual representation of the input proto.
// But when the proto-lite is used, it will contain the binary stream of the
// proto instead since it is not possible to build the textual representation in
// that case.
//
// The SatSolveProto() function will test if the proto-lite is used and expect a
// binary stream when it is the case. So in order for your code to be portable,
// you should always use this function to set the specific parameters.
std::string EncodeSatParametersAsString(const sat::SatParameters& parameters);
} // namespace operations_research
#endif // OR_TOOLS_LINEAR_SOLVER_SAT_PROTO_SOLVER_H_

113
ortools/linear_solver/sat_solver_utils.cc
View File

@@ -1,113 +0,0 @@
// Copyright 2010-2022 Google LLC
// Licensed under the Apache License, Version 2.0 (the "License");
// you may not use this file except in compliance with the License.
// You may obtain a copy of the License at
//
// http://www.apache.org/licenses/LICENSE-2.0
//
// Unless required by applicable law or agreed to in writing, software
// distributed under the License is distributed on an "AS IS" BASIS,
// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
// See the License for the specific language governing permissions and
// limitations under the License.
#include "ortools/linear_solver/sat_solver_utils.h"
#include <memory>
#include <string>
#include <utility>
#include <vector>
#include "absl/memory/memory.h"
#include "ortools/glop/parameters.pb.h"
#include "ortools/glop/preprocessor.h"
#include "ortools/lp_data/proto_utils.h"
namespace operations_research {
#define ADD_LP_PREPROCESSOR(name) \
names.push_back(#name); \
lp_preprocessors.push_back(absl::make_unique<name>(&glop_params));
glop::ProblemStatus ApplyMipPresolveSteps(
const glop::GlopParameters& glop_params, MPModelProto* model,
std::vector<std::unique_ptr<glop::Preprocessor>>* for_postsolve,
SolverLogger* logger) {
CHECK(model != nullptr);
// TODO(user): General constraints are currently not supported.
if (!model->general_constraint().empty()) {
return glop::ProblemStatus::INIT;
}
// We need to copy the hint because LinearProgramToMPModelProto() loose it.
const bool hint_is_present = model->has_solution_hint();
const auto copy_of_hint = model->solution_hint();
// TODO(user): Remove this back and forth conversion. We could convert
// the LinearProgram directly to a CpModelProto, or we could have a custom
// implementation of these presolve steps.
glop::LinearProgram lp;
glop::MPModelProtoToLinearProgram(*model, &lp);
// These presolve might change the problem size.
//
// TODO(user): transform the hint instead of disabling presolve.
if (!hint_is_present) {
const std::string header =
"Running basic LP presolve, initial problem dimensions: ";
SOLVER_LOG(logger, "");
SOLVER_LOG(logger, header, lp.GetDimensionString());
std::vector<std::string> names;
std::vector<std::unique_ptr<glop::Preprocessor>> lp_preprocessors;
ADD_LP_PREPROCESSOR(glop::FixedVariablePreprocessor);
ADD_LP_PREPROCESSOR(glop::SingletonPreprocessor);
ADD_LP_PREPROCESSOR(glop::ForcingAndImpliedFreeConstraintPreprocessor);
ADD_LP_PREPROCESSOR(glop::FreeConstraintPreprocessor);
// TODO(user): Usually it is good to run the ImpliedFreePreprocessor before
// this one. However this seems to cause problem on atm20-100.mps. Moreover,
// for the conversion, it is better to have tight bounds even if the bound
// propagator is supposed to undo what this presolve would have done.
ADD_LP_PREPROCESSOR(glop::UnconstrainedVariablePreprocessor);
for (int i = 0; i < lp_preprocessors.size(); ++i) {
auto& preprocessor = lp_preprocessors[i];
preprocessor->UseInMipContext();
const bool need_postsolve = preprocessor->Run(&lp);
names[i].resize(header.size(), ' '); // padding.
SOLVER_LOG(logger, names[i], lp.GetDimensionString());
const glop::ProblemStatus status = preprocessor->status();
if (status != glop::ProblemStatus::INIT) return status;
if (need_postsolve) for_postsolve->push_back(std::move(preprocessor));
}
}
// Finally, we make sure all domains contain zero.
if (!hint_is_present) {
auto shift_bounds =
std::make_unique<glop::ShiftVariableBoundsPreprocessor>(&glop_params);
shift_bounds->UseInMipContext();
const bool need_postsolve = shift_bounds->Run(&lp);
if (shift_bounds->status() != glop::ProblemStatus::INIT) {
return shift_bounds->status();
}
if (need_postsolve) {
for_postsolve->push_back(std::move(shift_bounds));
}
}
glop::LinearProgramToMPModelProto(lp, model);
// Restore the hint, note that none of the presolve steps we run here change
// the number of variables in the model.
if (hint_is_present) {
*model->mutable_solution_hint() = copy_of_hint;
}
return glop::ProblemStatus::INIT;
}
#undef ADD_LP_PREPROCESSOR
} // namespace operations_research

43
ortools/linear_solver/sat_solver_utils.h
View File

@@ -1,43 +0,0 @@
// Copyright 2010-2022 Google LLC
// Licensed under the Apache License, Version 2.0 (the "License");
// you may not use this file except in compliance with the License.
// You may obtain a copy of the License at
//
// http://www.apache.org/licenses/LICENSE-2.0
//
// Unless required by applicable law or agreed to in writing, software
// distributed under the License is distributed on an "AS IS" BASIS,
// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
// See the License for the specific language governing permissions and
// limitations under the License.
#ifndef OR_TOOLS_LINEAR_SOLVER_SAT_SOLVER_UTILS_H_
#define OR_TOOLS_LINEAR_SOLVER_SAT_SOLVER_UTILS_H_
#include <memory>
#include <vector>
#include "ortools/glop/preprocessor.h"
#include "ortools/linear_solver/linear_solver.pb.h"
#include "ortools/util/logging.h"
namespace operations_research {
// Applies presolve steps to improve the MIP -> IP imperfect conversion. The
// stricter the domain of the variables, the more room we have for scaling the
// constraint to integers and prevent overflow. Similarly if we can remove
// singleton continuous variables, it is just good to do so.
//
// Returns the presolve status which can currently be:
// - INIT for most cases were nothing was proven during this step.
// - PRIMAL_INFEASIBLE if the model was proven infeasible.
// - INFEASIBLE_OR_UNBOUNDED if the presolve couldn't distinguish between these
// two statuses.
// - ABNORMAL if an error occurred.
glop::ProblemStatus ApplyMipPresolveSteps(
const glop::GlopParameters& glop_params, MPModelProto* model,
std::vector<std::unique_ptr<glop::Preprocessor>>* for_postsolve,
SolverLogger* logger);
} // namespace operations_research
#endif // OR_TOOLS_LINEAR_SOLVER_SAT_SOLVER_UTILS_H_

945
ortools/linear_solver/scip_proto_solver.cc
View File

@@ -1,945 +0,0 @@
// Copyright 2010-2022 Google LLC
// Licensed under the Apache License, Version 2.0 (the "License");
// you may not use this file except in compliance with the License.
// You may obtain a copy of the License at
//
// http://www.apache.org/licenses/LICENSE-2.0
//
// Unless required by applicable law or agreed to in writing, software
// distributed under the License is distributed on an "AS IS" BASIS,
// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
// See the License for the specific language governing permissions and
// limitations under the License.
#if defined(USE_SCIP)
#include "ortools/linear_solver/scip_proto_solver.h"
#include <algorithm>
#include <cmath>
#include <limits>
#include <memory>
#include <numeric>
#include <set>
#include <string>
#include <vector>
#include "absl/status/status.h"
#include "absl/status/statusor.h"
#include "absl/strings/ascii.h"
#include "absl/strings/numbers.h"
#include "absl/strings/str_cat.h"
#include "absl/strings/str_format.h"
#include "absl/strings/str_split.h"
#include "absl/time/time.h"
#include "ortools/base/cleanup.h"
#include "ortools/base/commandlineflags.h"
#include "ortools/base/status_macros.h"
#include "ortools/base/timer.h"
#include "ortools/gscip/legacy_scip_params.h"
#include "ortools/linear_solver/linear_solver.pb.h"
#include "ortools/linear_solver/model_validator.h"
#include "ortools/linear_solver/scip_helper_macros.h"
#include "ortools/util/lazy_mutable_copy.h"
#include "scip/cons_disjunction.h"
#include "scip/cons_linear.h"
#include "scip/cons_quadratic.h"
#include "scip/pub_var.h"
#include "scip/scip.h"
#include "scip/scip_param.h"
#include "scip/scip_prob.h"
#include "scip/scip_var.h"
#include "scip/scipdefplugins.h"
#include "scip/set.h"
#include "scip/struct_paramset.h"
#include "scip/type_cons.h"
#include "scip/type_paramset.h"
#include "scip/type_var.h"
ABSL_FLAG(std::string, scip_proto_solver_output_cip_file, "",
"If given, saves the generated CIP file here. Useful for "
"reporting bugs to SCIP.");
namespace operations_research {
namespace {
// This function will create a new constraint if the indicator constraint has
// both a lower bound and an upper bound.
absl::Status AddIndicatorConstraint(const MPGeneralConstraintProto& gen_cst,
SCIP* scip, SCIP_CONS** scip_cst,
std::vector<SCIP_VAR*>* scip_variables,
std::vector<SCIP_CONS*>* scip_constraints,
std::vector<SCIP_VAR*>* tmp_variables,
std::vector<double>* tmp_coefficients) {
CHECK(scip != nullptr);
CHECK(scip_cst != nullptr);
CHECK(scip_variables != nullptr);
CHECK(scip_constraints != nullptr);
CHECK(tmp_variables != nullptr);
CHECK(tmp_coefficients != nullptr);
CHECK(gen_cst.has_indicator_constraint());
constexpr double kInfinity = std::numeric_limits<double>::infinity();
const auto& ind = gen_cst.indicator_constraint();
if (!ind.has_constraint()) return absl::OkStatus();
const MPConstraintProto& constraint = ind.constraint();
const int size = constraint.var_index_size();
tmp_variables->resize(size, nullptr);
tmp_coefficients->resize(size, 0);
for (int i = 0; i < size; ++i) {
(*tmp_variables)[i] = (*scip_variables)[constraint.var_index(i)];
(*tmp_coefficients)[i] = constraint.coefficient(i);
}
SCIP_VAR* ind_var = (*scip_variables)[ind.var_index()];
if (ind.var_value() == 0) {
RETURN_IF_SCIP_ERROR(
SCIPgetNegatedVar(scip, (*scip_variables)[ind.var_index()], &ind_var));
}
if (ind.constraint().upper_bound() < kInfinity) {
RETURN_IF_SCIP_ERROR(SCIPcreateConsIndicator(
scip, scip_cst, gen_cst.name().c_str(), ind_var, size,
tmp_variables->data(), tmp_coefficients->data(),
ind.constraint().upper_bound(),
/*initial=*/!ind.constraint().is_lazy(),
/*separate=*/true,
/*enforce=*/true,
/*check=*/true,
/*propagate=*/true,
/*local=*/false,
/*dynamic=*/false,
/*removable=*/ind.constraint().is_lazy(),
/*stickingatnode=*/false));
RETURN_IF_SCIP_ERROR(SCIPaddCons(scip, *scip_cst));
scip_constraints->push_back(nullptr);
scip_cst = &scip_constraints->back();
}
if (ind.constraint().lower_bound() > -kInfinity) {
for (int i = 0; i < size; ++i) {
(*tmp_coefficients)[i] *= -1;
}
RETURN_IF_SCIP_ERROR(SCIPcreateConsIndicator(
scip, scip_cst, gen_cst.name().c_str(), ind_var, size,
tmp_variables->data(), tmp_coefficients->data(),
-ind.constraint().lower_bound(),
/*initial=*/!ind.constraint().is_lazy(),
/*separate=*/true,
/*enforce=*/true,
/*check=*/true,
/*propagate=*/true,
/*local=*/false,
/*dynamic=*/false,
/*removable=*/ind.constraint().is_lazy(),
/*stickingatnode=*/false));
RETURN_IF_SCIP_ERROR(SCIPaddCons(scip, *scip_cst));
}
return absl::OkStatus();
}
absl::Status AddSosConstraint(const MPGeneralConstraintProto& gen_cst,
const std::vector<SCIP_VAR*>& scip_variables,
SCIP* scip, SCIP_CONS** scip_cst,
std::vector<SCIP_VAR*>* tmp_variables,
std::vector<double>* tmp_weights) {
CHECK(scip != nullptr);
CHECK(scip_cst != nullptr);
CHECK(tmp_variables != nullptr);
CHECK(tmp_weights != nullptr);
CHECK(gen_cst.has_sos_constraint());
const MPSosConstraint& sos_cst = gen_cst.sos_constraint();
// SOS constraints of type N indicate at most N variables are non-zero.
// Constraints with N variables or less are valid, but useless. They also
// crash SCIP, so we skip them.
if (sos_cst.var_index_size() <= 1) return absl::OkStatus();
if (sos_cst.type() == MPSosConstraint::SOS2 &&
sos_cst.var_index_size() <= 2) {
return absl::OkStatus();
}
tmp_variables->resize(sos_cst.var_index_size(), nullptr);
for (int v = 0; v < sos_cst.var_index_size(); ++v) {
(*tmp_variables)[v] = scip_variables[sos_cst.var_index(v)];
}
tmp_weights->resize(sos_cst.var_index_size(), 0);
if (sos_cst.weight_size() == sos_cst.var_index_size()) {
for (int w = 0; w < sos_cst.weight_size(); ++w) {
(*tmp_weights)[w] = sos_cst.weight(w);
}
} else {
// In theory, SCIP should accept empty weight arrays and use natural
// ordering, but in practice, this crashes their code.
std::iota(tmp_weights->begin(), tmp_weights->end(), 1);
}
switch (sos_cst.type()) {
case MPSosConstraint::SOS1_DEFAULT:
RETURN_IF_SCIP_ERROR(
SCIPcreateConsBasicSOS1(scip,
/*cons=*/scip_cst,
/*name=*/gen_cst.name().c_str(),
/*nvars=*/sos_cst.var_index_size(),
/*vars=*/tmp_variables->data(),
/*weights=*/tmp_weights->data()));
break;
case MPSosConstraint::SOS2:
RETURN_IF_SCIP_ERROR(
SCIPcreateConsBasicSOS2(scip,
/*cons=*/scip_cst,
/*name=*/gen_cst.name().c_str(),
/*nvars=*/sos_cst.var_index_size(),
/*vars=*/tmp_variables->data(),
/*weights=*/tmp_weights->data()));
break;
}
RETURN_IF_SCIP_ERROR(SCIPaddCons(scip, *scip_cst));
return absl::OkStatus();
}
absl::Status AddQuadraticConstraint(
const MPGeneralConstraintProto& gen_cst,
const std::vector<SCIP_VAR*>& scip_variables, SCIP* scip,
SCIP_CONS** scip_cst, std::vector<SCIP_VAR*>* tmp_variables,
std::vector<double>* tmp_coefficients,
std::vector<SCIP_VAR*>* tmp_qvariables1,
std::vector<SCIP_VAR*>* tmp_qvariables2,
std::vector<double>* tmp_qcoefficients) {
CHECK(scip != nullptr);
CHECK(scip_cst != nullptr);
CHECK(tmp_variables != nullptr);
CHECK(tmp_coefficients != nullptr);
CHECK(tmp_qvariables1 != nullptr);
CHECK(tmp_qvariables2 != nullptr);
CHECK(tmp_qcoefficients != nullptr);
CHECK(gen_cst.has_quadratic_constraint());
const MPQuadraticConstraint& quad_cst = gen_cst.quadratic_constraint();
// Process linear part of the constraint.
const int lsize = quad_cst.var_index_size();
CHECK_EQ(quad_cst.coefficient_size(), lsize);
tmp_variables->resize(lsize, nullptr);
tmp_coefficients->resize(lsize, 0.0);
for (int i = 0; i < lsize; ++i) {
(*tmp_variables)[i] = scip_variables[quad_cst.var_index(i)];
(*tmp_coefficients)[i] = quad_cst.coefficient(i);
}
// Process quadratic part of the constraint.
const int qsize = quad_cst.qvar1_index_size();
CHECK_EQ(quad_cst.qvar2_index_size(), qsize);
CHECK_EQ(quad_cst.qcoefficient_size(), qsize);
tmp_qvariables1->resize(qsize, nullptr);
tmp_qvariables2->resize(qsize, nullptr);
tmp_qcoefficients->resize(qsize, 0.0);
for (int i = 0; i < qsize; ++i) {
(*tmp_qvariables1)[i] = scip_variables[quad_cst.qvar1_index(i)];
(*tmp_qvariables2)[i] = scip_variables[quad_cst.qvar2_index(i)];
(*tmp_qcoefficients)[i] = quad_cst.qcoefficient(i);
}
RETURN_IF_SCIP_ERROR(
SCIPcreateConsBasicQuadratic(scip,
/*cons=*/scip_cst,
/*name=*/gen_cst.name().c_str(),
/*nlinvars=*/lsize,
/*linvars=*/tmp_variables->data(),
/*lincoefs=*/tmp_coefficients->data(),
/*nquadterms=*/qsize,
/*quadvars1=*/tmp_qvariables1->data(),
/*quadvars2=*/tmp_qvariables2->data(),
/*quadcoefs=*/tmp_qcoefficients->data(),
/*lhs=*/quad_cst.lower_bound(),
/*rhs=*/quad_cst.upper_bound()));
RETURN_IF_SCIP_ERROR(SCIPaddCons(scip, *scip_cst));
return absl::OkStatus();
}
// Models the constraint y = |x| as y >= 0 plus one disjunction constraint:
// y = x OR y = -x
absl::Status AddAbsConstraint(const MPGeneralConstraintProto& gen_cst,
const std::vector<SCIP_VAR*>& scip_variables,
SCIP* scip, SCIP_CONS** scip_cst) {
CHECK(scip != nullptr);
CHECK(scip_cst != nullptr);
CHECK(gen_cst.has_abs_constraint());
const auto& abs = gen_cst.abs_constraint();
SCIP_VAR* scip_var = scip_variables[abs.var_index()];
SCIP_VAR* scip_resultant_var = scip_variables[abs.resultant_var_index()];
// Set the resultant variable's lower bound to zero if it's negative.
if (SCIPvarGetLbLocal(scip_resultant_var) < 0.0) {
RETURN_IF_SCIP_ERROR(SCIPchgVarLb(scip, scip_resultant_var, 0.0));
}
std::vector<SCIP_VAR*> vars;
std::vector<double> vals;
std::vector<SCIP_CONS*> cons;
auto add_abs_constraint =
[&](const std::string& name_prefix) -> absl::Status {
SCIP_CONS* scip_cons = nullptr;
CHECK(vars.size() == vals.size());
const std::string name =
gen_cst.has_name() ? absl::StrCat(gen_cst.name(), name_prefix) : "";
RETURN_IF_SCIP_ERROR(SCIPcreateConsBasicLinear(
scip, /*cons=*/&scip_cons,
/*name=*/name.c_str(), /*nvars=*/vars.size(), /*vars=*/vars.data(),
/*vals=*/vals.data(), /*lhs=*/0.0, /*rhs=*/0.0));
// Note that the constraints are, by design, not added into the model using
// SCIPaddCons.
cons.push_back(scip_cons);
return absl::OkStatus();
};
// Create an intermediary constraint such that y = -x
vars = {scip_resultant_var, scip_var};
vals = {1, 1};
RETURN_IF_ERROR(add_abs_constraint("_neg"));
// Create an intermediary constraint such that y = x
vals = {1, -1};
RETURN_IF_ERROR(add_abs_constraint("_pos"));
// Activate at least one of the two above constraints.
const std::string name =
gen_cst.has_name() ? absl::StrCat(gen_cst.name(), "_disj") : "";
RETURN_IF_SCIP_ERROR(SCIPcreateConsBasicDisjunction(
scip, /*cons=*/scip_cst, /*name=*/name.c_str(),
/*nconss=*/cons.size(), /*conss=*/cons.data(), /*relaxcons=*/nullptr));
RETURN_IF_SCIP_ERROR(SCIPaddCons(scip, *scip_cst));
return absl::OkStatus();
}
absl::Status AddAndConstraint(const MPGeneralConstraintProto& gen_cst,
const std::vector<SCIP_VAR*>& scip_variables,
SCIP* scip, SCIP_CONS** scip_cst,
std::vector<SCIP_VAR*>* tmp_variables) {
CHECK(scip != nullptr);
CHECK(scip_cst != nullptr);
CHECK(tmp_variables != nullptr);
CHECK(gen_cst.has_and_constraint());
const auto& andcst = gen_cst.and_constraint();
tmp_variables->resize(andcst.var_index_size(), nullptr);
for (int i = 0; i < andcst.var_index_size(); ++i) {
(*tmp_variables)[i] = scip_variables[andcst.var_index(i)];
}
RETURN_IF_SCIP_ERROR(SCIPcreateConsBasicAnd(
scip, /*cons=*/scip_cst,
/*name=*/gen_cst.name().c_str(),
/*resvar=*/scip_variables[andcst.resultant_var_index()],
/*nvars=*/andcst.var_index_size(),
/*vars=*/tmp_variables->data()));
RETURN_IF_SCIP_ERROR(SCIPaddCons(scip, *scip_cst));
return absl::OkStatus();
}
absl::Status AddOrConstraint(const MPGeneralConstraintProto& gen_cst,
const std::vector<SCIP_VAR*>& scip_variables,
SCIP* scip, SCIP_CONS** scip_cst,
std::vector<SCIP_VAR*>* tmp_variables) {
CHECK(scip != nullptr);
CHECK(scip_cst != nullptr);
CHECK(tmp_variables != nullptr);
CHECK(gen_cst.has_or_constraint());
const auto& orcst = gen_cst.or_constraint();
tmp_variables->resize(orcst.var_index_size(), nullptr);
for (int i = 0; i < orcst.var_index_size(); ++i) {
(*tmp_variables)[i] = scip_variables[orcst.var_index(i)];
}
RETURN_IF_SCIP_ERROR(SCIPcreateConsBasicOr(
scip, /*cons=*/scip_cst,
/*name=*/gen_cst.name().c_str(),
/*resvar=*/scip_variables[orcst.resultant_var_index()],
/*nvars=*/orcst.var_index_size(),
/*vars=*/tmp_variables->data()));
RETURN_IF_SCIP_ERROR(SCIPaddCons(scip, *scip_cst));
return absl::OkStatus();
}
// Models the constraint y = min(x1, x2, ... xn, c) with c being a constant with
// - n + 1 constraints to ensure y <= min(x1, x2, ... xn, c)
// - one disjunction constraint among all of the possible y = x1, y = x2, ...
// y = xn, y = c constraints
// Does the equivalent thing for max (with y >= max(...) instead).
absl::Status AddMinMaxConstraint(const MPGeneralConstraintProto& gen_cst,
const std::vector<SCIP_VAR*>& scip_variables,
SCIP* scip, SCIP_CONS** scip_cst,
std::vector<SCIP_CONS*>* scip_constraints,
std::vector<SCIP_VAR*>* tmp_variables) {
CHECK(scip != nullptr);
CHECK(scip_cst != nullptr);
CHECK(tmp_variables != nullptr);
CHECK(gen_cst.has_min_constraint() || gen_cst.has_max_constraint());
const auto& minmax = gen_cst.has_min_constraint() ? gen_cst.min_constraint()
: gen_cst.max_constraint();
const std::set<int> unique_var_indices(minmax.var_index().begin(),
minmax.var_index().end());
SCIP_VAR* scip_resultant_var = scip_variables[minmax.resultant_var_index()];
std::vector<SCIP_VAR*> vars;
std::vector<double> vals;
std::vector<SCIP_CONS*> cons;
auto add_lin_constraint = [&](const std::string& name_prefix,
double lower_bound = 0.0,
double upper_bound = 0.0) -> absl::Status {
SCIP_CONS* scip_cons = nullptr;
CHECK(vars.size() == vals.size());
const std::string name =
gen_cst.has_name() ? absl::StrCat(gen_cst.name(), name_prefix) : "";
RETURN_IF_SCIP_ERROR(SCIPcreateConsBasicLinear(
scip, /*cons=*/&scip_cons,
/*name=*/name.c_str(), /*nvars=*/vars.size(), /*vars=*/vars.data(),
/*vals=*/vals.data(), /*lhs=*/lower_bound, /*rhs=*/upper_bound));
// Note that the constraints are, by design, not added into the model using
// SCIPaddCons.
cons.push_back(scip_cons);
return absl::OkStatus();
};
// Create intermediary constraints such that y = xi
for (const int var_index : unique_var_indices) {
vars = {scip_resultant_var, scip_variables[var_index]};
vals = {1, -1};
RETURN_IF_ERROR(add_lin_constraint(absl::StrCat("_", var_index)));
}
// Create an intermediary constraint such that y = c
if (minmax.has_constant()) {
vars = {scip_resultant_var};
vals = {1};
RETURN_IF_ERROR(
add_lin_constraint("_constant", minmax.constant(), minmax.constant()));
}
// Activate at least one of the above constraints.
const std::string name =
gen_cst.has_name() ? absl::StrCat(gen_cst.name(), "_disj") : "";
RETURN_IF_SCIP_ERROR(SCIPcreateConsBasicDisjunction(
scip, /*cons=*/scip_cst, /*name=*/name.c_str(),
/*nconss=*/cons.size(), /*conss=*/cons.data(), /*relaxcons=*/nullptr));
RETURN_IF_SCIP_ERROR(SCIPaddCons(scip, *scip_cst));
// Add all of the inequality constraints.
constexpr double kInfinity = std::numeric_limits<double>::infinity();
cons.clear();
for (const int var_index : unique_var_indices) {
vars = {scip_resultant_var, scip_variables[var_index]};
vals = {1, -1};
if (gen_cst.has_min_constraint()) {
RETURN_IF_ERROR(add_lin_constraint(absl::StrCat("_ineq_", var_index),
-kInfinity, 0.0));
} else {
RETURN_IF_ERROR(add_lin_constraint(absl::StrCat("_ineq_", var_index), 0.0,
kInfinity));
}
}
if (minmax.has_constant()) {
vars = {scip_resultant_var};
vals = {1};
if (gen_cst.has_min_constraint()) {
RETURN_IF_ERROR(add_lin_constraint(absl::StrCat("_ineq_constant"),
-kInfinity, minmax.constant()));
} else {
RETURN_IF_ERROR(add_lin_constraint(absl::StrCat("_ineq_constant"),
minmax.constant(), kInfinity));
}
}
for (SCIP_CONS* scip_cons : cons) {
scip_constraints->push_back(scip_cons);
RETURN_IF_SCIP_ERROR(SCIPaddCons(scip, scip_cons));
}
return absl::OkStatus();
}
absl::Status AddQuadraticObjective(const MPQuadraticObjective& quadobj,
SCIP* scip,
std::vector<SCIP_VAR*>* scip_variables,
std::vector<SCIP_CONS*>* scip_constraints) {
CHECK(scip != nullptr);
CHECK(scip_variables != nullptr);
CHECK(scip_constraints != nullptr);
constexpr double kInfinity = std::numeric_limits<double>::infinity();
const int size = quadobj.coefficient_size();
if (size == 0) return absl::OkStatus();
// SCIP supports quadratic objectives by adding a quadratic constraint. We
// need to create an extra variable to hold this quadratic objective.
scip_variables->push_back(nullptr);
RETURN_IF_SCIP_ERROR(SCIPcreateVarBasic(scip, /*var=*/&scip_variables->back(),
/*name=*/"quadobj",
/*lb=*/-kInfinity, /*ub=*/kInfinity,
/*obj=*/1,
/*vartype=*/SCIP_VARTYPE_CONTINUOUS));
RETURN_IF_SCIP_ERROR(SCIPaddVar(scip, scip_variables->back()));
scip_constraints->push_back(nullptr);
SCIP_VAR* linvars[1] = {scip_variables->back()};
double lincoefs[1] = {-1};
std::vector<SCIP_VAR*> quadvars1(size, nullptr);
std::vector<SCIP_VAR*> quadvars2(size, nullptr);
std::vector<double> quadcoefs(size, 0);
for (int i = 0; i < size; ++i) {
quadvars1[i] = scip_variables->at(quadobj.qvar1_index(i));
quadvars2[i] = scip_variables->at(quadobj.qvar2_index(i));
quadcoefs[i] = quadobj.coefficient(i);
}
RETURN_IF_SCIP_ERROR(SCIPcreateConsBasicQuadratic(
scip, /*cons=*/&scip_constraints->back(), /*name=*/"quadobj",
/*nlinvars=*/1, /*linvars=*/linvars, /*lincoefs=*/lincoefs,
/*nquadterms=*/size, /*quadvars1=*/quadvars1.data(),
/*quadvars2=*/quadvars2.data(), /*quadcoefs=*/quadcoefs.data(),
/*lhs=*/0, /*rhs=*/0));
RETURN_IF_SCIP_ERROR(SCIPaddCons(scip, scip_constraints->back()));
return absl::OkStatus();
}
absl::Status AddSolutionHint(const MPModelProto& model, SCIP* scip,
const std::vector<SCIP_VAR*>& scip_variables) {
CHECK(scip != nullptr);
if (!model.has_solution_hint()) return absl::OkStatus();
const PartialVariableAssignment& solution_hint = model.solution_hint();
SCIP_SOL* solution;
bool is_solution_partial =
solution_hint.var_index_size() != model.variable_size();
if (is_solution_partial) {
RETURN_IF_SCIP_ERROR(
SCIPcreatePartialSol(scip, /*sol=*/&solution, /*heur=*/nullptr));
} else {
RETURN_IF_SCIP_ERROR(
SCIPcreateSol(scip, /*sol=*/&solution, /*heur=*/nullptr));
}
for (int i = 0; i < solution_hint.var_index_size(); ++i) {
RETURN_IF_SCIP_ERROR(SCIPsetSolVal(
scip, solution, scip_variables[solution_hint.var_index(i)],
solution_hint.var_value(i)));
}
SCIP_Bool is_stored;
RETURN_IF_SCIP_ERROR(SCIPaddSolFree(scip, &solution, &is_stored));
return absl::OkStatus();
}
} // namespace
// Returns "" iff the model seems valid for SCIP, else returns a human-readable
// error message. Assumes that FindErrorInMPModelProto(model) found no error.
std::string FindErrorInMPModelForScip(const MPModelProto& model, SCIP* scip) {
CHECK(scip != nullptr);
const double infinity = SCIPinfinity(scip);
for (int v = 0; v < model.variable_size(); ++v) {
const MPVariableProto& variable = model.variable(v);
if (variable.lower_bound() >= infinity) {
return absl::StrFormat(
"Variable %i's lower bound is considered +infinity", v);
}
if (variable.upper_bound() <= -infinity) {
return absl::StrFormat(
"Variable %i's upper bound is considered -infinity", v);
}
const double coeff = variable.objective_coefficient();
if (coeff >= infinity || coeff <= -infinity) {
return absl::StrFormat(
"Variable %i's objective coefficient is considered infinite", v);
}
}
for (int c = 0; c < model.constraint_size(); ++c) {
const MPConstraintProto& cst = model.constraint(c);
if (cst.lower_bound() >= infinity) {
return absl::StrFormat(
"Constraint %d's lower_bound is considered +infinity", c);
}
if (cst.upper_bound() <= -infinity) {
return absl::StrFormat(
"Constraint %d's upper_bound is considered -infinity", c);
}
for (int i = 0; i < cst.coefficient_size(); ++i) {
if (std::abs(cst.coefficient(i)) >= infinity) {
return absl::StrFormat(
"Constraint %d's coefficient #%d is considered infinite", c, i);
}
}
}
for (int c = 0; c < model.general_constraint_size(); ++c) {
const MPGeneralConstraintProto& cst = model.general_constraint(c);
switch (cst.general_constraint_case()) {
case MPGeneralConstraintProto::kQuadraticConstraint:
if (cst.quadratic_constraint().lower_bound() >= infinity) {
return absl::StrFormat(
"Quadratic constraint %d's lower_bound is considered +infinity",
c);
}
if (cst.quadratic_constraint().upper_bound() <= -infinity) {
return absl::StrFormat(
"Quadratic constraint %d's upper_bound is considered -infinity",
c);
}
for (int i = 0; i < cst.quadratic_constraint().coefficient_size();
++i) {
const double coefficient = cst.quadratic_constraint().coefficient(i);
if (coefficient >= infinity || coefficient <= -infinity) {
return absl::StrFormat(
"Quadratic constraint %d's linear coefficient #%d considered "
"infinite",
c, i);
}
}
for (int i = 0; i < cst.quadratic_constraint().qcoefficient_size();
++i) {
const double qcoefficient =
cst.quadratic_constraint().qcoefficient(i);
if (qcoefficient >= infinity || qcoefficient <= -infinity) {
return absl::StrFormat(
"Quadratic constraint %d's quadratic coefficient #%d "
"considered infinite",
c, i);
}
}
break;
case MPGeneralConstraintProto::kMinConstraint:
if (cst.min_constraint().constant() >= infinity ||
cst.min_constraint().constant() <= -infinity) {
return absl::StrFormat(
"Min constraint %d's coefficient constant considered infinite",
c);
}
break;
case MPGeneralConstraintProto::kMaxConstraint:
if (cst.max_constraint().constant() >= infinity ||
cst.max_constraint().constant() <= -infinity) {
return absl::StrFormat(
"Max constraint %d's coefficient constant considered infinite",
c);
}
break;
default:
continue;
}
}
const MPQuadraticObjective& quad_obj = model.quadratic_objective();
for (int i = 0; i < quad_obj.coefficient_size(); ++i) {
if (std::abs(quad_obj.coefficient(i)) >= infinity) {
return absl::StrFormat(
"Quadratic objective term #%d's coefficient is considered infinite",
i);
}
}
if (model.has_solution_hint()) {
for (int i = 0; i < model.solution_hint().var_value_size(); ++i) {
const double value = model.solution_hint().var_value(i);
if (value >= infinity || value <= -infinity) {
return absl::StrFormat(
"Variable %i's solution hint is considered infinite",
model.solution_hint().var_index(i));
}
}
}
if (model.objective_offset() >= infinity ||
model.objective_offset() <= -infinity) {
return "Model's objective offset is considered infinite.";
}
return "";
}
absl::StatusOr<MPSolutionResponse> ScipSolveProto(
const MPModelRequest& request) {
MPSolutionResponse response;
const absl::optional<LazyMutableCopy<MPModelProto>> optional_model =
ExtractValidMPModelOrPopulateResponseStatus(request, &response);
if (!optional_model) return response;
const MPModelProto& model = optional_model->get();
SCIP* scip = nullptr;
std::vector<SCIP_VAR*> scip_variables(model.variable_size(), nullptr);
std::vector<SCIP_CONS*> scip_constraints(
model.constraint_size() + model.general_constraint_size(), nullptr);
auto delete_scip_objects = [&]() -> absl::Status {
// Release all created pointers.
if (scip == nullptr) return absl::OkStatus();
for (SCIP_VAR* variable : scip_variables) {
if (variable != nullptr) {
RETURN_IF_SCIP_ERROR(SCIPreleaseVar(scip, &variable));
}
}
for (SCIP_CONS* constraint : scip_constraints) {
if (constraint != nullptr) {
RETURN_IF_SCIP_ERROR(SCIPreleaseCons(scip, &constraint));
}
}
RETURN_IF_SCIP_ERROR(SCIPfree(&scip));
return absl::OkStatus();
};
auto scip_deleter = absl::MakeCleanup([delete_scip_objects]() {
const absl::Status deleter_status = delete_scip_objects();
LOG_IF(DFATAL, !deleter_status.ok()) << deleter_status;
});
RETURN_IF_SCIP_ERROR(SCIPcreate(&scip));
RETURN_IF_SCIP_ERROR(SCIPincludeDefaultPlugins(scip));
const std::string scip_model_invalid_error =
FindErrorInMPModelForScip(model, scip);
if (!scip_model_invalid_error.empty()) {
response.set_status(MPSOLVER_MODEL_INVALID);
response.set_status_str(scip_model_invalid_error);
return response;
}
const auto parameters_status = LegacyScipSetSolverSpecificParameters(
request.solver_specific_parameters(), scip);
if (!parameters_status.ok()) {
response.set_status(MPSOLVER_MODEL_INVALID_SOLVER_PARAMETERS);
response.set_status_str(
std::string(parameters_status.message())); // NOLINT
return response;
}
// Default clock type. We use wall clock time because getting CPU user seconds
// involves calling times() which is very expensive.
// NOTE(user): Also, time limit based on CPU user seconds is *NOT* thread
// safe. We observed that different instances of SCIP running concurrently
// in different threads consume the time limit *together*. E.g., 2 threads
// running SCIP with time limit 10s each will both terminate after ~5s.
RETURN_IF_SCIP_ERROR(
SCIPsetIntParam(scip, "timing/clocktype", SCIP_CLOCKTYPE_WALL));
if (request.solver_time_limit_seconds() > 0 &&
request.solver_time_limit_seconds() < 1e20) {
RETURN_IF_SCIP_ERROR(SCIPsetRealParam(scip, "limits/time",
request.solver_time_limit_seconds()));
}
SCIPsetMessagehdlrQuiet(scip, !request.enable_internal_solver_output());
RETURN_IF_SCIP_ERROR(SCIPcreateProbBasic(scip, model.name().c_str()));
if (model.maximize()) {
RETURN_IF_SCIP_ERROR(SCIPsetObjsense(scip, SCIP_OBJSENSE_MAXIMIZE));
}
for (int v = 0; v < model.variable_size(); ++v) {
const MPVariableProto& variable = model.variable(v);
RETURN_IF_SCIP_ERROR(SCIPcreateVarBasic(
scip, /*var=*/&scip_variables[v], /*name=*/variable.name().c_str(),
/*lb=*/variable.lower_bound(), /*ub=*/variable.upper_bound(),
/*obj=*/variable.objective_coefficient(),
/*vartype=*/variable.is_integer() ? SCIP_VARTYPE_INTEGER
: SCIP_VARTYPE_CONTINUOUS));
RETURN_IF_SCIP_ERROR(SCIPaddVar(scip, scip_variables[v]));
}
{
std::vector<SCIP_VAR*> ct_variables;
std::vector<double> ct_coefficients;
for (int c = 0; c < model.constraint_size(); ++c) {
const MPConstraintProto& constraint = model.constraint(c);
const int size = constraint.var_index_size();
ct_variables.resize(size, nullptr);
ct_coefficients.resize(size, 0);
for (int i = 0; i < size; ++i) {
ct_variables[i] = scip_variables[constraint.var_index(i)];
ct_coefficients[i] = constraint.coefficient(i);
}
RETURN_IF_SCIP_ERROR(SCIPcreateConsLinear(
scip, /*cons=*/&scip_constraints[c],
/*name=*/constraint.name().c_str(),
/*nvars=*/constraint.var_index_size(), /*vars=*/ct_variables.data(),
/*vals=*/ct_coefficients.data(),
/*lhs=*/constraint.lower_bound(), /*rhs=*/constraint.upper_bound(),
/*initial=*/!constraint.is_lazy(),
/*separate=*/true,
/*enforce=*/true,
/*check=*/true,
/*propagate=*/true,
/*local=*/false,
/*modifiable=*/false,
/*dynamic=*/false,
/*removable=*/constraint.is_lazy(),
/*stickingatnode=*/false));
RETURN_IF_SCIP_ERROR(SCIPaddCons(scip, scip_constraints[c]));
}
// These extra arrays are used by quadratic constraints.
std::vector<SCIP_VAR*> ct_qvariables1;
std::vector<SCIP_VAR*> ct_qvariables2;
std::vector<double> ct_qcoefficients;
const int lincst_size = model.constraint_size();
for (int c = 0; c < model.general_constraint_size(); ++c) {
const MPGeneralConstraintProto& gen_cst = model.general_constraint(c);
switch (gen_cst.general_constraint_case()) {
case MPGeneralConstraintProto::kIndicatorConstraint: {
RETURN_IF_ERROR(AddIndicatorConstraint(
gen_cst, scip, &scip_constraints[lincst_size + c],
&scip_variables, &scip_constraints, &ct_variables,
&ct_coefficients));
break;
}
case MPGeneralConstraintProto::kSosConstraint: {
RETURN_IF_ERROR(AddSosConstraint(gen_cst, scip_variables, scip,
&scip_constraints[lincst_size + c],
&ct_variables, &ct_coefficients));
break;
}
case MPGeneralConstraintProto::kQuadraticConstraint: {
RETURN_IF_ERROR(AddQuadraticConstraint(
gen_cst, scip_variables, scip, &scip_constraints[lincst_size + c],
&ct_variables, &ct_coefficients, &ct_qvariables1, &ct_qvariables2,
&ct_qcoefficients));
break;
}
case MPGeneralConstraintProto::kAbsConstraint: {
RETURN_IF_ERROR(AddAbsConstraint(gen_cst, scip_variables, scip,
&scip_constraints[lincst_size + c]));
break;
}
case MPGeneralConstraintProto::kAndConstraint: {
RETURN_IF_ERROR(AddAndConstraint(gen_cst, scip_variables, scip,
&scip_constraints[lincst_size + c],
&ct_variables));
break;
}
case MPGeneralConstraintProto::kOrConstraint: {
RETURN_IF_ERROR(AddOrConstraint(gen_cst, scip_variables, scip,
&scip_constraints[lincst_size + c],
&ct_variables));
break;
}
case MPGeneralConstraintProto::kMinConstraint:
case MPGeneralConstraintProto::kMaxConstraint: {
RETURN_IF_ERROR(AddMinMaxConstraint(
gen_cst, scip_variables, scip, &scip_constraints[lincst_size + c],
&scip_constraints, &ct_variables));
break;
}
default:
return absl::UnimplementedError(
absl::StrFormat("General constraints of type %i not supported.",
gen_cst.general_constraint_case()));
}
}
}
if (model.has_quadratic_objective()) {
RETURN_IF_ERROR(AddQuadraticObjective(model.quadratic_objective(), scip,
&scip_variables, &scip_constraints));
}
RETURN_IF_SCIP_ERROR(SCIPaddOrigObjoffset(scip, model.objective_offset()));
RETURN_IF_ERROR(AddSolutionHint(model, scip, scip_variables));
if (!absl::GetFlag(FLAGS_scip_proto_solver_output_cip_file).empty()) {
SCIPwriteOrigProblem(
scip, absl::GetFlag(FLAGS_scip_proto_solver_output_cip_file).c_str(),
nullptr, true);
}
const absl::Time time_before = absl::Now();
UserTimer user_timer;
user_timer.Start();
RETURN_IF_SCIP_ERROR(SCIPsolve(scip));
const absl::Duration solving_duration = absl::Now() - time_before;
user_timer.Stop();
VLOG(1) << "Finished solving in ScipSolveProto(), walltime = "
<< solving_duration << ", usertime = " << user_timer.GetDuration();
response.mutable_solve_info()->set_solve_wall_time_seconds(
absl::ToDoubleSeconds(solving_duration));
response.mutable_solve_info()->set_solve_user_time_seconds(
absl::ToDoubleSeconds(user_timer.GetDuration()));
const int solution_count =
std::min(SCIPgetNSols(scip),
std::min(request.populate_additional_solutions_up_to(),
std::numeric_limits<int32_t>::max() - 1) +
1);
if (solution_count > 0) {
// can't make 'scip_solution' const, as SCIPxxx does not offer const
// parameter functions.
auto scip_solution_to_repeated_field = [&](SCIP_SOL* scip_solution) {
google::protobuf::RepeatedField<double> variable_value;
variable_value.Reserve(model.variable_size());
for (int v = 0; v < model.variable_size(); ++v) {
double value = SCIPgetSolVal(scip, scip_solution, scip_variables[v]);
if (model.variable(v).is_integer()) {
value = std::round(value);
}
variable_value.AddAlreadyReserved(value);
}
return variable_value;
};
// NOTE(user): As of SCIP 8.0.1, getting the pointer to all
// solutions is as fast as getting the pointer to the best solution.
SCIP_SOL** const scip_solutions = SCIPgetSols(scip);
response.set_objective_value(SCIPgetSolOrigObj(scip, scip_solutions[0]));
response.set_best_objective_bound(SCIPgetDualbound(scip));
*response.mutable_variable_value() =
scip_solution_to_repeated_field(scip_solutions[0]);
for (int i = 1; i < solution_count; ++i) {
MPSolution* solution = response.add_additional_solutions();
solution->set_objective_value(SCIPgetSolOrigObj(scip, scip_solutions[i]));
*solution->mutable_variable_value() =
scip_solution_to_repeated_field(scip_solutions[i]);
}
}
const SCIP_STATUS scip_status = SCIPgetStatus(scip);
switch (scip_status) {
case SCIP_STATUS_OPTIMAL:
response.set_status(MPSOLVER_OPTIMAL);
break;
case SCIP_STATUS_GAPLIMIT:
// To be consistent with the other solvers.
response.set_status(MPSOLVER_OPTIMAL);
break;
case SCIP_STATUS_INFORUNBD:
// NOTE(user): After looking at the SCIP code on 2019-06-14, it seems
// that this will mostly happen for INFEASIBLE problems in practice.
// Since most (all?) users shouldn't have their application behave very
// differently upon INFEASIBLE or UNBOUNDED, the potential error that we
// are making here seems reasonable (and not worth a LOG, unless in
// debug mode).
DLOG(INFO) << "SCIP solve returned SCIP_STATUS_INFORUNBD, which we treat "
"as INFEASIBLE even though it may mean UNBOUNDED.";
response.set_status_str(
"The model may actually be unbounded: SCIP returned "
"SCIP_STATUS_INFORUNBD");
ABSL_FALLTHROUGH_INTENDED;
case SCIP_STATUS_INFEASIBLE:
response.set_status(MPSOLVER_INFEASIBLE);
break;
case SCIP_STATUS_UNBOUNDED:
response.set_status(MPSOLVER_UNBOUNDED);
break;
default:
if (solution_count > 0) {
response.set_status(MPSOLVER_FEASIBLE);
} else {
response.set_status(MPSOLVER_NOT_SOLVED);
response.set_status_str(absl::StrFormat("SCIP status code %d",
static_cast<int>(scip_status)));
}
break;
}
VLOG(1) << "ScipSolveProto() status="
<< MPSolverResponseStatus_Name(response.status()) << ".";
return response;
}
} // namespace operations_research
#endif // #if defined(USE_SCIP)

36
ortools/linear_solver/scip_proto_solver.h
View File

@@ -1,36 +0,0 @@
// Copyright 2010-2022 Google LLC
// Licensed under the Apache License, Version 2.0 (the "License");
// you may not use this file except in compliance with the License.
// You may obtain a copy of the License at
//
// http://www.apache.org/licenses/LICENSE-2.0
//
// Unless required by applicable law or agreed to in writing, software
// distributed under the License is distributed on an "AS IS" BASIS,
// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
// See the License for the specific language governing permissions and
// limitations under the License.
#ifndef OR_TOOLS_LINEAR_SOLVER_SCIP_PROTO_SOLVER_H_
#define OR_TOOLS_LINEAR_SOLVER_SCIP_PROTO_SOLVER_H_
#include <string>
#include "absl/status/statusor.h"
#include "ortools/linear_solver/linear_solver.pb.h"
#include "scip/type_scip.h"
namespace operations_research {
// Note, here we do not override any of SCIP default parameters. This behavior
// *differs* from `MPSolver::Solve()` which sets the feasibility tolerance to
// 1e-7, and the gap limit to 0.0001 (whereas SCIP defaults are 1e-6 and 0,
// respectively, and they are being used here).
absl::StatusOr<MPSolutionResponse> ScipSolveProto(
const MPModelRequest& request);
std::string FindErrorInMPModelForScip(const MPModelProto& model, SCIP* scip);
} // namespace operations_research
#endif // OR_TOOLS_LINEAR_SOLVER_SCIP_PROTO_SOLVER_H_

2
ortools/math_opt/solvers/cp_sat_solver.cc
View File

@@ -36,7 +36,7 @@
#include "ortools/base/protoutil.h"
#include "ortools/base/status_macros.h"
#include "ortools/linear_solver/linear_solver.pb.h"
#include "ortools/linear_solver/sat_proto_solver.h"
#include "ortools/linear_solver/proto_solver/sat_proto_solver.h"
#include "ortools/math_opt/callback.pb.h"
#include "ortools/math_opt/core/math_opt_proto_utils.h"
#include "ortools/math_opt/core/solve_interrupter.h"