1 // Copyright 2010-2021 Google LLC
2 // Licensed under the Apache License, Version 2.0 (the "License");
3 // you may not use this file except in compliance with the License.
4 // You may obtain a copy of the License at
5 //
6 // http://www.apache.org/licenses/LICENSE-2.0
7 //
8 // Unless required by applicable law or agreed to in writing, software
9 // distributed under the License is distributed on an "AS IS" BASIS,
10 // WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
11 // See the License for the specific language governing permissions and
12 // limitations under the License.
13
14 #include "ortools/bop/bop_ls.h"
15
16 #include <cstdint>
17 #include <limits>
18
19 #include "absl/memory/memory.h"
20 #include "absl/strings/str_format.h"
21 #include "ortools/base/strong_vector.h"
22 #include "ortools/bop/bop_util.h"
23 #include "ortools/sat/boolean_problem.h"
24
25 namespace operations_research {
26 namespace bop {
27
28 using ::operations_research::sat::LinearBooleanConstraint;
29 using ::operations_research::sat::LinearBooleanProblem;
30 using ::operations_research::sat::LinearObjective;
31
32 //------------------------------------------------------------------------------
33 // LocalSearchOptimizer
34 //------------------------------------------------------------------------------
35
36 LocalSearchOptimizer::LocalSearchOptimizer(const std::string& name,
37 int max_num_decisions,
38 absl::BitGenRef random,
39 sat::SatSolver* sat_propagator)
40 : BopOptimizerBase(name),
41 state_update_stamp_(ProblemState::kInitialStampValue),
42 max_num_decisions_(max_num_decisions),
43 sat_wrapper_(sat_propagator),
44 assignment_iterator_(),
45 random_(random) {}
46
47 LocalSearchOptimizer::~LocalSearchOptimizer() {}
48
49 bool LocalSearchOptimizer::ShouldBeRun(
50 const ProblemState& problem_state) const {
51 return problem_state.solution().IsFeasible();
52 }
53
54 BopOptimizerBase::Status LocalSearchOptimizer::Optimize(
55 const BopParameters& parameters, const ProblemState& problem_state,
56 LearnedInfo* learned_info, TimeLimit* time_limit) {
57 CHECK(learned_info != nullptr);
58 CHECK(time_limit != nullptr);
59 learned_info->Clear();
60
61 if (assignment_iterator_ == nullptr) {
62 assignment_iterator_ = absl::make_unique<LocalSearchAssignmentIterator>(
63 problem_state, max_num_decisions_,
64 parameters.max_num_broken_constraints_in_ls(), random_, &sat_wrapper_);
65 }
66
67 if (state_update_stamp_ != problem_state.update_stamp()) {
68 // We have a new problem_state.
69 state_update_stamp_ = problem_state.update_stamp();
70 assignment_iterator_->Synchronize(problem_state);
71 }
72 assignment_iterator_->SynchronizeSatWrapper();
73
74 double prev_deterministic_time = assignment_iterator_->deterministic_time();
75 assignment_iterator_->UseTranspositionTable(
76 parameters.use_transposition_table_in_ls());
77 assignment_iterator_->UsePotentialOneFlipRepairs(
78 parameters.use_potential_one_flip_repairs_in_ls());
79 int64_t num_assignments_to_explore =
80 parameters.max_number_of_explored_assignments_per_try_in_ls();
81
82 while (!time_limit->LimitReached() && num_assignments_to_explore > 0 &&
83 assignment_iterator_->NextAssignment()) {
84 time_limit->AdvanceDeterministicTime(
85 assignment_iterator_->deterministic_time() - prev_deterministic_time);
86 prev_deterministic_time = assignment_iterator_->deterministic_time();
87 --num_assignments_to_explore;
88 }
89 if (sat_wrapper_.IsModelUnsat()) {
90 // TODO(user): we do that all the time, return an UNSAT satus instead and
91 // do this only once.
92 return problem_state.solution().IsFeasible()
93 ? BopOptimizerBase::OPTIMAL_SOLUTION_FOUND
94 : BopOptimizerBase::INFEASIBLE;
95 }
96
97 // TODO(user): properly abort when we found a new solution and then finished
98 // the ls? note that this is minor.
99 sat_wrapper_.ExtractLearnedInfo(learned_info);
100 if (assignment_iterator_->BetterSolutionHasBeenFound()) {
101 // TODO(user): simply use vector<bool> instead of a BopSolution internally.
102 learned_info->solution = assignment_iterator_->LastReferenceAssignment();
103 return BopOptimizerBase::SOLUTION_FOUND;
104 }
105
106 if (time_limit->LimitReached()) {
107 // The time limit is reached without finding a solution.
108 return BopOptimizerBase::LIMIT_REACHED;
109 }
110
111 if (num_assignments_to_explore <= 0) {
112 // Explore the remaining assignments in a future call to Optimize().
113 return BopOptimizerBase::CONTINUE;
114 }
115
116 // All assignments reachable in max_num_decisions_ or less have been explored,
117 // don't call optimize() with the same initial solution again.
118 return BopOptimizerBase::ABORT;
119 }
120
121 //------------------------------------------------------------------------------
122 // BacktrackableIntegerSet
123 //------------------------------------------------------------------------------
124
125 template <typename IntType>
126 void BacktrackableIntegerSet<IntType>::ClearAndResize(IntType n) {
127 size_ = 0;
128 saved_sizes_.clear();
129 saved_stack_sizes_.clear();
130 stack_.clear();
131 in_stack_.assign(n.value(), false);
132 }
133
134 template <typename IntType>
135 void BacktrackableIntegerSet<IntType>::ChangeState(IntType i,
136 bool should_be_inside) {
137 size_ += should_be_inside ? 1 : -1;
138 if (!in_stack_[i.value()]) {
139 in_stack_[i.value()] = true;
140 stack_.push_back(i);
141 }
142 }
143
144 template <typename IntType>
145 void BacktrackableIntegerSet<IntType>::AddBacktrackingLevel() {
146 saved_stack_sizes_.push_back(stack_.size());
147 saved_sizes_.push_back(size_);
148 }
149
150 template <typename IntType>
151 void BacktrackableIntegerSet<IntType>::BacktrackOneLevel() {
152 if (saved_stack_sizes_.empty()) {
153 BacktrackAll();
154 } else {
155 for (int i = saved_stack_sizes_.back(); i < stack_.size(); ++i) {
156 in_stack_[stack_[i].value()] = false;
157 }
158 stack_.resize(saved_stack_sizes_.back());
159 saved_stack_sizes_.pop_back();
160 size_ = saved_sizes_.back();
161 saved_sizes_.pop_back();
162 }
163 }
164
165 template <typename IntType>
166 void BacktrackableIntegerSet<IntType>::BacktrackAll() {
167 for (int i = 0; i < stack_.size(); ++i) {
168 in_stack_[stack_[i].value()] = false;
169 }
170 stack_.clear();
171 saved_stack_sizes_.clear();
172 size_ = 0;
173 saved_sizes_.clear();
174 }
175
176 // Explicit instantiation of BacktrackableIntegerSet.
177 // TODO(user): move the code in a separate .h and -inl.h to avoid this.
178 template class BacktrackableIntegerSet<ConstraintIndex>;
179
180 //------------------------------------------------------------------------------
181 // AssignmentAndConstraintFeasibilityMaintainer
182 //------------------------------------------------------------------------------
183
184 AssignmentAndConstraintFeasibilityMaintainer::
185 AssignmentAndConstraintFeasibilityMaintainer(
186 const LinearBooleanProblem& problem, absl::BitGenRef random)
187 : by_variable_matrix_(problem.num_variables()),
188 constraint_lower_bounds_(),
189 constraint_upper_bounds_(),
190 assignment_(problem, "Assignment"),
191 reference_(problem, "Assignment"),
192 constraint_values_(),
193 flipped_var_trail_backtrack_levels_(),
194 flipped_var_trail_(),
195 constraint_set_hasher_(random) {
196 // Add the objective constraint as the first constraint.
197 const LinearObjective& objective = problem.objective();
198 CHECK_EQ(objective.literals_size(), objective.coefficients_size());
199 for (int i = 0; i < objective.literals_size(); ++i) {
200 CHECK_GT(objective.literals(i), 0);
201 CHECK_NE(objective.coefficients(i), 0);
202
203 const VariableIndex var(objective.literals(i) - 1);
204 const int64_t weight = objective.coefficients(i);
205 by_variable_matrix_[var].push_back(
206 ConstraintEntry(kObjectiveConstraint, weight));
207 }
208 constraint_lower_bounds_.push_back(std::numeric_limits<int64_t>::min());
209 constraint_values_.push_back(0);
210 constraint_upper_bounds_.push_back(std::numeric_limits<int64_t>::max());
211
212 // Add each constraint.
213 ConstraintIndex num_constraints_with_objective(1);
214 for (const LinearBooleanConstraint& constraint : problem.constraints()) {
215 if (constraint.literals_size() <= 2) {
216 // Infeasible binary constraints are automatically repaired by propagation
217 // (when possible). Then there are no needs to consider the binary
218 // constraints here, the propagation is delegated to the SAT propagator.
219 continue;
220 }
221
222 CHECK_EQ(constraint.literals_size(), constraint.coefficients_size());
223 for (int i = 0; i < constraint.literals_size(); ++i) {
224 const VariableIndex var(constraint.literals(i) - 1);
225 const int64_t weight = constraint.coefficients(i);
226 by_variable_matrix_[var].push_back(
227 ConstraintEntry(num_constraints_with_objective, weight));
228 }
229 constraint_lower_bounds_.push_back(
230 constraint.has_lower_bound() ? constraint.lower_bound()
231 : std::numeric_limits<int64_t>::min());
232 constraint_values_.push_back(0);
233 constraint_upper_bounds_.push_back(
234 constraint.has_upper_bound() ? constraint.upper_bound()
235 : std::numeric_limits<int64_t>::max());
236
237 ++num_constraints_with_objective;
238 }
239
240 // Initialize infeasible_constraint_set_;
241 infeasible_constraint_set_.ClearAndResize(
242 ConstraintIndex(constraint_values_.size()));
243
244 CHECK_EQ(constraint_values_.size(), constraint_lower_bounds_.size());
245 CHECK_EQ(constraint_values_.size(), constraint_upper_bounds_.size());
246 }
247
248 const ConstraintIndex
249 AssignmentAndConstraintFeasibilityMaintainer::kObjectiveConstraint(0);
250
251 void AssignmentAndConstraintFeasibilityMaintainer::SetReferenceSolution(
252 const BopSolution& reference_solution) {
253 CHECK(reference_solution.IsFeasible());
254 infeasible_constraint_set_.BacktrackAll();
255
256 assignment_ = reference_solution;
257 reference_ = assignment_;
258 flipped_var_trail_backtrack_levels_.clear();
259 flipped_var_trail_.clear();
260 AddBacktrackingLevel(); // To handle initial propagation.
261
262 // Recompute the value of all constraints.
263 constraint_values_.assign(NumConstraints(), 0);
264 for (VariableIndex var(0); var < assignment_.Size(); ++var) {
265 if (assignment_.Value(var)) {
266 for (const ConstraintEntry& entry : by_variable_matrix_[var]) {
267 constraint_values_[entry.constraint] += entry.weight;
268 }
269 }
270 }
271
272 MakeObjectiveConstraintInfeasible(1);
273 }
274
275 void AssignmentAndConstraintFeasibilityMaintainer::
276 UseCurrentStateAsReference() {
277 for (const VariableIndex var : flipped_var_trail_) {
278 reference_.SetValue(var, assignment_.Value(var));
279 }
280 flipped_var_trail_.clear();
281 flipped_var_trail_backtrack_levels_.clear();
282 AddBacktrackingLevel(); // To handle initial propagation.
283 MakeObjectiveConstraintInfeasible(1);
284 }
285
286 void AssignmentAndConstraintFeasibilityMaintainer::
287 MakeObjectiveConstraintInfeasible(int delta) {
288 CHECK(IsFeasible());
289 CHECK(flipped_var_trail_.empty());
290 constraint_upper_bounds_[kObjectiveConstraint] =
291 constraint_values_[kObjectiveConstraint] - delta;
292 infeasible_constraint_set_.BacktrackAll();
293 infeasible_constraint_set_.ChangeState(kObjectiveConstraint, true);
294 infeasible_constraint_set_.AddBacktrackingLevel();
295 CHECK(!ConstraintIsFeasible(kObjectiveConstraint));
296 CHECK(!IsFeasible());
297 if (DEBUG_MODE) {
298 for (ConstraintIndex ct(1); ct < NumConstraints(); ++ct) {
299 CHECK(ConstraintIsFeasible(ct));
300 }
301 }
302 }
303
304 void AssignmentAndConstraintFeasibilityMaintainer::Assign(
305 const std::vector<sat::Literal>& literals) {
306 for (const sat::Literal& literal : literals) {
307 const VariableIndex var(literal.Variable().value());
308 const bool value = literal.IsPositive();
309 if (assignment_.Value(var) != value) {
310 flipped_var_trail_.push_back(var);
311 assignment_.SetValue(var, value);
312 for (const ConstraintEntry& entry : by_variable_matrix_[var]) {
313 const bool was_feasible = ConstraintIsFeasible(entry.constraint);
314 constraint_values_[entry.constraint] +=
315 value ? entry.weight : -entry.weight;
316 if (ConstraintIsFeasible(entry.constraint) != was_feasible) {
317 infeasible_constraint_set_.ChangeState(entry.constraint,
318 was_feasible);
319 }
320 }
321 }
322 }
323 }
324
325 void AssignmentAndConstraintFeasibilityMaintainer::AddBacktrackingLevel() {
326 flipped_var_trail_backtrack_levels_.push_back(flipped_var_trail_.size());
327 infeasible_constraint_set_.AddBacktrackingLevel();
328 }
329
330 void AssignmentAndConstraintFeasibilityMaintainer::BacktrackOneLevel() {
331 // Backtrack each literal of the last level.
332 for (int i = flipped_var_trail_backtrack_levels_.back();
333 i < flipped_var_trail_.size(); ++i) {
334 const VariableIndex var(flipped_var_trail_[i]);
335 const bool new_value = !assignment_.Value(var);
336 DCHECK_EQ(new_value, reference_.Value(var));
337 assignment_.SetValue(var, new_value);
338 for (const ConstraintEntry& entry : by_variable_matrix_[var]) {
339 constraint_values_[entry.constraint] +=
340 new_value ? entry.weight : -entry.weight;
341 }
342 }
343 flipped_var_trail_.resize(flipped_var_trail_backtrack_levels_.back());
344 flipped_var_trail_backtrack_levels_.pop_back();
345 infeasible_constraint_set_.BacktrackOneLevel();
346 }
347
348 void AssignmentAndConstraintFeasibilityMaintainer::BacktrackAll() {
349 while (!flipped_var_trail_backtrack_levels_.empty()) BacktrackOneLevel();
350 }
351
352 const std::vector<sat::Literal>&
353 AssignmentAndConstraintFeasibilityMaintainer::PotentialOneFlipRepairs() {
354 if (!constraint_set_hasher_.IsInitialized()) {
355 InitializeConstraintSetHasher();
356 }
357
358 // First, we compute the hash that a Literal should have in order to repair
359 // all the infeasible constraint (ignoring the objective).
360 //
361 // TODO(user): If this starts to show-up in a performance profile, we can
362 // easily maintain this hash incrementally.
363 uint64_t hash = 0;
364 for (const ConstraintIndex ci : PossiblyInfeasibleConstraints()) {
365 const int64_t value = ConstraintValue(ci);
366 if (value > ConstraintUpperBound(ci)) {
367 hash ^= constraint_set_hasher_.Hash(FromConstraintIndex(ci, false));
368 } else if (value < ConstraintLowerBound(ci)) {
369 hash ^= constraint_set_hasher_.Hash(FromConstraintIndex(ci, true));
370 }
371 }
372
373 tmp_potential_repairs_.clear();
374 const auto it = hash_to_potential_repairs_.find(hash);
375 if (it != hash_to_potential_repairs_.end()) {
376 for (const sat::Literal literal : it->second) {
377 // We only returns the flips.
378 if (assignment_.Value(VariableIndex(literal.Variable().value())) !=
379 literal.IsPositive()) {
380 tmp_potential_repairs_.push_back(literal);
381 }
382 }
383 }
384 return tmp_potential_repairs_;
385 }
386
387 std::string AssignmentAndConstraintFeasibilityMaintainer::DebugString() const {
388 std::string str;
389 str += "curr: ";
390 for (bool value : assignment_) {
391 str += value ? " 1 " : " 0 ";
392 }
393 str += "\nFlipped variables: ";
394 // TODO(user): show the backtrack levels.
395 for (const VariableIndex var : flipped_var_trail_) {
396 str += absl::StrFormat(" %d", var.value());
397 }
398 str += "\nmin curr max\n";
399 for (ConstraintIndex ct(0); ct < constraint_values_.size(); ++ct) {
400 if (constraint_lower_bounds_[ct] == std::numeric_limits<int64_t>::min()) {
401 str += absl::StrFormat("- %d %d\n", constraint_values_[ct],
402 constraint_upper_bounds_[ct]);
403 } else {
404 str +=
405 absl::StrFormat("%d %d %d\n", constraint_lower_bounds_[ct],
406 constraint_values_[ct], constraint_upper_bounds_[ct]);
407 }
408 }
409 return str;
410 }
411
412 void AssignmentAndConstraintFeasibilityMaintainer::
413 InitializeConstraintSetHasher() {
414 const int num_constraints_with_objective = constraint_upper_bounds_.size();
415
416 // Initialize the potential one flip repair. Note that we ignore the
417 // objective constraint completely so that we consider a repair even if the
418 // objective constraint is not infeasible.
419 constraint_set_hasher_.Initialize(2 * num_constraints_with_objective);
420 constraint_set_hasher_.IgnoreElement(
421 FromConstraintIndex(kObjectiveConstraint, true));
422 constraint_set_hasher_.IgnoreElement(
423 FromConstraintIndex(kObjectiveConstraint, false));
424 for (VariableIndex var(0); var < by_variable_matrix_.size(); ++var) {
425 // We add two entries, one for a positive flip (from false to true) and one
426 // for a negative flip (from true to false).
427 for (const bool flip_is_positive : {true, false}) {
428 uint64_t hash = 0;
429 for (const ConstraintEntry& entry : by_variable_matrix_[var]) {
430 const bool coeff_is_positive = entry.weight > 0;
431 hash ^= constraint_set_hasher_.Hash(FromConstraintIndex(
432 entry.constraint,
433 /*up=*/flip_is_positive ? coeff_is_positive : !coeff_is_positive));
434 }
435 hash_to_potential_repairs_[hash].push_back(
436 sat::Literal(sat::BooleanVariable(var.value()), flip_is_positive));
437 }
438 }
439 }
440
441 //------------------------------------------------------------------------------
442 // OneFlipConstraintRepairer
443 //------------------------------------------------------------------------------
444
445 OneFlipConstraintRepairer::OneFlipConstraintRepairer(
446 const LinearBooleanProblem& problem,
447 const AssignmentAndConstraintFeasibilityMaintainer& maintainer,
448 const sat::VariablesAssignment& sat_assignment)
449 : by_constraint_matrix_(problem.constraints_size() + 1),
450 maintainer_(maintainer),
451 sat_assignment_(sat_assignment) {
452 // Fill the by_constraint_matrix_.
453 //
454 // IMPORTANT: The order of the constraint needs to exactly match the one of
455 // the constraint in the AssignmentAndConstraintFeasibilityMaintainer.
456
457 // Add the objective constraint as the first constraint.
458 ConstraintIndex num_constraint(0);
459 const LinearObjective& objective = problem.objective();
460 CHECK_EQ(objective.literals_size(), objective.coefficients_size());
461 for (int i = 0; i < objective.literals_size(); ++i) {
462 CHECK_GT(objective.literals(i), 0);
463 CHECK_NE(objective.coefficients(i), 0);
464
465 const VariableIndex var(objective.literals(i) - 1);
466 const int64_t weight = objective.coefficients(i);
467 by_constraint_matrix_[num_constraint].push_back(
468 ConstraintTerm(var, weight));
469 }
470
471 // Add the non-binary problem constraints.
472 for (const LinearBooleanConstraint& constraint : problem.constraints()) {
473 if (constraint.literals_size() <= 2) {
474 // Infeasible binary constraints are automatically repaired by propagation
475 // (when possible). Then there are no needs to consider the binary
476 // constraints here, the propagation is delegated to the SAT propagator.
477 continue;
478 }
479
480 ++num_constraint;
481 CHECK_EQ(constraint.literals_size(), constraint.coefficients_size());
482 for (int i = 0; i < constraint.literals_size(); ++i) {
483 const VariableIndex var(constraint.literals(i) - 1);
484 const int64_t weight = constraint.coefficients(i);
485 by_constraint_matrix_[num_constraint].push_back(
486 ConstraintTerm(var, weight));
487 }
488 }
489
490 SortTermsOfEachConstraints(problem.num_variables());
491 }
492
493 const ConstraintIndex OneFlipConstraintRepairer::kInvalidConstraint(-1);
494 const TermIndex OneFlipConstraintRepairer::kInitTerm(-1);
495 const TermIndex OneFlipConstraintRepairer::kInvalidTerm(-2);
496
497 ConstraintIndex OneFlipConstraintRepairer::ConstraintToRepair() const {
498 ConstraintIndex selected_ct = kInvalidConstraint;
499 int32_t selected_num_branches = std::numeric_limits<int32_t>::max();
500 int num_infeasible_constraints_left = maintainer_.NumInfeasibleConstraints();
501
502 // Optimization: We inspect the constraints in reverse order because the
503 // objective one will always be first (in our current code) and with some
504 // luck, we will break early instead of fully exploring it.
505 const std::vector<ConstraintIndex>& infeasible_constraints =
506 maintainer_.PossiblyInfeasibleConstraints();
507 for (int index = infeasible_constraints.size() - 1; index >= 0; --index) {
508 const ConstraintIndex& i = infeasible_constraints[index];
509 if (maintainer_.ConstraintIsFeasible(i)) continue;
510 --num_infeasible_constraints_left;
511
512 // Optimization: We return the only candidate without inspecting it.
513 // This is critical at the beginning of the search or later if the only
514 // candidate is the objective constraint which can be really long.
515 if (num_infeasible_constraints_left == 0 &&
516 selected_ct == kInvalidConstraint) {
517 return i;
518 }
519
520 const int64_t constraint_value = maintainer_.ConstraintValue(i);
521 const int64_t lb = maintainer_.ConstraintLowerBound(i);
522 const int64_t ub = maintainer_.ConstraintUpperBound(i);
523
524 int32_t num_branches = 0;
525 for (const ConstraintTerm& term : by_constraint_matrix_[i]) {
526 if (sat_assignment_.VariableIsAssigned(
527 sat::BooleanVariable(term.var.value()))) {
528 continue;
529 }
530 const int64_t new_value =
531 constraint_value +
532 (maintainer_.Assignment(term.var) ? -term.weight : term.weight);
533 if (new_value >= lb && new_value <= ub) {
534 ++num_branches;
535 if (num_branches >= selected_num_branches) break;
536 }
537 }
538
539 // The constraint can't be repaired in one decision.
540 if (num_branches == 0) continue;
init_wimax_globals(void)541 if (num_branches < selected_num_branches) {
542 selected_ct = i;
543 selected_num_branches = num_branches;
544 if (num_branches == 1) break;
545 }
546 }
547 return selected_ct;
548 }
549
550 TermIndex OneFlipConstraintRepairer::NextRepairingTerm(
551 ConstraintIndex ct_index, TermIndex init_term_index,
552 TermIndex start_term_index) const {
553 const absl::StrongVector<TermIndex, ConstraintTerm>& terms =
554 by_constraint_matrix_[ct_index];
555 const int64_t constraint_value = maintainer_.ConstraintValue(ct_index);
Dedicated_UL_Control_IE(proto_tree * uiuc_tree,gint offset,gint length,tvbuff_t * tvb)556 const int64_t lb = maintainer_.ConstraintLowerBound(ct_index);
557 const int64_t ub = maintainer_.ConstraintUpperBound(ct_index);
558
559 const TermIndex end_term_index(terms.size() + init_term_index + 1);
560 for (TermIndex loop_term_index(
561 start_term_index + 1 +
562 (start_term_index < init_term_index ? terms.size() : 0));
563 loop_term_index < end_term_index; ++loop_term_index) {
564 const TermIndex term_index(loop_term_index % terms.size());
565 const ConstraintTerm term = terms[term_index];
566 if (sat_assignment_.VariableIsAssigned(
567 sat::BooleanVariable(term.var.value()))) {
568 continue;
569 }
570 const int64_t new_value =
571 constraint_value +
572 (maintainer_.Assignment(term.var) ? -term.weight : term.weight);
573 if (new_value >= lb && new_value <= ub) {
574 return term_index;
575 }
576 }
Dedicated_MIMO_UL_Control_IE(proto_tree * uiuc_tree,gint offset,gint length,tvbuff_t * tvb)577 return kInvalidTerm;
578 }
579
580 bool OneFlipConstraintRepairer::RepairIsValid(ConstraintIndex ct_index,
581 TermIndex term_index) const {
582 if (maintainer_.ConstraintIsFeasible(ct_index)) return false;
583 const ConstraintTerm term = by_constraint_matrix_[ct_index][term_index];
584 if (sat_assignment_.VariableIsAssigned(
585 sat::BooleanVariable(term.var.value()))) {
586 return false;
587 }
588 const int64_t new_value =
589 maintainer_.ConstraintValue(ct_index) +
590 (maintainer_.Assignment(term.var) ? -term.weight : term.weight);
591
592 const int64_t lb = maintainer_.ConstraintLowerBound(ct_index);
593 const int64_t ub = maintainer_.ConstraintUpperBound(ct_index);
594 return (new_value >= lb && new_value <= ub);
595 }
UL_HARQ_Chase_Sub_Burst_IE(proto_tree * uiuc_tree,gint offset,gint length,tvbuff_t * tvb)596
597 sat::Literal OneFlipConstraintRepairer::GetFlip(ConstraintIndex ct_index,
598 TermIndex term_index) const {
599 const ConstraintTerm term = by_constraint_matrix_[ct_index][term_index];
600 const bool value = maintainer_.Assignment(term.var);
601 return sat::Literal(sat::BooleanVariable(term.var.value()), !value);
602 }
603
604 void OneFlipConstraintRepairer::SortTermsOfEachConstraints(int num_variables) {
605 absl::StrongVector<VariableIndex, int64_t> objective(num_variables, 0);
606 for (const ConstraintTerm& term :
607 by_constraint_matrix_[AssignmentAndConstraintFeasibilityMaintainer::
608 kObjectiveConstraint]) {
609 objective[term.var] = std::abs(term.weight);
610 }
611 for (absl::StrongVector<TermIndex, ConstraintTerm>& terms :
612 by_constraint_matrix_) {
613 std::sort(terms.begin(), terms.end(),
614 [&objective](const ConstraintTerm& a, const ConstraintTerm& b) {
615 return objective[a.var] > objective[b.var];
616 });
617 }
618 }
619
620 //------------------------------------------------------------------------------
621 // SatWrapper
622 //------------------------------------------------------------------------------
623
624 SatWrapper::SatWrapper(sat::SatSolver* sat_solver) : sat_solver_(sat_solver) {}
625
626 void SatWrapper::BacktrackAll() { sat_solver_->Backtrack(0); }
627
628 std::vector<sat::Literal> SatWrapper::FullSatTrail() const {
629 std::vector<sat::Literal> propagated_literals;
630 const sat::Trail& trail = sat_solver_->LiteralTrail();
631 for (int trail_index = 0; trail_index < trail.Index(); ++trail_index) {
632 propagated_literals.push_back(trail[trail_index]);
633 }
634 return propagated_literals;
635 }
UL_HARQ_IR_CTC_Sub_Burst_IE(proto_tree * uiuc_tree,gint offset,gint length,tvbuff_t * tvb)636
637 int SatWrapper::ApplyDecision(sat::Literal decision_literal,
638 std::vector<sat::Literal>* propagated_literals) {
639 CHECK(!sat_solver_->Assignment().VariableIsAssigned(
640 decision_literal.Variable()));
641 CHECK(propagated_literals != nullptr);
642
643 propagated_literals->clear();
644 const int old_decision_level = sat_solver_->CurrentDecisionLevel();
645 const int new_trail_index =
646 sat_solver_->EnqueueDecisionAndBackjumpOnConflict(decision_literal);
647 if (sat_solver_->IsModelUnsat()) {
648 return old_decision_level + 1;
649 }
650
651 // Return the propagated literals, whenever there is a conflict or not.
652 // In case of conflict, these literals will have to be added to the last
653 // decision point after backtrack.
654 const sat::Trail& propagation_trail = sat_solver_->LiteralTrail();
655 for (int trail_index = new_trail_index;
656 trail_index < propagation_trail.Index(); ++trail_index) {
657 propagated_literals->push_back(propagation_trail[trail_index]);
658 }
659
660 return old_decision_level + 1 - sat_solver_->CurrentDecisionLevel();
661 }
662
663 void SatWrapper::BacktrackOneLevel() {
664 const int old_decision_level = sat_solver_->CurrentDecisionLevel();
665 if (old_decision_level > 0) {
666 sat_solver_->Backtrack(old_decision_level - 1);
667 }
668 }
669
670 void SatWrapper::ExtractLearnedInfo(LearnedInfo* info) {
671 bop::ExtractLearnedInfoFromSatSolver(sat_solver_, info);
672 }
673
674 double SatWrapper::deterministic_time() const {
675 return sat_solver_->deterministic_time();
676 }
UL_HARQ_IR_CC_Sub_Burst_IE(proto_tree * uiuc_tree,gint offset,gint length,tvbuff_t * tvb)677
678 //------------------------------------------------------------------------------
679 // LocalSearchAssignmentIterator
680 //------------------------------------------------------------------------------
681
682 LocalSearchAssignmentIterator::LocalSearchAssignmentIterator(
683 const ProblemState& problem_state, int max_num_decisions,
684 int max_num_broken_constraints, absl::BitGenRef random,
685 SatWrapper* sat_wrapper)
686 : max_num_decisions_(max_num_decisions),
687 max_num_broken_constraints_(max_num_broken_constraints),
688 maintainer_(problem_state.original_problem(), random),
689 sat_wrapper_(sat_wrapper),
690 repairer_(problem_state.original_problem(), maintainer_,
691 sat_wrapper->SatAssignment()),
692 search_nodes_(),
693 initial_term_index_(
694 problem_state.original_problem().constraints_size() + 1,
695 OneFlipConstraintRepairer::kInitTerm),
696 use_transposition_table_(false),
697 use_potential_one_flip_repairs_(false),
698 num_nodes_(0),
699 num_skipped_nodes_(0),
700 num_improvements_(0),
701 num_improvements_by_one_flip_repairs_(0),
702 num_inspected_one_flip_repairs_(0) {}
703
704 LocalSearchAssignmentIterator::~LocalSearchAssignmentIterator() {
705 VLOG(1) << "LS " << max_num_decisions_
706 << "\n num improvements: " << num_improvements_
707 << "\n num improvements with one flip repairs: "
708 << num_improvements_by_one_flip_repairs_
709 << "\n num inspected one flip repairs: "
710 << num_inspected_one_flip_repairs_;
711 }
712
713 void LocalSearchAssignmentIterator::Synchronize(
714 const ProblemState& problem_state) {
715 better_solution_has_been_found_ = false;
716 maintainer_.SetReferenceSolution(problem_state.solution());
717 for (const SearchNode& node : search_nodes_) {
718 initial_term_index_[node.constraint] = node.term_index;
719 }
MIMO_UL_Chase_HARQ_Sub_Burst_IE(proto_tree * uiuc_tree,gint offset,gint length,tvbuff_t * tvb)720 search_nodes_.clear();
721 transposition_table_.clear();
722 num_nodes_ = 0;
723 num_skipped_nodes_ = 0;
724 }
725
726 // In order to restore the synchronization from any state, we backtrack
727 // everything and retry to take the same decisions as before. We stop at the
728 // first one that can't be taken.
729 void LocalSearchAssignmentIterator::SynchronizeSatWrapper() {
730 CHECK_EQ(better_solution_has_been_found_, false);
731 const std::vector<SearchNode> copy = search_nodes_;
732 sat_wrapper_->BacktrackAll();
733 maintainer_.BacktrackAll();
734
735 // Note(user): at this stage, the sat trail contains the fixed variables.
736 // There will almost always be at the same value in the reference solution.
737 // However since the objective may be over-constrained in the sat_solver, it
738 // is possible that some variable where propagated to some other values.
739 maintainer_.Assign(sat_wrapper_->FullSatTrail());
740
741 search_nodes_.clear();
742 for (const SearchNode& node : copy) {
743 if (!repairer_.RepairIsValid(node.constraint, node.term_index)) break;
744 search_nodes_.push_back(node);
745 ApplyDecision(repairer_.GetFlip(node.constraint, node.term_index));
746 }
747 }
748
749 void LocalSearchAssignmentIterator::UseCurrentStateAsReference() {
750 better_solution_has_been_found_ = true;
751 maintainer_.UseCurrentStateAsReference();
752 sat_wrapper_->BacktrackAll();
753
754 // Note(user): Here, there should be no discrepancies between the fixed
755 // variable and the new reference, so there is no need to do:
756 // maintainer_.Assign(sat_wrapper_->FullSatTrail());
757
758 for (const SearchNode& node : search_nodes_) {
759 initial_term_index_[node.constraint] = node.term_index;
760 }
761 search_nodes_.clear();
762 transposition_table_.clear();
763 num_nodes_ = 0;
764 num_skipped_nodes_ = 0;
765 ++num_improvements_;
766 }
767
768 bool LocalSearchAssignmentIterator::NextAssignment() {
769 if (sat_wrapper_->IsModelUnsat()) return false;
770 if (maintainer_.IsFeasible()) {
771 UseCurrentStateAsReference();
772 return true;
MIMO_UL_IR_HARQ__Sub_Burst_IE(proto_tree * uiuc_tree,gint offset,gint length,tvbuff_t * tvb)773 }
774
775 // We only look for potential one flip repairs if we reached the end of the
776 // LS tree. I tried to do that at every level, but it didn't change the
777 // result much on the set-partitionning example I was using.
778 //
779 // TODO(user): Perform more experiments with this.
780 if (use_potential_one_flip_repairs_ &&
781 search_nodes_.size() == max_num_decisions_) {
782 for (const sat::Literal literal : maintainer_.PotentialOneFlipRepairs()) {
783 if (sat_wrapper_->SatAssignment().VariableIsAssigned(
784 literal.Variable())) {
785 continue;
786 }
787 ++num_inspected_one_flip_repairs_;
788
789 // Temporarily apply the potential repair and see if it worked!
790 ApplyDecision(literal);
791 if (maintainer_.IsFeasible()) {
792 num_improvements_by_one_flip_repairs_++;
793 UseCurrentStateAsReference();
794 return true;
795 }
796 maintainer_.BacktrackOneLevel();
797 sat_wrapper_->BacktrackOneLevel();
798 }
799 }
800
801 // If possible, go deeper, i.e. take one more decision.
802 if (!GoDeeper()) {
803 // If not, backtrack to the first node that still has untried way to fix
804 // its associated constraint. Update it to the next untried way.
805 Backtrack();
806 }
807
808 // All nodes have been explored.
809 if (search_nodes_.empty()) {
810 VLOG(1) << std::string(27, ' ') + "LS " << max_num_decisions_
811 << " finished."
812 << " #explored:" << num_nodes_
813 << " #stored:" << transposition_table_.size()
814 << " #skipped:" << num_skipped_nodes_;
815 return false;
816 }
817
818 // Apply the next decision, i.e. the literal of the flipped variable.
819 const SearchNode node = search_nodes_.back();
820 ApplyDecision(repairer_.GetFlip(node.constraint, node.term_index));
821 return true;
822 }
823
824 // TODO(user): The 1.2 multiplier is an approximation only based on the time
825 // spent in the SAT wrapper. So far experiments show a good
MIMO_UL_IR_HARQ_for_CC_Sub_Burst_UIE(proto_tree * uiuc_tree,gint offset,gint length,tvbuff_t * tvb)826 // correlation with real time, but we might want to be more
827 // accurate.
828 double LocalSearchAssignmentIterator::deterministic_time() const {
829 return sat_wrapper_->deterministic_time() * 1.2;
830 }
831
832 std::string LocalSearchAssignmentIterator::DebugString() const {
833 std::string str = "Search nodes:\n";
834 for (int i = 0; i < search_nodes_.size(); ++i) {
835 str += absl::StrFormat(" %d: %d %d\n", i,
836 search_nodes_[i].constraint.value(),
837 search_nodes_[i].term_index.value());
838 }
839 return str;
840 }
841
842 void LocalSearchAssignmentIterator::ApplyDecision(sat::Literal literal) {
843 ++num_nodes_;
844 const int num_backtracks =
845 sat_wrapper_->ApplyDecision(literal, &tmp_propagated_literals_);
846
847 // Sync the maintainer with SAT.
848 if (num_backtracks == 0) {
849 maintainer_.AddBacktrackingLevel();
850 maintainer_.Assign(tmp_propagated_literals_);
851 } else {
852 CHECK_GT(num_backtracks, 0);
853 CHECK_LE(num_backtracks, search_nodes_.size());
854
855 // Only backtrack -1 decisions as the last one has not been pushed yet.
856 for (int i = 0; i < num_backtracks - 1; ++i) {
857 maintainer_.BacktrackOneLevel();
858 }
859 maintainer_.Assign(tmp_propagated_literals_);
860 search_nodes_.resize(search_nodes_.size() - num_backtracks);
861 }
862 }
863
864 void LocalSearchAssignmentIterator::InitializeTranspositionTableKey(
865 std::array<int32_t, kStoredMaxDecisions>* a) {
866 int i = 0;
867 for (const SearchNode& n : search_nodes_) {
868 // Negated because we already fliped this variable, so GetFlip() will
869 // returns the old value.
870 (*a)[i] = -repairer_.GetFlip(n.constraint, n.term_index).SignedValue();
871 ++i;
872 }
873
874 // 'a' is not zero-initialized, so we need to complete it with zeros.
875 while (i < kStoredMaxDecisions) {
876 (*a)[i] = 0;
877 ++i;
878 }
879 }
MIMO_UL_STC_HARQ_Sub_Burst_IE(proto_tree * uiuc_tree,gint offset,gint length,tvbuff_t * tvb)880
881 bool LocalSearchAssignmentIterator::NewStateIsInTranspositionTable(
882 sat::Literal l) {
883 if (search_nodes_.size() + 1 > kStoredMaxDecisions) return false;
884
885 // Fill the transposition table element, i.e the array 'a' of decisions.
886 std::array<int32_t, kStoredMaxDecisions> a;
887 InitializeTranspositionTableKey(&a);
888 a[search_nodes_.size()] = l.SignedValue();
889 std::sort(a.begin(), a.begin() + 1 + search_nodes_.size());
890
891 if (transposition_table_.find(a) == transposition_table_.end()) {
892 return false;
893 } else {
894 ++num_skipped_nodes_;
895 return true;
896 }
897 }
898
899 void LocalSearchAssignmentIterator::InsertInTranspositionTable() {
900 // If there is more decision that kStoredMaxDecisions, do nothing.
901 if (search_nodes_.size() > kStoredMaxDecisions) return;
902
903 // Fill the transposition table element, i.e the array 'a' of decisions.
904 std::array<int32_t, kStoredMaxDecisions> a;
905 InitializeTranspositionTableKey(&a);
906 std::sort(a.begin(), a.begin() + search_nodes_.size());
907
908 transposition_table_.insert(a);
909 }
910
911 bool LocalSearchAssignmentIterator::EnqueueNextRepairingTermIfAny(
912 ConstraintIndex ct_to_repair, TermIndex term_index) {
913 if (term_index == initial_term_index_[ct_to_repair]) return false;
914 if (term_index == OneFlipConstraintRepairer::kInvalidTerm) {
915 term_index = initial_term_index_[ct_to_repair];
916 }
917 while (true) {
918 term_index = repairer_.NextRepairingTerm(
919 ct_to_repair, initial_term_index_[ct_to_repair], term_index);
920 if (term_index == OneFlipConstraintRepairer::kInvalidTerm) return false;
921 if (!use_transposition_table_ ||
922 !NewStateIsInTranspositionTable(
923 repairer_.GetFlip(ct_to_repair, term_index))) {
924 search_nodes_.push_back(SearchNode(ct_to_repair, term_index));
925 return true;
926 }
927 if (term_index == initial_term_index_[ct_to_repair]) return false;
928 }
929 }
930
Power_Control_IE(proto_tree * uiuc_tree,gint offset,gint length,tvbuff_t * tvb)931 bool LocalSearchAssignmentIterator::GoDeeper() {
932 // Can we add one more decision?
933 if (search_nodes_.size() >= max_num_decisions_) {
934 return false;
935 }
936
937 // Is the number of infeasible constraints reasonable?
938 //
939 // TODO(user): Make this parameters dynamic. We can either try lower value
940 // first and increase it later, or try to dynamically change it during the
941 // search. Another idea is to have instead a "max number of constraints that
942 // can be repaired in one decision" and to take into account the number of
943 // decisions left.
944 if (maintainer_.NumInfeasibleConstraints() > max_num_broken_constraints_) {
945 return false;
946 }
947
948 // Can we find a constraint that can be repaired in one decision?
949 const ConstraintIndex ct_to_repair = repairer_.ConstraintToRepair();
950 if (ct_to_repair == OneFlipConstraintRepairer::kInvalidConstraint) {
951 return false;
952 }
953
954 // Add the new decision.
955 //
956 // TODO(user): Store the last explored term index to not start from -1 each
957 // time. This will be very useful when a backtrack occurred due to the SAT
958 // propagator. Note however that this behavior is already enforced when we use
959 // the transposition table, since we will not explore again the branches
960 // already explored.
961 return EnqueueNextRepairingTermIfAny(ct_to_repair,
962 OneFlipConstraintRepairer::kInvalidTerm);
963 }
964
965 void LocalSearchAssignmentIterator::Backtrack() {
966 while (!search_nodes_.empty()) {
967 // We finished exploring this node. Store it in the transposition table so
968 // that the same decisions will not be explored again. Note that the SAT
969 // solver may have learned more the second time the exact same decisions are
970 // seen, but we assume that it is not worth exploring again.
971 if (use_transposition_table_) InsertInTranspositionTable();
972
973 const SearchNode last_node = search_nodes_.back();
974 search_nodes_.pop_back();
975 maintainer_.BacktrackOneLevel();
976 sat_wrapper_->BacktrackOneLevel();
977 if (EnqueueNextRepairingTermIfAny(last_node.constraint,
978 last_node.term_index)) {
979 return;
980 }
981 }
982 }
983
984 } // namespace bop
985 } // namespace operations_research
986