From 8627cbf9b8adfa70f14a6b5192b09ae3634cd4b2 Mon Sep 17 00:00:00 2001 From: Laurent Perron Date: Fri, 16 Feb 2024 23:13:09 +0100 Subject: [PATCH] [CP-SAT] improve cumulative relaxation; more spans --- ortools/sat/BUILD.bazel | 3 + ortools/sat/clause.cc | 4 +- ortools/sat/cp_model_solver.cc | 4 + ortools/sat/cp_model_symmetries.cc | 3 +- ortools/sat/encoding.cc | 3 +- ortools/sat/encoding.h | 3 +- ortools/sat/linear_relaxation.cc | 130 ++++++++++++++++++++++++++--- 7 files changed, 134 insertions(+), 16 deletions(-) diff --git a/ortools/sat/BUILD.bazel b/ortools/sat/BUILD.bazel index 06d1525a84..64344b7953 100644 --- a/ortools/sat/BUILD.bazel +++ b/ortools/sat/BUILD.bazel @@ -1347,6 +1347,7 @@ cc_library( ":scheduling_cuts", ":util", "//ortools/base", + "//ortools/base:mathutil", "//ortools/base:stl_util", "//ortools/base:strong_vector", "//ortools/util:logging", @@ -1914,6 +1915,7 @@ cc_library( "//ortools/util:strong_integers", "@com_google_absl//absl/log:check", "@com_google_absl//absl/strings", + "@com_google_absl//absl/types:span", ], ) @@ -2154,6 +2156,7 @@ cc_library( "@com_google_absl//absl/meta:type_traits", "@com_google_absl//absl/status", "@com_google_absl//absl/strings", + "@com_google_absl//absl/types:span", "@com_google_protobuf//:protobuf", ], ) diff --git a/ortools/sat/clause.cc b/ortools/sat/clause.cc index 19f59cc570..ff05d46450 100644 --- a/ortools/sat/clause.cc +++ b/ortools/sat/clause.cc @@ -1627,8 +1627,8 @@ bool BinaryImplicationGraph::ComputeTransitiveReduction(bool log_info) { namespace { -int ElementInIntersectionOrMinusOne(const std::vector& a, - const std::vector& b) { +int ElementInIntersectionOrMinusOne(absl::Span a, + absl::Span b) { DCHECK(std::is_sorted(a.begin(), a.end())); DCHECK(std::is_sorted(b.begin(), b.end())); if (a.empty() || b.empty()) return -1; diff --git a/ortools/sat/cp_model_solver.cc b/ortools/sat/cp_model_solver.cc index fe7a87afc5..ae6140d709 100644 --- a/ortools/sat/cp_model_solver.cc +++ b/ortools/sat/cp_model_solver.cc @@ -940,6 +940,7 @@ IntegerVariable AddLPConstraints(bool objective_need_to_be_tight, const CpModelProto& model_proto, Model* m) { // Non const as we will std::move() stuff out of there. LinearRelaxation relaxation = ComputeLinearRelaxation(model_proto, m); + if (m->GetOrCreate()->ModelIsUnsat()) return kNoIntegerVariable; // The bipartite graph of LP constraints might be disconnected: // make a partition of the variables into connected components. @@ -1590,6 +1591,8 @@ void LoadFeasibilityPump(const CpModelProto& model_proto, Model* model) { // Add linear constraints to Feasibility Pump. const LinearRelaxation relaxation = ComputeLinearRelaxation(model_proto, model); + if (model->GetOrCreate()->ModelIsUnsat()) return; + const int num_lp_constraints = static_cast(relaxation.linear_constraints.size()); if (num_lp_constraints == 0) return; @@ -1701,6 +1704,7 @@ void LoadCpModel(const CpModelProto& model_proto, Model* model) { // Linearize some part of the problem and register LP constraint(s). objective_var = AddLPConstraints(objective_need_to_be_tight, model_proto, model); + if (sat_solver->ModelIsUnsat()) return unsat(); } else if (model_proto.has_objective()) { const CpObjectiveProto& obj = model_proto.objective(); std::vector> terms; diff --git a/ortools/sat/cp_model_symmetries.cc b/ortools/sat/cp_model_symmetries.cc index 339c9e59a5..9cd429cfa4 100644 --- a/ortools/sat/cp_model_symmetries.cc +++ b/ortools/sat/cp_model_symmetries.cc @@ -30,6 +30,7 @@ #include "absl/status/status.h" #include "absl/strings/str_cat.h" #include "absl/strings/str_join.h" +#include "absl/types/span.h" #include "google/protobuf/message.h" #include "ortools/algorithms/find_graph_symmetries.h" #include "ortools/algorithms/sparse_permutation.h" @@ -55,7 +56,7 @@ namespace sat { namespace { struct VectorHash { - std::size_t operator()(const std::vector& values) const { + std::size_t operator()(absl::Span values) const { size_t hash = 0; for (const int64_t value : values) { hash = util_hash::Hash(value, hash); diff --git a/ortools/sat/encoding.cc b/ortools/sat/encoding.cc index b6d24fa0aa..ca7efb2b8d 100644 --- a/ortools/sat/encoding.cc +++ b/ortools/sat/encoding.cc @@ -26,6 +26,7 @@ #include "absl/log/check.h" #include "absl/strings/str_cat.h" +#include "absl/types/span.h" #include "ortools/base/stl_util.h" #include "ortools/sat/boolean_problem.pb.h" #include "ortools/sat/pb_constraint.h" @@ -582,7 +583,7 @@ Coefficient MaxNodeWeightSmallerThan(const std::vector& nodes, return result; } -bool ObjectiveEncoder::ProcessCore(const std::vector& core, +bool ObjectiveEncoder::ProcessCore(absl::Span core, Coefficient min_weight, Coefficient gap, std::string* info) { // Backtrack to be able to add new constraints. diff --git a/ortools/sat/encoding.h b/ortools/sat/encoding.h index 744fae4d54..e6602c133a 100644 --- a/ortools/sat/encoding.h +++ b/ortools/sat/encoding.h @@ -28,6 +28,7 @@ #include #include "absl/log/check.h" +#include "absl/types/span.h" #include "ortools/base/logging.h" #include "ortools/base/macros.h" #include "ortools/base/types.h" @@ -238,7 +239,7 @@ class ObjectiveEncoder { // Updates the encoding using the given core. The literals in the core must // match the order in nodes. Returns false if the model become infeasible. - bool ProcessCore(const std::vector& core, Coefficient min_weight, + bool ProcessCore(absl::Span core, Coefficient min_weight, Coefficient gap, std::string* info); void AddBaseNode(EncodingNode node) { diff --git a/ortools/sat/linear_relaxation.cc b/ortools/sat/linear_relaxation.cc index 5fdf18a9db..d5a6382593 100644 --- a/ortools/sat/linear_relaxation.cc +++ b/ortools/sat/linear_relaxation.cc @@ -28,6 +28,7 @@ #include "absl/meta/type_traits.h" #include "absl/types/span.h" #include "ortools/base/logging.h" +#include "ortools/base/mathutil.h" #include "ortools/base/stl_util.h" #include "ortools/base/strong_vector.h" #include "ortools/sat/circuit.h" // for ReindexArcs. @@ -820,15 +821,39 @@ void AddCumulativeRelaxation(const AffineExpression& capacity, IntegerValue max_of_ends = kMinIntegerValue; int num_variable_energies = 0; int num_optionals = 0; + int64_t sizes_gcd = 0; + int64_t demands_gcd = 0; + int64_t num_active_intervals = 0; IntegerTrail* integer_trail = model->GetOrCreate(); for (int index = 0; index < num_intervals; ++index) { if (helper->IsAbsent(index)) continue; - if (helper->IsOptional(index) && demands_helper->EnergyMin(index) >= 0) { - num_optionals++; + if (helper->IsOptional(index)) { + if (demands_helper->EnergyMin(index) >= 0) { + num_optionals++; + } else { + continue; + } } + if (!helper->SizeIsFixed(index) || !demands_helper->DemandIsFixed(index)) { num_variable_energies++; } + if (sizes_gcd != 1) { + if (helper->SizeIsFixed(index)) { + sizes_gcd = MathUtil::GCD64(helper->SizeMin(index).value(), sizes_gcd); + } else { + sizes_gcd = 1; + } + } + if (demands_gcd != 1) { + if (demands_helper->DemandIsFixed(index)) { + demands_gcd = MathUtil::GCD64(demands_helper->DemandMin(index).value(), + demands_gcd); + } else { + demands_gcd = 1; + } + } + num_active_intervals++; min_of_starts = std::min(min_of_starts, helper->StartMin(index)); max_of_ends = std::max(max_of_ends, helper->EndMax(index)); } @@ -836,15 +861,88 @@ void AddCumulativeRelaxation(const AffineExpression& capacity, VLOG(2) << "Span [" << min_of_starts << ".." << max_of_ends << "] with " << num_optionals << " optional intervals, and " << num_variable_energies << " variable energy tasks out of " - << num_intervals << " intervals" - << (makespan.has_value() ? ", and 1 makespan" : ""); + << num_active_intervals << "/" << num_intervals << " intervals" + << (makespan.has_value() ? ", and 1 makespan" : "") + << " sizes_gcd: " << sizes_gcd << " demands_gcd: " << demands_gcd; - // If nothing is variable, the linear relaxation will already be enforced - // by the scheduling propagators. - if (num_variable_energies + num_optionals == 0) return; + // There are no active intervals, no need to add the relaxation. + if (num_active_intervals == 0) return; + + // If nothing is variable, and the coefficients cannot be reduced, the linear + // relaxation will already be enforced by the scheduling propagators. + if (num_variable_energies + num_optionals == 0 && sizes_gcd == 1 && + demands_gcd == 1) { + return; + } + + // Specialized case 1: sizes are fixed with a non 1 gcd and no makespan. + if (sizes_gcd != 1 && !makespan.has_value()) { + VLOG(2) << "Cumulative relaxation: sizes_gcd = " << sizes_gcd + << ", demands_gcd = " << demands_gcd + << ", no makespan, capacity is " << capacity.DebugString(); + // We can simplify the capacity only if it is fixed. + // TODO(user): We could use (capacity / demands_gcd) * demands_gcd. + const int64_t active_demand_gcd = + integer_trail->IsFixed(capacity) ? demands_gcd : 1; + LinearConstraintBuilder lc(model, kMinIntegerValue, IntegerValue(0)); + for (int index = 0; index < num_intervals; ++index) { + if (helper->IsAbsent(index)) continue; + if (helper->IsOptional(index)) { + const IntegerValue energy_min = demands_helper->EnergyMin(index); + if (energy_min == 0) continue; + DCHECK_EQ(energy_min % (sizes_gcd * active_demand_gcd), 0); + if (!lc.AddLiteralTerm(helper->PresenceLiteral(index), + energy_min / sizes_gcd / active_demand_gcd)) { + return; + } + } else { + // Copy the decomposed energy. + std::vector product = + demands_helper->DecomposedEnergies()[index]; + if (!product.empty()) { + // The energy is defined if the vector is not empty. + // Let's reduce the coefficients. + for (LiteralValueValue& entry : product) { + DCHECK_EQ(entry.left_value % sizes_gcd, 0); + entry.left_value /= sizes_gcd; + DCHECK_EQ(entry.right_value % demands_gcd, 0); + entry.right_value /= demands_gcd; + } + if (!lc.AddDecomposedProduct(product)) return; + } else { + // We know the size is fixed. + const IntegerValue local_size = + integer_trail->FixedValue(helper->Sizes()[index]); + DCHECK_EQ(local_size % sizes_gcd, 0); + if (active_demand_gcd == 1) { + lc.AddTerm(demands_helper->Demands()[index], + local_size / sizes_gcd); + } else { + const IntegerValue local_demand = + integer_trail->FixedValue(demands_helper->Demands()[index]); + DCHECK_EQ(local_demand % active_demand_gcd, 0); + lc.AddConstant(local_size * local_demand / sizes_gcd / + active_demand_gcd); + } + } + } + } + + // Add the available energy of the cumulative. + if (active_demand_gcd == 1) { + lc.AddTerm(capacity, -(max_of_ends - min_of_starts) / sizes_gcd); + } else { + const IntegerValue fixed_capacity = integer_trail->FixedValue(capacity); + lc.AddConstant(-fixed_capacity * (max_of_ends - min_of_starts) / + sizes_gcd / active_demand_gcd); + } + relaxation->linear_constraints.push_back(lc.Build()); + return; + } + + // TODO(user): Implement demands_gcd != 1 && capacity is fixed. LinearConstraintBuilder lc(model, kMinIntegerValue, IntegerValue(0)); - int num_intervals_added = 0; for (int index = 0; index < num_intervals; ++index) { if (helper->IsAbsent(index)) continue; if (helper->IsOptional(index)) { @@ -868,9 +966,7 @@ void AddCumulativeRelaxation(const AffineExpression& capacity, integer_trail); } } - ++num_intervals_added; } - if (num_intervals_added == 0) return; // Create and link span_start and span_end to the starts and ends of the // tasks. @@ -1839,7 +1935,19 @@ LinearRelaxation ComputeLinearRelaxation(const CpModelProto& model_proto, // so that we don't do extra work in the connected component computation. relaxation.at_most_ones.clear(); - // Remove size one LP constraints, they are not useful. + // Propagate unary constraints. + { + SatSolver* sat_solver = m->GetOrCreate(); + for (const LinearConstraint& lc : relaxation.linear_constraints) { + if (lc.num_terms > 1) continue; + LoadLinearConstraint(lc, m); + if (sat_solver->ModelIsUnsat()) return relaxation; + } + if (!sat_solver->FinishPropagation()) return relaxation; + } + + // Remove size one LP constraints from the main algorithms, they are not + // useful. relaxation.linear_constraints.erase( std::remove_if( relaxation.linear_constraints.begin(),