diff --git a/ortools/sat/BUILD b/ortools/sat/BUILD index 71eb1960a6..a311a8cf88 100644 --- a/ortools/sat/BUILD +++ b/ortools/sat/BUILD @@ -633,6 +633,7 @@ cc_library( ":sat_base", ":sat_parameters_cc_proto", ":sat_solver", + ":timetable", ":theta_tree", "//ortools/base", "//ortools/base:int_type", diff --git a/ortools/sat/cp_model_lns.cc b/ortools/sat/cp_model_lns.cc index 09002a4b2f..8da4699251 100644 --- a/ortools/sat/cp_model_lns.cc +++ b/ortools/sat/cp_model_lns.cc @@ -487,6 +487,18 @@ Neighborhood RelaxationInducedNeighborhoodGenerator::Generate( neighborhood.cp_model.mutable_variables(var)->add_domain(value); neighborhood.is_reduced = true; } + + for (const std::pair> + reduced_var : rins_neighborhood_opt.value().reduced_domain_vars) { + int var = reduced_var.first.model_var; + int64 lb = reduced_var.second.first; + int64 ub = reduced_var.second.second; + if (!helper_.IsActive(var)) continue; + Domain domain = ReadDomainFromProto(neighborhood.cp_model.variables(var)); + domain = domain.IntersectionWith(Domain(lb, ub)); + FillDomainInProto(domain, neighborhood.cp_model.mutable_variables(var)); + neighborhood.is_reduced = true; + } neighborhood.is_generated = true; return neighborhood; } diff --git a/ortools/sat/cp_model_lns.h b/ortools/sat/cp_model_lns.h index 4d976332f6..432404318a 100644 --- a/ortools/sat/cp_model_lns.h +++ b/ortools/sat/cp_model_lns.h @@ -140,6 +140,9 @@ class NeighborhoodGenerator { virtual Neighborhood Generate(const CpSolverResponse& initial_solution, int64 seed, double difficulty) const = 0; + // Returns true if the generator needs a solution to generate a neighborhood. + virtual bool NeedsFirstSolution() const { return true; } + // Returns a short description of the generator. std::string name() const { return name_; } @@ -242,6 +245,10 @@ class SchedulingTimeWindowNeighborhoodGenerator : public NeighborhoodGenerator { // as their linear relaxation. This was published in "Exploring relaxation // induced neighborhoods to improve MIP solutions" 2004 by E. Danna et. // +// If no solution is available, this generates a neighborhood using only the +// linear relaxation values. This was published in "RENS – The Relaxation +// Enforced Neighborhood" 2009 by Timo Berthold. +// // NOTE: The neighborhoods are generated outside of this generator and are // managed by SharedRINSNeighborhoodManager. class RelaxationInducedNeighborhoodGenerator : public NeighborhoodGenerator { @@ -254,6 +261,8 @@ class RelaxationInducedNeighborhoodGenerator : public NeighborhoodGenerator { Neighborhood Generate(const CpSolverResponse& initial_solution, int64 seed, double difficulty) const final; + bool NeedsFirstSolution() const override { return false; } + const Model* model_; }; diff --git a/ortools/sat/cp_model_solver.cc b/ortools/sat/cp_model_solver.cc index 8274afb99a..fc86b32a18 100644 --- a/ortools/sat/cp_model_solver.cc +++ b/ortools/sat/cp_model_solver.cc @@ -1733,13 +1733,6 @@ void SolveCpModelWithLNS(const CpModelProto& model_proto, int num_workers, SatParameters* parameters = model->GetOrCreate(); CpSolverResponse response = shared_response_manager->GetResponse(); - if (response.status() == CpSolverStatus::UNKNOWN) { - parameters->set_stop_after_first_solution(true); - LoadCpModel(model_proto, shared_response_manager, model); - SolveLoadedCpModel(model_proto, shared_response_manager, model); - response = shared_response_manager->GetResponse(); - } - if (response.status() != CpSolverStatus::FEASIBLE) return; const bool focus_on_decision_variables = parameters->lns_focus_on_decision_variables(); @@ -1772,6 +1765,28 @@ void SolveCpModelWithLNS(const CpModelProto& model_proto, int num_workers, absl::make_unique( &helper, model, "rins")); } + bool all_generators_require_solution = true; + for (int i = 0; i < generators.size(); ++i) { + if (!generators[i]->NeedsFirstSolution()) { + all_generators_require_solution = false; + break; + } + } + if (response.status() == CpSolverStatus::UNKNOWN) { + if (all_generators_require_solution) { + parameters->set_stop_after_first_solution(true); + LoadCpModel(model_proto, shared_response_manager, model); + SolveLoadedCpModel(model_proto, shared_response_manager, model); + response = shared_response_manager->GetResponse(); + if (response.status() != CpSolverStatus::FEASIBLE) { + return; + } + } + } + if (response.status() != CpSolverStatus::FEASIBLE && + response.status() != CpSolverStatus::UNKNOWN) { + return; + } // The "optimal" difficulties do not have to be the same for different // generators. @@ -1836,7 +1851,11 @@ void SolveCpModelWithLNS(const CpModelProto& model_proto, int num_workers, response.objective_value() == response.best_objective_bound(); }, [&](int64 seed) { - const int selected_generator = seed % generators.size(); + int selected_generator = seed % generators.size(); + while (response.status() == CpSolverStatus::UNKNOWN && + generators[selected_generator]->NeedsFirstSolution()) { + selected_generator = (selected_generator + 1) % generators.size(); + } AdaptiveParameterValue& difficulty = difficulties[selected_generator]; const double saved_difficulty = difficulty.value(); Neighborhood neighborhood = generators[selected_generator]->Generate( @@ -2130,7 +2149,7 @@ CpSolverResponse SolveCpModel(const CpModelProto& model_proto, Model* model) { { const std::string error = ValidateCpModel(model_proto); if (!error.empty()) { - VLOG(1) << error; + LOG_IF(INFO, log_search) << error; CpSolverResponse response; response.set_status(CpSolverStatus::MODEL_INVALID); LOG_IF(INFO, log_search) << CpSolverResponseStats(response); diff --git a/ortools/sat/doc/boolean_logic.md b/ortools/sat/doc/boolean_logic.md index 0e584a7cb8..1ac6600000 100644 --- a/ortools/sat/doc/boolean_logic.md +++ b/ortools/sat/doc/boolean_logic.md @@ -41,8 +41,8 @@ https://en.wikipedia.org/wiki/Boolean_satisfiability_problem#Basic_definitions_a ## Boolean variables and literals -We can create a Boolean variable 'x' and a literal 'not_x' equal to the logical -negation of 'x'. +We can create a Boolean variable `x` and a literal `not_x` equal to the logical +negation of `x`. ### Python code @@ -241,10 +241,10 @@ The CP-SAT solver supports *half-reified* constraints, also called x implies constraint -where the constraint must hold if *x* is true. +where the constraint must hold if `x` is true. Please note that this is not an equivalence relation. The constraint can still -be true if *x* is false. +be true if `x` is false. So we can write b => And(x, not y). That is, if b is true, then x is true and y is false. Note that in this particular example, there are multiple ways to @@ -401,13 +401,12 @@ A useful construct is the product `p` of two Boolean variables `x` and `y`. p == x * y - This is equivalent to the logical relation p <=> x and y -This is encoded using one bool_or constraint and two implications. The -following code samples output this truth table: +This is encoded using one bool_or constraint and two implications. The following +code samples output this truth table: x = 0 y = 0 p = 0 x = 1 y = 0 p = 0 diff --git a/ortools/sat/integer_search.cc b/ortools/sat/integer_search.cc index de59ff1c7e..8c3021936e 100644 --- a/ortools/sat/integer_search.cc +++ b/ortools/sat/integer_search.cc @@ -665,11 +665,29 @@ SatSolver::Status SolveIntegerProblem(Model* model) { } const int old_level = solver->CurrentDecisionLevel(); + // No decision means that we reached a leave of the search tree and that + // we have a feasible solution. if (decision == kNoLiteralIndex) { SolutionDetails* solution_details = model->Mutable(); if (solution_details != nullptr) { solution_details->LoadFromTrail(*integer_trail); } + + // Save the current polarity of all Booleans in the solution. It will be + // followed for the next SAT decisions. This is known to be a good policy + // for optimization problem. Note that for decision problem we don't care + // since we are just done as soon as a solution is found. + // + // This idea is kind of "well known", see for instance the "LinSBPS" + // submission to the maxSAT 2018 competition by Emir Demirovic and Peter + // Stuckey where they show it is a good idea and provide more references. + if (model->GetOrCreate()->use_optimization_hints()) { + auto* sat_decision = model->GetOrCreate(); + const auto& trail = *model->GetOrCreate(); + for (int i = 0; i < trail.Index(); ++i) { + sat_decision->SetAssignmentPreference(trail[i], 0.0); + } + } return SatSolver::FEASIBLE; } @@ -689,8 +707,10 @@ SatSolver::Status SolveIntegerProblem(Model* model) { solver->AdvanceDeterministicTime(time_limit); if (!solver->ReapplyAssumptionsIfNeeded()) return solver->UnsatStatus(); - if (model->GetOrCreate()->solution_count > 0 && - model->Get() != nullptr) { + if (model->Get() != nullptr) { + // If RINS is activated, we need to make sure the SolutionDetails is + // created. + model->GetOrCreate(); num_decisions_without_rins++; // TODO(user): Experiment more around dynamically changing the // threshold for trigerring RINS. Alternatively expose this as parameter diff --git a/ortools/sat/rins.cc b/ortools/sat/rins.cc index 67f0e2b350..4b4920a53c 100644 --- a/ortools/sat/rins.cc +++ b/ortools/sat/rins.cc @@ -24,15 +24,16 @@ bool SharedRINSNeighborhoodManager::AddNeighborhood( const RINSNeighborhood& rins_neighborhood) { // Don't store this neighborhood if the current storage is already too large. // TODO(user): Consider instead removing one of the older neighborhoods. - if (total_num_fixed_vars_ + rins_neighborhood.fixed_vars.size() > - max_fixed_vars()) { + const int64 neighborhood_size = rins_neighborhood.fixed_vars.size() + + rins_neighborhood.reduced_domain_vars.size(); + if (total_stored_vars_ + neighborhood_size > max_stored_vars()) { return false; } absl::MutexLock lock(&mutex_); - total_num_fixed_vars_ += rins_neighborhood.fixed_vars.size(); + total_stored_vars_ += neighborhood_size; neighborhoods_.push_back(std::move(rins_neighborhood)); - VLOG(1) << "total fixed vars: " << total_num_fixed_vars_; + VLOG(1) << "total stored vars: " << total_stored_vars_; return true; } @@ -47,8 +48,9 @@ SharedRINSNeighborhoodManager::GetUnexploredNeighborhood() { // return the last added neighborhood and remove it. const RINSNeighborhood neighborhood = std::move(neighborhoods_.back()); neighborhoods_.pop_back(); - total_num_fixed_vars_ -= neighborhood.fixed_vars.size(); - VLOG(1) << "total fixed vars: " << total_num_fixed_vars_; + total_stored_vars_ -= + neighborhood.fixed_vars.size() + neighborhood.reduced_domain_vars.size(); + VLOG(1) << "total stored vars: " << total_stored_vars_; return neighborhood; } @@ -57,36 +59,60 @@ void AddRINSNeighborhood(Model* model) { auto* solution_details = model->Get(); const RINSVariables& rins_vars = *model->GetOrCreate(); - if (solution_details == nullptr) return; - RINSNeighborhood rins_neighborhood; + const double tolerance = 1e-6; for (const RINSVariable& rins_var : rins_vars.vars) { const IntegerVariable positive_var = rins_var.positive_var; if (integer_trail->IsCurrentlyIgnored(positive_var)) continue; - // TODO(user): Perform caching to make this more efficient. LinearProgrammingConstraint* lp = rins_var.lp; if (lp == nullptr || !lp->HasSolution()) continue; - if (positive_var >= solution_details->best_solution.size()) continue; + const double lp_value = lp->GetSolutionValue(positive_var); + + if (solution_details == nullptr || solution_details->solution_count == 0) { + // The tolerance make sure that if the lp_values is close to an integer, + // then we fix the variable to this integer value. + const int64 domain_lb = + static_cast(std::floor(lp_value + tolerance)); + const int64 domain_ub = + static_cast(std::ceil(lp_value - tolerance)); + if (domain_lb >= integer_trail->LowerBound(positive_var) && + domain_ub <= integer_trail->UpperBound(positive_var)) { + if (domain_lb == domain_ub) { + rins_neighborhood.fixed_vars.push_back({rins_var, domain_lb}); + } else { + rins_neighborhood.reduced_domain_vars.push_back( + {rins_var, {domain_lb, domain_ub}}); + } - const IntegerValue best_solution_value = - solution_details->best_solution[positive_var]; - if (std::abs(best_solution_value.value() - - lp->GetSolutionValue(positive_var)) < 1e-4) { - if (best_solution_value >= integer_trail->LowerBound(positive_var) || - best_solution_value <= integer_trail->UpperBound(positive_var)) { - rins_neighborhood.fixed_vars.push_back( - {rins_var, best_solution_value.value()}); } else { // The last lp_solution might not always be up to date. - VLOG(2) << "RINS common value out of bounds: " << best_solution_value + VLOG(2) << "RENS lp value out of bounds: " << lp_value << " LB: " << integer_trail->LowerBound(positive_var) << " UB: " << integer_trail->UpperBound(positive_var); } + } else { + if (positive_var >= solution_details->best_solution.size()) continue; + + const IntegerValue best_solution_value = + solution_details->best_solution[positive_var]; + if (std::abs(best_solution_value.value() - lp_value) < 1e-4) { + if (best_solution_value >= integer_trail->LowerBound(positive_var) && + best_solution_value <= integer_trail->UpperBound(positive_var)) { + rins_neighborhood.fixed_vars.push_back( + {rins_var, best_solution_value.value()}); + } else { + // The last lp_solution might not always be up to date. + VLOG(2) << "RINS common value out of bounds: " << best_solution_value + << " LB: " << integer_trail->LowerBound(positive_var) + << " UB: " << integer_trail->UpperBound(positive_var); + } + } } } + model->Mutable()->AddNeighborhood( rins_neighborhood); } diff --git a/ortools/sat/rins.h b/ortools/sat/rins.h index f64682ca63..7f2188a576 100644 --- a/ortools/sat/rins.h +++ b/ortools/sat/rins.h @@ -53,7 +53,10 @@ struct RINSVariables { }; struct RINSNeighborhood { + // A variable will appear only once and not in both vectors. std::vector> fixed_vars; + std::vector>> + reduced_domain_vars; }; // Shared object to pass around generated RINS neighborhoods across workers. @@ -74,7 +77,7 @@ class SharedRINSNeighborhoodManager { // This is the memory limit we use to avoid storing too many neighborhoods. // The sum of the fixed variables of the stored neighborhood will always // stays smaller than this. - int64 max_fixed_vars() const { return 100 * num_model_vars_; } + int64 max_stored_vars() const { return 100 * num_model_vars_; } private: // Used while adding and removing neighborhoods. @@ -84,9 +87,9 @@ class SharedRINSNeighborhoodManager { // collection. std::vector neighborhoods_; - // This is the sum of number of fixed variables across all the shared - // neighborhoods. This is used for controlling the size of storage. - int64 total_num_fixed_vars_ = 0; + // This is the sum of number of fixed and reduced variables across all the + // shared neighborhoods. This is used for controlling the size of storage. + int64 total_stored_vars_ = 0; const int64 num_model_vars_; };