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

better presolve of int_max, make extraction of size 2 table constraint deterministic; internal cleanup of CP-SAT search internals

Browse Source
This commit is contained in:
Laurent Perron
2019-05-10 23:26:36 +02:00
parent 35a9dc3b8a
commit 0f2fe9e578
6 changed files with 215 additions and 148 deletions
Download Patch File Download Diff File Expand all files Collapse all files

65
ortools/sat/cp_model_presolve.cc
View File

@@ -739,8 +739,8 @@ bool PresolveIntMax(ConstraintProto* ct, PresolveContext* context) {
}
// Pass 1, compute the infered min of the target, and remove duplicates.
int64 target_min = context->MinOf(target_ref);
int64 target_max = kint64min;
int64 infered_min = kint64min;
int64 infered_max = kint64min;
bool contains_target_ref = false;
std::set<int> used_ref;
int new_size = 0;
@@ -749,12 +749,12 @@ bool PresolveIntMax(ConstraintProto* ct, PresolveContext* context) {
if (gtl::ContainsKey(used_ref, ref)) continue;
if (gtl::ContainsKey(used_ref, NegatedRef(ref)) ||
ref == NegatedRef(target_ref)) {
target_min = std::max(target_min, int64{0});
infered_min = std::max(infered_min, int64{0});
}
used_ref.insert(ref);
ct->mutable_int_max()->set_vars(new_size++, ref);
target_min = std::max(target_min, context->MinOf(ref));
target_max = std::max(target_max, context->MaxOf(ref));
infered_min = std::max(infered_min, context->MinOf(ref));
infered_max = std::max(infered_max, context->MaxOf(ref));
}
if (new_size < ct->int_max().vars_size()) {
context->UpdateRuleStats("int_max: removed dup");
@@ -777,24 +777,65 @@ bool PresolveIntMax(ConstraintProto* ct, PresolveContext* context) {
return RemoveConstraint(ct, context);
}
// Compute the infered target_domain.
Domain infered_domain;
for (const int ref : ct->int_max().vars()) {
infered_domain = infered_domain.UnionWith(
context->DomainOf(ref).IntersectionWith({infered_min, infered_max}));
}
// Update the target domain.
bool domain_reduced = false;
if (!HasEnforcementLiteral(*ct)) {
Domain infered_domain;
for (const int ref : ct->int_max().vars()) {
infered_domain = infered_domain.UnionWith(
context->DomainOf(ref).IntersectionWith({target_min, target_max}));
}
if (!context->IntersectDomainWith(target_ref, infered_domain,
&domain_reduced)) {
return true;
}
}
// If the target is only used here and if
// infered_domain ∩ [kint64min, target_ub] ⊂ target_domain
// then the constraint is really max(...) <= target_ub and we can simplify it.
if (context->VariableIsUniqueAndRemovable(target_ref)) {
const Domain target_domain = context->DomainOf(target_ref);
if (infered_domain.IntersectionWith(Domain(kint64min, target_domain.Max()))
.IsIncludedIn(context->DomainOf(target_ref))) {
if (infered_domain.Max() <= target_domain.Max()) {
// The constraint is always satisfiable.
context->UpdateRuleStats("int_max: always true");
} else if (ct->enforcement_literal().empty()) {
// The constraint just restrict the upper bound of its variable.
for (const int ref : ct->int_max().vars()) {
context->UpdateRuleStats("int_max: lower than constant");
if (!context->IntersectDomainWith(
ref, Domain(kint64min, target_domain.Max()))) {
return false;
}
}
} else {
// We simply transform this into n reified constraints
// enforcement => [var_i <= target_domain.Max()].
context->UpdateRuleStats("int_max: reified lower than constant");
for (const int ref : ct->int_max().vars()) {
ConstraintProto* new_ct = context->working_model->add_constraints();
*(new_ct->mutable_enforcement_literal()) = ct->enforcement_literal();
ct->mutable_linear()->add_vars(ref);
ct->mutable_linear()->add_coeffs(1);
ct->mutable_linear()->add_domain(kint64min);
ct->mutable_linear()->add_domain(target_domain.Max());
}
}
// In all cases we delete the original constraint.
*(context->mapping_model->add_constraints()) = *ct;
return RemoveConstraint(ct, context);
}
}
// Pass 2, update the argument domains. Filter them eventually.
new_size = 0;
const int size = ct->int_max().vars_size();
target_max = context->MaxOf(target_ref);
const int64 target_max = context->MaxOf(target_ref);
for (const int ref : ct->int_max().vars()) {
if (!HasEnforcementLiteral(*ct)) {
if (!context->IntersectDomainWith(ref, Domain(kint64min, target_max),
@@ -802,7 +843,7 @@ bool PresolveIntMax(ConstraintProto* ct, PresolveContext* context) {
return true;
}
}
if (context->MaxOf(ref) >= target_min) {
if (context->MaxOf(ref) >= infered_min) {
ct->mutable_int_max()->set_vars(new_size++, ref);
}
}

3
ortools/sat/csharp/sat.i
View File

@@ -28,9 +28,6 @@ using System.Collections;
#include "ortools/sat/swig_helper.h"
%}
typedef int64_t int64;
typedef uint64_t uint64;
%module(directors="1") operations_research_sat
PROTO_INPUT(operations_research::sat::CpModelProto,

52
ortools/sat/linear_relaxation.cc
View File

@@ -17,6 +17,7 @@
#include "ortools/base/iterator_adaptors.h"
#include "ortools/sat/cp_model_loader.h"
#include "ortools/sat/integer.h"
#include "ortools/sat/linear_constraint.h"
#include "ortools/sat/linear_programming_constraint.h"
namespace operations_research {
@@ -380,9 +381,39 @@ void TryToLinearizeConstraint(const CpModelProto& model_proto,
}
relaxation->linear_constraints.push_back(constraint.Build());
} else if (ct.constraint_case() ==
ConstraintProto::ConstraintCase::kInterval) {
if (linearization_level < 3) return;
if (HasEnforcementLiteral(ct)) return;
const IntegerVariable start = mapping->Integer(ct.interval().start());
const IntegerVariable size = mapping->Integer(ct.interval().size());
const IntegerVariable end = mapping->Integer(ct.interval().end());
LinearConstraintBuilder lc(model, IntegerValue(0), IntegerValue(0));
lc.AddTerm(start, IntegerValue(1));
lc.AddTerm(size, IntegerValue(1));
lc.AddTerm(end, IntegerValue(-1));
relaxation->linear_constraints.push_back(lc.Build());
}
}
namespace {
// Adds enforcing_var => target <= bounding_var to relaxation.
void AppendEnforcedUpperBound(const IntegerVariable enforcing_var,
const IntegerVariable target,
const IntegerVariable bounding_var, Model* model,
LinearRelaxation* relaxation) {
IntegerTrail* integer_trail = model->GetOrCreate<IntegerTrail>();
const IntegerValue max_target_value = integer_trail->UpperBound(target);
const IntegerValue min_var_value = integer_trail->LowerBound(bounding_var);
const IntegerValue max_term_value = max_target_value - min_var_value;
LinearConstraintBuilder lc(model, kMinIntegerValue, max_term_value);
lc.AddTerm(target, IntegerValue(1));
lc.AddTerm(bounding_var, IntegerValue(-1));
lc.AddTerm(enforcing_var, max_term_value);
relaxation->linear_constraints.push_back(lc.Build());
}
} // namespace
void AppendMaxRelaxation(IntegerVariable target,
const std::vector<IntegerVariable>& vars,
int linearization_level, Model* model,
@@ -400,11 +431,14 @@ void AppendMaxRelaxation(IntegerVariable target,
}
// Part 2: Encode upper bound on X.
// NOTE(user) for size = 2, we can do this with 1 less variable.
if (vars.size() == 2 || linearization_level >= 2) {
IntegerTrail* integer_trail = model->GetOrCreate<IntegerTrail>();
const IntegerValue max_target_value = integer_trail->UpperBound(target);
// For size = 2, we do this with 1 less variable.
if (vars.size() == 2) {
IntegerVariable y = model->Add(NewIntegerVariable(0, 1));
AppendEnforcedUpperBound(y, target, vars[0], model, relaxation);
AppendEnforcedUpperBound(NegationOf(y), target, vars[1], model, relaxation);
return;
}
if (linearization_level >= 2) {
// For each X_i, we encode y_i => X <= X_i. And at least one of the y_i
// is true. Note that the correct y_i will be chosen because of the first
// part in linearlization (X >= X_i).
@@ -417,13 +451,7 @@ void AppendMaxRelaxation(IntegerVariable target,
// <=> max_term_value * y + X - X_i <= max_term_value.
// where max_tern_value is X_ub - X_i_lb.
IntegerVariable y = model->Add(NewIntegerVariable(0, 1));
const IntegerValue min_var_value = integer_trail->LowerBound(var);
const IntegerValue max_term_value = max_target_value - min_var_value;
LinearConstraintBuilder lc(model, kMinIntegerValue, max_term_value);
lc.AddTerm(target, IntegerValue(1));
lc.AddTerm(var, IntegerValue(-1));
lc.AddTerm(y, max_term_value);
relaxation->linear_constraints.push_back(lc.Build());
AppendEnforcedUpperBound(y, target, var, model, relaxation);
lc_exactly_one.AddTerm(y, IntegerValue(1));
}

21
ortools/sat/optimization.cc
View File

@@ -1090,21 +1090,17 @@ SatSolver::Status MinimizeIntegerVariableWithLinearScanAndLazyEncoding(
const std::function<void()>& feasible_solution_observer, Model* model) {
SatSolver* sat_solver = model->GetOrCreate<SatSolver>();
IntegerTrail* integer_trail = model->GetOrCreate<IntegerTrail>();
const SatParameters& parameters = *(model->GetOrCreate<SatParameters>());
// Simple linear scan algorithm to find the optimal.
SatSolver::Status result;
bool model_is_feasible = false;
while (true) {
result = ResetAndSolveIntegerProblem(/*assumptions=*/{}, model);
if (result != SatSolver::FEASIBLE) break;
const SatSolver::Status result = SolveIntegerProblem(model);
if (result != SatSolver::FEASIBLE) return result;
// The objective is the current lower bound of the objective_var.
const IntegerValue objective = integer_trail->LowerBound(objective_var);
// We have a solution!
model_is_feasible = true;
if (feasible_solution_observer != nullptr) {
feasible_solution_observer();
}
@@ -1117,20 +1113,9 @@ SatSolver::Status MinimizeIntegerVariableWithLinearScanAndLazyEncoding(
if (!integer_trail->Enqueue(
IntegerLiteral::LowerOrEqual(objective_var, objective - 1), {},
{})) {
result = SatSolver::INFEASIBLE;
break;
return SatSolver::INFEASIBLE;
}
}
CHECK_NE(result, SatSolver::FEASIBLE);
if (result == SatSolver::INFEASIBLE && model_is_feasible) {
// We proved the optimal and use the FEASIBLE value for this.
result = SatSolver::FEASIBLE;
} else {
sat_solver->Backtrack(0);
}
return result;
}
void RestrictObjectiveDomainWithBinarySearch(

6
ortools/sat/optimization.h
View File

@@ -127,10 +127,14 @@ SatSolver::Status SolveWithCardinalityEncodingAndCore(
// Model-based API, for now we just provide a basic algorithm that minimizes a
// given IntegerVariable by solving a sequence of decision problem by using
// SolveIntegerProblem().
// SolveIntegerProblem(). Returns the status of the last solved decision
// problem.
//
// The feasible_solution_observer function will be called each time a new
// feasible solution is found.
//
// Note that this function will resume the search from the current state of the
// solver, and it is up to the client to backtrack to the root node if needed.
SatSolver::Status MinimizeIntegerVariableWithLinearScanAndLazyEncoding(
IntegerVariable objective_var,
const std::function<void()>& feasible_solution_observer, Model* model);

216
ortools/sat/table.cc
View File

@@ -83,12 +83,14 @@ void ProcessOneColumn(
}
}
// Simpler encoding for table constraints with 2 variables.
void AddSizeTwoTable(
absl::Span<const IntegerVariable> vars,
const std::vector<std::vector<int64>>& tuples,
const std::vector<absl::flat_hash_set<int64>>& values_per_var,
Model* model) {
const int n = vars.size();
CHECK_EQ(n, 2);
IntegerTrail* const integer_trail = model->GetOrCreate<IntegerTrail>();
std::vector<absl::flat_hash_map<IntegerValue, Literal>> encodings(n);
@@ -106,10 +108,8 @@ void AddSizeTwoTable(
// One variable is fixed. Propagation is complete.
if (values_per_var[0].size() == 1 || values_per_var[1].size() == 1) return;
absl::flat_hash_map<LiteralIndex, absl::flat_hash_set<LiteralIndex>>
left_to_right;
absl::flat_hash_map<LiteralIndex, absl::flat_hash_set<LiteralIndex>>
right_to_left;
std::map<LiteralIndex, std::vector<Literal>> left_to_right;
std::map<LiteralIndex, std::vector<Literal>> right_to_left;
for (const auto& tuple : tuples) {
const IntegerValue left_value(tuple[0]);
@@ -119,105 +119,113 @@ void AddSizeTwoTable(
continue;
}
Literal left = gtl::FindOrDie(encodings[0], left_value);
Literal right = gtl::FindOrDie(encodings[1], right_value);
left_to_right[left.Index()].insert(right.Index());
right_to_left[right.Index()].insert(left.Index());
const Literal left = gtl::FindOrDie(encodings[0], left_value);
const Literal right = gtl::FindOrDie(encodings[1], right_value);
left_to_right[left.Index()].push_back(right);
right_to_left[right.Index()].push_back(left);
}
int implications = 0;
int clause_added = 0;
int large_clause_added = 0;
std::vector<Literal> clauses;
auto add_support =
[model, &clause_added, &large_clause_added, &implications, &clauses](
LiteralIndex lit, const absl::flat_hash_set<LiteralIndex>& supports,
int max_support_size) {
int num_implications = 0;
int num_clause_added = 0;
int num_large_clause_added = 0;
std::vector<Literal> clause;
auto add_support_constraint =
[model, &num_clause_added, &num_large_clause_added, &num_implications,
&clause](LiteralIndex lit, const std::vector<Literal>& supports,
int max_support_size) {
if (supports.size() == max_support_size) return;
if (supports.size() == 1) {
model->Add(Implication(Literal(lit), Literal(*supports.begin())));
implications++;
model->Add(Implication(Literal(lit), supports.front()));
num_implications++;
} else {
clauses.clear();
for (const LiteralIndex index : supports) {
clauses.push_back(Literal(index));
}
clauses.push_back(Literal(lit).Negated());
model->Add(ClauseConstraint(clauses));
clause_added++;
clause.assign(supports.begin(), supports.end());
clause.push_back(Literal(lit).Negated());
model->Add(ClauseConstraint(clause));
num_clause_added++;
if (supports.size() > max_support_size / 2) {
large_clause_added++;
num_large_clause_added++;
}
}
};
for (const auto& it : left_to_right) {
add_support(it.first, it.second, values_per_var[1].size());
add_support_constraint(it.first, it.second, values_per_var[1].size());
}
for (const auto& it : right_to_left) {
add_support(it.first, it.second, values_per_var[0].size());
add_support_constraint(it.first, it.second, values_per_var[0].size());
}
VLOG(2) << "Table: 2 variables, " << tuples.size() << " tuples encoded using "
<< clause_added << " clauses, " << large_clause_added
<< " large clauses, " << implications << " implications";
<< num_clause_added << " clauses, " << num_large_clause_added
<< " large clauses, " << num_implications << " implications";
}
void ExplorePrefixes(const std::vector<std::vector<int64>>& tuples,
const std::vector<std::vector<int64>>& var_domains,
absl::Span<const IntegerVariable> vars, Model* model) {
auto explore_prefix_span = [&](int start, int end) {
// Compute the maximum number of such prefix tuples.
int64 max_num_prefix_tuples = 1;
for (int i = start; i <= end; ++i) {
max_num_prefix_tuples =
CapProd(max_num_prefix_tuples, var_domains[i].size());
}
// This method heuristically explores subsets of variables and decide if the
// projection of all tuples nearly fills all the possible combination of
// projected variables domains.
//
// In that case, it creates the complement of the projected tuples and add that
// as a forbidden assignment constraint.
void ExploreSubsetOfVariablesAndAddNegatedTables(
const std::vector<std::vector<int64>>& tuples,
const std::vector<std::vector<int64>>& var_domains,
absl::Span<const IntegerVariable> vars, Model* model) {
const int num_vars = var_domains.size();
for (int start = 0; start < num_vars; ++start) {
const int limit = start == 0 ? num_vars : std::min(num_vars, start + 3);
for (int end = start + 1; end < limit; ++end) {
// TODO(user,user): If we add negated table for more than one value of
// end, because the set of variables will be included in each other, we
// could reduce the number of clauses added. I.e if we excluded
// (x=2, y=3) there is no need to exclude any of the tuples
// (x=2, y=3, z=*).
// Abort early.
if (max_num_prefix_tuples > 2 * tuples.size()) return;
// Compute the maximum number of such prefix tuples.
int64 max_num_prefix_tuples = 1;
for (int i = start; i <= end; ++i) {
max_num_prefix_tuples =
CapProd(max_num_prefix_tuples, var_domains[i].size());
}
absl::flat_hash_set<absl::Span<const int64>> prefixes;
for (const std::vector<int64>& tuple : tuples) {
prefixes.insert(absl::MakeSpan(&tuple[start], end - start + 1));
if (prefixes.size() == max_num_prefix_tuples) return;
}
const int num_prefix_tuples = prefixes.size();
// Abort early.
if (max_num_prefix_tuples > 2 * tuples.size()) break;
std::vector<std::vector<int64>> negated_tuples;
int created = 0;
if (num_prefix_tuples < max_num_prefix_tuples &&
max_num_prefix_tuples < num_prefix_tuples * 2) {
std::vector<int64> tmp_tuple;
for (int i = 0; i < max_num_prefix_tuples; ++i) {
tmp_tuple.clear();
int index = i;
for (int j = start; j <= end; ++j) {
tmp_tuple.push_back(var_domains[j][index % var_domains[j].size()]);
index /= var_domains[j].size();
}
if (!prefixes.contains(tmp_tuple)) {
negated_tuples.push_back(tmp_tuple);
created++;
absl::flat_hash_set<absl::Span<const int64>> prefixes;
bool skip = false;
for (const std::vector<int64>& tuple : tuples) {
prefixes.insert(absl::MakeSpan(&tuple[start], end - start + 1));
if (prefixes.size() == max_num_prefix_tuples) {
// Nothing to add with this range [start..end].
skip = true;
break;
}
}
AddNegatedTableConstraint(vars.subspan(start, end - start + 1),
negated_tuples, model);
VLOG(2) << " created = " << created << " for " << start << " .. " << end;
}
};
if (skip) continue;
const int num_prefix_tuples = prefixes.size();
for (int end = 1; end < var_domains.size(); ++end) {
explore_prefix_span(0, end);
}
for (int start = 1; start + 1 < var_domains.size(); ++start) {
explore_prefix_span(start, start + 1);
}
for (int start = 1; start + 2 < var_domains.size(); ++start) {
explore_prefix_span(start, start + 2);
}
for (int start = 1; start + 3 < var_domains.size(); ++start) {
explore_prefix_span(start, start + 3);
std::vector<std::vector<int64>> negated_tuples;
int created = 0;
if (num_prefix_tuples < max_num_prefix_tuples &&
max_num_prefix_tuples < num_prefix_tuples * 2) {
std::vector<int64> tmp_tuple;
for (int i = 0; i < max_num_prefix_tuples; ++i) {
tmp_tuple.clear();
int index = i;
for (int j = start; j <= end; ++j) {
tmp_tuple.push_back(var_domains[j][index % var_domains[j].size()]);
index /= var_domains[j].size();
}
if (!prefixes.contains(tmp_tuple)) {
negated_tuples.push_back(tmp_tuple);
created++;
}
}
AddNegatedTableConstraint(vars.subspan(start, end - start + 1),
negated_tuples, model);
VLOG(2) << " add negated tables with " << created
<< " tuples on the range [" << start << "," << end << "]";
}
}
}
}
@@ -303,35 +311,29 @@ void AddTableConstraint(absl::Span<const IntegerVariable> vars,
return;
}
// Detect the case when the first n-1 columns are all different.
// This encodes the implication table (tuple of size n - 1) implies value.
absl::flat_hash_set<absl::Span<const int64>> prefixes;
for (const std::vector<int64>& tuple : tuples) {
prefixes.insert(absl::MakeSpan(tuple.data(), n - 1));
}
const int num_prefix_tuples = prefixes.size();
// Compute the maximum number of such prefix tuples.
int64 max_num_prefix_tuples = 1;
for (int i = 0; i + 1 < n; ++i) {
max_num_prefix_tuples =
CapProd(max_num_prefix_tuples, values_per_var[i].size());
}
if (n == 2) {
AddSizeTwoTable(vars, tuples, values_per_var, model);
return;
}
// It is easier to compute this before compression, as compression will merge
// tuples.
int num_prefix_tuples = 0;
{
absl::flat_hash_set<absl::Span<const int64>> prefixes;
for (const std::vector<int64>& tuple : tuples) {
prefixes.insert(absl::MakeSpan(tuple.data(), n - 1));
}
num_prefix_tuples = prefixes.size();
}
std::vector<std::vector<int64>> var_domains(n);
for (int j = 0; j < n; ++j) {
var_domains[j].assign(values_per_var[j].begin(), values_per_var[j].end());
std::sort(var_domains[j].begin(), var_domains[j].end());
}
if (vars.size() > 2) {
ExplorePrefixes(tuples, var_domains, vars, model);
}
// Detect if prefix tuples are all different.
const bool prefixes_are_all_different = num_prefix_tuples == num_valid_tuples;
CHECK_GT(vars.size(), 2);
ExploreSubsetOfVariablesAndAddNegatedTables(tuples, var_domains, vars, model);
// The variable domains have been computed. Fully encode variables.
// Note that in some corner cases (like duplicate vars), as we call
@@ -359,7 +361,16 @@ void AddTableConstraint(absl::Span<const IntegerVariable> vars,
CompressTuples(domain_sizes, any_value, &tuples);
const int num_compressed_tuples = tuples.size();
// Detect if prefix tuples are all different.
const bool prefixes_are_all_different = num_prefix_tuples == num_valid_tuples;
if (VLOG_IS_ON(2)) {
// Compute the maximum number of prefix tuples.
int64 max_num_prefix_tuples = 1;
for (int i = 0; i + 1 < n; ++i) {
max_num_prefix_tuples =
CapProd(max_num_prefix_tuples, values_per_var[i].size());
}
std::string message = absl::StrCat(
"Table: ", n, " variables, original tuples = ", num_original_tuples);
if (num_valid_tuples != num_original_tuples) {
@@ -373,7 +384,7 @@ void AddTableConstraint(absl::Span<const IntegerVariable> vars,
absl::StrAppend(&message, ", full prefix = true");
}
} else {
absl::StrAppend(&message, ", num prefix tuples = ", prefixes.size());
absl::StrAppend(&message, ", num prefix tuples = ", num_prefix_tuples);
}
if (num_compressed_tuples != num_valid_tuples) {
absl::StrAppend(&message,
@@ -437,8 +448,9 @@ void AddTableConstraint(absl::Span<const IntegerVariable> vars,
}
if (prefixes_are_all_different) {
// This is optional propagation wise, but it should lead
// to better explanation.
// The first n-1 columns are all different, this encodes the implication
// table (tuple of size n - 1) implies value. We can add an optional
// propagation that should lead to better explanation.
// For each tuple, we add a clause prefix => last value.
std::vector<Literal> clause;
for (int j = 0; j < tuples.size(); ++j) {