From a841258da74eab4ef70285bb86d3a709f28f8320 Mon Sep 17 00:00:00 2001 From: Laurent Perron Date: Tue, 16 Sep 2025 18:01:41 +0200 Subject: [PATCH] [CP-SAT] misc improvements to the internals; more work on general constraints with enforcements --- ortools/sat/BUILD.bazel | 6 +- ortools/sat/circuit.cc | 73 +++-- ortools/sat/cp_model_expand.cc | 52 ---- ortools/sat/cp_model_loader.cc | 13 +- ortools/sat/cp_model_presolve.cc | 2 +- ortools/sat/cp_model_search.cc | 56 ++-- ortools/sat/cp_model_search.h | 2 +- ortools/sat/cp_model_solver_helpers.cc | 22 +- ortools/sat/diffn.cc | 38 ++- ortools/sat/enforcement.cc | 14 +- ortools/sat/enforcement.h | 4 + ortools/sat/integer_search.cc | 25 +- ortools/sat/linear_propagation.cc | 12 +- ortools/sat/pb_constraint.cc | 364 +++++++++++++++++++------ ortools/sat/pb_constraint.h | 70 ++++- ortools/sat/pb_constraint_test.cc | 155 ++++++++++- ortools/sat/sat_parameters.proto | 2 +- ortools/sat/sat_solver.cc | 65 ++++- ortools/sat/sat_solver.h | 15 +- ortools/sat/util.cc | 4 +- ortools/sat/work_assignment.cc | 3 +- 21 files changed, 712 insertions(+), 285 deletions(-) diff --git a/ortools/sat/BUILD.bazel b/ortools/sat/BUILD.bazel index 8333f1ebee..ae05db8de0 100644 --- a/ortools/sat/BUILD.bazel +++ b/ortools/sat/BUILD.bazel @@ -1696,6 +1696,7 @@ cc_library( srcs = ["pb_constraint.cc"], hdrs = ["pb_constraint.h"], deps = [ + ":enforcement", ":model", ":sat_base", ":sat_parameters_cc_proto", @@ -1705,10 +1706,10 @@ cc_library( "//ortools/util:saturated_arithmetic", "//ortools/util:stats", "//ortools/util:strong_integers", + "@abseil-cpp//absl/cleanup", "@abseil-cpp//absl/container:flat_hash_map", "@abseil-cpp//absl/hash", "@abseil-cpp//absl/log:check", - "@abseil-cpp//absl/strings", "@abseil-cpp//absl/strings:str_format", "@abseil-cpp//absl/types:span", ], @@ -1719,6 +1720,7 @@ cc_test( size = "small", srcs = ["pb_constraint_test.cc"], deps = [ + ":enforcement", ":model", ":pb_constraint", ":sat_base", @@ -3745,6 +3747,7 @@ cc_library( ":intervals", ":model", ":no_overlap_2d_helper", + ":pb_constraint", ":sat_base", ":sat_parameters_cc_proto", ":sat_solver", @@ -3801,6 +3804,7 @@ cc_library( ":enforcement_helper", ":integer", ":model", + ":pb_constraint", ":sat_base", ":sat_solver", ":util", diff --git a/ortools/sat/circuit.cc b/ortools/sat/circuit.cc index a21b41efc6..f18867cdc9 100644 --- a/ortools/sat/circuit.cc +++ b/ortools/sat/circuit.cc @@ -28,6 +28,7 @@ #include "ortools/sat/enforcement.h" #include "ortools/sat/integer.h" #include "ortools/sat/model.h" +#include "ortools/sat/pb_constraint.h" #include "ortools/sat/sat_base.h" #include "ortools/sat/sat_solver.h" #include "ortools/util/strong_integers.h" @@ -696,6 +697,25 @@ std::function ExactlyOnePerRowAndPerColumn( }; } +namespace { +bool AddAtMostOne(absl::Span enforcement_literals, + absl::Span literals, Model* model) { + if (enforcement_literals.empty()) { + return model->GetOrCreate()->AddAtMostOne(literals); + } + std::vector enforcement(enforcement_literals.begin(), + enforcement_literals.end()); + std::vector cst; + cst.reserve(literals.size()); + for (const Literal l : literals) { + cst.emplace_back(l, Coefficient(1)); + } + return model->GetOrCreate()->AddLinearConstraint( + /*use_lower_bound=*/false, Coefficient(0), + /*use_upper_bound=*/true, Coefficient(1), &enforcement, &cst); +} +} // namespace + void LoadSubcircuitConstraint(int num_nodes, absl::Span tails, absl::Span heads, absl::Span enforcement_literals, @@ -709,39 +729,32 @@ void LoadSubcircuitConstraint(int num_nodes, absl::Span tails, // If a node has no outgoing or no incoming arc, the model will be unsat // as soon as we add the corresponding ExactlyOneConstraint(). auto sat_solver = model->GetOrCreate(); - auto implications = model->GetOrCreate(); - if (enforcement_literals.empty()) { - // TODO(user): how to generalize this to support enforcement literals? - // (we can easily add the negated enforcement literals to AddProblemClause() - // calls, but the AddAtMostOne() calls are trickier). For now this is done - // with additional constraints added by ExpandCircuit() during expansion. - std::vector> exactly_one_incoming(num_nodes); - std::vector> exactly_one_outgoing(num_nodes); - for (int arc = 0; arc < num_arcs; arc++) { - const int tail = tails[arc]; - const int head = heads[arc]; - exactly_one_outgoing[tail].push_back(literals[arc]); - exactly_one_incoming[head].push_back(literals[arc]); + std::vector> exactly_one_incoming(num_nodes); + std::vector> exactly_one_outgoing(num_nodes); + for (int arc = 0; arc < num_arcs; arc++) { + const int tail = tails[arc]; + const int head = heads[arc]; + exactly_one_outgoing[tail].push_back(literals[arc]); + exactly_one_incoming[head].push_back(literals[arc]); + } + for (int i = 0; i < exactly_one_incoming.size(); ++i) { + if (i == 0 && multiple_subcircuit_through_zero) continue; + if (!AddAtMostOne(enforcement_literals, exactly_one_incoming[i], model)) { + sat_solver->NotifyThatModelIsUnsat(); + return; } - for (int i = 0; i < exactly_one_incoming.size(); ++i) { - if (i == 0 && multiple_subcircuit_through_zero) continue; - if (!implications->AddAtMostOne(exactly_one_incoming[i])) { - sat_solver->NotifyThatModelIsUnsat(); - return; - } - sat_solver->AddProblemClause(exactly_one_incoming[i]); - if (sat_solver->ModelIsUnsat()) return; - } - for (int i = 0; i < exactly_one_outgoing.size(); ++i) { - if (i == 0 && multiple_subcircuit_through_zero) continue; - if (!implications->AddAtMostOne(exactly_one_outgoing[i])) { - sat_solver->NotifyThatModelIsUnsat(); - return; - } - sat_solver->AddProblemClause(exactly_one_outgoing[i]); - if (sat_solver->ModelIsUnsat()) return; + model->Add(EnforcedClause(enforcement_literals, exactly_one_incoming[i])); + if (sat_solver->ModelIsUnsat()) return; + } + for (int i = 0; i < exactly_one_outgoing.size(); ++i) { + if (i == 0 && multiple_subcircuit_through_zero) continue; + if (!AddAtMostOne(enforcement_literals, exactly_one_outgoing[i], model)) { + sat_solver->NotifyThatModelIsUnsat(); + return; } + model->Add(EnforcedClause(enforcement_literals, exactly_one_outgoing[i])); + if (sat_solver->ModelIsUnsat()) return; } CircuitPropagator::Options options; diff --git a/ortools/sat/cp_model_expand.cc b/ortools/sat/cp_model_expand.cc index bc6e4d9ee8..8d345ef85f 100644 --- a/ortools/sat/cp_model_expand.cc +++ b/ortools/sat/cp_model_expand.cc @@ -2703,50 +2703,6 @@ void MaybeExpandAllDiff(ConstraintProto* ct, PresolveContext* context, if (!keep_after_expansion) ct->Clear(); } -template -void ExpandCircuitOrRoutes(absl::Span enforcement_literals, - const T& graph_proto, PresolveContext* context, - bool multiple_subcircuit_through_zero) { - // The constraints added below are added by LoadSubcircuitConstraint() when - // enforcement_literals is empty. - // TODO(user): ideally we don't want to do that here, but only create - // these constraints at loading time so that the presolve does not have to - // deal with redundant constraints. - if (enforcement_literals.empty()) return; - const int num_arcs = graph_proto.tails_size(); - // At this point the node indices can have arbitrary values. - absl::flat_hash_map dense_node_index; - for (int arc = 0; arc < num_arcs; arc++) { - dense_node_index.insert({graph_proto.tails(arc), dense_node_index.size()}); - dense_node_index.insert({graph_proto.heads(arc), dense_node_index.size()}); - } - const int num_nodes = dense_node_index.size(); - std::vector> exactly_one_incoming(num_nodes); - std::vector> exactly_one_outgoing(num_nodes); - for (int arc = 0; arc < num_arcs; arc++) { - const int tail = dense_node_index[graph_proto.tails(arc)]; - const int head = dense_node_index[graph_proto.heads(arc)]; - exactly_one_outgoing[tail].push_back(graph_proto.literals(arc)); - exactly_one_incoming[head].push_back(graph_proto.literals(arc)); - } - for (int i = 0; i < exactly_one_incoming.size(); ++i) { - if (i == 0 && multiple_subcircuit_through_zero) continue; - ConstraintProto* const new_ct = context->working_model->add_constraints(); - *new_ct->mutable_enforcement_literal() = {enforcement_literals.begin(), - enforcement_literals.end()}; - LiteralsToLinear(exactly_one_incoming[i], /*lb=*/1, /*ub=*/1, - new_ct->mutable_linear()); - } - for (int i = 0; i < exactly_one_outgoing.size(); ++i) { - if (i == 0 && multiple_subcircuit_through_zero) continue; - ConstraintProto* const new_ct = context->working_model->add_constraints(); - *new_ct->mutable_enforcement_literal() = {enforcement_literals.begin(), - enforcement_literals.end()}; - LiteralsToLinear(exactly_one_outgoing[i], /*lb=*/1, /*ub=*/1, - new_ct->mutable_linear()); - } -} - } // namespace void ExpandCpModel(PresolveContext* context) { @@ -2824,14 +2780,6 @@ void ExpandCpModel(PresolveContext* context) { has_all_diffs = true; skip = true; break; - case ConstraintProto::kCircuit: - ExpandCircuitOrRoutes(ct->enforcement_literal(), ct->circuit(), context, - /*multiple_subcircuit_through_zero=*/false); - break; - case ConstraintProto::kRoutes: - ExpandCircuitOrRoutes(ct->enforcement_literal(), ct->routes(), context, - /*multiple_subcircuit_through_zero=*/true); - break; default: skip = true; break; diff --git a/ortools/sat/cp_model_loader.cc b/ortools/sat/cp_model_loader.cc index 7fe304f98a..ffe9284e90 100644 --- a/ortools/sat/cp_model_loader.cc +++ b/ortools/sat/cp_model_loader.cc @@ -29,11 +29,9 @@ #include "absl/log/check.h" #include "absl/log/log.h" #include "absl/log/vlog_is_on.h" -#include "absl/meta/type_traits.h" #include "absl/strings/str_cat.h" #include "absl/types/span.h" #include "ortools/algorithms/sparse_permutation.h" -#include "ortools/base/mathutil.h" #include "ortools/base/stl_util.h" #include "ortools/base/strong_vector.h" #include "ortools/sat/all_different.h" @@ -1380,8 +1378,7 @@ void LoadLinearConstraint(const ConstraintProto& ct, Model* m) { // Note that the domain/enforcement of the main constraint do not change. // Same for the min/sum and max_sum. The intermediate variables are always // equal to the intermediate sum, independently of the enforcement. - const bool pseudo_boolean = !HasEnforcementLiteral(ct) && - ct.linear().domain_size() == 2 && all_booleans; + const bool pseudo_boolean = ct.linear().domain_size() == 2 && all_booleans; if (!pseudo_boolean && ct.linear().vars().size() > params.linear_split_size()) { const auto& domain = ct.linear().domain(); @@ -1393,11 +1390,9 @@ void LoadLinearConstraint(const ConstraintProto& ct, Model* m) { if (ct.linear().domain_size() == 2) { const int64_t lb = ct.linear().domain(0); const int64_t ub = ct.linear().domain(1); - const std::vector enforcement_literals = + std::vector enforcement_literals = mapping->Literals(ct.enforcement_literal()); - if (all_booleans && enforcement_literals.empty()) { - // TODO(user): we should probably also implement an - // half-reified version of this constraint. + if (all_booleans) { std::vector cst; for (int i = 0; i < vars.size(); ++i) { const int ref = ct.linear().vars(i); @@ -1405,7 +1400,7 @@ void LoadLinearConstraint(const ConstraintProto& ct, Model* m) { } m->GetOrCreate()->AddLinearConstraint( /*use_lower_bound=*/(min_sum < lb), lb, - /*use_upper_bound=*/(max_sum > ub), ub, &cst); + /*use_upper_bound=*/(max_sum > ub), ub, &enforcement_literals, &cst); } else { if (min_sum < lb) { AddWeightedSumGreaterOrEqual(enforcement_literals, vars, coeffs, lb, m); diff --git a/ortools/sat/cp_model_presolve.cc b/ortools/sat/cp_model_presolve.cc index b157dd8433..611ef32088 100644 --- a/ortools/sat/cp_model_presolve.cc +++ b/ortools/sat/cp_model_presolve.cc @@ -3184,7 +3184,7 @@ bool CpModelPresolver::PresolveDiophantine(ConstraintProto* ct) { // TODO(user): Make sure the newly generated linear constraint // satisfy our no-overflow precondition on the min/max activity. // We should check that the model still satisfy conditions in - // https://source.corp.google.com/piper///depot/ortools/sat/cp_model_checker.cc;l=165;bpv=0 + // `PossibleIntegerOverflow` (sat/cp_model_checker.cc) // Create new variables. std::vector new_variables(num_new_variables); diff --git a/ortools/sat/cp_model_search.cc b/ortools/sat/cp_model_search.cc index 9d8011422b..2876df4c9a 100644 --- a/ortools/sat/cp_model_search.cc +++ b/ortools/sat/cp_model_search.cc @@ -168,18 +168,10 @@ struct VarValue { namespace { // TODO(user): Save this somewhere instead of recomputing it. -bool ModelHasSchedulingConstraints(const CpModelProto& cp_model_proto, - const SatParameters& params) { - const bool use_cumulative_relaxation_for_no_overlap_2d = - params.use_timetabling_in_no_overlap_2d() || - params.use_energetic_reasoning_in_no_overlap_2d(); - +bool ModelHasSchedulingConstraints(const CpModelProto& cp_model_proto) { for (const ConstraintProto& ct : cp_model_proto.constraints()) { if (ct.constraint_case() == ConstraintProto::kNoOverlap) return true; if (ct.constraint_case() == ConstraintProto::kCumulative) return true; - if (ct.constraint_case() == ConstraintProto::kNoOverlap2D && - use_cumulative_relaxation_for_no_overlap_2d) - return true; } return false; } @@ -197,8 +189,6 @@ void AddExtraSchedulingPropagators(SatParameters& new_params) { new_params.set_use_area_energetic_reasoning_in_no_overlap_2d(true); new_params.set_use_try_edge_reasoning_in_no_overlap_2d(true); new_params.set_no_overlap_2d_boolean_relations_limit(100); - new_params.set_use_dynamic_precedence_in_cumulative(true); - new_params.set_use_dynamic_precedence_in_disjunctive(true); } // We want a random tie breaking among variables with equivalent values. @@ -361,12 +351,12 @@ std::function ConstructUserSearchStrategy( }; } -// TODO(user): Implement a routing search. +// TODO(user): Implement a routing search strategy. std::function ConstructHeuristicSearchStrategy( const CpModelProto& cp_model_proto, Model* model) { - const auto& params = *model->GetOrCreate(); - if (ModelHasSchedulingConstraints(cp_model_proto, params)) { + if (ModelHasSchedulingConstraints(cp_model_proto)) { std::vector> heuristics; + const auto& params = *model->GetOrCreate(); bool possible_new_constraints = false; if (params.use_dynamic_precedence_in_disjunctive()) { possible_new_constraints = true; @@ -387,9 +377,10 @@ std::function ConstructHeuristicSearchStrategy( } heuristics.push_back(SchedulingSearchHeuristic(model)); + CHECK(!heuristics.empty()); return SequentialSearch(std::move(heuristics)); } - return PseudoCost(model); + return nullptr; } std::function @@ -440,7 +431,7 @@ std::function ConstructHintSearchStrategy( std::function ConstructFixedSearchStrategy( std::function user_search, std::function heuristic_search, - std::function integer_completion) { + std::function integer_completion, Model* model) { // We start by the user specified heuristic. std::vector> heuristics; if (user_search != nullptr) { @@ -449,11 +440,14 @@ std::function ConstructFixedSearchStrategy( if (heuristic_search != nullptr) { heuristics.push_back(heuristic_search); } + if (heuristics.empty()) { + heuristics.push_back(PseudoCost(model)); + } if (integer_completion != nullptr) { heuristics.push_back(integer_completion); } - return SequentialSearch(heuristics); + return SequentialSearch(std::move(heuristics)); } std::function InstrumentSearchStrategy( @@ -688,6 +682,29 @@ absl::flat_hash_map GetNamedParameters( new_params.set_use_dynamic_precedence_in_disjunctive(false); new_params.set_use_dynamic_precedence_in_cumulative(false); strategies["fixed"] = new_params; + + new_params.set_linearization_level(0); + strategies["fixed_no_lp"] = new_params; + + new_params.set_linearization_level(2); + new_params.set_add_lp_constraints_lazily(false); + new_params.set_root_lp_iterations(100'000); + strategies["fixed_max_lp"] = new_params; + } + + // Portfolio search. + { + SatParameters new_params = base_params; + new_params.set_search_branching(SatParameters::PORTFOLIO_SEARCH); + strategies["portfolio"] = new_params; + + new_params.set_linearization_level(0); + strategies["portfolio_no_lp"] = new_params; + + new_params.set_linearization_level(2); + new_params.set_add_lp_constraints_lazily(false); + new_params.set_root_lp_iterations(100'000); + strategies["portfolio_max_lp"] = new_params; } // Quick restart. @@ -841,9 +858,8 @@ std::vector GetFullWorkerParameters( // TODO(user): For scheduling, this is important to find good first solution // but afterwards it is not really great and should probably be replaced by a // LNS worker. - const bool use_fixed_strategy = - !cp_model.search_strategy().empty() || - ModelHasSchedulingConstraints(cp_model, base_params); + const bool use_fixed_strategy = !cp_model.search_strategy().empty() || + ModelHasSchedulingConstraints(cp_model); // Our current set of strategies // diff --git a/ortools/sat/cp_model_search.h b/ortools/sat/cp_model_search.h index 9c2dc3cfc4..a14f192619 100644 --- a/ortools/sat/cp_model_search.h +++ b/ortools/sat/cp_model_search.h @@ -94,7 +94,7 @@ std::function ConstructHintSearchStrategy( std::function ConstructFixedSearchStrategy( std::function user_search, std::function heuristic_search, - std::function integer_completion); + std::function integer_completion, Model* model); // For debugging fixed-search: display information about the named variables // domain before taking each decision. Note that we copy the instrumented diff --git a/ortools/sat/cp_model_solver_helpers.cc b/ortools/sat/cp_model_solver_helpers.cc index 1e8f696ffc..d2aaaecb67 100644 --- a/ortools/sat/cp_model_solver_helpers.cc +++ b/ortools/sat/cp_model_solver_helpers.cc @@ -180,7 +180,7 @@ void InitializeDebugSolution(const CpModelProto& model_proto, Model* model) { // We also register a DEBUG callback to check our reasons. auto* encoder = model->GetOrCreate(); - const auto checker = [mapping, encoder, debug_sol, model]( + const auto checker = [mapping = mapping, encoder, debug_sol, model]( absl::Span clause, absl::Span integers) { bool is_satisfied = false; @@ -257,7 +257,8 @@ void InitializeDebugSolution(const CpModelProto& model_proto, Model* model) { } return is_satisfied; }; - const auto lit_checker = [checker](absl::Span clause) { + const auto lit_checker = [checker = + checker](absl::Span clause) { return checker(clause, {}); }; @@ -779,8 +780,8 @@ void RegisterVariableBoundsLevelZeroImport( const int id = shared_bounds_manager->RegisterNewId(); const auto& import_level_zero_bounds = [&model_proto, shared_bounds_manager, - name, sat_solver, integer_trail, - trail, id, mapping]() { + name = name, sat_solver, + integer_trail, trail, id, mapping]() { std::vector model_variables; std::vector new_lower_bounds; std::vector new_upper_bounds; @@ -941,8 +942,8 @@ void RegisterObjectiveBoundsImport( auto* integer_trail = model->GetOrCreate(); auto* objective = model->GetOrCreate(); const std::string name = model->Name(); - const auto import_objective_bounds = [name, solver, integer_trail, objective, - shared_response_manager]() { + const auto import_objective_bounds = [name = name, solver, integer_trail, + objective, shared_response_manager]() { if (solver->AssumptionLevel() != 0) return true; bool tighter_bounds = false; @@ -1102,6 +1103,7 @@ int RegisterClausesLevelZeroImport(int id, } } clause_manager->SetAddClauseCallback(std::move(callback)); + if (new_clauses > 0 && !sat_solver->FinishPropagation()) return false; if (minimize_shared_clauses && new_clauses > 0) { // The new clauses may be subsumed, so try to minimize them to reduce // overhead of sharing. @@ -1596,7 +1598,7 @@ void LoadCpModel(const CpModelProto& model_proto, Model* model) { objective_var, model); search_heuristics->fixed_search = ConstructFixedSearchStrategy( search_heuristics->user_search, search_heuristics->heuristic_search, - search_heuristics->integer_completion_search); + search_heuristics->integer_completion_search, model); if (VLOG_IS_ON(3)) { search_heuristics->fixed_search = InstrumentSearchStrategy(model_proto, mapping->GetVariableMapping(), @@ -1803,8 +1805,10 @@ void QuickSolveWithHint(const CpModelProto& model_proto, Model* model) { parameters->set_search_branching(SatParameters::HINT_SEARCH); parameters->set_optimize_with_core(false); parameters->set_use_sat_inprocessing(false); - auto cleanup = ::absl::MakeCleanup( - [parameters, saved_params]() { *parameters = saved_params; }); + auto cleanup = + ::absl::MakeCleanup([parameters, saved_params = saved_params]() { + *parameters = saved_params; + }); // Solve decision problem. ConfigureSearchHeuristics(model); diff --git a/ortools/sat/diffn.cc b/ortools/sat/diffn.cc index 03faa7f93e..6b1ef7d2aa 100644 --- a/ortools/sat/diffn.cc +++ b/ortools/sat/diffn.cc @@ -46,6 +46,7 @@ #include "ortools/sat/intervals.h" #include "ortools/sat/model.h" #include "ortools/sat/no_overlap_2d_helper.h" +#include "ortools/sat/pb_constraint.h" #include "ortools/sat/sat_base.h" #include "ortools/sat/sat_parameters.pb.h" #include "ortools/sat/sat_solver.h" @@ -158,9 +159,6 @@ void AddDiffnCumulativeRelationOnX( model->GetOrCreate()->GetOrCreateDemandHelper( x, y->Sizes()); - model->GetOrCreate()->RegisterCumulative( - {.capacity = capacity, .task_helper = x, .demand_helper = demands}); - // Propagator responsible for applying Timetabling filtering rule. It // increases the minimum of the start variables, decrease the maximum of the // end variables, and increase the minimum of the capacity variable. @@ -179,6 +177,23 @@ void AddDiffnCumulativeRelationOnX( } } +bool AddAtMostOne(absl::Span enforcement_literals, + absl::Span literals, Model* model) { + if (enforcement_literals.empty()) { + return model->GetOrCreate()->AddAtMostOne(literals); + } + std::vector enforcement(enforcement_literals.begin(), + enforcement_literals.end()); + std::vector cst; + cst.reserve(literals.size()); + for (const Literal l : literals) { + cst.emplace_back(l, Coefficient(1)); + } + return model->GetOrCreate()->AddLinearConstraint( + /*use_lower_bound=*/false, Coefficient(0), + /*use_upper_bound=*/true, Coefficient(1), &enforcement, &cst); +} + } // namespace void AddNonOverlappingRectangles( @@ -283,14 +298,8 @@ void AddNonOverlappingRectangles( // of the general repository, however, as we add constraints, this // propagates which might swap/change the underlying x_helper and // y_helper... - // - // TODO(user): add support for enforcement_literals. The - // AddProblemClause below is easy to extend, but the AddAtMostOne calls - // are trickier... const int num_boxes = x.size(); - if (num_boxes < params.no_overlap_2d_boolean_relations_limit() && - enforcement_literals.empty()) { - auto* implications = model->GetOrCreate(); + if (num_boxes < params.no_overlap_2d_boolean_relations_limit()) { auto* sat_solver = model->GetOrCreate(); auto* integer_trail = model->GetOrCreate(); DCHECK_EQ(sat_solver->CurrentDecisionLevel(), 0); @@ -320,7 +329,7 @@ void AddNonOverlappingRectangles( repository->End(x[j]), repository->Start(x[i])); if ((integer_trail->LowerBound(repository->Size(x[i])) > 0 || integer_trail->LowerBound(repository->Size(x[j])) > 0) && - !implications->AddAtMostOne({x_ij, x_ji})) { + !AddAtMostOne(enforcement_literals, {x_ij, x_ji}, model)) { sat_solver->NotifyThatModelIsUnsat(); return; } @@ -333,7 +342,7 @@ void AddNonOverlappingRectangles( repository->End(y[j]), repository->Start(y[i])); if ((integer_trail->LowerBound(repository->Size(y[i])) > 0 || integer_trail->LowerBound(repository->Size(y[j])) > 0) && - !implications->AddAtMostOne({y_ij, y_ji})) { + !AddAtMostOne(enforcement_literals, {y_ij, y_ji}, model)) { sat_solver->NotifyThatModelIsUnsat(); return; } @@ -352,9 +361,8 @@ void AddNonOverlappingRectangles( if (repository->IsOptional(y[j])) { clause.push_back(repository->PresenceLiteral(y[j]).Negated()); } - if (!sat_solver->AddProblemClause(clause)) { - return; - } + model->Add(EnforcedClause(enforcement_literals, clause)); + if (sat_solver->ModelIsUnsat()) return; } } } diff --git a/ortools/sat/enforcement.cc b/ortools/sat/enforcement.cc index 44f551e73a..55ff0264ee 100644 --- a/ortools/sat/enforcement.cc +++ b/ortools/sat/enforcement.cc @@ -293,21 +293,25 @@ void EnforcementPropagator::ChangeStatus(EnforcementId id, if (callbacks_[id] != nullptr) callbacks_[id](id, new_status); } -EnforcementStatus EnforcementPropagator::DebugStatus(EnforcementId id) { - if (id < 0) return EnforcementStatus::IS_ENFORCED; - +EnforcementStatus EnforcementPropagator::Status( + absl::Span enforcement) const { int num_true = 0; - for (const Literal l : GetSpan(id)) { + for (const Literal l : enforcement) { if (assignment_.LiteralIsFalse(l)) { return EnforcementStatus::IS_FALSE; } if (assignment_.LiteralIsTrue(l)) ++num_true; } - const int size = GetSpan(id).size(); + const int size = enforcement.size(); if (num_true == size) return EnforcementStatus::IS_ENFORCED; if (num_true + 1 == size) return EnforcementStatus::CAN_PROPAGATE_ENFORCEMENT; return EnforcementStatus::CANNOT_PROPAGATE; } +EnforcementStatus EnforcementPropagator::DebugStatus(EnforcementId id) { + if (id < 0) return EnforcementStatus::IS_ENFORCED; + return Status(GetSpan(id)); +} + } // namespace sat } // namespace operations_research diff --git a/ortools/sat/enforcement.h b/ortools/sat/enforcement.h index 1c2624f22d..8805b925be 100644 --- a/ortools/sat/enforcement.h +++ b/ortools/sat/enforcement.h @@ -77,6 +77,10 @@ class EnforcementPropagator : public SatPropagator { return statuses_[id]; } + // This recomputes the current status by scanning the given list. + // It thus has a linear complexity and should not be used in hot loops. + EnforcementStatus Status(absl::Span enforcement) const; + // Recompute the status from the current assignment. // This should only used in DCHECK(). EnforcementStatus DebugStatus(EnforcementId id); diff --git a/ortools/sat/integer_search.cc b/ortools/sat/integer_search.cc index eef82f00b7..65b0396422 100644 --- a/ortools/sat/integer_search.cc +++ b/ortools/sat/integer_search.cc @@ -18,6 +18,7 @@ #include #include #include +#include #include #include @@ -259,7 +260,7 @@ std::function LpPseudoCostHeuristic(Model* model) { if (info.score > best_score) { best_score = info.score; - // This direction works better than the inverse in the benchs. But + // This direction works better than the inverse in the benchmarks. But // always branching up seems even better. TODO(user): investigate. if (info.down_score > info.up_score) { decision = BooleanOrIntegerLiteral(info.down_branch); @@ -296,6 +297,9 @@ UnassignedVarWithLowestMinAtItsMinHeuristic( std::function SequentialSearch( std::vector> heuristics) { + for (const auto& h : heuristics) { + CHECK(h != nullptr); + } return [heuristics]() { for (const auto& h : heuristics) { const BooleanOrIntegerLiteral decision = h(); @@ -990,10 +994,21 @@ std::function RandomizeOnRestartHeuristic( weights.push_back(1); } - // Always add heuristic search. - policies.push_back(SequentialSearch({heuristics.heuristic_search, sat_policy, - heuristics.integer_completion_search})); - weights.push_back(1); + // Add heuristic search if present. + if (heuristics.heuristic_search != nullptr) { + policies.push_back( + SequentialSearch({heuristics.heuristic_search, sat_policy, + heuristics.integer_completion_search})); + weights.push_back(1); + } + + if (policies.size() == 1) { + CHECK(heuristics.fixed_search != nullptr); + policies.push_back( + SequentialSearch({heuristics.fixed_search, sat_policy, + heuristics.integer_completion_search})); + weights.push_back(1); + } // The higher weight for the sat policy is because this policy actually // contains a lot of variation as we randomize the sat parameters. diff --git a/ortools/sat/linear_propagation.cc b/ortools/sat/linear_propagation.cc index b2e30467f1..8e23d9087f 100644 --- a/ortools/sat/linear_propagation.cc +++ b/ortools/sat/linear_propagation.cc @@ -458,13 +458,11 @@ bool LinearPropagator::AddConstraint( // Propagate this new constraint. // TODO(user): Do we want to do that? - num_terms_for_dtime_update_ = 0; - const auto cleanup = ::absl::MakeCleanup([this]() { - time_limit_->AdvanceDeterministicTime( - static_cast(num_terms_for_dtime_update_) * 1e-9); - }); - if (!PropagateOneConstraint(id)) return false; - return true; + // + // Tricky: we cannot just call PropagateOneConstraint(id) if some variable + // bounds where modified since the last time we propagated, otherwise we + // might wrongly detect cycles for instance. + return Propagate(); } absl::Span LinearPropagator::GetCoeffs( diff --git a/ortools/sat/pb_constraint.cc b/ortools/sat/pb_constraint.cc index 26cfa62dd8..0e8616f703 100644 --- a/ortools/sat/pb_constraint.cc +++ b/ortools/sat/pb_constraint.cc @@ -15,11 +15,14 @@ #include #include +#include #include #include +#include #include #include +#include "absl/cleanup/cleanup.h" #include "absl/container/flat_hash_map.h" #include "absl/hash/hash.h" #include "absl/log/check.h" @@ -27,6 +30,7 @@ #include "absl/types/span.h" #include "ortools/base/logging.h" #include "ortools/base/strong_vector.h" +#include "ortools/sat/enforcement.h" #include "ortools/sat/sat_base.h" #include "ortools/sat/sat_parameters.pb.h" #include "ortools/util/bitset.h" @@ -61,7 +65,7 @@ bool ComputeBooleanLinearExpressionCanonicalForm( *max_value = 0; // First, sort by literal to remove duplicate literals. - // This also remove term with a zero coefficient. + // This also removes terms with a zero coefficient. std::sort(cst->begin(), cst->end(), LiteralComparator); int index = 0; LiteralWithCoeff* representative = nullptr; @@ -148,7 +152,12 @@ bool ApplyLiteralMapping( // TODO(user): Also check for no duplicates literals + unit tests. bool BooleanLinearExpressionIsCanonical( + absl::Span enforcement_literals, absl::Span cst) { + if (!std::is_sorted(enforcement_literals.begin(), + enforcement_literals.end())) { + return false; + } Coefficient previous(1); for (LiteralWithCoeff term : cst) { if (term.coefficient < previous) return false; @@ -200,12 +209,12 @@ Coefficient ComputeNegatedCanonicalRhs(Coefficient lower_bound, // Positive overflow. The constraint is infeasible. return Coefficient(-1); } else { - // Negative overflow. The constraint is trivialy satisfiable. + // Negative overflow. The constraint is trivially satisfiable. return max_value; } } if (shifted_lb <= 0) { - // If shifted_lb <= 0 then the constraint is trivialy satisfiable. We test + // If shifted_lb <= 0 then the constraint is trivially satisfiable. We test // this so we are sure that max_value - shifted_lb doesn't overflow below. return max_value; } @@ -353,7 +362,7 @@ void MutableUpperBoundedLinearConstraint::ReduceSlackTo( CHECK_LE(target, slack); CHECK_GE(target, 0); - // This is not stricly needed, but true in our use case: + // This is not strictly needed, but true in our use case: // The variable assigned at trail_index was causing a conflict. const Coefficient coeff = GetCoefficient(trail[trail_index].Variable()); CHECK_LT(slack, coeff); @@ -402,12 +411,16 @@ Coefficient MutableUpperBoundedLinearConstraint::ComputeMaxSum() const { } UpperBoundedLinearConstraint::UpperBoundedLinearConstraint( + const std::vector& enforcement_literals, const std::vector& cst) : is_marked_for_deletion_(false), is_learned_(false), first_reason_trail_index_(-1), - activity_(0.0) { + activity_(0.0), + enforcement_id_(-1) { DCHECK(!cst.empty()); + DCHECK( + std::is_sorted(enforcement_literals.begin(), enforcement_literals.end())); DCHECK(std::is_sorted(cst.begin(), cst.end(), CoeffComparator)); literals_.reserve(cst.size()); @@ -438,11 +451,14 @@ UpperBoundedLinearConstraint::UpperBoundedLinearConstraint( // Sentinel. starts_.push_back(literals_.size()); - hash_ = absl::Hash>()(cst); + hash_ = absl::Hash< + std::pair, std::vector>>()( + {enforcement_literals, cst}); } void UpperBoundedLinearConstraint::AddToConflict( MutableUpperBoundedLinearConstraint* conflict) { + CHECK_LT(enforcement_id_, 0) << "Enforcement literals are not supported"; int literal_index = 0; int coeff_index = 0; for (Literal literal : literals_) { @@ -453,8 +469,14 @@ void UpperBoundedLinearConstraint::AddToConflict( conflict->AddToRhs(rhs_); } -bool UpperBoundedLinearConstraint::HasIdenticalTerms( - absl::Span cst) { +bool UpperBoundedLinearConstraint::HasIdenticalTermsAndEnforcement( + absl::Span enforcement_literals, + absl::Span cst, + EnforcementPropagator* enforcement_propagator) { + if (enforcement_literals != + enforcement_propagator->GetEnforcementLiterals(enforcement_id_)) { + return false; + } if (cst.size() != literals_.size()) return false; int literal_index = 0; int coeff_index = 0; @@ -470,7 +492,9 @@ bool UpperBoundedLinearConstraint::HasIdenticalTerms( } bool UpperBoundedLinearConstraint::InitializeRhs( - Coefficient rhs, int trail_index, Coefficient* threshold, Trail* trail, + EnforcementStatus enforcement_status, + absl::Span enforcement_literals, Coefficient rhs, + int trail_index, Coefficient* threshold, Trail* trail, PbConstraintsEnqueueHelper* helper) { // Compute the slack from the assigned variables with a trail index // smaller than the given trail_index. The variable at trail_index has not @@ -504,7 +528,9 @@ bool UpperBoundedLinearConstraint::InitializeRhs( } // The constraint is infeasible provided the current propagated trail. - if (slack < 0) return false; + if (slack < 0 && enforcement_literals.empty()) { + return false; + } // Cumulative sum. for (int i = 1; i < sum_at_previous_level.size(); ++i) { @@ -538,18 +564,50 @@ bool UpperBoundedLinearConstraint::InitializeRhs( index_ = coeffs_.size() - 1; already_propagated_end_ = literals_.size(); Update(slack, threshold); - return *threshold < 0 - ? Propagate(max_relevant_trail_index, threshold, trail, helper) - : true; + if (enforcement_status == EnforcementStatus::IS_FALSE || + enforcement_status == EnforcementStatus::CANNOT_PROPAGATE || + *threshold >= 0) { + return true; + } + return Propagate(max_relevant_trail_index, threshold, trail, + enforcement_status, enforcement_literals, helper); } bool UpperBoundedLinearConstraint::Propagate( int trail_index, Coefficient* threshold, Trail* trail, - PbConstraintsEnqueueHelper* helper) { + EnforcementStatus enforcement_status, + absl::Span enforcement_literals, + PbConstraintsEnqueueHelper* helper, bool* need_untrail_inspection) { DCHECK_LT(*threshold, 0); + DCHECK_NE(enforcement_status, EnforcementStatus::IS_FALSE); + DCHECK_NE(enforcement_status, EnforcementStatus::CANNOT_PROPAGATE); const Coefficient slack = GetSlackFromThreshold(*threshold); - DCHECK_GE(slack, 0) << "The constraint is already a conflict!"; + + if (slack < 0) { + // This can happen if the slack became negative while the status was + // IS_FALSE or CANNOT_PROPAGATE. + DCHECK(!enforcement_literals.empty()); + if (enforcement_status == EnforcementStatus::CAN_PROPAGATE_ENFORCEMENT) { + // Propagate the unique unassigned enforcement literal. + for (const Literal literal : enforcement_literals) { + if (!trail->Assignment().LiteralIsAssigned(literal)) { + helper->Enqueue(literal.Negated(), trail_index, this, trail); + break; + } + } + return true; + } else { + // Conflict. + FillReason(*trail, trail_index, enforcement_literals, + /*propagated_variable=*/kNoBooleanVariable, &helper->conflict); + return false; + } + } + while (index_ >= 0 && coeffs_[index_] > slack) --index_; + if (need_untrail_inspection != nullptr) { + *need_untrail_inspection = true; + } // Check propagation. BooleanVariable first_propagated_variable(-1); @@ -557,14 +615,26 @@ bool UpperBoundedLinearConstraint::Propagate( if (trail->Assignment().LiteralIsFalse(literals_[i])) continue; if (trail->Assignment().LiteralIsTrue(literals_[i])) { if (trail->Info(literals_[i].Variable()).trail_index > trail_index) { - // Conflict. - FillReason(*trail, trail_index, literals_[i].Variable(), - &helper->conflict); - helper->conflict.push_back(literals_[i].Negated()); - Update(slack, threshold); - return false; + if (enforcement_status == EnforcementStatus::IS_ENFORCED) { + // Conflict. + FillReason(*trail, trail_index, enforcement_literals, + literals_[i].Variable(), &helper->conflict); + helper->conflict.push_back(literals_[i].Negated()); + Update(slack, threshold); + return false; + } else { + // Propagate the unique unassigned enforcement literal. + for (const Literal literal : enforcement_literals) { + if (!trail->Assignment().LiteralIsAssigned(literal)) { + helper->Enqueue(literal.Negated(), trail_index, this, trail); + break; + } + } + Update(slack, threshold); + return true; + } } - } else { + } else if (enforcement_status == EnforcementStatus::IS_ENFORCED) { // Propagation. if (first_propagated_variable < 0) { if (first_reason_trail_index_ == -1) { @@ -588,72 +658,143 @@ bool UpperBoundedLinearConstraint::Propagate( void UpperBoundedLinearConstraint::FillReason( const Trail& trail, int source_trail_index, + absl::Span enforcement_literals, BooleanVariable propagated_variable, std::vector* reason) { + bool enforcement_propagation = false; reason->clear(); - - // Optimization for an "at most one" constraint. - // Note that the source_trail_index set by InitializeRhs() is ok in this case. - if (rhs_ == 1) { - reason->push_back(trail[source_trail_index].Negated()); - return; + for (const Literal literal : enforcement_literals) { + if (trail.Assignment().LiteralIsTrue(literal)) { + reason->push_back(literal.Negated()); + } else if (literal.Variable() == propagated_variable) { + enforcement_propagation = true; + } } + const int enforcement_reason_size = reason->size(); + Coefficient slack = rhs_; + Coefficient propagated_variable_coefficient(0); + Literal extra_literal_reason; // Optimization: This will be set to the index of the last literal in the - // reason. + // reason (they are sorted in decreasing coefficient order). int last_i = 0; int last_coeff_index = 0; - // Compute the initial reason which is formed by all the literals of the - // constraint that were assigned to true at the time of the propagation. - // We remove literals with a level of 0 since they are not needed. - // We also compute the slack at the time. - Coefficient slack = rhs_; - Coefficient propagated_variable_coefficient(0); - int coeff_index = coeffs_.size() - 1; - for (int i = literals_.size() - 1; i >= 0; --i) { - const Literal literal = literals_[i]; - if (literal.Variable() == propagated_variable) { - propagated_variable_coefficient = coeffs_[coeff_index]; - } else { - if (trail.Assignment().LiteralIsTrue(literal) && - trail.Info(literal.Variable()).trail_index <= source_trail_index) { - if (trail.Info(literal.Variable()).level > 0) { - reason->push_back(literal.Negated()); - last_i = i; - last_coeff_index = coeff_index; - } - slack -= coeffs_[coeff_index]; + // propagated_variable is set to kNoBooleanVariable when the constraint + // becomes enforced when the slack is already negative. In this case, or when + // the enforcement can be propagated, the reason must include all the terms + // explaining the negative slack. The code below does that (while the 'else' + // part computes a reason for propagated_variable). extra_literal_reason is + // the literal which makes the slack become negative. + if (enforcement_propagation || propagated_variable == kNoBooleanVariable) { + int literal_index = 0; + int coeff_index = 0; + // Vector of (trail_index, literal_index, coeff_index) tuples. + std::vector> true_literals; + for (Literal literal : literals_) { + if (trail.Assignment().LiteralIsTrue(literal)) { + true_literals.push_back({trail.Info(literal.Variable()).trail_index, + literal_index, coeff_index}); } + ++literal_index; + if (literal_index == starts_[coeff_index + 1]) ++coeff_index; } - if (i == starts_[coeff_index]) { - --coeff_index; + std::sort(true_literals.begin(), true_literals.end()); + // Vector of (literal_index, coeff_index) pairs. + std::vector> reason_indices; + for (const auto& [trail_index, literal_index, coeff_index] : + true_literals) { + const Literal literal = literals_[literal_index]; + const Coefficient coeff = coeffs_[coeff_index]; + if (coeff > slack) { + propagated_variable_coefficient = coeff; + // This literal is added to the reason at the very end (see the cleanup + // below) because the code minimizing the reason assumes that it is not + // part of it. Another solution would be insert it at the beginning (and + // to increment enforcement_reason_size), but this is less efficient. + extra_literal_reason = literal.Negated(); + break; + } + if (trail.Info(literal.Variable()).level > 0) { + reason_indices.push_back({literal_index, coeff_index}); + } + slack -= coeff.value(); + } + std::sort(reason_indices.begin(), reason_indices.end(), std::greater<>()); + for (const auto& [literal_index, coeff_index] : reason_indices) { + reason->push_back(literals_[literal_index].Negated()); + } + if (!reason_indices.empty()) { + last_i = reason_indices.back().first; + last_coeff_index = reason_indices.back().second; + } + } else { + // Optimization for an "at most one" constraint. Note that the + // source_trail_index set by InitializeRhs() is ok in this case. + if (rhs_ == 1) { + reason->push_back(trail[source_trail_index].Negated()); + return; + } + + // Compute the initial reason which is formed by all the literals of the + // constraint that were assigned to true at the time of the propagation. + // We remove literals with a level of 0 since they are not needed. + // We also compute the slack at the time. + int coeff_index = coeffs_.size() - 1; + for (int i = literals_.size() - 1; i >= 0; --i) { + const Literal literal = literals_[i]; + if (literal.Variable() == propagated_variable) { + propagated_variable_coefficient = coeffs_[coeff_index]; + } else { + if (trail.Assignment().LiteralIsTrue(literal) && + trail.Info(literal.Variable()).trail_index <= source_trail_index) { + if (trail.Info(literal.Variable()).level > 0) { + reason->push_back(literal.Negated()); + last_i = i; + last_coeff_index = coeff_index; + } + slack -= coeffs_[coeff_index]; + } + } + if (i == starts_[coeff_index]) { + --coeff_index; + } } } DCHECK_GT(propagated_variable_coefficient, slack); DCHECK_GE(propagated_variable_coefficient, 0); + auto cleanup = absl::MakeCleanup([&] { + if (enforcement_propagation || propagated_variable == kNoBooleanVariable) { + reason->push_back(extra_literal_reason); + } + }); // In both cases, we can't minimize the reason further. - if (reason->size() <= 1 || coeffs_.size() == 1) return; + if (reason->size() <= enforcement_reason_size + 1 || coeffs_.size() == 1) { + return; + } Coefficient limit = propagated_variable_coefficient - slack; DCHECK_GE(limit, 1); // Remove literals with small coefficients from the reason as long as the - // limit is still stricly positive. - coeff_index = last_coeff_index; - if (coeffs_[coeff_index] >= limit) return; + // limit is still strictly positive. + int coeff_index = last_coeff_index; + if (coeffs_[coeff_index] >= limit) { + return; + } for (int i = last_i; i < literals_.size(); ++i) { const Literal literal = literals_[i]; if (i == starts_[coeff_index + 1]) { ++coeff_index; if (coeffs_[coeff_index] >= limit) break; } + DCHECK_GT(reason->size(), enforcement_reason_size); if (literal.Negated() != reason->back()) continue; limit -= coeffs_[coeff_index]; reason->pop_back(); + if (reason->size() == enforcement_reason_size) break; if (coeffs_[coeff_index] >= limit) break; } - DCHECK(!reason->empty()); DCHECK_GE(limit, 1); } @@ -678,6 +819,7 @@ void UpperBoundedLinearConstraint::ResolvePBConflict( const Trail& trail, BooleanVariable var, MutableUpperBoundedLinearConstraint* conflict, Coefficient* conflict_slack) { + CHECK_LT(enforcement_id_, 0) << "Enforcement literals are not supported"; const int limit_trail_index = trail.Info(var).trail_index; // Compute the constraint activity at the time and the coefficient of the @@ -834,10 +976,13 @@ void UpperBoundedLinearConstraint::Untrail(Coefficient* threshold, // TODO(user): This is relatively slow. Take the "transpose" all at once, and // maybe put small constraints first on the to_update_ lists. -bool PbConstraints::AddConstraint(const std::vector& cst, - Coefficient rhs, Trail* trail) { +bool PbConstraints::AddConstraint( + const std::vector& enforcement_literals, + const std::vector& cst, Coefficient rhs, Trail* trail) { SCOPED_TIME_STAT(&stats_); DCHECK(!cst.empty()); + DCHECK( + std::is_sorted(enforcement_literals.begin(), enforcement_literals.end())); DCHECK(std::is_sorted(cst.begin(), cst.end(), CoeffComparator)); // Special case if this is the first constraint. @@ -849,13 +994,14 @@ bool PbConstraints::AddConstraint(const std::vector& cst, } std::unique_ptr c( - new UpperBoundedLinearConstraint(cst)); + new UpperBoundedLinearConstraint(enforcement_literals, cst)); std::vector& duplicate_candidates = possible_duplicates_[c->hash()]; // Optimization if the constraint terms are duplicates. for (UpperBoundedLinearConstraint* candidate : duplicate_candidates) { - if (candidate->HasIdenticalTerms(cst)) { + if (candidate->HasIdenticalTermsAndEnforcement(enforcement_literals, cst, + enforcement_propagator_)) { if (rhs < candidate->Rhs()) { // TODO(user): the index is needed to give the correct thresholds_ entry // to InitializeRhs() below, but this linear scan is not super @@ -866,9 +1012,10 @@ bool PbConstraints::AddConstraint(const std::vector& cst, ++i; } CHECK_LT(i, constraints_.size()); - return candidate->InitializeRhs(rhs, propagation_trail_index_, - &thresholds_[i], trail, - &enqueue_helper_); + return candidate->InitializeRhs( + enforcement_propagator_->Status(candidate->enforcement_id()), + enforcement_literals, rhs, propagation_trail_index_, + &thresholds_[i], trail, &enqueue_helper_); } else { // The constraint is redundant, so there is nothing to do. return true; @@ -877,13 +1024,27 @@ bool PbConstraints::AddConstraint(const std::vector& cst, } thresholds_.push_back(Coefficient(0)); - if (!c->InitializeRhs(rhs, propagation_trail_index_, &thresholds_.back(), - trail, &enqueue_helper_)) { + const EnforcementStatus enforcement_status = + enforcement_propagator_->Status(enforcement_literals); + if (!c->InitializeRhs(enforcement_status, enforcement_literals, rhs, + propagation_trail_index_, &thresholds_.back(), trail, + &enqueue_helper_)) { thresholds_.pop_back(); return false; } const ConstraintIndex cst_index(constraints_.size()); + enforcement_status_changed_.Resize(cst_index + 1); + if (!enforcement_literals.empty()) { + c->set_enforcement_id(enforcement_propagator_->Register( + enforcement_literals, + [&, cst_index](EnforcementId, EnforcementStatus status) { + if (status == EnforcementStatus::IS_ENFORCED || + status == EnforcementStatus::CAN_PROPAGATE_ENFORCEMENT) { + enforcement_status_changed_.Set(cst_index); + } + })); + } duplicate_candidates.push_back(c.get()); constraints_.emplace_back(c.release()); for (LiteralWithCoeff term : cst) { @@ -899,7 +1060,8 @@ bool PbConstraints::AddLearnedConstraint( const std::vector& cst, Coefficient rhs, Trail* trail) { DeleteSomeLearnedConstraintIfNeeded(); const int old_num_constraints = constraints_.size(); - const bool result = AddConstraint(cst, rhs, trail); + const bool result = + AddConstraint(/*enforcement_literals=*/{}, cst, rhs, trail); // The second test is to avoid marking a problem constraint as learned because // of the "reuse last constraint" optimization. @@ -909,6 +1071,36 @@ bool PbConstraints::AddLearnedConstraint( return result; } +bool PbConstraints::PropagateConstraint(ConstraintIndex index, Trail* trail, + int source_trail_index, + bool* need_untrail_inspection) { + UpperBoundedLinearConstraint* const cst = constraints_[index.value()].get(); + const EnforcementStatus enforcement_status = + enforcement_propagator_->Status(cst->enforcement_id()); + if (enforcement_status == EnforcementStatus::IS_FALSE || + enforcement_status == EnforcementStatus::CANNOT_PROPAGATE) { + return true; + } + bool conflict = false; + ++num_constraint_lookups_; + const int old_value = cst->already_propagated_end(); + if (!cst->Propagate(source_trail_index, &thresholds_[index], trail, + enforcement_status, + enforcement_propagator_->GetEnforcementLiterals( + cst->enforcement_id()), + &enqueue_helper_, need_untrail_inspection)) { + trail->MutableConflict()->swap(enqueue_helper_.conflict); + conflicting_constraint_index_ = index; + conflict = true; + + // We bump the activity of the conflict. + BumpActivity(constraints_[index.value()].get()); + } + num_inspected_constraint_literals_ += + old_value - cst->already_propagated_end(); + return !conflict; +} + bool PbConstraints::PropagateNext(Trail* trail) { SCOPED_TIME_STAT(&stats_); const int source_trail_index = propagation_trail_index_; @@ -924,28 +1116,24 @@ bool PbConstraints::PropagateNext(Trail* trail) { thresholds_[update.index] - update.coefficient; thresholds_[update.index] = threshold; if (threshold < 0 && !conflict) { - UpperBoundedLinearConstraint* const cst = - constraints_[update.index.value()].get(); - update.need_untrail_inspection = true; - ++num_constraint_lookups_; - const int old_value = cst->already_propagated_end(); - if (!cst->Propagate(source_trail_index, &thresholds_[update.index], trail, - &enqueue_helper_)) { - trail->MutableConflict()->swap(enqueue_helper_.conflict); - conflicting_constraint_index_ = update.index; - conflict = true; - - // We bump the activity of the conflict. - BumpActivity(constraints_[update.index.value()].get()); - } - num_inspected_constraint_literals_ += - old_value - cst->already_propagated_end(); + conflict = !PropagateConstraint(update.index, trail, source_trail_index, + &update.need_untrail_inspection); } } return !conflict; } bool PbConstraints::Propagate(Trail* trail) { + for (const ConstraintIndex index : + enforcement_status_changed_.PositionsSetAtLeastOnce()) { + if (thresholds_[index] < 0 && + !PropagateConstraint(index, trail, propagation_trail_index_)) { + enforcement_status_changed_.ResetAllToFalse(); + return false; + } + } + enforcement_status_changed_.ResetAllToFalse(); + const int old_index = trail->Index(); while (trail->Index() == old_index && propagation_trail_index_ < old_index) { if (!PropagateNext(trail)) return false; @@ -983,8 +1171,11 @@ absl::Span PbConstraints::Reason(const Trail& trail, const PbConstraintsEnqueueHelper::ReasonInfo& reason_info = enqueue_helper_.reasons[trail_index]; std::vector* reason = trail.GetEmptyVectorToStoreReason(trail_index); - reason_info.pb_constraint->FillReason(trail, reason_info.source_trail_index, - trail[trail_index].Variable(), reason); + reason_info.pb_constraint->FillReason( + trail, reason_info.source_trail_index, + enforcement_propagator_->GetEnforcementLiterals( + reason_info.pb_constraint->enforcement_id()), + trail[trail_index].Variable(), reason); return *reason; } @@ -1023,6 +1214,8 @@ void PbConstraints::DeleteSomeLearnedConstraintIfNeeded() { std::vector activities; for (int i = 0; i < constraints_.size(); ++i) { const UpperBoundedLinearConstraint& constraint = *(constraints_[i].get()); + CHECK_LT(constraint.enforcement_id(), 0) + << "Enforcement literals are not supported"; if (constraint.is_learned() && !constraint.is_used_as_a_reason()) { activities.push_back(constraint.activity()); } @@ -1115,6 +1308,9 @@ void PbConstraints::DeleteConstraintMarkedForDeletion() { } else { // Remove it from possible_duplicates_. UpperBoundedLinearConstraint* c = constraints_[i.value()].get(); + // We only delete learned constraints, and learned constraints never have + // enforcement literals. + CHECK_LT(c->enforcement_id(), 0); std::vector& ref = possible_duplicates_[c->hash()]; for (int i = 0; i < ref.size(); ++i) { diff --git a/ortools/sat/pb_constraint.h b/ortools/sat/pb_constraint.h index 68b208093f..d70b6e4f72 100644 --- a/ortools/sat/pb_constraint.h +++ b/ortools/sat/pb_constraint.h @@ -24,10 +24,10 @@ #include "absl/container/flat_hash_map.h" #include "absl/log/check.h" -#include "absl/strings/string_view.h" #include "absl/types/span.h" #include "ortools/base/logging.h" #include "ortools/base/strong_vector.h" +#include "ortools/sat/enforcement.h" #include "ortools/sat/model.h" #include "ortools/sat/sat_base.h" #include "ortools/sat/sat_parameters.pb.h" @@ -125,8 +125,13 @@ Coefficient ComputeNegatedCanonicalRhs(Coefficient lower_bound, Coefficient bound_shift, Coefficient max_value); -// Returns true iff the Boolean linear expression is in canonical form. -bool BooleanLinearExpressionIsCanonical(absl::Span cst); +// Returns true iff the enforced Boolean linear expression is in canonical form. +// The enforcement literals must be sorted and unique, and cst must be in the +// form returned by ComputeBooleanLinearExpressionCanonicalForm(). Moreover the +// enforcement literals should not appear in cst. +bool BooleanLinearExpressionIsCanonical( + absl::Span enforcement_literals, + absl::Span cst); // Given a Boolean linear constraint in canonical form, simplify its // coefficients using simple heuristics. @@ -211,7 +216,7 @@ class MutableUpperBoundedLinearConstraint { // // If we take a constraint sum ci.xi <= rhs, take its negation and add max_sum // on both side, we have sum ci.(1 - xi) >= max_sum - rhs - // So every ci > (max_sum - rhs) can be replacend by (max_sum - rhs). + // So every ci > (max_sum - rhs) can be replaced by (max_sum - rhs). // Not that this operation also change the original rhs of the constraint. void ReduceCoefficients(); @@ -240,7 +245,7 @@ class MutableUpperBoundedLinearConstraint { // This allow to DCHECK some assumptions on what coefficients can be reduced // or not. // - // TODO(user): Ideally the slack should be maitainable incrementally. + // TODO(user): Ideally the slack should be maintainable incrementally. Coefficient ReduceCoefficientsAndComputeSlackForTrailPrefix( const Trail& trail, int trail_index); @@ -385,11 +390,20 @@ struct PbConstraintsEnqueueHelper { class UpperBoundedLinearConstraint { public: // Takes a pseudo-Boolean formula in canonical form. - explicit UpperBoundedLinearConstraint( - const std::vector& cst); + UpperBoundedLinearConstraint(const std::vector& enforcement_literals, + const std::vector& cst); - // Returns true if the given terms are the same as the one in this constraint. - bool HasIdenticalTerms(absl::Span cst); + EnforcementId enforcement_id() const { return enforcement_id_; }; + void set_enforcement_id(EnforcementId enforcement_id) { + enforcement_id_ = enforcement_id; + } + + // Returns true if the given terms and enforcement literals are the same as + // the one in this constraint. + bool HasIdenticalTermsAndEnforcement( + absl::Span enforcement_literals, + absl::Span cst, + EnforcementPropagator* enforcement_propagator); Coefficient Rhs() const { return rhs_; } // Sets the rhs of this constraint. Compute the initial threshold value using @@ -398,7 +412,9 @@ class UpperBoundedLinearConstraint { // // Returns false if the preconditions described in // PbConstraints::AddConstraint() are not meet. - bool InitializeRhs(Coefficient rhs, int trail_index, Coefficient* threshold, + bool InitializeRhs(EnforcementStatus enforcement_status, + absl::Span enforcement_literals, + Coefficient rhs, int trail_index, Coefficient* threshold, Trail* trail, PbConstraintsEnqueueHelper* helper); // Tests for propagation and enqueues propagated literals on the trail. @@ -415,7 +431,10 @@ class UpperBoundedLinearConstraint { // // The threshold is updated to its new value. bool Propagate(int trail_index, Coefficient* threshold, Trail* trail, - PbConstraintsEnqueueHelper* helper); + EnforcementStatus enforcement_status, + absl::Span enforcement_literals, + PbConstraintsEnqueueHelper* helper, + bool* need_untrail_inspection = nullptr); // Updates the given threshold and the internal state. This is the opposite of // Propagate(). Each time a literal in unassigned, the threshold value must @@ -424,7 +443,7 @@ class UpperBoundedLinearConstraint { void Untrail(Coefficient* threshold, int trail_index); // Provided that the literal with given source_trail_index was the one that - // propagated the conflict or the literal we wants to explain, then this will + // propagated the conflict or the literal we want to explain, then this will // compute the reason. // // Some properties of the reason: @@ -439,6 +458,7 @@ class UpperBoundedLinearConstraint { // better to use during conflict minimization (namely the one already in the // 1-UIP conflict). void FillReason(const Trail& trail, int source_trail_index, + absl::Span enforcement_literals, BooleanVariable propagated_variable, std::vector* reason); @@ -470,7 +490,10 @@ class UpperBoundedLinearConstraint { // Only learned constraints are considered for deletion during the constraint // cleanup phase. We also can't delete variables used as a reason. - void set_is_learned(bool is_learned) { is_learned_ = is_learned; } + void set_is_learned(bool is_learned) { + CHECK(!is_learned || enforcement_id_ < 0); + is_learned_ = is_learned; + } bool is_learned() const { return is_learned_; } bool is_used_as_a_reason() const { return first_reason_trail_index_ != -1; } @@ -488,7 +511,7 @@ class UpperBoundedLinearConstraint { int already_propagated_end() const { return already_propagated_end_; } private: - Coefficient GetSlackFromThreshold(Coefficient threshold) { + Coefficient GetSlackFromThreshold(Coefficient threshold) const { return (index_ < 0) ? threshold : coeffs_[index_] + threshold; } void Update(Coefficient slack, Coefficient* threshold) { @@ -514,6 +537,7 @@ class UpperBoundedLinearConstraint { // - coeffs_ contains unique increasing coefficients. // - starts_[i] is the index in literals_ of the first literal with // coefficient coeffs_[i]. + EnforcementId enforcement_id_; std::vector coeffs_; std::vector starts_; std::vector literals_; @@ -528,6 +552,7 @@ class PbConstraints : public SatPropagator { public: explicit PbConstraints(Model* model) : SatPropagator("PbConstraints"), + enforcement_propagator_(model->GetOrCreate()), conflicting_constraint_index_(-1), num_learned_constraint_before_cleanup_(0), constraint_activity_increment_(1.0), @@ -577,6 +602,11 @@ class PbConstraints : public SatPropagator { // - The constraint cannot be conflicting. // - The constraint cannot have propagated at an earlier decision level. bool AddConstraint(const std::vector& cst, Coefficient rhs, + Trail* trail) { + return AddConstraint(/*enforcement_literals=*/{}, cst, rhs, trail); + } + bool AddConstraint(const std::vector& enforcement_literals, + const std::vector& cst, Coefficient rhs, Trail* trail); // Same as AddConstraint(), but also marks the added constraint as learned @@ -623,7 +653,12 @@ class PbConstraints : public SatPropagator { int64_t num_threshold_updates() const { return num_threshold_updates_; } private: + DEFINE_STRONG_INDEX_TYPE(ConstraintIndex); + bool PropagateNext(Trail* trail); + bool PropagateConstraint(ConstraintIndex index, Trail* trail, + int source_trail_index, + bool* need_untrail_inspection = nullptr); // Same function as the clause related one is SatSolver(). // TODO(user): Remove duplication. @@ -642,7 +677,6 @@ class PbConstraints : public SatPropagator { // about two times faster with this implementation than one with direct // pointer to an UpperBoundedLinearConstraint. The main reason for this is // probably that the thresholds_ vector is a lot more efficient cache-wise. - DEFINE_STRONG_INDEX_TYPE(ConstraintIndex); struct ConstraintIndexWithCoeff { ConstraintIndexWithCoeff() = default; // Needed for vector.resize() ConstraintIndexWithCoeff(bool n, ConstraintIndex i, Coefficient c) @@ -663,6 +697,10 @@ class PbConstraints : public SatPropagator { util_intops::StrongVector> to_update_; + // The indices of the constraints that need to be updated because of an + // enforcement status change. + SparseBitset enforcement_status_changed_; + // Bitset used to optimize the Untrail() function. SparseBitset to_untrail_; @@ -674,6 +712,8 @@ class PbConstraints : public SatPropagator { // Helper to enqueue propagated literals on the trail and store their reasons. PbConstraintsEnqueueHelper enqueue_helper_; + EnforcementPropagator* enforcement_propagator_; + // Last conflicting PB constraint index. This is reset to -1 when // ClearConflictingConstraint() is called. ConstraintIndex conflicting_constraint_index_; diff --git a/ortools/sat/pb_constraint_test.cc b/ortools/sat/pb_constraint_test.cc index 94225a3a12..482adc92d9 100644 --- a/ortools/sat/pb_constraint_test.cc +++ b/ortools/sat/pb_constraint_test.cc @@ -23,6 +23,7 @@ #include "gtest/gtest.h" #include "ortools/base/gmock.h" #include "ortools/base/strong_vector.h" +#include "ortools/sat/enforcement.h" #include "ortools/sat/model.h" #include "ortools/sat/sat_base.h" #include "ortools/util/strong_integers.h" @@ -128,7 +129,7 @@ TEST(ComputeBooleanLinearExpressionCanonicalForm, BigIntCase) { TEST(ApplyLiteralMappingTest, BasicTest) { Coefficient bound_shift, max_value; - // This is needed to initizalize the ITIVector below. + // This is needed to initialize the ITIVector below. std::vector temp{ kTrueLiteralIndex, kFalseLiteralIndex, // var1 fixed to true. Literal(-1).Index(), Literal(+1).Index(), // var2 mapped to not(var1) @@ -208,7 +209,7 @@ TEST(CanonicalBooleanLinearProblem, OverflowCases) { reference = MakePb({{+1, 10}, {-1, 10}}); } - // All These constraint are trivially satisfiables, so no new constraints + // All these constraints are trivially satisfiable, so no new constraints // should be added. cst = reference; EXPECT_TRUE(problem.AddLinearConstraint(true, -kCoefficientMax, true, @@ -252,15 +253,18 @@ TEST(UpperBoundedLinearConstraintTest, ConstructionAndBasicPropagation) { trail.Resize(10); UpperBoundedLinearConstraint cst( + /*enforcement_literals=*/{}, MakePb({{+1, 4}, {+2, 4}, {-3, 5}, {+4, 10}})); - cst.InitializeRhs(Coefficient(7), 0, &threshold, &trail, &helper); + cst.InitializeRhs(EnforcementStatus::IS_ENFORCED, /*enforcement_literals=*/{}, + Coefficient(7), 0, &threshold, &trail, &helper); EXPECT_EQ(threshold, 2); EXPECT_THAT(TrailToVector(trail), LiteralsAre(-4)); trail.Enqueue(Literal(-3), AssignmentType::kSearchDecision); threshold -= 5; // The coeff of -3 in cst. EXPECT_TRUE(cst.Propagate(trail.Info(Literal(-3).Variable()).trail_index, - &threshold, &trail, &helper)); + &threshold, &trail, EnforcementStatus::IS_ENFORCED, + /*enforcement_literals=*/{}, &helper)); EXPECT_EQ(threshold, 2); EXPECT_THAT(TrailToVector(trail), LiteralsAre(-4, -3, -1, -2)); @@ -280,8 +284,10 @@ TEST(UpperBoundedLinearConstraintTest, Conflict) { // At most one constraint. UpperBoundedLinearConstraint cst( + /*enforcement_literals=*/{}, MakePb({{+1, 1}, {+2, 1}, {+3, 1}, {+4, 1}})); - cst.InitializeRhs(Coefficient(1), 0, &threshold, &trail, &helper); + cst.InitializeRhs(EnforcementStatus::IS_ENFORCED, /*enforcement_literals=*/{}, + Coefficient(1), 0, &threshold, &trail, &helper); EXPECT_EQ(threshold, 0); // Two assignment from other part of the solver. @@ -293,7 +299,8 @@ TEST(UpperBoundedLinearConstraintTest, Conflict) { // We propagate only +1. threshold -= 1; EXPECT_FALSE(cst.Propagate(trail.Info(Literal(+1).Variable()).trail_index, - &threshold, &trail, &helper)); + &threshold, &trail, EnforcementStatus::IS_ENFORCED, + /*enforcement_literals=*/{}, &helper)); EXPECT_THAT(helper.conflict, LiteralsAre(-1, -2)); } @@ -304,10 +311,11 @@ TEST(UpperBoundedLinearConstraintTest, CompactReason) { PbConstraintsEnqueueHelper helper; helper.reasons.resize(10); - // At most one constraint. UpperBoundedLinearConstraint cst( + /*enforcement_literals=*/{}, MakePb({{+1, 1}, {+2, 2}, {+3, 3}, {+4, 4}})); - cst.InitializeRhs(Coefficient(7), 0, &threshold, &trail, &helper); + cst.InitializeRhs(EnforcementStatus::IS_ENFORCED, /*enforcement_literals=*/{}, + Coefficient(7), 0, &threshold, &trail, &helper); EXPECT_EQ(threshold, 3); // Two assignment from other part of the solver. @@ -321,17 +329,138 @@ TEST(UpperBoundedLinearConstraintTest, CompactReason) { // We propagate when +3 is processed. threshold = -3; const int source_trail_index = trail.Info(Literal(+3).Variable()).trail_index; - EXPECT_TRUE(cst.Propagate(source_trail_index, &threshold, &trail, &helper)); + EXPECT_TRUE(cst.Propagate(source_trail_index, &threshold, &trail, + EnforcementStatus::IS_ENFORCED, + /*enforcement_literals=*/{}, &helper)); EXPECT_EQ(trail.Index(), 4); EXPECT_EQ(trail[3], Literal(-4)); // -1 do not need to be in the reason since {-3, -2} propagates exactly // the same way. - cst.FillReason(trail, source_trail_index, Literal(-4).Variable(), - &helper.conflict); + cst.FillReason(trail, source_trail_index, /*enforcement_literals=*/{}, + Literal(-4).Variable(), &helper.conflict); EXPECT_THAT(helper.conflict, LiteralsAre(-3, -2)); } +TEST(UpperBoundedLinearConstraintTest, ConflictAfterEnforcementStatusChange) { + Coefficient threshold; + Trail trail; + trail.Resize(10); + PbConstraintsEnqueueHelper helper; + helper.reasons.resize(10); + + std::vector enforcement_literals = {Literal(+9)}; + UpperBoundedLinearConstraint cst( + enforcement_literals, MakePb({{+1, 1}, {+2, 2}, {+3, 3}, {+4, 4}})); + cst.InitializeRhs(EnforcementStatus::IS_FALSE, enforcement_literals, + Coefficient(7), 0, &threshold, &trail, &helper); + EXPECT_EQ(threshold, 3); + + // Some assignments from other parts of the solver. + trail.SetDecisionLevel(1); + trail.Enqueue(Literal(+1), AssignmentType::kSearchDecision); + trail.SetDecisionLevel(2); + trail.Enqueue(Literal(+2), AssignmentType::kSearchDecision); + trail.SetDecisionLevel(3); + trail.Enqueue(Literal(+3), AssignmentType::kSearchDecision); + trail.SetDecisionLevel(4); + trail.Enqueue(Literal(+4), AssignmentType::kSearchDecision); + trail.SetDecisionLevel(5); + trail.Enqueue(Literal(+9), AssignmentType::kSearchDecision); + + // We detect a conflict when +9 is processed. + threshold = -7; + const int source_trail_index = trail.Info(Literal(+9).Variable()).trail_index; + EXPECT_FALSE(cst.Propagate(source_trail_index, &threshold, &trail, + EnforcementStatus::IS_ENFORCED, + enforcement_literals, &helper)); + + // -1 do not need to be in the reason since {-4, -3, -2} propagates exactly + // the same way. + EXPECT_THAT(helper.conflict, LiteralsAre(-9, -3, -2, -4)); +} + +TEST(UpperBoundedLinearConstraintTest, PropagateEnforcementAfterStatusChange) { + Coefficient threshold; + Trail trail; + trail.Resize(10); + PbConstraintsEnqueueHelper helper; + helper.reasons.resize(10); + + std::vector enforcement_literals = {Literal(+8), Literal(+9)}; + UpperBoundedLinearConstraint cst( + enforcement_literals, MakePb({{+1, 1}, {+2, 2}, {+3, 3}, {+4, 4}})); + cst.InitializeRhs(EnforcementStatus::IS_FALSE, enforcement_literals, + Coefficient(7), 0, &threshold, &trail, &helper); + EXPECT_EQ(threshold, 3); + + // Some assignments from other parts of the solver. + trail.SetDecisionLevel(1); + trail.Enqueue(Literal(+1), AssignmentType::kSearchDecision); + trail.SetDecisionLevel(2); + trail.Enqueue(Literal(+2), AssignmentType::kSearchDecision); + trail.SetDecisionLevel(3); + trail.Enqueue(Literal(+3), AssignmentType::kSearchDecision); + trail.SetDecisionLevel(4); + trail.Enqueue(Literal(+4), AssignmentType::kSearchDecision); + trail.SetDecisionLevel(5); + trail.Enqueue(Literal(+9), AssignmentType::kSearchDecision); + + // We should propagate -8 when +9 is processed. + threshold = -7; + const int source_trail_index = trail.Info(Literal(+9).Variable()).trail_index; + EXPECT_TRUE(cst.Propagate(source_trail_index, &threshold, &trail, + EnforcementStatus::CAN_PROPAGATE_ENFORCEMENT, + enforcement_literals, &helper)); + EXPECT_EQ(trail.Index(), 6); + EXPECT_EQ(trail[5], Literal(-8)); + + // -1 do not need to be in the reason since {-4, -3, -2} propagates exactly + // the same way. + const PbConstraintsEnqueueHelper::ReasonInfo& reason = helper.reasons[5]; + cst.FillReason(trail, reason.source_trail_index, enforcement_literals, + Literal(-8).Variable(), &helper.conflict); + EXPECT_THAT(helper.conflict, LiteralsAre(-9, -3, -2, -4)); +} + +TEST(UpperBoundedLinearConstraintTest, + PropagateEnforcementAfterTermAssignment) { + Coefficient threshold; + Trail trail; + trail.Resize(10); + PbConstraintsEnqueueHelper helper; + helper.reasons.resize(10); + + // At most one constraint with enforcement literal. + std::vector enforcement_literals = {Literal(+9)}; + UpperBoundedLinearConstraint cst( + enforcement_literals, MakePb({{+1, 1}, {+2, 1}, {+3, 1}, {+4, 1}})); + cst.InitializeRhs(EnforcementStatus::CAN_PROPAGATE_ENFORCEMENT, + enforcement_literals, Coefficient(1), 0, &threshold, &trail, + &helper); + EXPECT_EQ(threshold, 0); + + // Some assignments from other parts of the solver. + trail.SetDecisionLevel(1); + trail.Enqueue(Literal(+1), AssignmentType::kSearchDecision); + trail.SetDecisionLevel(2); + trail.Enqueue(Literal(+2), AssignmentType::kSearchDecision); + + // We should propagate -9 when +2 is processed. + const int source_trail_index = trail.Info(Literal(+1).Variable()).trail_index; + threshold = -1; + EXPECT_TRUE(cst.Propagate(source_trail_index, &threshold, &trail, + EnforcementStatus::CAN_PROPAGATE_ENFORCEMENT, + enforcement_literals, &helper)); + EXPECT_EQ(trail.Index(), 3); + EXPECT_EQ(trail[2], Literal(-9)); + + const PbConstraintsEnqueueHelper::ReasonInfo& reason = helper.reasons[2]; + cst.FillReason(trail, reason.source_trail_index, enforcement_literals, + Literal(-9).Variable(), &helper.conflict); + EXPECT_THAT(helper.conflict, LiteralsAre(-1, -2)); +} + TEST(PbConstraintsTest, Duplicates) { Model model; PbConstraints& csts = *(model.GetOrCreate()); @@ -419,7 +548,7 @@ TEST(PbConstraintsTest, BasicDeletion) { while (!csts.PropagationIsDone(trail)) EXPECT_TRUE(csts.Propagate(&trail)); EXPECT_EQ("-1", trail.DebugString()); - // But also enqueing -2 should. + // But also enqueuing -2 should. trail.Enqueue(Literal(-2), AssignmentType::kSearchDecision); while (!csts.PropagationIsDone(trail)) EXPECT_TRUE(csts.Propagate(&trail)); EXPECT_EQ("-1 -2 -3 -4", trail.DebugString()); @@ -491,7 +620,7 @@ TEST(PbConstraintsTest, AddConstraintWithLevel0Propagation) { TEST(PbConstraintsTest, AddConstraintUMR) { const auto cst = MakePb({{+3, 7}}); - UpperBoundedLinearConstraint c(cst); + UpperBoundedLinearConstraint c(/*enforcement_literals=*/{}, cst); // Calling hashing on c generates an UMR that is triggered during the hash_map // lookup below. const uint64_t ct_hash = c.hash(); diff --git a/ortools/sat/sat_parameters.proto b/ortools/sat/sat_parameters.proto index 7bbdc61a3d..7e441b53cb 100644 --- a/ortools/sat/sat_parameters.proto +++ b/ortools/sat/sat_parameters.proto @@ -1219,7 +1219,7 @@ message SatParameters { // the same subtrees as one another. // Specifying a negative number uses a heuristic to select an appropriate // number of shared tree workeres based on the total number of workers. - optional int32 shared_tree_num_workers = 235 [default = 0]; + optional int32 shared_tree_num_workers = 235 [default = -1]; // Set on shared subtree workers. Users should not set this directly. optional bool use_shared_tree_search = 236 [default = false]; diff --git a/ortools/sat/sat_solver.cc b/ortools/sat/sat_solver.cc index 9c1bd8f5d5..efa9f9e614 100644 --- a/ortools/sat/sat_solver.cc +++ b/ortools/sat/sat_solver.cc @@ -281,12 +281,24 @@ bool SatSolver::AddProblemClauseInternal(absl::Span literals) { } bool SatSolver::AddLinearConstraintInternal( + const std::vector& enforcement_literals, const std::vector& cst, Coefficient rhs, Coefficient max_value) { SCOPED_TIME_STAT(&stats_); - DCHECK(BooleanLinearExpressionIsCanonical(cst)); - if (rhs < 0) return SetModelUnsat(); // Unsatisfiable constraint. - if (rhs >= max_value) return true; // Always satisfied constraint. + DCHECK(BooleanLinearExpressionIsCanonical(enforcement_literals, cst)); + if (rhs < 0) { + // Unsatisfiable constraint if enforced. + if (enforcement_literals.empty()) { + return SetModelUnsat(); + } else { + literals_scratchpad_.clear(); + for (const Literal& literal : enforcement_literals) { + literals_scratchpad_.push_back(literal.Negated()); + } + return AddProblemClauseInternal(literals_scratchpad_); + } + } + if (rhs >= max_value) return true; // Always satisfied constraint. // Since the constraint is in canonical form, the coefficients are sorted. const Coefficient min_coeff = cst.front().coefficient; @@ -296,17 +308,28 @@ bool SatSolver::AddLinearConstraintInternal( // assignment is the one where all the literals are true. if (max_value - min_coeff <= rhs) { // This constraint is actually a clause. It is faster to treat it as one. - literals_scratchpad_.clear(); - for (const LiteralWithCoeff& term : cst) { - literals_scratchpad_.push_back(term.literal.Negated()); + if (enforcement_literals.empty()) { + literals_scratchpad_.clear(); + for (const LiteralWithCoeff& term : cst) { + literals_scratchpad_.push_back(term.literal.Negated()); + } + return AddProblemClauseInternal(literals_scratchpad_); + } else { + std::vector literals; + for (const Literal& literal : enforcement_literals) { + literals.push_back(literal.Negated()); + } + for (const LiteralWithCoeff& term : cst) { + literals.push_back(term.literal.Negated()); + } + return AddProblemClause(literals); } - return AddProblemClauseInternal(literals_scratchpad_); } // Detect at most one constraints. Note that this use the fact that the // coefficient are sorted. if (!parameters_->use_pb_resolution() && max_coeff <= rhs && - 2 * min_coeff > rhs) { + 2 * min_coeff > rhs && enforcement_literals.empty()) { literals_scratchpad_.clear(); for (const LiteralWithCoeff& term : cst) { literals_scratchpad_.push_back(term.literal); @@ -317,9 +340,13 @@ bool SatSolver::AddLinearConstraintInternal( return true; } + // TODO(user): fix literals with coefficient larger than rhs to false, or + // add implication enforcement => not(literal) (and remove them from the + // constraint)? + // TODO(user): If this constraint forces all its literal to false (when rhs is // zero for instance), we still add it. Optimize this? - return pb_constraints_->AddConstraint(cst, rhs, trail_); + return pb_constraints_->AddConstraint(enforcement_literals, cst, rhs, trail_); } void SatSolver::CanonicalizeLinear(std::vector* cst, @@ -355,11 +382,25 @@ bool SatSolver::AddLinearConstraint(bool use_lower_bound, Coefficient lower_bound, bool use_upper_bound, Coefficient upper_bound, + std::vector* enforcement_literals, std::vector* cst) { SCOPED_TIME_STAT(&stats_); CHECK_EQ(CurrentDecisionLevel(), 0); if (model_is_unsat_) return false; + gtl::STLSortAndRemoveDuplicates(enforcement_literals); + int num_enforcement_literals = 0; + for (int i = 0; i < enforcement_literals->size(); ++i) { + const Literal literal = (*enforcement_literals)[i]; + if (trail_->Assignment().LiteralIsFalse(literal)) { + return true; + } + if (!trail_->Assignment().LiteralIsTrue(literal)) { + (*enforcement_literals)[num_enforcement_literals++] = literal; + } + } + enforcement_literals->resize(num_enforcement_literals); + Coefficient bound_shift(0); if (use_upper_bound) { @@ -367,7 +408,8 @@ bool SatSolver::AddLinearConstraint(bool use_lower_bound, CanonicalizeLinear(cst, &bound_shift, &max_value); const Coefficient rhs = ComputeCanonicalRhs(upper_bound, bound_shift, max_value); - if (!AddLinearConstraintInternal(*cst, rhs, max_value)) { + if (!AddLinearConstraintInternal(*enforcement_literals, *cst, rhs, + max_value)) { return SetModelUnsat(); } } @@ -384,7 +426,8 @@ bool SatSolver::AddLinearConstraint(bool use_lower_bound, } const Coefficient rhs = ComputeNegatedCanonicalRhs(lower_bound, bound_shift, max_value); - if (!AddLinearConstraintInternal(*cst, rhs, max_value)) { + if (!AddLinearConstraintInternal(*enforcement_literals, *cst, rhs, + max_value)) { return SetModelUnsat(); } } diff --git a/ortools/sat/sat_solver.h b/ortools/sat/sat_solver.h index 28559af0b6..63bfab81b4 100644 --- a/ortools/sat/sat_solver.h +++ b/ortools/sat/sat_solver.h @@ -138,8 +138,17 @@ class SatSolver { // TODO(user): Instead of failing, implement an error handling code. bool AddLinearConstraint(bool use_lower_bound, Coefficient lower_bound, bool use_upper_bound, Coefficient upper_bound, + std::vector* enforcement_literals, std::vector* cst); + bool AddLinearConstraint(bool use_lower_bound, Coefficient lower_bound, + bool use_upper_bound, Coefficient upper_bound, + std::vector* cst) { + std::vector enforcement_literals; + return AddLinearConstraint(use_lower_bound, lower_bound, use_upper_bound, + upper_bound, &enforcement_literals, cst); + } + // Returns true if the model is UNSAT. Note that currently the status is // "sticky" and once this happen, nothing else can be done with the solver. // @@ -600,8 +609,10 @@ class SatSolver { // This is used by all the Add*LinearConstraint() functions. It detects // infeasible/trivial constraints or clause constraints and takes the proper // action. - bool AddLinearConstraintInternal(const std::vector& cst, - Coefficient rhs, Coefficient max_value); + bool AddLinearConstraintInternal( + const std::vector& enforcement_literals, + const std::vector& cst, Coefficient rhs, + Coefficient max_value); // Makes sure a pseudo boolean constraint is in canonical form. void CanonicalizeLinear(std::vector* cst, diff --git a/ortools/sat/util.cc b/ortools/sat/util.cc index 8e29e93780..a58d356de7 100644 --- a/ortools/sat/util.cc +++ b/ortools/sat/util.cc @@ -128,8 +128,8 @@ void RandomizeDecisionHeuristic(absl::BitGenRef random, namespace { -// This will be optimized into one division. I tested that in other places: -// https://source.corp.google.com/piper///depot/ortools/sat/integer_test.cc;l=1223-1228;bpv=0 +// This will be optimized into one division. I tested that in +// `BM_DivisionAndRemainderAlternative` (sat/integer_test.cc) // // Note that I am not 100% sure we need the indirection for the optimization // to kick in though, but this seemed safer given our weird r[i ^ 1] inputs. diff --git a/ortools/sat/work_assignment.cc b/ortools/sat/work_assignment.cc index 010e359323..a9a48286e4 100644 --- a/ortools/sat/work_assignment.cc +++ b/ortools/sat/work_assignment.cc @@ -232,7 +232,7 @@ absl::Span ProtoTrail::Implications(int level) const { SharedTreeManager::SharedTreeManager(Model* model) : params_(*model->GetOrCreate()), - num_workers_(params_.shared_tree_num_workers()), + num_workers_(std::max(0, params_.shared_tree_num_workers())), max_path_depth_(MaxPossibleLeafDepth(params_)), shared_response_manager_(model->GetOrCreate()), num_splits_wanted_( @@ -242,7 +242,6 @@ SharedTreeManager::SharedTreeManager(Model* model) std::numeric_limits::max() / std::max(num_workers_, 1) ? std::numeric_limits::max() : num_workers_ * params_.shared_tree_max_nodes_per_worker()) { - CHECK_GE(num_workers_, 0); // Create the root node with a fake literal. nodes_.push_back( {.literal = ProtoLiteral(),