25#include "absl/log/check.h"
26#include "absl/random/bit_gen_ref.h"
27#include "absl/random/distributions.h"
28#include "absl/strings/str_format.h"
29#include "absl/strings/string_view.h"
34#include "ortools/bop/bop_parameters.pb.h"
39#include "ortools/glop/parameters.pb.h"
43#include "ortools/sat/boolean_problem.pb.h"
48#include "ortools/sat/sat_parameters.pb.h"
59using ::operations_research::glop::ColIndex;
61using ::operations_research::glop::GlopParameters;
62using ::operations_research::glop::RowIndex;
65 int64_t lower_bound) {
74bool AllIntegralValues(
const DenseRow& values,
double tolerance) {
78 if (value >= tolerance && value + tolerance < 1.0) {
87 CHECK_EQ(
solution->Size(), values.size());
88 for (VariableIndex var(0); var <
solution->Size(); ++var) {
89 solution->SetValue(var, round(values[ColIndex(var.value())]));
99 absl::string_view name, Policy policy)
110 if (state_update_stamp_ == problem_state.
update_stamp()) {
117 sat_solver_ = std::make_unique<sat::SatSolver>();
121 .exploit_symmetry_in_sat_first_solution()) {
122 std::vector<std::unique_ptr<SparsePermutation>> generators;
125 std::unique_ptr<sat::SymmetryPropagator> propagator(
127 for (
int i = 0; i < generators.size(); ++i) {
128 propagator->AddSymmetry(std::move(generators[i]));
130 sat_solver_->AddPropagator(propagator.get());
131 sat_solver_->TakePropagatorOwnership(std::move(propagator));
143 for (ColIndex col(0); col < problem_state.
lp_values().size(); ++col) {
144 const double value = problem_state.
lp_values()[col];
145 sat_solver_->SetAssignmentPreference(
146 sat::Literal(sat::BooleanVariable(col.value()), round(value) == 1),
147 1 - fabs(value - round(value)));
156 sat_solver_->SetAssignmentPreference(
157 sat::Literal(sat::BooleanVariable(i),
168 if (abort_)
return false;
173 problem_state.assignment_preference().empty()) {
180 const BopParameters& parameters,
const ProblemState& problem_state,
182 CHECK(learned_info !=
nullptr);
184 learned_info->
Clear();
187 SynchronizeIfNeeded(problem_state);
190 sat::SatParameters sat_params;
191 sat_params.set_max_time_in_seconds(
time_limit->GetTimeLeft());
192 sat_params.set_max_deterministic_time(
time_limit->GetDeterministicTimeLeft());
193 sat_params.set_random_seed(parameters.random_seed());
199 sat_params.set_max_number_of_conflicts(
200 parameters.guided_sat_conflicts_chunk());
201 sat_solver_->SetParameters(sat_params);
203 const double initial_deterministic_time = sat_solver_->deterministic_time();
205 time_limit->AdvanceDeterministicTime(sat_solver_->deterministic_time() -
206 initial_deterministic_time);
210 if (problem_state.upper_bound() != std::numeric_limits<int64_t>::max()) {
212 learned_info->
lower_bound = problem_state.upper_bound();
223 return SolutionStatus(learned_info->
solution, problem_state.lower_bound());
233 absl::string_view name,
const BopParameters& parameters,
237 sat_propagator_(sat_propagator) {}
248 const BopParameters& parameters,
const ProblemState& problem_state,
250 CHECK(learned_info !=
nullptr);
252 learned_info->
Clear();
255 const sat::SatParameters saved_params = sat_propagator_->parameters();
256 const std::vector<std::pair<sat::Literal, float>> saved_prefs =
257 sat_propagator_->AllPreferences();
259 const int kMaxNumConflicts = 10;
262 : std::numeric_limits<int64_t>::max();
263 int64_t remaining_num_conflicts =
264 parameters.max_number_of_conflicts_in_random_solution_generation();
265 int64_t old_num_failures = 0;
269 bool objective_need_to_be_overconstrained =
270 (best_cost != std::numeric_limits<int64_t>::max());
272 bool solution_found =
false;
273 while (remaining_num_conflicts > 0 && !
time_limit->LimitReached()) {
274 sat_propagator_->Backtrack(0);
275 old_num_failures = sat_propagator_->num_failures();
277 sat::SatParameters sat_params = saved_params;
279 sat_params.set_max_number_of_conflicts(kMaxNumConflicts);
280 sat_propagator_->SetParameters(sat_params);
281 sat_propagator_->ResetDecisionHeuristic();
283 if (objective_need_to_be_overconstrained) {
284 if (!AddObjectiveConstraint(
286 true, sat::Coefficient(best_cost) - 1, sat_propagator_)) {
289 return best_cost == std::numeric_limits<int64_t>::max()
293 objective_need_to_be_overconstrained =
false;
297 const int preference = absl::Uniform(random_, 0, 4);
298 if (preference == 0) {
301 }
else if (preference == 1 && !problem_state.
lp_values().empty()) {
303 for (ColIndex col(0); col < problem_state.
lp_values().size(); ++col) {
304 const double value = problem_state.
lp_values()[col];
306 sat::Literal(sat::BooleanVariable(col.value()), round(value) == 1),
307 1 - fabs(value - round(value)));
312 sat_propagator_->SolveWithTimeLimit(
time_limit);
314 objective_need_to_be_overconstrained =
true;
315 solution_found =
true;
323 return best_cost == std::numeric_limits<int64_t>::max()
330 remaining_num_conflicts -=
331 sat_propagator_->num_failures() - old_num_failures;
339 CHECK_EQ(0, sat_propagator_->AssumptionLevel());
340 sat_propagator_->RestoreSolverToAssumptionLevel();
341 sat_propagator_->SetParameters(saved_params);
342 sat_propagator_->ResetDecisionHeuristic();
343 for (
const auto [literal, weight] : saved_prefs) {
344 sat_propagator_->SetAssignmentPreference(literal, weight);
348 if (sat_propagator_->IsModelUnsat()) {
351 return best_cost == std::numeric_limits<int64_t>::max()
366 absl::string_view name)
368 parameters_(parameters),
370 lp_model_loaded_(
false),
376 num_fixed_variables_(-1),
377 problem_already_solved_(false),
384 if (state_update_stamp_ == problem_state.
update_stamp()) {
392 parameters_.max_lp_solve_for_feasibility_problems() >= 0 &&
393 num_full_solves_ >= parameters_.max_lp_solve_for_feasibility_problems()) {
399 int num_fixed_variables = 0;
400 for (
const bool is_fixed : problem_state.
is_fixed()) {
402 ++num_fixed_variables;
405 problem_already_solved_ =
406 problem_already_solved_ && num_fixed_variables_ >= num_fixed_variables;
410 num_fixed_variables_ = num_fixed_variables;
411 if (!lp_model_loaded_) {
415 lp_model_loaded_ =
true;
417 for (VariableIndex var(0); var < problem_state.
is_fixed().size(); ++var) {
421 lp_model_.SetVariableBounds(ColIndex(var.value()), value, value);
426 if (parameters_.use_learned_binary_clauses_in_lp()) {
427 for (
const sat::BinaryClause& clause :
429 const RowIndex constraint_index = lp_model_.CreateNewConstraint();
430 const int64_t coefficient_a = clause.a.IsPositive() ? 1 : -1;
431 const int64_t coefficient_b = clause.b.IsPositive() ? 1 : -1;
432 const int64_t rhs = 1 + (clause.a.IsPositive() ? 0 : -1) +
433 (clause.b.IsPositive() ? 0 : -1);
434 const ColIndex col_a(clause.a.Variable().value());
435 const ColIndex col_b(clause.b.Variable().value());
436 const std::string name_a = lp_model_.GetVariableName(col_a);
437 const std::string name_b = lp_model_.GetVariableName(col_b);
439 lp_model_.SetConstraintName(
441 (clause.a.IsPositive() ? name_a :
"not(" + name_a +
")") +
" or " +
442 (clause.b.IsPositive() ? name_b :
"not(" + name_b +
")"));
443 lp_model_.SetCoefficient(constraint_index, col_a, coefficient_a);
444 lp_model_.SetCoefficient(constraint_index, col_b, coefficient_b);
451 scaled_solution_cost_ =
465 return problem_state.original_problem().objective().literals_size() > 0 ||
466 parameters_.max_lp_solve_for_feasibility_problems() != 0;
472 CHECK(learned_info !=
nullptr);
474 learned_info->
Clear();
477 SynchronizeIfNeeded(problem_state);
485 <<
" status: " << GetProblemStatusString(lp_status);
490 problem_already_solved_ =
true;
501 learned_info->
lp_values = lp_solver_.variable_values();
506 double lower_bound = lp_solver_.GetObjectiveValue();
507 if (parameters_.use_lp_strong_branching()) {
509 ComputeLowerBoundUsingStrongBranching(learned_info,
time_limit);
511 << absl::StrFormat(
"%.6f", lower_bound)
512 <<
" using strong branching.";
515 const int tolerance_sign = scaling_ < 0 ? 1 : -1;
516 const double unscaled_cost =
519 lp_solver_.GetParameters().solution_feasibility_tolerance()) /
522 learned_info->
lower_bound =
static_cast<int64_t
>(ceil(unscaled_cost));
524 if (AllIntegralValues(
526 lp_solver_.GetParameters().primal_feasibility_tolerance())) {
542 GlopParameters glop_params;
543 if (incremental_solve) {
544 glop_params.set_use_dual_simplex(
true);
545 glop_params.set_allow_simplex_algorithm_change(
true);
546 glop_params.set_use_preprocessing(
false);
547 lp_solver_.SetParameters(glop_params);
550 parameters_.lp_max_deterministic_time());
552 lp_model_, nested_time_limit.GetTimeLimit());
556double LinearRelaxation::ComputeLowerBoundUsingStrongBranching(
558 const glop::DenseRow initial_lp_values = lp_solver_.variable_values();
559 const double tolerance =
560 lp_solver_.GetParameters().primal_feasibility_tolerance();
561 double best_lp_objective = lp_solver_.GetObjectiveValue();
562 for (glop::ColIndex col(0); col < initial_lp_values.size(); ++col) {
570 if (lp_model_.variable_lower_bounds()[col] ==
571 lp_model_.variable_upper_bounds()[col]) {
574 CHECK_EQ(0.0, lp_model_.variable_lower_bounds()[col]);
575 CHECK_EQ(1.0, lp_model_.variable_upper_bounds()[col]);
586 (initial_lp_values[col] < tolerance ||
587 initial_lp_values[col] + tolerance > 1)) {
591 double objective_true = best_lp_objective;
592 double objective_false = best_lp_objective;
595 lp_model_.SetVariableBounds(col, 1.0, 1.0);
603 objective_true = lp_solver_.GetObjectiveValue();
606 lp_model_.SetVariableBounds(col, 0.0, 0.0);
610 objective_false = lp_solver_.GetObjectiveValue();
614 lp_model_.IsMaximizationProblem()
615 ? std::min(best_lp_objective,
616 std::max(objective_true, objective_false))
617 : std::max(best_lp_objective,
618 std::min(objective_true, objective_false));
622 if (CostIsWorseThanSolution(objective_true, tolerance)) {
625 lp_model_.SetVariableBounds(col, 0.0, 0.0);
626 learned_info->fixed_literals.push_back(
627 sat::Literal(sat::BooleanVariable(col.value()),
false));
628 }
else if (CostIsWorseThanSolution(objective_false, tolerance)) {
631 lp_model_.SetVariableBounds(col, 1.0, 1.0);
632 learned_info->fixed_literals.push_back(
633 sat::Literal(sat::BooleanVariable(col.value()),
true));
636 lp_model_.SetVariableBounds(col, 0.0, 1.0);
639 return best_lp_objective;
642bool LinearRelaxation::CostIsWorseThanSolution(
double scaled_cost,
643 double tolerance)
const {
644 return lp_model_.IsMaximizationProblem()
645 ? scaled_cost + tolerance < scaled_solution_cost_
646 : scaled_cost > scaled_solution_cost_ + tolerance;
BopOptimizerBase(absl::string_view name)
@ ABORT
There is no need to call this optimizer again on the same problem state.
const std::string & name() const
Returns the name given at construction.
~BopRandomFirstSolutionGenerator() override
BopRandomFirstSolutionGenerator(absl::string_view name, const BopParameters ¶meters, sat::SatSolver *sat_propagator, absl::BitGenRef random)
Status Optimize(const BopParameters ¶meters, const ProblemState &problem_state, LearnedInfo *learned_info, TimeLimit *time_limit) override
bool ShouldBeRun(const ProblemState &problem_state) const override
Only run the RandomFirstSolution when there is an objective to minimize.
double GetScaledCost() const
Status Optimize(const BopParameters ¶meters, const ProblemState &problem_state, LearnedInfo *learned_info, TimeLimit *time_limit) override
bool ShouldBeRun(const ProblemState &problem_state) const override
~GuidedSatFirstSolutionGenerator() override
GuidedSatFirstSolutionGenerator(absl::string_view name, Policy policy)
~LinearRelaxation() override
Status Optimize(const BopParameters ¶meters, const ProblemState &problem_state, LearnedInfo *learned_info, TimeLimit *time_limit) override
LinearRelaxation(const BopParameters ¶meters, absl::string_view name)
bool ShouldBeRun(const ProblemState &problem_state) const override
bool IsVariableFixed(VariableIndex var) const
int64_t update_stamp() const
const BopSolution & solution() const
const std::vector< sat::BinaryClause > & NewlyAddedBinaryClauses() const
Returns the newly added binary clause since the last SynchronizationDone().
std::vector< bool > assignment_preference() const
const glop::DenseRow & lp_values() const
const BopParameters & GetParameters() const
const util_intops::StrongVector< VariableIndex, bool > & is_fixed() const
const sat::LinearBooleanProblem & original_problem() const
bool GetVariableFixedValue(VariableIndex var) const
Fractional GetObjectiveValue() const
Returns the objective value of the solution with its offset and scaling.
void SetAssignmentPreference(Literal literal, float weight)
constexpr double kInfinity
Infinity for type Fractional.
BopOptimizerBase::Status LoadStateProblemToSatSolver(const ProblemState &problem_state, sat::SatSolver *sat_solver)
void SatAssignmentToBopSolution(const sat::VariablesAssignment &assignment, BopSolution *solution)
void ExtractLearnedInfoFromSatSolver(sat::SatSolver *solver, LearnedInfo *info)
constexpr double kInfinity
Infinity for type Fractional.
StrictITIVector< ColIndex, Fractional > DenseRow
Row-vector types. Row-vector types are indexed by a column index.
ProblemStatus
Different statuses for a given problem.
@ INIT
The solver didn't had a chance to prove anything.
For infeasible and unbounded see Not checked if options check_solutions_if_inf_or_unbounded and the If options first_solution_only is false
problem is infeasible or unbounded (default).
SolutionStatus
Feasibility of a primal or dual solution as claimed by the solver.
void RandomizeDecisionHeuristic(absl::BitGenRef random, SatParameters *parameters)
Randomizes the decision heuristic of the given SatParameters.
void ConvertBooleanProblemToLinearProgram(const LinearBooleanProblem &problem, glop::LinearProgram *lp)
Converts a Boolean optimization problem to its lp formulation.
void FindLinearBooleanProblemSymmetries(const LinearBooleanProblem &problem, std::vector< std::unique_ptr< SparsePermutation > > *generators)
void UseObjectiveForSatAssignmentPreference(const LinearBooleanProblem &problem, SatSolver *solver)
In SWIG mode, we don't want anything besides these top-level includes.
Select next search node to expand Select next item_i to add this new search node to the search Generate a new search node where item_i is not in the knapsack Check validity of this new partial solution(using propagators) - If valid
int64_t lower_bound
A lower bound (for multi-threading purpose).
BopSolution solution
New solution. Note that the solution might be infeasible.