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

improve lp cut and constraint management

Browse Source
This commit is contained in:
Laurent Perron
2019-11-14 12:42:35 -08:00
parent 09063ae0c2
commit 28033c996f
5 changed files with 132 additions and 93 deletions
Download Patch File Download Diff File Expand all files Collapse all files

2
ortools/sat/cp_model_expand.cc
View File

@@ -43,7 +43,7 @@ void ExpandReservoir(ConstraintProto* ct, PresolveContext* context) {
auto is_optional = [&context, &reservoir](int index) {
if (reservoir.actives_size() == 0) return false;
const int literal = reservoir.actives(index);
return !context->DomainOf(PositiveRef(literal)).IsFixed();
return !context->IsFixed(literal);
};
const int true_literal = context->GetOrCreateConstantVar(1);
auto active = [&reservoir, true_literal](int index) {

94
ortools/sat/linear_constraint_manager.cc
View File

@@ -82,22 +82,48 @@ LinearConstraintManager::~LinearConstraintManager() {
}
}
// TODO(user,user): Also update the revised simplex starting basis for the
// next solve.
void LinearConstraintManager::RemoveMarkedConstraints() {
int new_index = 0;
for (int i = 0; i < lp_constraints_.size(); ++i) {
const ConstraintIndex constraint = lp_constraints_[i];
if (constraints_removal_list_.contains(constraint)) {
constraint_infos_[constraint].is_in_lp = false;
continue;
bool LinearConstraintManager::MaybeRemoveSomeInactiveConstraints(
glop::BasisState* solution_state) {
if (solution_state->IsEmpty()) return false; // Mainly to simplify tests.
const glop::RowIndex num_rows(lp_constraints_.size());
const glop::ColIndex num_cols =
solution_state->statuses.size() - RowToColIndex(num_rows);
int new_size = 0;
for (int i = 0; i < num_rows; ++i) {
const ConstraintIndex constraint_index = lp_constraints_[i];
// Constraints that are not tight in the current solution have a basic
// status. We remove the ones that have been inactive in the last recent
// solves.
//
// TODO(user): More advanced heuristics might perform better, I didn't do
// a lot of tuning experiments yet.
const glop::VariableStatus row_status =
solution_state->statuses[num_cols + glop::ColIndex(i)];
if (row_status == glop::VariableStatus::BASIC) {
constraint_infos_[constraint_index].inactive_count++;
if (constraint_infos_[constraint_index].inactive_count >
sat_parameters_.max_consecutive_inactive_count()) {
constraint_infos_[constraint_index].is_in_lp = false;
continue; // Remove it.
}
} else {
// Only count consecutive inactivities.
constraint_infos_[constraint_index].inactive_count = 0;
}
lp_constraints_[new_index] = constraint;
new_index++;
lp_constraints_[new_size] = constraint_index;
solution_state->statuses[num_cols + glop::ColIndex(new_size)] = row_status;
new_size++;
}
lp_constraints_.resize(new_index);
VLOG(2) << "Removed " << constraints_removal_list_.size() << " constraints.";
constraints_removal_list_.clear();
const int num_removed_constraints = lp_constraints_.size() - new_size;
lp_constraints_.resize(new_size);
solution_state->statuses.resize(num_cols + glop::ColIndex(new_size));
if (num_removed_constraints > 0) {
VLOG(2) << "Removed " << num_removed_constraints << " constraints";
}
return num_removed_constraints > 0;
}
// Because sometimes we split a == constraint in two (>= and <=), it makes sense
@@ -343,7 +369,8 @@ bool LinearConstraintManager::SimplifyConstraint(LinearConstraint* ct) {
}
bool LinearConstraintManager::ChangeLp(
const gtl::ITIVector<IntegerVariable, double>& lp_solution) {
const gtl::ITIVector<IntegerVariable, double>& lp_solution,
glop::BasisState* solution_state) {
VLOG(3) << "Enter ChangeLP, scan " << constraint_infos_.size()
<< " constraints";
std::vector<ConstraintIndex> new_constraints;
@@ -375,25 +402,16 @@ bool LinearConstraintManager::ChangeLp(
equiv_constraints_[constraint_infos_[i].hash] = i;
}
// TODO(user,user): Use constraint status from GLOP for constraints
// already in LP.
if (constraint_infos_[i].is_in_lp) continue;
const double activity =
ComputeActivity(constraint_infos_[i].constraint, lp_solution);
const double lb_violation =
ToDouble(constraint_infos_[i].constraint.lb) - activity;
const double ub_violation =
activity - ToDouble(constraint_infos_[i].constraint.ub);
const double violation = std::max(lb_violation, ub_violation);
if (constraint_infos_[i].is_in_lp && violation < tolerance) {
constraint_infos_[i].inactive_count++;
if (constraint_infos_[i].is_cut &&
constraint_infos_[i].inactive_count >
sat_parameters_.max_inactive_count()) {
// Mark cut for removal.
constraints_removal_list_.insert(i);
}
} else if (!constraint_infos_[i].is_in_lp && violation >= tolerance) {
if (violation >= tolerance) {
constraint_infos_[i].inactive_count = 0;
new_constraints.push_back(i);
new_constraints_efficacies.push_back(violation /
@@ -413,10 +431,11 @@ bool LinearConstraintManager::ChangeLp(
}
}
// Remove constraints in a batch to avoid changing the LP too frequently.
if (constraints_removal_list_.size() >=
sat_parameters_.constraint_removal_batch_size()) {
RemoveMarkedConstraints();
// Remove constraints from the current LP that have been inactive for a while.
// We do that after we computed new_constraints so we do not need to iterate
// over the just deleted constraints.
if (MaybeRemoveSomeInactiveConstraints(solution_state)) {
current_lp_is_changed_ = true;
}
// Note that the algo below is in O(limit * new_constraint). In order to
@@ -430,7 +449,8 @@ bool LinearConstraintManager::ChangeLp(
// TODO(user): This blowup factor could be adaptative w.r.t. the constraint
// limit.
const int kBlowupFactor = 4;
int constraint_limit = std::min(50, static_cast<int>(new_constraints.size()));
int constraint_limit = std::min(sat_parameters_.new_constraints_batch_size(),
static_cast<int>(new_constraints.size()));
if (lp_constraints_.empty()) {
constraint_limit = std::min(1000, static_cast<int>(new_constraints.size()));
}
@@ -448,6 +468,7 @@ bool LinearConstraintManager::ChangeLp(
new_constraints.resize(kBlowupFactor * constraint_limit);
}
int num_added = 0;
int num_skipped_checks = 0;
const int kCheckFrequency = 100;
ConstraintIndex last_added_candidate = kInvalidConstraintIndex;
@@ -507,13 +528,20 @@ bool LinearConstraintManager::ChangeLp(
// Note that it is important for LP incremental solving that the old
// constraints stays at the same position in this list (and thus in the
// returned GetLp()).
++num_added;
current_lp_is_changed_ = true;
lp_constraints_.push_back(best_candidate);
last_added_candidate = best_candidate;
}
}
VLOG(3) << " - Exit ChangeLP";
if (num_added > 0) {
// We update the solution sate to match the new LP size.
VLOG(2) << "Added " << num_added << " constraints.";
solution_state->statuses.resize(solution_state->statuses.size() + num_added,
glop::VariableStatus::BASIC);
}
// The LP changed only if we added new constraints or if some constraints
// already inside changed (simplification or tighter bounds).
if (current_lp_is_changed_) {

24
ortools/sat/linear_constraint_manager.h
View File

@@ -19,6 +19,7 @@
#include "absl/container/flat_hash_map.h"
#include "absl/container/flat_hash_set.h"
#include "ortools/glop/revised_simplex.h"
#include "ortools/sat/linear_constraint.h"
#include "ortools/sat/model.h"
#include "ortools/sat/sat_parameters.pb.h"
@@ -81,8 +82,13 @@ class LinearConstraintManager {
// should usually be the optimal solution of the LP returned by GetLp() before
// this call, but is just used as an heuristic.
//
// The current solution state is used for detecting inactive constraints. It
// is also updated correctly on constraint deletion/addition so that the
// simplex can be fully iterative on restart by loading this modified state.
//
// Returns true iff LpConstraints() will return a different LP than before.
bool ChangeLp(const gtl::ITIVector<IntegerVariable, double>& lp_solution);
bool ChangeLp(const gtl::ITIVector<IntegerVariable, double>& lp_solution,
glop::BasisState* solution_state);
// This can be called initially to add all the current constraint to the LP
// returned by GetLp().
@@ -105,8 +111,16 @@ class LinearConstraintManager {
int64 num_coeff_strenghtening() const { return num_coeff_strenghtening_; }
private:
// Removes the marked constraints from the LP.
void RemoveMarkedConstraints();
// Heuristic that decide which constraints we should remove from the current
// LP. Note that such constraints can be added back later by the heuristic
// responsible for adding new constraints from the pool.
//
// Returns true iff one or more constraints where removed.
//
// If the solutions_state is empty, then this function does nothing and
// returns false (this is used for tests). Otherwise, the solutions_state is
// assumed to correspond to the current LP and to be of the correct size.
bool MaybeRemoveSomeInactiveConstraints(glop::BasisState* solution_state);
// Apply basic inprocessing simplification rules:
// - remove fixed variable
@@ -133,10 +147,6 @@ class LinearConstraintManager {
gtl::ITIVector<ConstraintIndex, ConstraintInfo> constraint_infos_;
// Temporary list of constraints marked for removal. Note that we remove
// constraints in batch to avoid changing LP too frequently.
absl::flat_hash_set<ConstraintIndex> constraints_removal_list_;
// The subset of constraints currently in the lp.
std::vector<ConstraintIndex> lp_constraints_;

92
ortools/sat/linear_programming_constraint.cc
View File

@@ -886,15 +886,10 @@ bool LinearProgrammingConstraint::Propagate() {
// this beeing called again on the next IncrementalPropagate() call, but that
// might not always happen at level zero.
if (simplex_.GetProblemStatus() == glop::ProblemStatus::OPTIMAL) {
// First add any new lazy constraints or cuts that where previsouly
// generated and are now cutting the current solution.
if (constraint_manager_.ChangeLp(expanded_lp_solution_)) {
CreateLpFromConstraintManager();
if (!SolveLp()) return true;
} else if (constraint_manager_.num_cuts() <
sat_parameters_.max_num_cuts()) {
const int old_num_cuts = constraint_manager_.num_cuts();
// We wait for the first batch of problem constraints to be added before we
// begin to generate cuts.
if (!integer_lp_.empty() &&
constraint_manager_.num_cuts() < sat_parameters_.max_num_cuts()) {
// The "generic" cuts are currently part of this class as they are using
// data from the current LP.
//
@@ -912,22 +907,23 @@ bool LinearProgrammingConstraint::Propagate() {
generator.generate_cuts(expanded_lp_solution_, &constraint_manager_);
}
}
}
if (constraint_manager_.num_cuts() > old_num_cuts &&
constraint_manager_.ChangeLp(expanded_lp_solution_)) {
CreateLpFromConstraintManager();
const double old_obj = simplex_.GetObjectiveValue();
if (!SolveLp()) return true;
if (simplex_.GetProblemStatus() == glop::ProblemStatus::OPTIMAL) {
VLOG(1) << "Cuts relaxation improvement " << old_obj << " -> "
<< simplex_.GetObjectiveValue()
<< " diff: " << simplex_.GetObjectiveValue() - old_obj
<< " level: " << trail_->CurrentDecisionLevel();
}
} else {
if (trail_->CurrentDecisionLevel() == 0) {
lp_at_level_zero_is_final_ = true;
}
glop::BasisState state = simplex_.GetState();
if (constraint_manager_.ChangeLp(expanded_lp_solution_, &state)) {
simplex_.LoadStateForNextSolve(state);
CreateLpFromConstraintManager();
const double old_obj = simplex_.GetObjectiveValue();
if (!SolveLp()) return true;
if (simplex_.GetProblemStatus() == glop::ProblemStatus::OPTIMAL) {
VLOG(1) << "Relaxation improvement " << old_obj << " -> "
<< simplex_.GetObjectiveValue()
<< " diff: " << simplex_.GetObjectiveValue() - old_obj
<< " level: " << trail_->CurrentDecisionLevel();
}
} else {
if (trail_->CurrentDecisionLevel() == 0) {
lp_at_level_zero_is_final_ = true;
}
}
}
@@ -1408,41 +1404,41 @@ void LinearProgrammingConstraint::AdjustNewLinearConstraint(
for (const auto entry : integer_lp_[row].terms) {
const ColIndex col = entry.first;
const IntegerValue coeff = entry.second;
const IntegerValue abs_coef = IntTypeAbs(coeff);
CHECK_NE(coeff, 0);
const IntegerVariable var = integer_variables_[col.value()];
const IntegerValue lb = integer_trail_->LowerBound(var);
const IntegerValue ub = integer_trail_->UpperBound(var);
// Moving a variable away from zero seems to improve the bound even
// if it reduces the number of non-zero. Note that this is because of
// this that positive_diff and negative_diff are not the same.
const IntegerValue current = (*dense_terms)[col];
if (current == 0) {
const IntegerValue overflow_limit(
FloorRatio(kMaxWantedCoeff, IntTypeAbs(coeff)));
FloorRatio(kMaxWantedCoeff, abs_coef));
positive_limit = std::min(positive_limit, overflow_limit);
negative_limit = std::min(negative_limit, overflow_limit);
const IntegerVariable var = integer_variables_[col.value()];
if (coeff > 0) {
positive_diff -= coeff * integer_trail_->LowerBound(var);
negative_diff -= coeff * integer_trail_->UpperBound(var);
positive_diff -= coeff * lb;
negative_diff -= coeff * ub;
} else {
positive_diff -= coeff * integer_trail_->UpperBound(var);
negative_diff -= coeff * integer_trail_->LowerBound(var);
positive_diff -= coeff * ub;
negative_diff -= coeff * lb;
}
continue;
}
// We don't want to change the sign of current or to have an overflow.
IntegerValue before_sign_change(
FloorRatio(IntTypeAbs(current), IntTypeAbs(coeff)));
// If the variable is fixed, we don't actually care about changing the
// sign.
const IntegerVariable var = integer_variables_[col.value()];
if (integer_trail_->LowerBound(var) == integer_trail_->UpperBound(var)) {
before_sign_change = kMaxWantedCoeff;
// We don't want to change the sign of current (except if the variable is
// fixed) or to have an overflow.
IntegerValue before_sign_change = kMaxWantedCoeff;
if (lb != ub) {
before_sign_change = FloorRatio(IntTypeAbs(current), abs_coef);
}
const IntegerValue overflow_limit(
FloorRatio(kMaxWantedCoeff - IntTypeAbs(current), IntTypeAbs(coeff)));
FloorRatio(kMaxWantedCoeff - IntTypeAbs(current), abs_coef));
if (current > 0 == coeff > 0) { // Same sign.
negative_limit = std::min(negative_limit, before_sign_change);
positive_limit = std::min(positive_limit, overflow_limit);
@@ -1452,12 +1448,11 @@ void LinearProgrammingConstraint::AdjustNewLinearConstraint(
}
// This is how diff change.
const IntegerValue implied = current > 0
? integer_trail_->LowerBound(var)
: integer_trail_->UpperBound(var);
positive_diff -= coeff * implied;
negative_diff -= coeff * implied;
const IntegerValue implied = current > 0 ? lb : ub;
if (implied != 0) {
positive_diff -= coeff * implied;
negative_diff -= coeff * implied;
}
}
// Only add a multiple of this row if it tighten the final constraint.
@@ -1553,6 +1548,11 @@ bool LinearProgrammingConstraint::ExactLpReasonning() {
IntegerSumLE* cp_constraint =
new IntegerSumLE({}, new_constraint.vars, new_constraint.coeffs,
new_constraint.ub, model_);
if (trail_->CurrentDecisionLevel() == 0) {
// Since we will never ask the reason for a constraint at level 0, we just
// keep the last one.
optimal_constraints_.clear();
}
optimal_constraints_.emplace_back(cp_constraint);
rev_optimal_constraints_size_ = optimal_constraints_.size();
return cp_constraint->Propagate();

13
ortools/sat/sat_parameters.proto
View File

@@ -590,13 +590,14 @@ message SatParameters {
// feature.
optional double min_orthogonality_for_lp_constraints = 115 [default = 0.0];
// If a constraint/cut in LP is not active for the 'max_inactive_count',
// remove it from the LP.
optional int64 max_inactive_count = 121 [default = 1000];
// If a constraint/cut in LP is not active for that many consecutive OPTIMAL
// solves, remove it from the LP. Note that it might be added again later if
// it become violated by the current LP solution.
optional int32 max_consecutive_inactive_count = 121 [default = 100];
// Remove constraints only if at least 'constraint_removal_batch_size'
// constraints are marked for removal.
optional int64 constraint_removal_batch_size = 122 [default = 100];
// Add that many lazy contraints (or cuts) at once in the LP. Note that at the
// beginning of the solve, we do add more than this.
optional int32 new_constraints_batch_size = 122 [default = 50];
// The search branching will be used to decide how to branch on unfixed nodes.
enum SearchBranching {