Skip to content

Commit

Permalink
[CP-SAT] improve memory usage for LNS
Browse files Browse the repository at this point in the history
  • Loading branch information
@lperron
lperron committed Feb 10, 2025
1 parent fc237ea commit d621036
Show file tree
Hide file tree
Showing 4 changed files with 74 additions and 40 deletions.
41 changes: 28 additions & 13 deletions ortools/sat/cp_model_lns.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -19,6 +19,7 @@
#include <deque>
#include <functional>
#include <limits>
#include <memory>
#include <random>
#include <string>
#include <tuple>
Expand Down Expand Up @@ -77,16 +78,19 @@ NeighborhoodGeneratorHelper::NeighborhoodGeneratorHelper(
shared_bounds_(shared_bounds),
shared_response_(shared_response) {
// Initialize proto memory.
local_arena_storage_.assign(Neighborhood::kDefaultArenaSizePerVariable *
model_proto_.variables_size(),
0);
local_arena_ = std::make_unique<google::protobuf::Arena>(
local_arena_storage_.data(), local_arena_storage_.size());
simplified_model_proto_ =
google::protobuf::Arena::Create<CpModelProto>(&local_arena_);
model_proto_with_only_variables_ =
google::protobuf::Arena::Create<CpModelProto>(&local_arena_);
google::protobuf::Arena::Create<CpModelProto>(local_arena_.get());

CHECK(shared_response_ != nullptr);
if (shared_bounds_ != nullptr) {
shared_bounds_id_ = shared_bounds_->RegisterNewId();
}
*model_proto_with_only_variables_->mutable_variables() =
*model_proto_with_only_variables_.mutable_variables() =
model_proto_.variables();
InitializeHelperData();
RecomputeHelperData();
Expand All @@ -112,15 +116,15 @@ void NeighborhoodGeneratorHelper::Synchronize() {
const int64_t new_ub = new_upper_bounds[i];
if (VLOG_IS_ON(3)) {
const auto& domain =
model_proto_with_only_variables_->variables(var).domain();
model_proto_with_only_variables_.variables(var).domain();
const int64_t old_lb = domain.Get(0);
const int64_t old_ub = domain.Get(domain.size() - 1);
VLOG(3) << "Variable: " << var << " old domain: [" << old_lb << ", "
<< old_ub << "] new domain: [" << new_lb << ", " << new_ub
<< "]";
}
const Domain old_domain = ReadDomainFromProto(
model_proto_with_only_variables_->variables(var));
model_proto_with_only_variables_.variables(var));
const Domain new_domain =
old_domain.IntersectionWith(Domain(new_lb, new_ub));
if (new_domain.IsEmpty()) {
Expand All @@ -141,7 +145,7 @@ void NeighborhoodGeneratorHelper::Synchronize() {
}
FillDomainInProto(
new_domain,
model_proto_with_only_variables_->mutable_variables(var));
model_proto_with_only_variables_.mutable_variables(var));
new_variables_have_been_fixed |= new_domain.IsFixed();
}
}
Expand All @@ -164,7 +168,7 @@ bool NeighborhoodGeneratorHelper::ObjectiveDomainIsConstraining() const {
const int var = PositiveRef(model_proto_.objective().vars(i));
const int64_t coeff = model_proto_.objective().coeffs(i);
const auto& var_domain =
model_proto_with_only_variables_->variables(var).domain();
model_proto_with_only_variables_.variables(var).domain();
const int64_t v1 = coeff * var_domain[0];
const int64_t v2 = coeff * var_domain[var_domain.size() - 1];
min_activity += std::min(v1, v2);
Expand Down Expand Up @@ -218,9 +222,20 @@ void NeighborhoodGeneratorHelper::RecomputeHelperData() {
{
Model local_model;
CpModelProto mapping_proto;
// We want to replace the simplified_model_proto_ by a new one. Since
// deleting an object in the arena doesn't free the memory, we also delete
// and recreate the arena, but reusing the same storage.
int64_t new_size = local_arena_->SpaceUsed();
new_size += new_size / 2;
simplified_model_proto_->Clear();
local_arena_.reset();
local_arena_storage_.resize(new_size);
local_arena_ = std::make_unique<google::protobuf::Arena>(
local_arena_storage_.data(), local_arena_storage_.size());
simplified_model_proto_ =
google::protobuf::Arena::Create<CpModelProto>(local_arena_.get());
*simplified_model_proto_->mutable_variables() =
model_proto_with_only_variables_->variables();
model_proto_with_only_variables_.variables();
PresolveContext context(&local_model, simplified_model_proto_,
&mapping_proto);
ModelCopy copier(&context);
Expand Down Expand Up @@ -383,7 +398,7 @@ bool NeighborhoodGeneratorHelper::IsActive(int var) const {
}

bool NeighborhoodGeneratorHelper::IsConstant(int var) const {
const auto& var_proto = model_proto_with_only_variables_->variables(var);
const auto& var_proto = model_proto_with_only_variables_.variables(var);
return var_proto.domain_size() == 2 &&
var_proto.domain(0) == var_proto.domain(1);
}
Expand All @@ -395,7 +410,7 @@ Neighborhood NeighborhoodGeneratorHelper::FullNeighborhood() const {
{
absl::ReaderMutexLock lock(&domain_mutex_);
*neighborhood.delta.mutable_variables() =
model_proto_with_only_variables_->variables();
model_proto_with_only_variables_.variables();
}
return neighborhood;
}
Expand Down Expand Up @@ -1057,7 +1072,7 @@ Neighborhood NeighborhoodGeneratorHelper::FixGivenVariables(
absl::ReaderMutexLock domain_lock(&domain_mutex_);
for (int i = 0; i < num_variables; ++i) {
const IntegerVariableProto& current_var =
model_proto_with_only_variables_->variables(i);
model_proto_with_only_variables_.variables(i);
IntegerVariableProto* new_var = neighborhood.delta.add_variables();

// We only copy the name in debug mode.
Expand Down Expand Up @@ -1203,7 +1218,7 @@ CpModelProto NeighborhoodGeneratorHelper::UpdatedModelProtoCopy() const {
{
absl::MutexLock domain_lock(&domain_mutex_);
*updated_model.mutable_variables() =
model_proto_with_only_variables_->variables();
model_proto_with_only_variables_.variables();
}
return updated_model;
}
Expand Down
5 changes: 3 additions & 2 deletions ortools/sat/cp_model_lns.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -319,14 +319,15 @@ class NeighborhoodGeneratorHelper : public SubSolver {
// Arena holding the memory of the CpModelProto* of this class. This saves the
// destruction cost that can take time on problem with millions of
// variables/constraints.
google::protobuf::Arena local_arena_;
std::vector<char> local_arena_storage_;
std::unique_ptr<google::protobuf::Arena> local_arena_;

// This proto will only contain the field variables() with an updated version
// of the domains compared to model_proto_.variables(). We do it like this to
// reduce the memory footprint of the helper when the model is large.
//
// TODO(user): Use custom domain repository rather than a proto?
CpModelProto* model_proto_with_only_variables_ ABSL_GUARDED_BY(domain_mutex_);
CpModelProto model_proto_with_only_variables_ ABSL_GUARDED_BY(domain_mutex_);

// Constraints by types. This never changes.
std::vector<std::vector<int>> type_to_constraints_;
Expand Down
67 changes: 42 additions & 25 deletions ortools/sat/diffn_cuts.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -391,6 +391,28 @@ std::string DiffnCtEvent::DebugString() const {
// The original cut is:
// sum(end_min_i * duration_min_i) >=
// (sum(duration_min_i^2) + sum(duration_min_i)^2) / 2
//
// Let's build a figure where each horizontal rectangle represent a task. It
// ends at the end of the task, and its height is the duration of the task.
// For a given order, we pack each rectangle to the left while not overlapping,
// that is one rectangle starts when the previous one ends.
//
// e1
// -----
// :\ | s1
// : \| e2
// -------------
// :\ |
// : \ | s2
// : \| e3
// ----------------
// : \| s3
// ----------------
//
// We can notice that the total area is independent of the order of tasks.
// The first term of the rhs is the area above the diagonal.
// The second term of the rhs is the area below the diagonal.
//
// We apply the following changes (see the code for cumulative constraints):
// - we strengthen this cuts by noticing that if all tasks starts after S,
// then replacing end_min_i by (end_min_i - S) is still valid.
Expand Down Expand Up @@ -461,12 +483,11 @@ void GenerateNoOvelap2dCompletionTimeCuts(absl::string_view cut_name,
// Best cut so far for this loop.
int best_end = -1;
double best_efficacy = 0.01;
IntegerValue best_min_contrib = 0;
IntegerValue best_capacity = 0;
IntegerValue best_y_range = 0;
IntegerValue best_min_rhs = 0;
bool best_use_subset_sum = false;

// Used in the first term of the rhs of the equation.
IntegerValue sum_event_contributions = 0;
IntegerValue sum_event_areas = 0;
// Used in the second term of the rhs of the equation.
IntegerValue sum_energy = 0;
// For normalization.
Expand All @@ -486,8 +507,7 @@ void GenerateNoOvelap2dCompletionTimeCuts(absl::string_view cut_name,
DCHECK_GE(event.x_start_min, sequence_start_min);
// Make sure we do not overflow.
if (!AddTo(event.energy_min, &sum_energy)) break;
if (!AddProductTo(event.energy_min, event.x_size_min,
&sum_event_contributions)) {
if (!AddProductTo(event.energy_min, event.x_size_min, &sum_event_areas)) {
break;
}
if (!AddSquareTo(event.energy_min, &sum_square_energy)) break;
Expand All @@ -505,7 +525,8 @@ void GenerateNoOvelap2dCompletionTimeCuts(absl::string_view cut_name,
if (i == 0) {
dp.Reset((y_max_of_subset - y_min_of_subset).value());
} else {
if (y_max_of_subset - y_min_of_subset != dp.Bound()) {
// TODO(user): Can we increase the bound dynamically ?
if (y_max_of_subset - y_min_of_subset > dp.Bound()) {
use_dp = false;
}
}
Expand All @@ -522,23 +543,21 @@ void GenerateNoOvelap2dCompletionTimeCuts(absl::string_view cut_name,
if (sum_of_y_size_min <= reachable_capacity) continue;

// Do we have a violated cut ?
const IntegerValue large_rectangle_contrib =
CapProdI(sum_energy, sum_energy);
if (AtMinOrMaxInt64I(large_rectangle_contrib)) break;
const IntegerValue mean_large_rectangle_contrib =
CeilRatio(large_rectangle_contrib, reachable_capacity);
const IntegerValue square_sum_energy = CapProdI(sum_energy, sum_energy);
if (AtMinOrMaxInt64I(square_sum_energy)) break;
const IntegerValue rhs_second_term =
CeilRatio(square_sum_energy, reachable_capacity);

IntegerValue min_contrib =
CapAddI(sum_event_contributions, mean_large_rectangle_contrib);
if (AtMinOrMaxInt64I(min_contrib)) break;
min_contrib = CeilRatio(min_contrib, 2);
IntegerValue min_rhs = CapAddI(sum_event_areas, rhs_second_term);
if (AtMinOrMaxInt64I(min_rhs)) break;
min_rhs = CeilRatio(min_rhs, 2);

// shift contribution by current_start_min.
if (!AddProductTo(sum_energy, current_start_min, &min_contrib)) break;
if (!AddProductTo(sum_energy, current_start_min, &min_rhs)) break;

// The efficacy of the cut is the normalized violation of the above
// equation. We will normalize by the sqrt of the sum of squared energies.
const double efficacy = (ToDouble(min_contrib) - lp_contrib) /
const double efficacy = (ToDouble(min_rhs) - lp_contrib) /
std::sqrt(ToDouble(sum_square_energy));

// For a given start time, we only keep the best cut.
Expand All @@ -550,13 +569,13 @@ void GenerateNoOvelap2dCompletionTimeCuts(absl::string_view cut_name,
if (efficacy > best_efficacy) {
best_efficacy = efficacy;
best_end = i;
best_min_contrib = min_contrib;
best_capacity = reachable_capacity;
best_y_range = y_max_of_subset - y_min_of_subset;
best_min_rhs = min_rhs;
best_use_subset_sum =
reachable_capacity < y_max_of_subset - y_min_of_subset;
}
}
if (best_end != -1) {
LinearConstraintBuilder cut(model, best_min_contrib, kMaxIntegerValue);
LinearConstraintBuilder cut(model, best_min_rhs, kMaxIntegerValue);
bool is_lifted = false;
bool add_energy_to_name = false;
for (int i = 0; i <= best_end; ++i) {
Expand All @@ -568,9 +587,7 @@ void GenerateNoOvelap2dCompletionTimeCuts(absl::string_view cut_name,
std::string full_name(cut_name);
if (is_lifted) full_name.append("_lifted");
if (add_energy_to_name) full_name.append("_energy");
if (best_capacity < best_y_range) {
full_name.append("_subsetsum");
}
if (best_use_subset_sum) full_name.append("_subsetsum");
top_n_cuts.AddCut(cut.Build(), full_name, manager->LpValues());
}
}
Expand Down
1 change: 1 addition & 0 deletions ortools/util/fixed_shape_binary_tree.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -218,6 +218,7 @@ class FixedShapeBinaryTree {
template <typename TypeWithPushBack>
void PartitionIntervalIntoNodes(LeafIndex first_leaf, LeafIndex last_leaf,
TypeWithPushBack* result) const {
DCHECK_LE(first_leaf, last_leaf);
TreeNodeIndex prev(0);
TreeNodeIndex current = GetNodeStartOfRange(first_leaf, last_leaf);
if (current == Root()) {
Expand Down

0 comments on commit d621036

Please sign in to comment.