From 5f52c87b50e945ef45e4ee62c40d654a26bcce18 Mon Sep 17 00:00:00 2001 From: Laurent Perron Date: Fri, 8 Dec 2023 14:48:57 +0100 Subject: [PATCH] update gscip code --- ortools/gscip/gscip.cc | 5 +- ortools/gscip/gscip.h | 48 +++++- ortools/gscip/gscip_constraint_handler.cc | 191 ++++++++++++---------- ortools/gscip/gscip_constraint_handler.h | 83 ++++++++-- 4 files changed, 218 insertions(+), 109 deletions(-) diff --git a/ortools/gscip/gscip.cc b/ortools/gscip/gscip.cc index 789dc13d62..8b3c87edb1 100644 --- a/ortools/gscip/gscip.cc +++ b/ortools/gscip/gscip.cc @@ -937,8 +937,9 @@ absl::StatusOr GScip::Solve( "GScip is in an error state due to a previous GScip::Solve()"); } interrupt_event_handler_.set_interrupter(interrupter); - const absl::Cleanup interrupt_cleanup( - [this]() { interrupt_event_handler_.set_interrupter(nullptr); }); + const absl::Cleanup interrupt_cleanup = [this]() { + interrupt_event_handler_.set_interrupter(nullptr); + }; // A four step process: // 1. Apply parameters. diff --git a/ortools/gscip/gscip.h b/ortools/gscip/gscip.h index 6b16655b72..0bf6684ade 100644 --- a/ortools/gscip/gscip.h +++ b/ortools/gscip/gscip.h @@ -72,9 +72,11 @@ #include "absl/strings/string_view.h" #include "absl/synchronization/mutex.h" #include "absl/types/span.h" +#include "ortools/base/status_macros.h" #include "ortools/gscip/gscip.pb.h" #include "ortools/gscip/gscip_event_handler.h" #include "ortools/gscip/gscip_message_handler.h" // IWYU pragma: export +#include "ortools/linear_solver/scip_helper_macros.h" #include "scip/scip.h" #include "scip/scip_prob.h" #include "scip/type_cons.h" @@ -166,6 +168,7 @@ class GScip { private: std::atomic interrupted_{false}; }; + // Create a new GScip (the constructor is private). The default objective // direction is minimization. static absl::StatusOr> Create( @@ -413,6 +416,23 @@ class GScip { // CHECK fails if status is OK. void InterruptSolveFromCallbackOnCallbackError(absl::Status error_status); + // Internal use only. Users should instead call + // GScipConstraintHandler::AddCallbackConstraint(). + template + inline absl::StatusOr AddConstraintForHandler( + ConsHandler* handler, ConsData* data, const std::string& name = "", + const GScipConstraintOptions& options = DefaultGScipConstraintOptions()); + + // Internal use only. + // + // Replaces +/- inf by +/- ScipInf(), fails when |d| is in [ScipInf(), inf). + absl::StatusOr ScipInfClamp(double d); + + // Internal use only. + // + // Returns +/- inf if |d| >= ScipInf(), otherwise returns d. + double ScipInfUnclamp(double d); + private: // Event handler that it used to call SCIPinterruptSolve() is a safe manner. // @@ -463,12 +483,6 @@ class GScip { absl::string_view legacy_params); absl::Status FreeTransform(); - // Replaces +/- inf by +/- ScipInf(), fails when |d| is in [ScipInf(), inf). - absl::StatusOr ScipInfClamp(double d); - - // Returns +/- inf if |d| >= ScipInf(), otherwise returns d. - double ScipInfUnclamp(double d); - // Returns an error if |d| >= ScipInf(). absl::Status CheckScipFinite(double d); @@ -654,6 +668,28 @@ struct GScipConstraintOptions { bool keep_alive = true; }; +//////////////////////////////////////////////////////////////////////////////// +// Template function implementations +//////////////////////////////////////////////////////////////////////////////// + +template +absl::StatusOr GScip::AddConstraintForHandler( + ConsHandler* handler, ConsData* data, const std::string& name, + const GScipConstraintOptions& options) { + SCIP_CONS* constraint = nullptr; + RETURN_IF_SCIP_ERROR(SCIPcreateCons( + scip_, &constraint, name.data(), handler, data, options.initial, + options.separate, options.enforce, options.check, options.propagate, + options.local, options.modifiable, options.dynamic, options.removable, + options.sticking_at_node)); + if (constraint == nullptr) { + return absl::InternalError("SCIP failed to create constraint"); + } + RETURN_IF_SCIP_ERROR(SCIPaddCons(scip_, constraint)); + RETURN_IF_ERROR(MaybeKeepConstraintAlive(constraint, options)); + return constraint; +} + } // namespace operations_research #endif // OR_TOOLS_GSCIP_GSCIP_H_ diff --git a/ortools/gscip/gscip_constraint_handler.cc b/ortools/gscip/gscip_constraint_handler.cc index f545249041..2bdb068455 100644 --- a/ortools/gscip/gscip_constraint_handler.cc +++ b/ortools/gscip/gscip_constraint_handler.cc @@ -20,7 +20,6 @@ #include "absl/log/check.h" #include "absl/status/status.h" -#include "absl/strings/string_view.h" #include "absl/types/span.h" #include "ortools/base/logging.h" #include "ortools/base/status_macros.h" @@ -35,7 +34,6 @@ #include "scip/scip_cut.h" #include "scip/scip_general.h" #include "scip/scip_lp.h" -#include "scip/scip_prob.h" #include "scip/scip_sol.h" #include "scip/scip_solvingstats.h" #include "scip/scip_tree.h" @@ -66,22 +64,6 @@ using GScipHandler = namespace operations_research { namespace { -// Enum with supported user-implementable callback functions in the SCIP -// constraint handler. Non-user-implementable functions are not included here -// (e.g. CONSFREE). Same order as in type_cons.h. -enum class ConstraintHandlerCallbackType { - kSepaLp, // CONSSEPALP - kSepaSol, // CONSSEPASOL - kEnfoLp, // CONSENFOLP - // Unsupported: kEnfoRelax, // CONSENFORELAX - kEnfoPs, // CONSENFOPS - kConsCheck, // CONSCHECK - // Unsupported: kConsProp, // CONSPROP - // Unsupported: kConsPresol, // CONSPRESOL - // Unsupported: kConsResProp, // CONSRESPROP - kConsLock // CONSLOCK -}; - // Returns the "do nothing" callback result. Used to handle the edge case where // Enfo* callbacks need to return kFeasible instead of kDidNotRun. GScipCallbackResult DidNotRunCallbackResult( @@ -94,55 +76,6 @@ GScipCallbackResult DidNotRunCallbackResult( return GScipCallbackResult::kDidNotRun; } -// In callbacks, SCIP requires the first SCIP_RESULT in a priority list be -// returned when multiple results are applicable. This is a unified order of the -// priorities extracted from type_cons.h. The higher the result, the higher -// priority it is. -int ConstraintHandlerResultPriority( - const GScipCallbackResult result, - const ConstraintHandlerCallbackType callback_type) { - // In type_cons.h, callback results are consistently ordered across all - // constraint handler callback methods except that SCIP_SOLVELP (kSolveLp) - // takes higher priority than SCIP_BRANCHED (kBranched) in CONSENFOLP, and the - // reverse is true for CONSENFORELAX and CONSENFOLP. - switch (result) { - case GScipCallbackResult::kUnbounded: - return 14; - case GScipCallbackResult::kCutOff: - return 13; - case GScipCallbackResult::kSuccess: - return 12; - case GScipCallbackResult::kConstraintAdded: - return 11; - case GScipCallbackResult::kReducedDomain: - return 10; - case GScipCallbackResult::kSeparated: - return 9; - case GScipCallbackResult::kBranched: - return callback_type == ConstraintHandlerCallbackType::kEnfoLp ? 7 : 8; - case GScipCallbackResult::kSolveLp: - return callback_type == ConstraintHandlerCallbackType::kEnfoLp ? 8 : 7; - case GScipCallbackResult::kInfeasible: - return 6; - case GScipCallbackResult::kFeasible: - return 5; - case GScipCallbackResult::kNewRound: - return 4; - case GScipCallbackResult::kDidNotFind: - return 3; - case GScipCallbackResult::kDidNotRun: - return 2; - case GScipCallbackResult::kDelayed: - return 1; - case GScipCallbackResult::kDelayNode: - return 0; - default: - // kConstraintChanged, kFoundSolution, and kSuspend are not used in - // constraint handlers. - return -1; - } -} - constexpr GScipCallbackResult kMinPriority = GScipCallbackResult::kDelayNode; // Calls the callback function over a span of constraints, returning the highest @@ -207,7 +140,8 @@ absl::StatusOr ApplyCallback( return result; } -GScipCallbackStats GetCallbackStats(SCIP* scip) { +GScipCallbackStats GetCallbackStats(GScip* gscip) { + SCIP* scip = gscip->scip(); GScipCallbackStats stats; stats.num_processed_nodes = SCIPgetNNodes(scip); stats.num_processed_nodes_total = SCIPgetNTotalNodes(scip); @@ -223,6 +157,38 @@ GScipCallbackStats GetCallbackStats(SCIP* scip) { default: stats.current_node_id = stats.num_processed_nodes; } + stats.primal_bound = gscip->ScipInfUnclamp(SCIPgetPrimalbound(scip)); + stats.dual_bound = gscip->ScipInfUnclamp(SCIPgetDualbound(scip)); + const SCIP_STAGE stage = SCIPgetStage(scip); + switch (stage) { + case SCIP_STAGE_PRESOLVED: + case SCIP_STAGE_SOLVING: + case SCIP_STAGE_SOLVED: + stats.primal_simplex_iterations = SCIPgetNPrimalLPIterations(scip); + stats.dual_simplex_iterations = SCIPgetNDualLPIterations(scip); + stats.num_nodes_left = SCIPgetNNodesLeft(scip); + break; + default: + break; + } + // SCIP counts the focus node (the current node) as explored, but to be + // consistent with gurobi, we want to count it as open instead. In particular, + // for callbacks at the root, we want num_processed_nodes=0 + if (stats.num_processed_nodes > 0) { + stats.num_processed_nodes--; + stats.num_processed_nodes_total--; + stats.num_nodes_left++; + } + switch (stage) { + case SCIP_STAGE_SOLVING: + case SCIP_STAGE_SOLVED: + case SCIP_STAGE_EXITSOLVE: + stats.num_cuts_in_lp = SCIPgetNPoolCuts(scip); + break; + default: + break; + } + stats.num_solutions_found = SCIPgetNLimSolsFound(scip); return stats; } @@ -245,6 +211,62 @@ GScipConstraintOptions CallbackLazyConstraintOptions(const bool local, } // namespace +int ConstraintHandlerResultPriority( + const GScipCallbackResult result, + const ConstraintHandlerCallbackType callback_type) { + // In type_cons.h, callback results are consistently ordered across all + // constraint handler callback methods except that SCIP_SOLVELP (kSolveLp) + // takes higher priority than SCIP_BRANCHED (kBranched) in CONSENFOLP, and the + // reverse is true for CONSENFORELAX and CONSENFOLP. + switch (result) { + case GScipCallbackResult::kUnbounded: + return 14; + case GScipCallbackResult::kCutOff: + return 13; + case GScipCallbackResult::kSuccess: + return 12; + case GScipCallbackResult::kConstraintAdded: + return 11; + case GScipCallbackResult::kReducedDomain: + return 10; + case GScipCallbackResult::kSeparated: + return 9; + case GScipCallbackResult::kBranched: + return callback_type == ConstraintHandlerCallbackType::kEnfoLp ? 7 : 8; + case GScipCallbackResult::kSolveLp: + return callback_type == ConstraintHandlerCallbackType::kEnfoLp ? 8 : 7; + case GScipCallbackResult::kInfeasible: + return 6; + case GScipCallbackResult::kFeasible: + return 5; + case GScipCallbackResult::kNewRound: + return 4; + case GScipCallbackResult::kDidNotFind: + return 3; + case GScipCallbackResult::kDidNotRun: + return 2; + case GScipCallbackResult::kDelayed: + return 1; + case GScipCallbackResult::kDelayNode: + return 0; + default: + // kConstraintChanged, kFoundSolution, and kSuspend are not used in + // constraint handlers. + return -1; + } +} + +GScipCallbackResult MergeConstraintHandlerResults( + const GScipCallbackResult result1, const GScipCallbackResult result2, + const ConstraintHandlerCallbackType callback_type) { + const int priority1 = ConstraintHandlerResultPriority(result1, callback_type); + const int priority2 = ConstraintHandlerResultPriority(result2, callback_type); + if (priority2 > priority1) { + return result2; + } + return result1; +} + double GScipConstraintHandlerContext::VariableValue(SCIP_VAR* variable) const { return SCIPgetSolVal(gscip_->scip(), current_solution_, variable); } @@ -355,7 +377,7 @@ static SCIP_DECL_CONSENFOLP(EnforceLpC) { SCIP_CONSHDLRDATA* scip_handler_data = SCIPconshdlrGetData(conshdlr); operations_research::GScip* gscip = scip_handler_data->gscip; operations_research::GScipCallbackStats stats = - operations_research::GetCallbackStats(gscip->scip()); + operations_research::GetCallbackStats(gscip); operations_research::GScipConstraintHandlerContext context(gscip, &stats, conshdlr, nullptr); const bool solution_known_infeasible = static_cast(solinfeasible); @@ -380,7 +402,7 @@ static SCIP_DECL_CONSENFOPS(EnforcePseudoSolutionC) { SCIP_CONSHDLRDATA* scip_handler_data = SCIPconshdlrGetData(conshdlr); operations_research::GScip* gscip = scip_handler_data->gscip; operations_research::GScipCallbackStats stats = - operations_research::GetCallbackStats(gscip->scip()); + operations_research::GetCallbackStats(gscip); operations_research::GScipConstraintHandlerContext context(gscip, &stats, conshdlr, nullptr); const bool solution_known_infeasible = static_cast(solinfeasible); @@ -408,7 +430,7 @@ static SCIP_DECL_CONSCHECK(CheckFeasibilityC) { SCIP_CONSHDLRDATA* scip_handler_data = SCIPconshdlrGetData(conshdlr); operations_research::GScip* gscip = scip_handler_data->gscip; operations_research::GScipCallbackStats stats = - operations_research::GetCallbackStats(gscip->scip()); + operations_research::GetCallbackStats(gscip); operations_research::GScipConstraintHandlerContext context(gscip, &stats, conshdlr, sol); const bool check_integrality = static_cast(checkintegrality); @@ -437,7 +459,7 @@ static SCIP_DECL_CONSSEPALP(SeparateLpC) { SCIP_CONSHDLRDATA* scip_handler_data = SCIPconshdlrGetData(conshdlr); operations_research::GScip* gscip = scip_handler_data->gscip; operations_research::GScipCallbackStats stats = - operations_research::GetCallbackStats(gscip->scip()); + operations_research::GetCallbackStats(gscip); operations_research::GScipConstraintHandlerContext context(gscip, &stats, conshdlr, nullptr); GScipHandler* gscip_handler = scip_handler_data->gscip_handler.get(); @@ -460,7 +482,7 @@ static SCIP_DECL_CONSSEPASOL(SeparatePrimalSolutionC) { SCIP_CONSHDLRDATA* scip_handler_data = SCIPconshdlrGetData(conshdlr); operations_research::GScip* gscip = scip_handler_data->gscip; operations_research::GScipCallbackStats stats = - operations_research::GetCallbackStats(gscip->scip()); + operations_research::GetCallbackStats(gscip); operations_research::GScipConstraintHandlerContext context(gscip, &stats, conshdlr, sol); GScipHandler* gscip_handler = scip_handler_data->gscip_handler.get(); @@ -549,10 +571,10 @@ absl::Status RegisterConstraintHandler( return absl::OkStatus(); } -absl::Status AddCallbackConstraint(GScip* gscip, absl::string_view handler_name, - absl::string_view constraint_name, - void* constraint_data, - const GScipConstraintOptions& options) { +absl::StatusOr AddCallbackConstraint( + GScip* gscip, const std::string& handler_name, + const std::string& constraint_name, void* constraint_data, + const GScipConstraintOptions& options) { if (constraint_data == nullptr) { return absl::InvalidArgumentError( "Constraint data missing when adding a constraint handler callback"); @@ -567,18 +589,9 @@ absl::Status AddCallbackConstraint(GScip* gscip, absl::string_view handler_name, } SCIP_CONSDATA* consdata = new SCIP_CONSDATA; consdata->data = constraint_data; - SCIP_CONS* constraint = nullptr; - RETURN_IF_SCIP_ERROR(SCIPcreateCons( - scip, &constraint, constraint_name.data(), conshdlr, consdata, - options.initial, options.separate, options.enforce, options.check, - options.propagate, options.local, options.modifiable, options.dynamic, - options.removable, options.sticking_at_node)); - if (constraint == nullptr) { - return absl::InternalError("SCIP failed to create constraint"); - } - RETURN_IF_SCIP_ERROR(SCIPaddCons(scip, constraint)); - RETURN_IF_SCIP_ERROR(SCIPreleaseCons(scip, &constraint)); - return absl::OkStatus(); + + return gscip->AddConstraintForHandler(conshdlr, consdata, constraint_name, + options); } } // namespace internal diff --git a/ortools/gscip/gscip_constraint_handler.h b/ortools/gscip/gscip_constraint_handler.h index 695f27d431..cdfccedf4a 100644 --- a/ortools/gscip/gscip_constraint_handler.h +++ b/ortools/gscip/gscip_constraint_handler.h @@ -21,13 +21,11 @@ #include #include #include -#include #include #include #include "absl/status/status.h" #include "absl/status/statusor.h" -#include "absl/strings/string_view.h" #include "ortools/gscip/gscip.h" #include "ortools/gscip/gscip_callback_result.h" #include "scip/type_cons.h" @@ -140,10 +138,35 @@ struct GScipCallbackStats { // nodes before a restart), including the focus node. See SCIPgetNNodes(). int64_t num_processed_nodes = 0; + // The number of open nodes left in the current run. See SCIPgetNNodesLeft(). + int64_t num_nodes_left = 0; + // The total number of processed nodes in all runs, including the focus node. // If the solver restarts > 1 time, will be larger than // num_processed_nodes, otherwise is equal. See SCIPgetNTotalNodes(). int64_t num_processed_nodes_total = 0; + + // The global primal bound on the problem. See SCIPgetPrimalBound(). + double primal_bound = 0.0; + + // The global dual bound on the problem. See SCIPgetDualBound(). + double dual_bound = 0.0; + + // The number of pivots taken by the primal simplex method. See + // SCIPgetNPrimalLPIterations(). + int64_t primal_simplex_iterations = 0; + + // The number of pivots taken by the dual simplex method. See + // SCIPgetNDualLPIterations(). + int64_t dual_simplex_iterations = 0; + + // The number of feasible solutions found that were at least as good as the + // cutoff limit. See SCIPgetNLimSolsFound(). + int32_t num_solutions_found = 0; + + // The number of rows in the current LP that were added as cuts. See + // SCIPgetNPoolCuts(). + int32_t num_cuts_in_lp = 0; }; // Interface to the callback context and underlying problem. Supports adding @@ -281,10 +304,10 @@ class GScipConstraintHandler { // Adds a callback constraint to the model. That is, it attaches to the // constraint handler a constraint for the given constraint data. - absl::Status AddCallbackConstraint(GScip* gscip, - std::string_view constraint_name, - const ConstraintData* constraint_data, - const GScipConstraintOptions& options); + absl::StatusOr AddCallbackConstraint( + GScip* gscip, const std::string& constraint_name, + const ConstraintData* constraint_data, + const GScipConstraintOptions& options = DefaultGScipConstraintOptions()); // Callback function called at SCIP's CONSENFOLP. Must check if an LP solution // at a node is feasible, and if not, resolve the infeasibility if possible by @@ -418,6 +441,41 @@ class GScipConstraintHandler { GScipConstraintHandlerProperties properties_; }; +//////////////////////////////////////////////////////////////////////////////// +// Helpers for implementing callbacks +//////////////////////////////////////////////////////////////////////////////// + +// Enum with supported user-implementable callback functions in the SCIP +// constraint handler. Non-user-implementable functions are not included here +// (e.g. CONSFREE). Same order as in type_cons.h. +enum class ConstraintHandlerCallbackType { + kSepaLp, // CONSSEPALP + kSepaSol, // CONSSEPASOL + kEnfoLp, // CONSENFOLP + // Unsupported: kEnfoRelax, // CONSENFORELAX + kEnfoPs, // CONSENFOPS + kConsCheck, // CONSCHECK + // Unsupported: kConsProp, // CONSPROP + // Unsupported: kConsPresol, // CONSPRESOL + // Unsupported: kConsResProp, // CONSRESPROP + kConsLock // CONSLOCK +}; + +// In callbacks, SCIP requires the first SCIP_RESULT in a priority list be +// returned when multiple results are applicable. This is a unified order of the +// priorities extracted from type_cons.h. The higher the result, the higher +// priority it is. +// +// Larger number indicates higher priority in what should be returned. +int ConstraintHandlerResultPriority( + GScipCallbackResult result, ConstraintHandlerCallbackType callback_type); + +// Returns the GScipCallbackResult with larger priority, see +// ConstraintHandlerResultPriority(). +GScipCallbackResult MergeConstraintHandlerResults( + GScipCallbackResult result1, GScipCallbackResult result2, + ConstraintHandlerCallbackType callback_type); + namespace internal { // These classes implement a void pointer version of GScipConstraintHandler that @@ -495,10 +553,10 @@ absl::Status RegisterConstraintHandler( GScip* gscip, std::unique_ptr constraint_handler); -absl::Status AddCallbackConstraint(GScip* gscip, std::string_view handler_name, - std::string_view constraint_name, - void* constraint_data, - const GScipConstraintOptions& options); +absl::StatusOr AddCallbackConstraint( + GScip* gscip, const std::string& handler_name, + const std::string& constraint_name, void* constraint_data, + const GScipConstraintOptions& options); } // namespace internal @@ -513,8 +571,9 @@ absl::Status GScipConstraintHandler::Register(GScip* gscip) { } template -absl::Status GScipConstraintHandler::AddCallbackConstraint( - GScip* gscip, std::string_view constraint_name, +absl::StatusOr +GScipConstraintHandler::AddCallbackConstraint( + GScip* gscip, const std::string& constraint_name, const ConstraintData* constraint_data, const GScipConstraintOptions& options) { return internal::AddCallbackConstraint(