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

dynamic tightening of the LP relaxation in the CP-SAT solver

Browse Source
This commit is contained in:
Laurent Perron
2019-04-05 17:00:40 +02:00
parent e685805fc8
commit 8ce29056b5
7 changed files with 228 additions and 16 deletions
Download Patch File Download Diff File Expand all files Collapse all files

5
ortools/sat/cp_model_presolve.cc
View File

@@ -1287,6 +1287,11 @@ bool PresolveLinear(ConstraintProto* ct, PresolveContext* context) {
// list.
//
// This operation is similar to coefficient strengthening in the MIP world.
//
// TODO(user): If the constraint is boxed, split it in two if enforcement
// literals can be extracted this way. This would just be better for the CP
// propagation, and for the LP it will improve the relaxation at the cost of
// an extra row in the matrix, but this should still be better.
void ExtractEnforcementLiteralFromLinearConstraint(ConstraintProto* ct,
PresolveContext* context) {
if (context->ModelIsUnsat()) return;

1
ortools/sat/integer.cc
View File

@@ -1166,6 +1166,7 @@ bool IntegerTrail::EnqueueInternal(
// Special case for level zero.
if (integer_search_levels_.empty()) {
++num_level_zero_enqueues_;
vars_[i_lit.var].current_bound = i_lit.bound;
integer_trail_[i_lit.var.value()].bound = i_lit.bound;

4
ortools/sat/integer.h
View File

@@ -692,6 +692,9 @@ class IntegerTrail : public SatPropagator {
// doesn't change since we directly update the "fixed" bounds.
int64 num_enqueues() const { return num_enqueues_; }
// Same as num_enqueues but only count the level zero changes.
int64 num_level_zero_enqueues() const { return num_level_zero_enqueues_; }
// All the registered bitsets will be set to one each time a LbVar is
// modified. It is up to the client to clear it if it wants to be notified
// with the newly modified variables.
@@ -907,6 +910,7 @@ class IntegerTrail : public SatPropagator {
std::vector<int> boolean_trail_index_to_integer_one_;
int64 num_enqueues_ = 0;
int64 num_level_zero_enqueues_ = 0;
std::vector<SparseBitset<IntegerVariable>*> watchers_;
std::vector<ReversibleInterface*> reversible_classes_;

5
ortools/sat/linear_constraint.h
View File

@@ -58,10 +58,11 @@ struct LinearConstraint {
for (int i = 0; i < vars.size(); ++i) {
const IntegerValue coeff =
VariableIsPositive(vars[i]) ? coeffs[i] : -coeffs[i];
absl::StrAppend(&result, coeff.value(), "*X", vars[i].value() / 2, " ");
absl::StrAppend(&result, i > 0 ? " " : "", coeff.value(), "*X",
vars[i].value() / 2);
}
if (ub.value() < kMaxIntegerValue) {
absl::StrAppend(&result, "<= ", ub.value());
absl::StrAppend(&result, " <= ", ub.value());
}
return result;
}

169
ortools/sat/linear_constraint_manager.cc
View File

@@ -57,6 +57,15 @@ LinearConstraintManager::~LinearConstraintManager() {
if (num_merged_constraints_ > 0) {
VLOG(2) << "num_merged_constraints: " << num_merged_constraints_;
}
if (num_shortened_constraints_ > 0) {
VLOG(2) << "num_shortened_constraints: " << num_shortened_constraints_;
}
if (num_splitted_constraints_ > 0) {
VLOG(2) << "num_splitted_constraints: " << num_splitted_constraints_;
}
if (num_coeff_strenghtening_ > 0) {
VLOG(2) << "num_coeff_strenghtening: " << num_coeff_strenghtening_;
}
for (const auto entry : type_to_num_cuts_) {
VLOG(1) << "Added " << entry.second << " cuts of type '" << entry.first
<< "'.";
@@ -86,6 +95,7 @@ void LinearConstraintManager::RemoveMarkedConstraints() {
// we regenerate identical cuts for some reason.
//
// TODO(user): Also compute GCD and divide by it.
// TODO(user): Call SimplifyConstraint() too.
void LinearConstraintManager::Add(const LinearConstraint& ct) {
LinearConstraint canonicalized = ct;
const Terms terms = CanonicalizeConstraintAndGetTerms(&canonicalized);
@@ -109,7 +119,7 @@ void LinearConstraintManager::Add(const LinearConstraint& ct) {
++num_merged_constraints_;
} else {
constraint_l2_norms_.push_back(ComputeL2Norm(canonicalized));
equiv_constraints_[terms] = constraints_.size();
equiv_constraints_[terms] = ConstraintIndex(constraints_.size());
constraint_is_in_lp_.push_back(false);
constraint_is_cut_.push_back(false);
constraint_inactive_count_.push_back(0);
@@ -118,6 +128,21 @@ void LinearConstraintManager::Add(const LinearConstraint& ct) {
}
}
void LinearConstraintManager::RemoveConstraintFromEquivTable(
ConstraintIndex ct_index) {
const Terms terms =
CanonicalizeConstraintAndGetTerms(&constraints_[ct_index]);
equiv_constraints_.erase(terms);
}
void LinearConstraintManager::NotifyConstraintChanged(
ConstraintIndex ct_index) {
const Terms terms =
CanonicalizeConstraintAndGetTerms(&constraints_[ct_index]);
equiv_constraints_[terms] = ct_index;
constraint_l2_norms_[ct_index] = ComputeL2Norm(constraints_[ct_index]);
}
// Same as Add(), but logs some information about the newly added constraint.
// Cuts are also handled slightly differently than normal constraints.
void LinearConstraintManager::AddCut(
@@ -144,6 +169,140 @@ void LinearConstraintManager::AddCut(
<< " violation=" << violation << " eff=" << violation / l2_norm;
}
void LinearConstraintManager::SimplifyConstraint(ConstraintIndex ct_index) {
LinearConstraint& ct = constraints_[ct_index];
IntegerValue min_sum(0);
IntegerValue max_sum(0);
IntegerValue max_magnitude(0);
int new_size = 0;
const int num_terms = ct.vars.size();
for (int i = 0; i < num_terms; ++i) {
const IntegerVariable var = ct.vars[i];
const IntegerValue coeff = ct.coeffs[i];
const IntegerValue lb = integer_trail_.LevelZeroLowerBound(var);
const IntegerValue ub = integer_trail_.LevelZeroUpperBound(var);
// For now we do not change ct, but just compute its new_size if we where
// to remove a fixed term.
if (lb == ub) continue;
++new_size;
max_magnitude = std::max(max_magnitude, IntTypeAbs(coeff));
if (coeff > 0.0) {
min_sum += coeff * lb;
max_sum += coeff * ub;
} else {
min_sum += coeff * ub;
max_sum += coeff * lb;
}
}
// Shorten the constraint if needed.
if (new_size < num_terms) {
++num_shortened_constraints_;
RemoveConstraintFromEquivTable(ct_index);
new_size = 0;
for (int i = 0; i < num_terms; ++i) {
const IntegerVariable var = ct.vars[i];
const IntegerValue coeff = ct.coeffs[i];
const IntegerValue lb = integer_trail_.LevelZeroLowerBound(var);
const IntegerValue ub = integer_trail_.LevelZeroUpperBound(var);
if (lb == ub) {
const IntegerValue rhs_adjust = lb * coeff;
if (ct.lb > kMinIntegerValue) ct.lb -= rhs_adjust;
if (ct.ub < kMaxIntegerValue) ct.ub -= rhs_adjust;
continue;
}
ct.vars[new_size] = var;
ct.coeffs[new_size] = coeff;
++new_size;
}
ct.vars.resize(new_size);
ct.coeffs.resize(new_size);
NotifyConstraintChanged(ct_index);
}
// Relax the bound if needed, note that this doesn't require a change to
// the equiv map.
if (min_sum >= ct.lb) ct.lb = kMinIntegerValue;
if (max_sum <= ct.ub) ct.ub = kMaxIntegerValue;
// Clear constraints that are always true.
// We rely on the deletion code to remove them eventually.
if (ct.lb == kMinIntegerValue && ct.ub == kMaxIntegerValue) {
RemoveConstraintFromEquivTable(ct_index);
ct.vars.clear();
ct.coeffs.clear();
constraint_l2_norms_[ct_index] = 0.0;
return;
}
// TODO(user): Split constraint in two if it is boxed and there is possible
// reduction?
//
// TODO(user): Make sure there cannot be any overflow. They should't, but
// I am not sure all the generated cuts are safe regarding min/max sum
// computation. We should check this.
if (ct.ub != kMaxIntegerValue && max_magnitude > max_sum - ct.ub) {
if (ct.lb != kMinIntegerValue) {
++num_splitted_constraints_;
} else {
++num_coeff_strenghtening_;
RemoveConstraintFromEquivTable(ct_index);
const int num_terms = ct.vars.size();
const IntegerValue target = max_sum - ct.ub;
for (int i = 0; i < num_terms; ++i) {
const IntegerValue coeff = ct.coeffs[i];
if (coeff > target) {
const IntegerVariable var = ct.vars[i];
const IntegerValue ub = integer_trail_.LevelZeroUpperBound(var);
ct.coeffs[i] = target;
ct.ub -= (coeff - target) * ub;
} else if (coeff < -target) {
const IntegerVariable var = ct.vars[i];
const IntegerValue lb = integer_trail_.LevelZeroLowerBound(var);
ct.coeffs[i] = -target;
ct.ub += (-target - coeff) * lb;
}
}
NotifyConstraintChanged(ct_index);
}
}
if (ct.lb != kMinIntegerValue && max_magnitude > ct.lb - min_sum) {
if (ct.ub != kMaxIntegerValue) {
++num_splitted_constraints_;
} else {
++num_coeff_strenghtening_;
RemoveConstraintFromEquivTable(ct_index);
const int num_terms = ct.vars.size();
const IntegerValue target = ct.lb - min_sum;
for (int i = 0; i < num_terms; ++i) {
const IntegerValue coeff = ct.coeffs[i];
if (coeff > target) {
const IntegerVariable var = ct.vars[i];
const IntegerValue lb = integer_trail_.LevelZeroLowerBound(var);
ct.coeffs[i] = target;
ct.lb -= (coeff - target) * lb;
} else if (coeff < -target) {
const IntegerVariable var = ct.vars[i];
const IntegerValue ub = integer_trail_.LevelZeroUpperBound(var);
ct.coeffs[i] = -target;
ct.lb += (-target - coeff) * ub;
}
}
NotifyConstraintChanged(ct_index);
}
}
}
bool LinearConstraintManager::ChangeLp(
const gtl::ITIVector<IntegerVariable, double>& lp_solution) {
const int old_num_constraints = lp_constraints_.size();
@@ -152,12 +311,19 @@ bool LinearConstraintManager::ChangeLp(
std::vector<double> new_constraints_efficacies;
std::vector<double> new_constraints_orthogonalities;
const bool simplify_constraints =
integer_trail_.num_level_zero_enqueues() > last_simplification_timestamp_;
last_simplification_timestamp_ = integer_trail_.num_level_zero_enqueues();
// We keep any constraints that is already present, and otherwise, we add the
// ones that are currently not satisfied by at least "tolerance".
const double tolerance = 1e-6;
for (ConstraintIndex i(0); i < constraints_.size(); ++i) {
if (constraint_permanently_removed_[i]) continue;
// Inprocessing of the constraint.
if (simplify_constraints) SimplifyConstraint(i);
// TODO(user,user): Use constraint status from GLOP for constraints
// already in LP.
const double activity = ComputeActivity(constraints_[i], lp_solution);
@@ -179,6 +345,7 @@ bool LinearConstraintManager::ChangeLp(
new_constraints_orthogonalities.push_back(1.0);
}
}
// Remove constraints in a batch to avoid changing the LP too frequently.
if (constraints_removal_list_.size() >=
sat_parameters_.constraint_removal_batch_size()) {

38
ortools/sat/linear_constraint_manager.h
View File

@@ -39,7 +39,9 @@ namespace sat {
// which constraint should go into the current LP.
class LinearConstraintManager {
public:
LinearConstraintManager() {}
explicit LinearConstraintManager(Model* model)
: sat_parameters_(*model->GetOrCreate<SatParameters>()),
integer_trail_(*model->GetOrCreate<IntegerTrail>()) {}
~LinearConstraintManager();
// Add a new constraint to the manager. Note that we canonicalize constraints
@@ -75,19 +77,41 @@ class LinearConstraintManager {
return lp_constraints_;
}
void SetParameters(const SatParameters& params) { sat_parameters_ = params; }
int num_cuts() const { return num_cuts_; }
int64 num_shortened_constraints() const { return num_shortened_constraints_; }
int64 num_coeff_strenghtening() const { return num_coeff_strenghtening_; }
// Visible for testing.
//
// Apply basic inprocessing simplification rules:
// - remove fixed variable
// - reduce large coefficient (i.e. coeff strenghtenning or big-M reduction).
// This uses level-zero bounds.
void SimplifyConstraint(ConstraintIndex index);
private:
// Removes the marked constraints from the LP.
void RemoveMarkedConstraints();
SatParameters sat_parameters_;
// Important: before modifying a constraint in equiv_constraints_, one
// must call RemoveConstraintFromEquivTable() and when done, call
// NotifyConstraintChanged().
//
// TODO(user): That might be a bit slow and brittle. Redesign the equiv map
// so that we compare actual constraints (so we never merge different ones),
// and we don't need to keep duplicate of the constraints in memory.
void RemoveConstraintFromEquivTable(ConstraintIndex ct_index);
void NotifyConstraintChanged(ConstraintIndex ct_index);
const SatParameters& sat_parameters_;
const IntegerTrail& integer_trail_;
// Set at true by Add() and at false by ChangeLp().
bool some_lp_constraint_bounds_changed_ = false;
// Optimization to avoid calling SimplifyConstraint() when not needed.
int64 last_simplification_timestamp_ = -1;
// TODO(user): Merge all the constraint related info in a struct and store
// a vector of struct instead. The global list of constraint.
gtl::ITIVector<ConstraintIndex, LinearConstraint> constraints_;
@@ -123,8 +147,12 @@ class LinearConstraintManager {
return hash;
}
};
absl::flat_hash_map<Terms, int, TermsHash> equiv_constraints_;
absl::flat_hash_map<Terms, ConstraintIndex, TermsHash> equiv_constraints_;
int64 num_merged_constraints_ = 0;
int64 num_shortened_constraints_ = 0;
int64 num_splitted_constraints_ = 0;
int64 num_coeff_strenghtening_ = 0;
int num_cuts_ = 0;
std::map<std::string, int> type_to_num_cuts_;

22
ortools/sat/linear_programming_constraint.cc
View File

@@ -42,7 +42,8 @@ const double LinearProgrammingConstraint::kLpEpsilon = 1e-6;
// TODO(user): make SatParameters singleton too, otherwise changing them after
// a constraint was added will have no effect on this class.
LinearProgrammingConstraint::LinearProgrammingConstraint(Model* model)
: sat_parameters_(*(model->GetOrCreate<SatParameters>())),
: constraint_manager_(model),
sat_parameters_(*(model->GetOrCreate<SatParameters>())),
model_(model),
time_limit_(model->GetOrCreate<TimeLimit>()),
integer_trail_(model->GetOrCreate<IntegerTrail>()),
@@ -116,7 +117,6 @@ void LinearProgrammingConstraint::CreateLpFromConstraintManager() {
// Fill integer_lp_.
integer_lp_.clear();
infinity_norms_.clear();
constraint_manager_.SetParameters(sat_parameters_);
const auto& all_constraints = constraint_manager_.AllConstraints();
for (const auto index : constraint_manager_.LpConstraints()) {
const LinearConstraint& ct = all_constraints[index];
@@ -1168,9 +1168,12 @@ bool LinearProgrammingConstraint::ExactLpReasonning() {
gtl::ITIVector<ColIndex, IntegerValue> reduced_costs;
IntegerValue rc_ub;
CHECK(ComputeNewLinearConstraint(
/*use_constraint_status=*/false, integer_multipliers, &reduced_costs,
&rc_ub));
if (!ComputeNewLinearConstraint(
/*use_constraint_status=*/false, integer_multipliers, &reduced_costs,
&rc_ub)) {
VLOG(1) << "Issue while computing the exact LP reason. Aborting.";
return true;
}
// The "objective constraint" behave like if the unscaled cp multiplier was
// 1.0, so we will multiply it by this number and add it to reduced_costs.
@@ -1215,9 +1218,12 @@ bool LinearProgrammingConstraint::FillExactDualRayReason() {
gtl::ITIVector<ColIndex, IntegerValue> dense_new_constraint;
IntegerValue new_constraint_ub;
CHECK(ComputeNewLinearConstraint(
/*use_constraint_status=*/false, integer_multipliers,
&dense_new_constraint, &new_constraint_ub));
if (!ComputeNewLinearConstraint(
/*use_constraint_status=*/false, integer_multipliers,
&dense_new_constraint, &new_constraint_ub)) {
VLOG(1) << "Isse while computing the exact dual ray reason. Aborting.";
return false;
}
AdjustNewLinearConstraint(&integer_multipliers, &dense_new_constraint,
&new_constraint_ub);