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

[CP-SAT] fix best bound callback non C++ code; more work on pseudo-costs

Browse Source
This commit is contained in:
Laurent Perron
2024-05-02 14:04:47 +02:00
parent 2758b6cb5f
commit e5d09e1e12
11 changed files with 194 additions and 110 deletions
Download Patch File Download Diff File Expand all files Collapse all files

21
ortools/java/com/google/ortools/sat/CpSolver.java
View File

@@ -36,7 +36,7 @@ public final class CpSolver {
/** Solves the given model, and returns the solve status. */
public CpSolverStatus solve(CpModel model) {
return solveWithSolutionCallback(model, null);
return solve(model, null);
}
/**
@@ -97,9 +97,9 @@ public final class CpSolver {
public CpSolverStatus searchAllSolutions(CpModel model, CpSolverSolutionCallback cb) {
boolean oldValue = solveParameters.getEnumerateAllSolutions();
solveParameters.setEnumerateAllSolutions(true);
solve(model, cb);
CpSolverStatus status = solve(model, cb);
solveParameters.setEnumerateAllSolutions(oldValue);
return solveResponse.getStatus();
return status;
}
private synchronized void createSolveWrapper() {
@@ -189,6 +189,21 @@ public final class CpSolver {
this.logCallback = cb;
}
/** Clears the log callback. */
public void clearLogCallback() {
this.logCallback = null;
}
/** Sets the best bound callback for the solver. */
public void setBestBoundCallback(Consumer<Double> cb) {
this.bestBoundCallback = cb;
}
/** Clears the best bound callback. */
public void clearBestBoundCallback() {
this.bestBoundCallback = null;
}
/** Returns some statistics on the solution found as a string. */
public String responseStats() {
return CpSatHelper.solverResponseStats(solveResponse);

3
ortools/sat/cp_model_solver.h
View File

@@ -106,8 +106,7 @@ std::function<void(Model*)> NewFeasibleSolutionLogCallback(
/** Creates a callbacks that will be called on each new best objective bound
* found.
*
* Note that adding a callback is not free since some computation will be done
* before this is called, and it will only be called on optimization models.
* Note that this function is called before the update takes place.
*/
std::function<void(Model*)> NewBestBoundCallback(
const std::function<void(double)>& callback);

10
ortools/sat/csharp/CpSolver.cs
View File

@@ -145,11 +145,21 @@ public class CpSolver
log_callback_ = new LogCallbackDelegate(del);
}
public void ClearLogCallback()
{
log_callback_ = null;
}
public void SetBestBoundCallback(DoubleToVoidDelegate del)
{
best_bound_callback_ = new BestBoundCallbackDelegate(del);
}
public void ClearBestBoundCallback()
{
best_bound_callback_ = null;
}
public CpSolverResponse Response
{
get {

92
ortools/sat/integer_search.cc
View File

@@ -233,20 +233,14 @@ std::function<BooleanOrIntegerLiteral()> BoolPseudoCostHeuristic(Model* model) {
std::function<BooleanOrIntegerLiteral()> LpPseudoCostHeuristic(Model* model) {
auto* lp_values = model->GetOrCreate<ModelLpValues>();
auto* integer_trail = model->GetOrCreate<IntegerTrail>();
auto* pseudo_costs = model->GetOrCreate<PseudoCosts>();
auto* encoder = model->GetOrCreate<IntegerEncoder>();
return [lp_values, pseudo_costs, integer_trail, encoder, model]() {
return [lp_values, pseudo_costs]() {
double best_score = 0.0;
BooleanOrIntegerLiteral decision;
for (IntegerVariable var(0); var < lp_values->size(); var += 2) {
const IntegerValue lb = integer_trail->LowerBound(var);
const IntegerValue ub = integer_trail->UpperBound(var);
if (lb == ub) continue;
const double lp_value = (*lp_values)[var];
const bool is_reliable = pseudo_costs->LpReliability(var) >= 4;
const bool is_integer = std::abs(lp_value - std::round(lp_value)) < 1e-6;
const PseudoCosts::BranchingInfo info =
pseudo_costs->EvaluateVar(var, *lp_values);
if (info.is_fixed) continue;
// When not reliable, we skip integer.
//
@@ -254,35 +248,18 @@ std::function<BooleanOrIntegerLiteral()> LpPseudoCostHeuristic(Model* model) {
// TODO(user): do not branch on integer lp? however it seems better to
// do that !? Maybe this is because if it has a high pseudo cost
// average, it is good anyway?
if (!is_reliable && is_integer) continue;
// There are some corner cases if we are at the bound. Note that it is
// important to be in sync with the SplitAroundLpValue() below.
double down_fractionality = lp_value - std::floor(lp_value);
IntegerValue down_target = IntegerValue(std::floor(lp_value));
if (lp_value >= ToDouble(ub)) {
down_fractionality = 1.0;
down_target = ub - 1;
} else if (lp_value <= ToDouble(lb)) {
down_fractionality = 0.0;
down_target = lb;
}
const auto [down_score, up_score] =
pseudo_costs->LpPseudoCost(var, down_fractionality);
const double score = pseudo_costs->CombineScores(down_score, up_score);
if (!info.is_reliable && info.is_integer) continue;
// We delay to subsequent heuristic if the score is 0.0.
if (score > best_score) {
best_score = score;
if (info.score > best_score) {
best_score = info.score;
// This direction works better than the inverse in the benchs. But
// always branching up seems even better. TODO(user): investigate.
if (down_score > up_score) {
decision = BooleanOrIntegerLiteral(
IntegerLiteral::LowerOrEqual(var, down_target));
if (info.down_score > info.up_score) {
decision = BooleanOrIntegerLiteral(info.down_branch);
} else {
decision = BooleanOrIntegerLiteral(
IntegerLiteral::GreaterOrEqual(var, down_target + 1));
decision = BooleanOrIntegerLiteral(info.down_branch.Negated());
}
}
}
@@ -1380,6 +1357,31 @@ bool IntegerSearchHelper::BeforeTakingDecision() {
return true;
}
LiteralIndex IntegerSearchHelper::GetDecisionLiteral(
const BooleanOrIntegerLiteral& decision) {
DCHECK(decision.HasValue());
// Convert integer decision to literal one if needed.
//
// TODO(user): Ideally it would be cool to delay the creation even more
// until we have a conflict with these decisions, but it is currently
// hard to do so.
const Literal literal =
decision.boolean_literal_index != kNoLiteralIndex
? Literal(decision.boolean_literal_index)
: encoder_->GetOrCreateAssociatedLiteral(decision.integer_literal);
if (sat_solver_->Assignment().LiteralIsAssigned(literal)) {
// TODO(user): It would be nicer if this can never happen. For now, it
// does because of the Propagate() not reaching the fixed point as
// mentioned in a TODO above. As a work-around, we display a message
// but do not crash and recall the decision heuristic.
VLOG(1) << "Trying to take a decision that is already assigned!"
<< " Fix this. Continuing for now...";
return kNoLiteralIndex;
}
return literal.Index();
}
bool IntegerSearchHelper::GetDecision(
const std::function<BooleanOrIntegerLiteral()>& f, LiteralIndex* decision) {
*decision = kNoLiteralIndex;
@@ -1411,28 +1413,8 @@ bool IntegerSearchHelper::GetDecision(
}
if (!new_decision.HasValue()) break;
// Convert integer decision to literal one if needed.
//
// TODO(user): Ideally it would be cool to delay the creation even more
// until we have a conflict with these decisions, but it is currently
// hard to do so.
if (new_decision.boolean_literal_index != kNoLiteralIndex) {
*decision = new_decision.boolean_literal_index;
} else {
*decision =
encoder_->GetOrCreateAssociatedLiteral(new_decision.integer_literal)
.Index();
}
if (sat_solver_->Assignment().LiteralIsAssigned(Literal(*decision))) {
// TODO(user): It would be nicer if this can never happen. For now, it
// does because of the Propagate() not reaching the fixed point as
// mentioned in a TODO above. As a work-around, we display a message
// but do not crash and recall the decision heuristic.
VLOG(1) << "Trying to take a decision that is already assigned!"
<< " Fix this. Continuing for now...";
continue;
}
break;
*decision = GetDecisionLiteral(new_decision);
if (*decision != kNoLiteralIndex) break;
}
return true;
}

4
ortools/sat/integer_search.h
View File

@@ -281,6 +281,10 @@ class IntegerSearchHelper {
bool GetDecision(const std::function<BooleanOrIntegerLiteral()>& f,
LiteralIndex* decision);
// Inner function used by GetDecision().
// It will create a new associated literal if needed.
LiteralIndex GetDecisionLiteral(const BooleanOrIntegerLiteral& decision);
// Functions passed to GetDecision() might call this to notify a conflict
// was detected.
void NotifyThatConflictWasFoundDuringGetDecision() {

12
ortools/sat/linear_programming_constraint.cc
View File

@@ -557,6 +557,11 @@ void LinearProgrammingConstraint::SetLevel(int level) {
if (lp_solution_is_set_ && level < lp_solution_level_) {
lp_solution_is_set_ = false;
}
if (level < previous_level_) {
lp_at_optimal_ = false;
lp_objective_lower_bound_ = -std::numeric_limits<double>::infinity();
}
previous_level_ = level;
// Special case for level zero, we "reload" any previously known optimal
// solution from that level.
@@ -697,6 +702,12 @@ bool LinearProgrammingConstraint::SolveLp() {
<< objective_definition_->ScaleObjective(
simplex_.GetObjectiveValue());
if (simplex_.GetProblemStatus() == glop::ProblemStatus::OPTIMAL ||
simplex_.GetProblemStatus() == glop::ProblemStatus::DUAL_FEASIBLE) {
lp_objective_lower_bound_ = simplex_.GetObjectiveValue();
}
lp_at_optimal_ = simplex_.GetProblemStatus() == glop::ProblemStatus::OPTIMAL;
if (simplex_.GetProblemStatus() == glop::ProblemStatus::OPTIMAL) {
lp_solution_is_set_ = true;
lp_solution_level_ = trail_->CurrentDecisionLevel();
@@ -799,7 +810,6 @@ bool LinearProgrammingConstraint::AnalyzeLp() {
if (simplex_.GetProblemStatus() == glop::ProblemStatus::OPTIMAL) {
CHECK(lp_solution_is_set_);
lp_objective_ = simplex_.GetObjectiveValue();
lp_solution_is_integer_ = true;
const int num_vars = integer_variables_.size();
for (int i = 0; i < num_vars; i++) {

15
ortools/sat/linear_programming_constraint.h
View File

@@ -164,11 +164,17 @@ class LinearProgrammingConstraint : public PropagatorInterface,
// Note that this solution is always an OPTIMAL solution of an LP above or
// at the current decision level. We "erase" it when we backtrack over it.
bool HasSolution() const { return lp_solution_is_set_; }
double SolutionObjectiveValue() const { return lp_objective_; }
double GetSolutionValue(IntegerVariable variable) const;
double GetSolutionReducedCost(IntegerVariable variable) const;
bool SolutionIsInteger() const { return lp_solution_is_integer_; }
// Returns a valid lp lower bound for the current branch, and indicates if
// the latest LP solve in that branch was solved to optimality or not.
// During normal operation, we will always propagate the LP, so this should
// not refer to an optimality status lower in that branch.
bool AtOptimal() const { return lp_at_optimal_; }
double ObjectiveLpLowerBound() const { return lp_objective_lower_bound_; }
// PropagatorInterface API.
bool Propagate() override;
bool IncrementalPropagate(const std::vector<int>& watch_indices) override;
@@ -526,10 +532,15 @@ class LinearProgrammingConstraint : public PropagatorInterface,
int lp_solution_level_ = 0;
bool lp_solution_is_set_ = false;
bool lp_solution_is_integer_ = false;
double lp_objective_;
std::vector<double> lp_solution_;
std::vector<double> lp_reduced_cost_;
// Last objective lower bound found by the LP Solver.
// We erase this on backtrack.
int previous_level_ = 0;
bool lp_at_optimal_ = false;
double lp_objective_lower_bound_;
// If non-empty, this is the last known optimal lp solution at root-node. If
// the variable bounds changed, or cuts where added, it is possible that this
// solution is no longer optimal though.

1
ortools/sat/linear_propagation.h
View File

@@ -18,6 +18,7 @@
#include <deque>
#include <functional>
#include <limits>
#include <ostream>
#include <string>
#include <utility>

111
ortools/sat/pseudo_costs.cc
View File

@@ -14,7 +14,9 @@
#include "ortools/sat/pseudo_costs.h"
#include <algorithm>
#include <cmath>
#include <cstdint>
#include <cstdlib>
#include <limits>
#include <string>
#include <tuple>
@@ -85,30 +87,70 @@ PseudoCosts::ObjectiveInfo PseudoCosts::GetCurrentObjectiveInfo() {
result.lp_bound = 0.0;
result.lp_at_optimal = true;
for (const auto* lp : *lps_) {
if (!lp->HasSolution()) result.lp_at_optimal = false;
result.lp_bound += lp->SolutionObjectiveValue();
if (!lp->AtOptimal()) result.lp_at_optimal = false;
result.lp_bound += lp->ObjectiveLpLowerBound();
}
return result;
}
void PseudoCosts::BeforeTakingDecision(Literal decision) {
if (objective_var_ == kNoIntegerVariable) return;
bool PseudoCosts::SaveLpInfo() {
saved_info_ = GetCurrentObjectiveInfo();
bound_changes_ = GetBoundChanges(decision);
return saved_info_.lp_at_optimal;
}
std::pair<double, double> PseudoCosts::LpPseudoCost(
IntegerVariable var, double down_fractionality) const {
const int max_index = std::max(var.value(), NegationOf(var).value());
if (max_index >= average_unit_objective_increase_.size()) return {0.0, 0.0};
void PseudoCosts::SaveBoundChanges(Literal decision,
absl::Span<const double> lp_values) {
bound_changes_ = GetBoundChanges(decision, lp_values);
}
const double up_fractionality = 1.0 - down_fractionality;
const double up_branch =
up_fractionality * average_unit_objective_increase_[var].CurrentAverage();
const double down_branch =
down_fractionality *
average_unit_objective_increase_[NegationOf(var)].CurrentAverage();
return {down_branch, up_branch};
void PseudoCosts::BeforeTakingDecision(Literal decision) {
if (objective_var_ == kNoIntegerVariable) return;
SaveLpInfo();
SaveBoundChanges(decision, *lp_values_);
}
PseudoCosts::BranchingInfo PseudoCosts::EvaluateVar(
IntegerVariable var, absl::Span<const double> lp_values) {
DCHECK_NE(var, kNoIntegerVariable);
BranchingInfo result;
const IntegerValue lb = integer_trail_->LowerBound(var);
const IntegerValue ub = integer_trail_->UpperBound(var);
if (lb == ub) {
result.is_fixed = true;
return result;
}
const double lp_value = lp_values[var.value()];
double down_fractionality = lp_value - std::floor(lp_value);
IntegerValue down_target = IntegerValue(std::floor(lp_value));
if (lp_value >= ToDouble(ub)) {
down_fractionality = 1.0;
down_target = ub - 1;
} else if (lp_value <= ToDouble(lb)) {
down_fractionality = 0.0;
down_target = lb;
}
result.is_integer = std::abs(lp_value - std::round(lp_value)) < 1e-6;
result.down_fractionality = down_fractionality;
result.down_branch = IntegerLiteral::LowerOrEqual(var, down_target);
const int max_index = std::max(var.value(), NegationOf(var).value());
if (max_index < average_unit_objective_increase_.size()) {
result.down_score =
down_fractionality *
average_unit_objective_increase_[NegationOf(var)].CurrentAverage();
result.up_score = (1.0 - down_fractionality) *
average_unit_objective_increase_[var].CurrentAverage();
result.score = CombineScores(result.down_score, result.up_score);
const int reliablitity = std::min(
average_unit_objective_increase_[var].NumRecords(),
average_unit_objective_increase_[NegationOf(var)].NumRecords());
result.is_reliable = reliablitity >= 4;
}
return result;
}
void PseudoCosts::UpdateBoolPseudoCosts(absl::Span<const Literal> reason,
@@ -136,13 +178,19 @@ double PseudoCosts::BoolPseudoCost(Literal lit, double lp_value) const {
return CombineScores(down_branch, up_branch);
}
int PseudoCosts::LpReliability(IntegerVariable var) const {
const int max_index = std::max(var.value(), NegationOf(var).value());
if (max_index >= average_unit_objective_increase_.size()) return 0;
double PseudoCosts::ObjectiveIncrease(bool conflict) {
const ObjectiveInfo new_info = GetCurrentObjectiveInfo();
const double obj_lp_diff =
std::max(0.0, new_info.lp_bound - saved_info_.lp_bound);
const IntegerValue obj_int_diff = new_info.lb - saved_info_.lb;
return std::min(
average_unit_objective_increase_[var].NumRecords(),
average_unit_objective_increase_[NegationOf(var)].NumRecords());
double obj_increase =
obj_lp_diff > 0.0 ? obj_lp_diff : ToDouble(obj_int_diff);
if (conflict) {
// We count a conflict as a max increase + 1.0
obj_increase = ToDouble(saved_info_.ub) - ToDouble(saved_info_.lb) + 1.0;
}
return obj_increase;
}
void PseudoCosts::AfterTakingDecision(bool conflict) {
@@ -156,23 +204,14 @@ void PseudoCosts::AfterTakingDecision(bool conflict) {
// just be the "artificial" continuing of the current lp solve that create the
// increase.
if (saved_info_.lp_at_optimal) {
// Compute the increase in objective.
const double obj_lp_diff =
std::max(0.0, new_info.lp_bound - saved_info_.lp_bound);
const IntegerValue obj_int_diff = new_info.lb - saved_info_.lb;
double obj_diff = obj_lp_diff > 0.0 ? obj_lp_diff : ToDouble(obj_int_diff);
if (conflict) {
// We count a conflict as a max increase + 1.0
obj_diff = ToDouble(saved_info_.ub) - ToDouble(saved_info_.lb) + 1.0;
}
// Update the average unit increases.
const double obj_increase = ObjectiveIncrease(conflict);
for (const auto [var, lb_change, lp_increase] : bound_changes_) {
if (lp_increase < 1e-6) continue;
if (var >= average_unit_objective_increase_.size()) {
average_unit_objective_increase_.resize(var + 1);
}
average_unit_objective_increase_[var].AddData(obj_diff / lp_increase);
average_unit_objective_increase_[var].AddData(obj_increase / lp_increase);
}
}
@@ -243,16 +282,16 @@ IntegerVariable PseudoCosts::GetBestDecisionVar() {
}
std::vector<PseudoCosts::VariableBoundChange> PseudoCosts::GetBoundChanges(
Literal decision) {
Literal decision, absl::Span<const double> lp_values) {
std::vector<PseudoCosts::VariableBoundChange> bound_changes;
for (const IntegerLiteral l : encoder_->GetIntegerLiterals(decision)) {
PseudoCosts::VariableBoundChange entry;
entry.var = l.var;
entry.lower_bound_change = l.bound - integer_trail_->LowerBound(l.var);
if (l.var < lp_values_->size()) {
if (l.var < lp_values.size()) {
entry.lp_increase =
std::max(0.0, ToDouble(l.bound) - (*lp_values_)[l.var]);
std::max(0.0, ToDouble(l.bound) - lp_values[l.var.value()]);
}
bound_changes.push_back(entry);
}

32
ortools/sat/pseudo_costs.h
View File

@@ -44,6 +44,12 @@ class PseudoCosts {
// BeforeTakingDecision().
void AfterTakingDecision(bool conflict = false);
// Advanced usage. Internal functions used by Before/AfterTakingDecision(),
// that are exposed for strong branching.
double ObjectiveIncrease(bool conflict);
bool SaveLpInfo();
void SaveBoundChanges(Literal decision, absl::Span<const double> lp_values);
// Returns the variable with best reliable pseudo cost that is not fixed.
IntegerVariable GetBestDecisionVar();
@@ -60,15 +66,24 @@ class PseudoCosts {
return pseudo_costs_[var].NumRecords();
}
// Alternative pseudo-costs. This relies on the LP more heavily and is more
// in line with what a MIP solver would do. Return the (down, up) costs which
// can be combined with CombineScores();
// Combines the score of the two branch into one score.
double CombineScores(double down_branch, double up_branch) const;
std::pair<double, double> LpPseudoCost(IntegerVariable var,
double down_fractionality) const;
// Returns the pseudo cost "reliability".
int LpReliability(IntegerVariable var) const;
// Alternative pseudo-cost. This relies on the LP more heavily and is more in
// line with what a MIP solver would do. Returns all the info about taking a
// branch around the current lp_value of var.
struct BranchingInfo {
bool is_fixed = false;
bool is_reliable = false;
bool is_integer = false;
double down_fractionality = 0.0;
double score = 0.0;
double down_score = 0.0;
double up_score = 0.0;
IntegerLiteral down_branch;
};
BranchingInfo EvaluateVar(IntegerVariable var,
absl::Span<const double> lp_values);
// Experimental alternative pseudo cost based on the explanation for bound
// increases.
@@ -83,7 +98,8 @@ class PseudoCosts {
IntegerValue lower_bound_change = IntegerValue(0);
double lp_increase = 0.0;
};
std::vector<VariableBoundChange> GetBoundChanges(Literal decision);
std::vector<VariableBoundChange> GetBoundChanges(
Literal decision, absl::Span<const double> lp_values);
private:
// Returns the current objective info.

3
ortools/sat/synchronization.h
View File

@@ -272,9 +272,6 @@ class SharedResponseManager {
// Adds a callback that will be called on each new best objective bound
// found. Returns its id so it can be unregistered if needed.
//
// Note that adding a callback is not free since some computation will be done
// before this is called.
//
// Note that currently the class is waiting for the callback to finish before
// accepting any new updates. That could be changed if needed.
int AddBestBoundCallback(std::function<void(double)> callback);