OR-Tools  8.2
cp_model_presolve.cc
Go to the documentation of this file.
1// Copyright 2010-2018 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
15
16#include <sys/stat.h>
17
18#include <algorithm>
19#include <cstdlib>
20#include <deque>
21#include <map>
22#include <memory>
23#include <numeric>
24#include <set>
25#include <string>
26#include <utility>
27#include <vector>
28
29#include "absl/container/flat_hash_map.h"
30#include "absl/container/flat_hash_set.h"
31#include "absl/random/random.h"
32#include "absl/strings/str_join.h"
33#include "ortools/base/hash.h"
40#include "ortools/sat/circuit.h"
41#include "ortools/sat/cp_model.pb.h"
49#include "ortools/sat/probing.h"
51#include "ortools/sat/sat_parameters.pb.h"
54
55namespace operations_research {
56namespace sat {
57
58bool CpModelPresolver::RemoveConstraint(ConstraintProto* ct) {
59 ct->Clear();
60 return true;
61}
62
64 // Remove all empty constraints. Note that we need to remap the interval
65 // references.
66 std::vector<int> interval_mapping(context_->working_model->constraints_size(),
67 -1);
68 int new_num_constraints = 0;
69 const int old_num_non_empty_constraints =
70 context_->working_model->constraints_size();
71 for (int c = 0; c < old_num_non_empty_constraints; ++c) {
72 const auto type = context_->working_model->constraints(c).constraint_case();
73 if (type == ConstraintProto::ConstraintCase::CONSTRAINT_NOT_SET) continue;
74 if (type == ConstraintProto::ConstraintCase::kInterval) {
75 interval_mapping[c] = new_num_constraints;
76 }
77 context_->working_model->mutable_constraints(new_num_constraints++)
78 ->Swap(context_->working_model->mutable_constraints(c));
79 }
80 context_->working_model->mutable_constraints()->DeleteSubrange(
81 new_num_constraints, old_num_non_empty_constraints - new_num_constraints);
82 for (ConstraintProto& ct_ref :
83 *context_->working_model->mutable_constraints()) {
85 [&interval_mapping](int* ref) {
86 *ref = interval_mapping[*ref];
87 CHECK_NE(-1, *ref);
88 },
89 &ct_ref);
90 }
91}
92
93bool CpModelPresolver::PresolveEnforcementLiteral(ConstraintProto* ct) {
94 if (context_->ModelIsUnsat()) return false;
95 if (!HasEnforcementLiteral(*ct)) return false;
96
97 int new_size = 0;
98 const int old_size = ct->enforcement_literal().size();
99 for (const int literal : ct->enforcement_literal()) {
100 if (context_->LiteralIsTrue(literal)) {
101 // We can remove a literal at true.
102 context_->UpdateRuleStats("true enforcement literal");
103 continue;
104 }
105
106 if (context_->LiteralIsFalse(literal)) {
107 context_->UpdateRuleStats("false enforcement literal");
108 return RemoveConstraint(ct);
109 }
110
111 if (context_->VariableIsUniqueAndRemovable(literal)) {
112 // We can simply set it to false and ignore the constraint in this case.
113 context_->UpdateRuleStats("enforcement literal not used");
114 CHECK(context_->SetLiteralToFalse(literal));
115 return RemoveConstraint(ct);
116 }
117
118 // If the literal only appear in the objective, we might be able to fix it
119 // to false. TODO(user): generalize if the literal always appear with the
120 // same polarity.
122 const int64 obj_coeff =
124 if (RefIsPositive(literal) == (obj_coeff > 0)) {
125 // It is just more advantageous to set it to false!
126 context_->UpdateRuleStats("enforcement literal with unique direction");
127 CHECK(context_->SetLiteralToFalse(literal));
128 return RemoveConstraint(ct);
129 }
130 }
131
132 ct->set_enforcement_literal(new_size++, literal);
133 }
134 ct->mutable_enforcement_literal()->Truncate(new_size);
135 return new_size != old_size;
136}
137
138bool CpModelPresolver::PresolveBoolXor(ConstraintProto* ct) {
139 if (context_->ModelIsUnsat()) return false;
140 if (HasEnforcementLiteral(*ct)) return false;
141
142 int new_size = 0;
143 bool changed = false;
144 int num_true_literals = 0;
145 int true_literal = kint32min;
146 for (const int literal : ct->bool_xor().literals()) {
147 // TODO(user): More generally, if a variable appear in only bool xor
148 // constraints, we can simply eliminate it using linear algebra on Z/2Z.
149 // This should solve in polynomial time the parity-learning*.fzn problems
150 // for instance. This seems low priority, but it is also easy to do. Even
151 // better would be to have a dedicated propagator with all bool_xor
152 // constraints that do the necessary linear algebra.
153 if (context_->VariableIsUniqueAndRemovable(literal)) {
154 context_->UpdateRuleStats("TODO bool_xor: remove constraint");
155 }
156
157 if (context_->LiteralIsFalse(literal)) {
158 context_->UpdateRuleStats("bool_xor: remove false literal");
159 changed = true;
160 continue;
161 } else if (context_->LiteralIsTrue(literal)) {
162 true_literal = literal; // Keep if we need to put one back.
163 num_true_literals++;
164 continue;
165 }
166
167 ct->mutable_bool_xor()->set_literals(new_size++, literal);
168 }
169 if (new_size == 1) {
170 context_->UpdateRuleStats("TODO bool_xor: one active literal");
171 } else if (new_size == 2) {
172 context_->UpdateRuleStats("TODO bool_xor: two active literals");
173 }
174 if (num_true_literals % 2 == 1) {
175 CHECK_NE(true_literal, kint32min);
176 ct->mutable_bool_xor()->set_literals(new_size++, true_literal);
177 }
178 if (num_true_literals > 1) {
179 context_->UpdateRuleStats("bool_xor: remove even number of true literals");
180 changed = true;
181 }
182 ct->mutable_bool_xor()->mutable_literals()->Truncate(new_size);
183 return changed;
184}
185
186bool CpModelPresolver::PresolveBoolOr(ConstraintProto* ct) {
187 if (context_->ModelIsUnsat()) return false;
188
189 // Move the enforcement literal inside the clause if any. Note that we do not
190 // mark this as a change since the literal in the constraint are the same.
192 context_->UpdateRuleStats("bool_or: removed enforcement literal");
193 for (const int literal : ct->enforcement_literal()) {
194 ct->mutable_bool_or()->add_literals(NegatedRef(literal));
195 }
196 ct->clear_enforcement_literal();
197 }
198
199 // Inspects the literals and deal with fixed ones.
200 bool changed = false;
201 context_->tmp_literals.clear();
202 context_->tmp_literal_set.clear();
203 for (const int literal : ct->bool_or().literals()) {
204 if (context_->LiteralIsFalse(literal)) {
205 changed = true;
206 continue;
207 }
208 if (context_->LiteralIsTrue(literal)) {
209 context_->UpdateRuleStats("bool_or: always true");
210 return RemoveConstraint(ct);
211 }
212 // We can just set the variable to true in this case since it is not
213 // used in any other constraint (note that we artificially bump the
214 // objective var usage by 1).
215 if (context_->VariableIsUniqueAndRemovable(literal)) {
216 context_->UpdateRuleStats("bool_or: singleton");
217 if (!context_->SetLiteralToTrue(literal)) return true;
218 return RemoveConstraint(ct);
219 }
220 if (context_->tmp_literal_set.contains(NegatedRef(literal))) {
221 context_->UpdateRuleStats("bool_or: always true");
222 return RemoveConstraint(ct);
223 }
224
225 if (context_->tmp_literal_set.contains(literal)) {
226 changed = true;
227 } else {
228 context_->tmp_literal_set.insert(literal);
229 context_->tmp_literals.push_back(literal);
230 }
231 }
232 context_->tmp_literal_set.clear();
233
234 if (context_->tmp_literals.empty()) {
235 context_->UpdateRuleStats("bool_or: empty");
236 return context_->NotifyThatModelIsUnsat();
237 }
238 if (context_->tmp_literals.size() == 1) {
239 context_->UpdateRuleStats("bool_or: only one literal");
240 if (!context_->SetLiteralToTrue(context_->tmp_literals[0])) return true;
241 return RemoveConstraint(ct);
242 }
243 if (context_->tmp_literals.size() == 2) {
244 // For consistency, we move all "implication" into half-reified bool_and.
245 // TODO(user): merge by enforcement literal and detect implication cycles.
246 context_->UpdateRuleStats("bool_or: implications");
247 ct->add_enforcement_literal(NegatedRef(context_->tmp_literals[0]));
248 ct->mutable_bool_and()->add_literals(context_->tmp_literals[1]);
249 return changed;
250 }
251
252 if (changed) {
253 context_->UpdateRuleStats("bool_or: fixed literals");
254 ct->mutable_bool_or()->mutable_literals()->Clear();
255 for (const int lit : context_->tmp_literals) {
256 ct->mutable_bool_or()->add_literals(lit);
257 }
258 }
259 return changed;
260}
261
262ABSL_MUST_USE_RESULT bool CpModelPresolver::MarkConstraintAsFalse(
263 ConstraintProto* ct) {
265 // Change the constraint to a bool_or.
266 ct->mutable_bool_or()->clear_literals();
267 for (const int lit : ct->enforcement_literal()) {
268 ct->mutable_bool_or()->add_literals(NegatedRef(lit));
269 }
270 ct->clear_enforcement_literal();
271 PresolveBoolOr(ct);
272 return true;
273 } else {
274 return context_->NotifyThatModelIsUnsat();
275 }
276}
277
278bool CpModelPresolver::PresolveBoolAnd(ConstraintProto* ct) {
279 if (context_->ModelIsUnsat()) return false;
280
281 if (!HasEnforcementLiteral(*ct)) {
282 context_->UpdateRuleStats("bool_and: non-reified.");
283 for (const int literal : ct->bool_and().literals()) {
284 if (!context_->SetLiteralToTrue(literal)) return true;
285 }
286 return RemoveConstraint(ct);
287 }
288
289 bool changed = false;
290 context_->tmp_literals.clear();
291 for (const int literal : ct->bool_and().literals()) {
292 if (context_->LiteralIsFalse(literal)) {
293 context_->UpdateRuleStats("bool_and: always false");
294 return MarkConstraintAsFalse(ct);
295 }
296 if (context_->LiteralIsTrue(literal)) {
297 changed = true;
298 continue;
299 }
300 if (context_->VariableIsUniqueAndRemovable(literal)) {
301 changed = true;
302 if (!context_->SetLiteralToTrue(literal)) return true;
303 continue;
304 }
305 context_->tmp_literals.push_back(literal);
306 }
307
308 // Note that this is not the same behavior as a bool_or:
309 // - bool_or means "at least one", so it is false if empty.
310 // - bool_and means "all literals inside true", so it is true if empty.
311 if (context_->tmp_literals.empty()) return RemoveConstraint(ct);
312
313 if (changed) {
314 ct->mutable_bool_and()->mutable_literals()->Clear();
315 for (const int lit : context_->tmp_literals) {
316 ct->mutable_bool_and()->add_literals(lit);
317 }
318 context_->UpdateRuleStats("bool_and: fixed literals");
319 }
320
321 // If a variable can move freely in one direction except for this constraint,
322 // we can make it an equality.
323 //
324 // TODO(user): also consider literal on the other side of the =>.
325 if (ct->enforcement_literal().size() == 1 &&
326 ct->bool_and().literals().size() == 1) {
327 const int enforcement = ct->enforcement_literal(0);
328 if (context_->VariableWithCostIsUniqueAndRemovable(enforcement)) {
329 int var = PositiveRef(enforcement);
330 int64 obj_coeff = gtl::FindOrDie(context_->ObjectiveMap(), var);
331 if (!RefIsPositive(enforcement)) obj_coeff = -obj_coeff;
332
333 // The other case where the constraint is redundant is treated elsewhere.
334 if (obj_coeff < 0) {
335 context_->UpdateRuleStats("bool_and: dual equality.");
336 context_->StoreBooleanEqualityRelation(enforcement,
337 ct->bool_and().literals(0));
338 }
339 }
340 }
341
342 return changed;
343}
344
345bool CpModelPresolver::PresolveAtMostOrExactlyOne(ConstraintProto* ct) {
346 bool is_at_most_one = ct->constraint_case() == ConstraintProto::kAtMostOne;
347 const std::string name = is_at_most_one ? "at_most_one: " : "exactly_one: ";
348 auto* literals = is_at_most_one
349 ? ct->mutable_at_most_one()->mutable_literals()
350 : ct->mutable_exactly_one()->mutable_literals();
351
352 // Deal with duplicate variable reference.
353 context_->tmp_literal_set.clear();
354 for (const int literal : *literals) {
355 if (context_->tmp_literal_set.contains(literal)) {
356 if (!context_->SetLiteralToFalse(literal)) return true;
357 context_->UpdateRuleStats(absl::StrCat(name, "duplicate literals"));
358 }
359 if (context_->tmp_literal_set.contains(NegatedRef(literal))) {
360 int num_positive = 0;
361 int num_negative = 0;
362 for (const int other : *literals) {
363 if (PositiveRef(other) != PositiveRef(literal)) {
364 if (!context_->SetLiteralToFalse(other)) return true;
365 context_->UpdateRuleStats(absl::StrCat(name, "x and not(x)"));
366 } else {
367 if (other == literal) {
368 ++num_positive;
369 } else {
370 ++num_negative;
371 }
372 }
373 }
374
375 // This is tricky for the case where the at most one reduce to (lit,
376 // not(lit), not(lit)) for instance.
377 if (num_positive > 1 && !context_->SetLiteralToFalse(literal)) {
378 return true;
379 }
380 if (num_negative > 1 && !context_->SetLiteralToTrue(literal)) {
381 return true;
382 }
383 return RemoveConstraint(ct);
384 }
385 context_->tmp_literal_set.insert(literal);
386 }
387
388 // Remove fixed variables.
389 bool changed = false;
390 bool transform_to_at_most_one = false;
391 context_->tmp_literals.clear();
392 for (const int literal : *literals) {
393 if (context_->LiteralIsTrue(literal)) {
394 context_->UpdateRuleStats(absl::StrCat(name, "satisfied"));
395 for (const int other : *literals) {
396 if (other != literal) {
397 if (!context_->SetLiteralToFalse(other)) return true;
398 }
399 }
400 return RemoveConstraint(ct);
401 }
402
403 if (context_->LiteralIsFalse(literal)) {
404 changed = true;
405 continue;
406 }
407
408 // A singleton variable in an at most one can just be set to zero.
409 //
410 // In an exactly one, it can be left to the postsolve to decide, and the
411 // rest of the constraint can be transformed to an at most one.
412 bool is_removable = context_->VariableIsUniqueAndRemovable(literal);
413 if (is_at_most_one && !is_removable &&
415 const auto it = context_->ObjectiveMap().find(PositiveRef(literal));
416 CHECK(it != context_->ObjectiveMap().end());
417 const int64 coeff = it->second;
418
419 // Fixing it to zero need to go in the correct direction.
420 is_removable = (coeff > 0) == RefIsPositive(literal);
421 }
422
423 if (is_removable) {
424 if (is_at_most_one) {
425 context_->UpdateRuleStats("at_most_one: singleton");
426 if (!context_->SetLiteralToFalse(literal)) return false;
427 changed = true;
428 continue;
429 } else {
430 changed = true;
431 is_at_most_one = true;
432 transform_to_at_most_one = true;
433 *(context_->mapping_model->add_constraints()) = *ct;
434 context_->UpdateRuleStats("exactly_one: singleton");
435 continue;
436 }
437 }
438
439 context_->tmp_literals.push_back(literal);
440 }
441
442 if (transform_to_at_most_one) {
443 CHECK(changed);
444 ct->Clear();
445 literals = ct->mutable_at_most_one()->mutable_literals();
446 }
447 if (changed) {
448 literals->Clear();
449 for (const int lit : context_->tmp_literals) {
450 literals->Add(lit);
451 }
452 context_->UpdateRuleStats(absl::StrCat(name, "removed literals"));
453 }
454 return changed;
455}
456
457bool CpModelPresolver::PresolveAtMostOne(ConstraintProto* ct) {
458 if (context_->ModelIsUnsat()) return false;
460 const bool changed = PresolveAtMostOrExactlyOne(ct);
461 if (ct->constraint_case() != ConstraintProto::kAtMostOne) return changed;
462
463 // Size zero: ok.
464 const auto& literals = ct->at_most_one().literals();
465 if (literals.empty()) {
466 context_->UpdateRuleStats("at_most_one: empty or all false");
467 return RemoveConstraint(ct);
468 }
469
470 // Size one: always satisfied.
471 if (literals.size() == 1) {
472 context_->UpdateRuleStats("at_most_one: size one");
473 return RemoveConstraint(ct);
474 }
475
476 return changed;
477}
478
479bool CpModelPresolver::PresolveExactlyOne(ConstraintProto* ct) {
480 if (context_->ModelIsUnsat()) return false;
482 const bool changed = PresolveAtMostOrExactlyOne(ct);
483 if (ct->constraint_case() != ConstraintProto::kExactlyOne) return changed;
484
485 // Size zero: UNSAT.
486 const auto& literals = ct->exactly_one().literals();
487 if (literals.empty()) {
488 return context_->NotifyThatModelIsUnsat("exactly_one: empty or all false");
489 }
490
491 // Size one: fix variable.
492 if (literals.size() == 1) {
493 context_->UpdateRuleStats("exactly_one: size one");
494 if (!context_->SetLiteralToTrue(literals[0])) return false;
495 return RemoveConstraint(ct);
496 }
497
498 // Size two: Equivalence.
499 if (literals.size() == 2) {
500 context_->UpdateRuleStats("exactly_one: size two");
501 context_->StoreBooleanEqualityRelation(literals[0],
502 NegatedRef(literals[1]));
503 return RemoveConstraint(ct);
504 }
505
506 return changed;
507}
508
509bool CpModelPresolver::PresolveIntMax(ConstraintProto* ct) {
510 if (context_->ModelIsUnsat()) return false;
511 if (ct->int_max().vars().empty()) {
512 context_->UpdateRuleStats("int_max: no variables!");
513 return MarkConstraintAsFalse(ct);
514 }
515 const int target_ref = ct->int_max().target();
516
517 // Pass 1, compute the infered min of the target, and remove duplicates.
518 int64 infered_min = kint64min;
519 int64 infered_max = kint64min;
520 bool contains_target_ref = false;
521 std::set<int> used_ref;
522 int new_size = 0;
523 for (const int ref : ct->int_max().vars()) {
524 if (ref == target_ref) contains_target_ref = true;
525 if (gtl::ContainsKey(used_ref, ref)) continue;
526 if (gtl::ContainsKey(used_ref, NegatedRef(ref)) ||
527 ref == NegatedRef(target_ref)) {
528 infered_min = std::max(infered_min, int64{0});
529 }
530 used_ref.insert(ref);
531 ct->mutable_int_max()->set_vars(new_size++, ref);
532 infered_min = std::max(infered_min, context_->MinOf(ref));
533 infered_max = std::max(infered_max, context_->MaxOf(ref));
534 }
535 if (new_size < ct->int_max().vars_size()) {
536 context_->UpdateRuleStats("int_max: removed dup");
537 }
538 ct->mutable_int_max()->mutable_vars()->Truncate(new_size);
539 if (contains_target_ref) {
540 context_->UpdateRuleStats("int_max: x = max(x, ...)");
541 for (const int ref : ct->int_max().vars()) {
542 if (ref == target_ref) continue;
543 ConstraintProto* new_ct = context_->working_model->add_constraints();
544 *new_ct->mutable_enforcement_literal() = ct->enforcement_literal();
545 auto* arg = new_ct->mutable_linear();
546 arg->add_vars(target_ref);
547 arg->add_coeffs(1);
548 arg->add_vars(ref);
549 arg->add_coeffs(-1);
550 arg->add_domain(0);
551 arg->add_domain(kint64max);
552 }
553 return RemoveConstraint(ct);
554 }
555
556 // Compute the infered target_domain.
557 Domain infered_domain;
558 for (const int ref : ct->int_max().vars()) {
559 infered_domain = infered_domain.UnionWith(
560 context_->DomainOf(ref).IntersectionWith({infered_min, infered_max}));
561 }
562
563 // Update the target domain.
564 bool domain_reduced = false;
565 if (!HasEnforcementLiteral(*ct)) {
566 if (!context_->IntersectDomainWith(target_ref, infered_domain,
567 &domain_reduced)) {
568 return true;
569 }
570 }
571
572 // If the target is only used here and if
573 // infered_domain ∩ [kint64min, target_ub] ⊂ target_domain
574 // then the constraint is really max(...) <= target_ub and we can simplify it.
575 if (context_->VariableIsUniqueAndRemovable(target_ref)) {
576 const Domain& target_domain = context_->DomainOf(target_ref);
577 if (infered_domain.IntersectionWith(Domain(kint64min, target_domain.Max()))
578 .IsIncludedIn(target_domain)) {
579 if (infered_domain.Max() <= target_domain.Max()) {
580 // The constraint is always satisfiable.
581 context_->UpdateRuleStats("int_max: always true");
582 } else if (ct->enforcement_literal().empty()) {
583 // The constraint just restrict the upper bound of its variable.
584 for (const int ref : ct->int_max().vars()) {
585 context_->UpdateRuleStats("int_max: lower than constant");
586 if (!context_->IntersectDomainWith(
587 ref, Domain(kint64min, target_domain.Max()))) {
588 return false;
589 }
590 }
591 } else {
592 // We simply transform this into n reified constraints
593 // enforcement => [var_i <= target_domain.Max()].
594 context_->UpdateRuleStats("int_max: reified lower than constant");
595 for (const int ref : ct->int_max().vars()) {
596 ConstraintProto* new_ct = context_->working_model->add_constraints();
597 *(new_ct->mutable_enforcement_literal()) = ct->enforcement_literal();
598 ct->mutable_linear()->add_vars(ref);
599 ct->mutable_linear()->add_coeffs(1);
600 ct->mutable_linear()->add_domain(kint64min);
601 ct->mutable_linear()->add_domain(target_domain.Max());
602 }
603 }
604
605 // In all cases we delete the original constraint.
606 context_->MarkVariableAsRemoved(target_ref);
607 *(context_->mapping_model->add_constraints()) = *ct;
608 return RemoveConstraint(ct);
609 }
610 }
611
612 // Pass 2, update the argument domains. Filter them eventually.
613 new_size = 0;
614 const int size = ct->int_max().vars_size();
615 const int64 target_max = context_->MaxOf(target_ref);
616 for (const int ref : ct->int_max().vars()) {
617 if (!HasEnforcementLiteral(*ct)) {
618 if (!context_->IntersectDomainWith(ref, Domain(kint64min, target_max),
619 &domain_reduced)) {
620 return true;
621 }
622 }
623 if (context_->MaxOf(ref) >= infered_min) {
624 ct->mutable_int_max()->set_vars(new_size++, ref);
625 }
626 }
627 if (domain_reduced) {
628 context_->UpdateRuleStats("int_max: reduced domains");
629 }
630
631 bool modified = false;
632 if (new_size < size) {
633 context_->UpdateRuleStats("int_max: removed variables");
634 ct->mutable_int_max()->mutable_vars()->Truncate(new_size);
635 modified = true;
636 }
637
638 if (new_size == 0) {
639 context_->UpdateRuleStats("int_max: no variables!");
640 return MarkConstraintAsFalse(ct);
641 }
642 if (new_size == 1) {
643 // Convert to an equality. Note that we create a new constraint otherwise it
644 // might not be processed again.
645 context_->UpdateRuleStats("int_max: converted to equality");
646 ConstraintProto* new_ct = context_->working_model->add_constraints();
647 *new_ct = *ct; // copy name and potential reification.
648 auto* arg = new_ct->mutable_linear();
649 arg->add_vars(target_ref);
650 arg->add_coeffs(1);
651 arg->add_vars(ct->int_max().vars(0));
652 arg->add_coeffs(-1);
653 arg->add_domain(0);
654 arg->add_domain(0);
656 return RemoveConstraint(ct);
657 }
658 return modified;
659}
660
661bool CpModelPresolver::PresolveLinMin(ConstraintProto* ct) {
662 if (context_->ModelIsUnsat()) return false;
663 // Convert to lin_max and presolve lin_max.
664 const auto copy = ct->lin_min();
665 SetToNegatedLinearExpression(copy.target(),
666 ct->mutable_lin_max()->mutable_target());
667 for (const LinearExpressionProto& expr : copy.exprs()) {
668 LinearExpressionProto* const new_expr = ct->mutable_lin_max()->add_exprs();
669 SetToNegatedLinearExpression(expr, new_expr);
670 }
671 return PresolveLinMax(ct);
672}
673
674bool CpModelPresolver::PresolveLinMax(ConstraintProto* ct) {
675 if (context_->ModelIsUnsat()) return false;
676 if (ct->lin_max().exprs().empty()) {
677 context_->UpdateRuleStats("lin_max: no exprs");
678 return MarkConstraintAsFalse(ct);
679 }
680
681 // Canonicalize all involved expression.
682 //
683 // TODO(user): If we start to have many constraints like this, we should
684 // use reflexion (see cp_model_util) to do that generically.
685 bool changed = CanonicalizeLinearExpression(
686 *ct, ct->mutable_lin_max()->mutable_target());
687 for (LinearExpressionProto& exp : *(ct->mutable_lin_max()->mutable_exprs())) {
688 changed |= CanonicalizeLinearExpression(*ct, &exp);
689 }
690
691 // TODO(user): Remove duplicate expressions. This might be expensive.
692
693 // Pass 1, Compute the infered min of the target.
694 int64 infered_min = context_->MinOf(ct->lin_max().target());
695 for (const LinearExpressionProto& expr : ct->lin_max().exprs()) {
696 // TODO(user): Check if the expressions contain target.
697
698 // TODO(user): Check if the negated expression is already present and
699 // reduce inferred domain if so.
700
701 infered_min = std::max(infered_min, context_->MinOf(expr));
702 }
703
704 // Pass 2, Filter the expressions which are smaller than inferred min.
705 int new_size = 0;
706 for (int i = 0; i < ct->lin_max().exprs_size(); ++i) {
707 const LinearExpressionProto& expr = ct->lin_max().exprs(i);
708 if (context_->MaxOf(expr) >= infered_min) {
709 *ct->mutable_lin_max()->mutable_exprs(new_size) = expr;
710 new_size++;
711 }
712 }
713
714 if (new_size < ct->lin_max().exprs_size()) {
715 context_->UpdateRuleStats("lin_max: Removed exprs");
716 ct->mutable_lin_max()->mutable_exprs()->DeleteSubrange(
717 new_size, ct->lin_max().exprs_size() - new_size);
718 return true;
719 }
720
721 return changed;
722}
723
724bool CpModelPresolver::PresolveIntAbs(ConstraintProto* ct) {
725 CHECK_EQ(ct->enforcement_literal_size(), 0);
726 if (context_->ModelIsUnsat()) return false;
727 const int target_ref = ct->int_max().target();
728 const int var = PositiveRef(ct->int_max().vars(0));
729
730 // Propagate from the variable domain to the target variable.
731 const Domain var_domain = context_->DomainOf(var);
732 const Domain new_target_domain = var_domain.UnionWith(var_domain.Negation())
734 if (!context_->DomainOf(target_ref).IsIncludedIn(new_target_domain)) {
735 if (!context_->IntersectDomainWith(target_ref, new_target_domain)) {
736 return true;
737 }
738 context_->UpdateRuleStats("int_abs: propagate domain x to abs(x)");
739 }
740
741 // Propagate from target domain to variable.
742 const Domain target_domain = context_->DomainOf(target_ref);
743 const Domain new_var_domain =
744 target_domain.UnionWith(target_domain.Negation());
745 if (!context_->DomainOf(var).IsIncludedIn(new_var_domain)) {
746 if (!context_->IntersectDomainWith(var, new_var_domain)) {
747 return true;
748 }
749 context_->UpdateRuleStats("int_abs: propagate domain abs(x) to x");
750 }
751
752 if (context_->MinOf(var) >= 0 && !context_->IsFixed(var)) {
753 context_->UpdateRuleStats("int_abs: converted to equality");
754 ConstraintProto* new_ct = context_->working_model->add_constraints();
755 new_ct->set_name(ct->name());
756 auto* arg = new_ct->mutable_linear();
757 arg->add_vars(target_ref);
758 arg->add_coeffs(1);
759 arg->add_vars(var);
760 arg->add_coeffs(-1);
761 arg->add_domain(0);
762 arg->add_domain(0);
764 return RemoveConstraint(ct);
765 }
766
767 if (context_->MaxOf(var) <= 0 && !context_->IsFixed(var)) {
768 context_->UpdateRuleStats("int_abs: converted to equality");
769 ConstraintProto* new_ct = context_->working_model->add_constraints();
770 new_ct->set_name(ct->name());
771 auto* arg = new_ct->mutable_linear();
772 arg->add_vars(target_ref);
773 arg->add_coeffs(1);
774 arg->add_vars(var);
775 arg->add_coeffs(1);
776 arg->add_domain(0);
777 arg->add_domain(0);
779 return RemoveConstraint(ct);
780 }
781
782 // Remove the abs constraint if the target is removable or fixed, as domains
783 // have been propagated.
784 if (context_->VariableIsUniqueAndRemovable(target_ref) ||
785 context_->IsFixed(target_ref)) {
786 if (!context_->IsFixed(target_ref)) {
787 context_->MarkVariableAsRemoved(target_ref);
788 *context_->mapping_model->add_constraints() = *ct;
789 }
790 context_->UpdateRuleStats("int_abs: remove constraint");
791 return RemoveConstraint(ct);
792 }
793
794 if (context_->StoreAbsRelation(target_ref, var)) {
795 context_->UpdateRuleStats("int_abs: store abs(x) == y");
796 }
797
798 return false;
799}
800
801bool CpModelPresolver::PresolveIntMin(ConstraintProto* ct) {
802 if (context_->ModelIsUnsat()) return false;
803
804 const auto copy = ct->int_min();
805 ct->mutable_int_max()->set_target(NegatedRef(copy.target()));
806 for (const int ref : copy.vars()) {
807 ct->mutable_int_max()->add_vars(NegatedRef(ref));
808 }
809 return PresolveIntMax(ct);
810}
811
812bool CpModelPresolver::PresolveIntProd(ConstraintProto* ct) {
813 if (context_->ModelIsUnsat()) return false;
814 if (HasEnforcementLiteral(*ct)) return false;
815
816 if (ct->int_prod().vars().empty()) {
817 if (!context_->IntersectDomainWith(ct->int_prod().target(), Domain(1))) {
818 return false;
819 }
820 context_->UpdateRuleStats("int_prod: empty_product");
821 return RemoveConstraint(ct);
822 }
823 bool changed = false;
824
825 // Replace any affine relation without offset.
826 // TODO(user): Also remove constant rhs variables.
827 int64 constant = 1;
828 for (int i = 0; i < ct->int_prod().vars().size(); ++i) {
829 const int ref = ct->int_prod().vars(i);
830 const AffineRelation::Relation& r = context_->GetAffineRelation(ref);
831 if (r.representative != ref && r.offset == 0) {
832 changed = true;
833 ct->mutable_int_prod()->set_vars(i, r.representative);
834 constant *= r.coeff;
835 }
836 }
837
838 // TODO(user): Probably better to add a fixed variable to the product
839 // instead in this case. But we do need to support product with more than
840 // two variables properly for that.
841 //
842 // TODO(user): We might do that too early since the other presolve step below
843 // might simplify the constraint in such a way that there is no need to create
844 // a new variable!
845 if (constant != 1) {
846 context_->UpdateRuleStats("int_prod: extracted product by constant.");
847
848 const int old_target = ct->int_prod().target();
849 const int new_target = context_->working_model->variables_size();
850
851 IntegerVariableProto* var_proto = context_->working_model->add_variables();
853 context_->DomainOf(old_target).InverseMultiplicationBy(constant),
854 var_proto);
855 context_->InitializeNewDomains();
856 if (context_->ModelIsUnsat()) return false;
857
858 ct->mutable_int_prod()->set_target(new_target);
859 if (context_->IsFixed(new_target)) {
860 // We need to fix old_target too.
861 if (!context_->IntersectDomainWith(
862 old_target,
863 context_->DomainOf(new_target).MultiplicationBy(constant))) {
864 return false;
865 }
866 } else {
867 if (!context_->StoreAffineRelation(old_target, new_target, constant, 0)) {
868 // We cannot store the affine relation because the old target seems
869 // to already be in affine relation with another variable. This is rare
870 // and we need to add a new constraint in that case.
871 ConstraintProto* new_ct = context_->working_model->add_constraints();
872 LinearConstraintProto* lin = new_ct->mutable_linear();
873 lin->add_vars(old_target);
874 lin->add_coeffs(1);
875 lin->add_vars(new_target);
876 lin->add_coeffs(-constant);
877 lin->add_domain(0);
878 lin->add_domain(0);
880 }
881 }
882 }
883
884 // Restrict the target domain if possible.
885 Domain implied(1);
886 for (const int ref : ct->int_prod().vars()) {
887 implied = implied.ContinuousMultiplicationBy(context_->DomainOf(ref));
888 }
889 bool modified = false;
890 if (!context_->IntersectDomainWith(ct->int_prod().target(), implied,
891 &modified)) {
892 return false;
893 }
894 if (modified) {
895 context_->UpdateRuleStats("int_prod: reduced target domain.");
896 }
897
898 if (ct->int_prod().vars_size() == 2) {
899 int a = ct->int_prod().vars(0);
900 int b = ct->int_prod().vars(1);
901 const int product = ct->int_prod().target();
902
903 if (context_->IsFixed(b)) std::swap(a, b);
904 if (context_->IsFixed(a)) {
905 const int64 value_a = context_->MinOf(a);
906 if (value_a == 0) { // Fix target to 0.
907 if (!context_->IntersectDomainWith(product, Domain(0, 0))) {
908 return false;
909 }
910 context_->UpdateRuleStats("int_prod: fix target to zero.");
911 return RemoveConstraint(ct);
912 } else if (b != product) {
913 ConstraintProto* const lin = context_->working_model->add_constraints();
914 lin->mutable_linear()->add_vars(b);
915 lin->mutable_linear()->add_coeffs(value_a);
916 lin->mutable_linear()->add_vars(product);
917 lin->mutable_linear()->add_coeffs(-1);
918 lin->mutable_linear()->add_domain(0);
919 lin->mutable_linear()->add_domain(0);
921 context_->UpdateRuleStats("int_prod: linearize product by constant.");
922 return RemoveConstraint(ct);
923 } else if (value_a != 1) {
924 if (!context_->IntersectDomainWith(product, Domain(0, 0))) {
925 return false;
926 }
927 context_->UpdateRuleStats("int_prod: fix variable to zero.");
928 return RemoveConstraint(ct);
929 } else {
930 context_->UpdateRuleStats("int_prod: remove identity.");
931 return RemoveConstraint(ct);
932 }
933 } else if (a == b && a == product) { // x = x * x, only true for {0, 1}.
934 if (!context_->IntersectDomainWith(product, Domain(0, 1))) {
935 return false;
936 }
937 context_->UpdateRuleStats("int_prod: fix variable to zero or one.");
938 return RemoveConstraint(ct);
939 }
940 }
941
942 // For now, we only presolve the case where all variables are Booleans.
943 const int target_ref = ct->int_prod().target();
944 if (!RefIsPositive(target_ref)) return changed;
945 for (const int var : ct->int_prod().vars()) {
946 if (!RefIsPositive(var)) return changed;
947 if (context_->MinOf(var) < 0) return changed;
948 if (context_->MaxOf(var) > 1) return changed;
949 }
950
951 // This is a bool constraint!
952 if (!context_->IntersectDomainWith(target_ref, Domain(0, 1))) {
953 return false;
954 }
955 context_->UpdateRuleStats("int_prod: all Boolean.");
956 {
957 ConstraintProto* new_ct = context_->working_model->add_constraints();
958 new_ct->add_enforcement_literal(target_ref);
959 auto* arg = new_ct->mutable_bool_and();
960 for (const int var : ct->int_prod().vars()) {
961 arg->add_literals(var);
962 }
963 }
964 {
965 ConstraintProto* new_ct = context_->working_model->add_constraints();
966 auto* arg = new_ct->mutable_bool_or();
967 arg->add_literals(target_ref);
968 for (const int var : ct->int_prod().vars()) {
969 arg->add_literals(NegatedRef(var));
970 }
971 }
973 return RemoveConstraint(ct);
974}
975
976bool CpModelPresolver::PresolveIntDiv(ConstraintProto* ct) {
977 if (context_->ModelIsUnsat()) return false;
978
979 // For now, we only presolve the case where the divisor is constant.
980 const int target = ct->int_div().target();
981 const int ref_x = ct->int_div().vars(0);
982 const int ref_div = ct->int_div().vars(1);
983 if (!RefIsPositive(target) || !RefIsPositive(ref_x) ||
984 !RefIsPositive(ref_div) || context_->DomainIsEmpty(ref_div) ||
985 !context_->IsFixed(ref_div)) {
986 return false;
987 }
988
989 const int64 divisor = context_->MinOf(ref_div);
990 if (divisor == 1) {
991 LinearConstraintProto* const lin =
992 context_->working_model->add_constraints()->mutable_linear();
993 lin->add_vars(ref_x);
994 lin->add_coeffs(1);
995 lin->add_vars(target);
996 lin->add_coeffs(-1);
997 lin->add_domain(0);
998 lin->add_domain(0);
1000 context_->UpdateRuleStats("int_div: rewrite to equality");
1001 return RemoveConstraint(ct);
1002 }
1003 bool domain_modified = false;
1004 if (context_->IntersectDomainWith(
1005 target, context_->DomainOf(ref_x).DivisionBy(divisor),
1006 &domain_modified)) {
1007 if (domain_modified) {
1008 context_->UpdateRuleStats(
1009 "int_div: updated domain of target in target = X / cte");
1010 }
1011 } else {
1012 // Model is unsat.
1013 return false;
1014 }
1015
1016 // Linearize if everything is positive.
1017 // TODO(user,user): Deal with other cases where there is no change of
1018 // sign.We can also deal with target = cte, div variable.
1019
1020 if (context_->MinOf(target) >= 0 && context_->MinOf(ref_x) >= 0 &&
1021 divisor > 1) {
1022 LinearConstraintProto* const lin =
1023 context_->working_model->add_constraints()->mutable_linear();
1024 lin->add_vars(ref_x);
1025 lin->add_coeffs(1);
1026 lin->add_vars(target);
1027 lin->add_coeffs(-divisor);
1028 lin->add_domain(0);
1029 lin->add_domain(divisor - 1);
1031 context_->UpdateRuleStats(
1032 "int_div: linearize positive division with a constant divisor");
1033 return RemoveConstraint(ct);
1034 }
1035
1036 // TODO(user): reduce the domain of X by introducing an
1037 // InverseDivisionOfSortedDisjointIntervals().
1038 return false;
1039}
1040
1041bool CpModelPresolver::ExploitEquivalenceRelations(int c, ConstraintProto* ct) {
1042 bool changed = false;
1043
1044 // Optim: Special case for the linear constraint. We just remap the
1045 // enforcement literals, the normal variables will be replaced by their
1046 // representative in CanonicalizeLinear().
1047 if (ct->constraint_case() == ConstraintProto::ConstraintCase::kLinear) {
1048 for (int& ref : *ct->mutable_enforcement_literal()) {
1049 const int rep = this->context_->GetLiteralRepresentative(ref);
1050 if (rep != ref) {
1051 changed = true;
1052 ref = rep;
1053 }
1054 }
1055 return changed;
1056 }
1057
1058 // Optim: This extra loop is a lot faster than reparsing the variable from the
1059 // proto when there is nothing to do, which is quite often.
1060 bool work_to_do = false;
1061 for (const int var : context_->ConstraintToVars(c)) {
1062 const AffineRelation::Relation r = context_->GetAffineRelation(var);
1063 if (r.representative != var) {
1064 work_to_do = true;
1065 break;
1066 }
1067 }
1068 if (!work_to_do) return false;
1069
1070 // Remap equal and negated variables to their representative.
1072 [&changed, this](int* ref) {
1073 const int rep = context_->GetVariableRepresentative(*ref);
1074 if (rep != *ref) {
1075 changed = true;
1076 *ref = rep;
1077 }
1078 },
1079 ct);
1080
1081 // Remap literal and negated literal to their representative.
1083 [&changed, this](int* ref) {
1084 const int rep = this->context_->GetLiteralRepresentative(*ref);
1085 if (rep != *ref) {
1086 changed = true;
1087 *ref = rep;
1088 }
1089 },
1090 ct);
1091 return changed;
1092}
1093
1094void CpModelPresolver::DivideLinearByGcd(ConstraintProto* ct) {
1095 if (context_->ModelIsUnsat()) return;
1096
1097 // Compute the GCD of all coefficients.
1098 int64 gcd = 0;
1099 const int num_vars = ct->linear().vars().size();
1100 for (int i = 0; i < num_vars; ++i) {
1101 const int64 magnitude = std::abs(ct->linear().coeffs(i));
1102 gcd = MathUtil::GCD64(gcd, magnitude);
1103 if (gcd == 1) break;
1104 }
1105 if (gcd > 1) {
1106 context_->UpdateRuleStats("linear: divide by GCD");
1107 for (int i = 0; i < num_vars; ++i) {
1108 ct->mutable_linear()->set_coeffs(i, ct->linear().coeffs(i) / gcd);
1109 }
1110 const Domain rhs = ReadDomainFromProto(ct->linear());
1111 FillDomainInProto(rhs.InverseMultiplicationBy(gcd), ct->mutable_linear());
1112 if (ct->linear().domain_size() == 0) {
1113 return (void)MarkConstraintAsFalse(ct);
1114 }
1115 }
1116}
1117
1118void CpModelPresolver::PresolveLinearEqualityModuloTwo(ConstraintProto* ct) {
1119 if (!ct->enforcement_literal().empty()) return;
1120 if (ct->linear().domain().size() != 2) return;
1121 if (ct->linear().domain(0) != ct->linear().domain(1)) return;
1122 if (context_->ModelIsUnsat()) return;
1123
1124 // Any equality must be true modulo n.
1125 // The case modulo 2 is interesting if the non-zero terms are Booleans.
1126 std::vector<int> literals;
1127 for (int i = 0; i < ct->linear().vars().size(); ++i) {
1128 const int64 coeff = ct->linear().coeffs(i);
1129 const int ref = ct->linear().vars(i);
1130 if (coeff % 2 == 0) continue;
1131 if (!context_->CanBeUsedAsLiteral(ref)) return;
1132 literals.push_back(PositiveRef(ref));
1133 if (literals.size() > 2) return;
1134 }
1135 if (literals.size() == 1) {
1136 const int64 rhs = std::abs(ct->linear().domain(0));
1137 context_->UpdateRuleStats("linear: only one odd Boolean in equality");
1138 if (!context_->IntersectDomainWith(literals[0], Domain(rhs % 2))) return;
1139 } else if (literals.size() == 2) {
1140 const int64 rhs = std::abs(ct->linear().domain(0));
1141 context_->UpdateRuleStats("linear: only two odd Booleans in equality");
1142 if (rhs % 2) {
1143 context_->StoreBooleanEqualityRelation(literals[0],
1144 NegatedRef(literals[1]));
1145 } else {
1146 context_->StoreBooleanEqualityRelation(literals[0], literals[1]);
1147 }
1148 }
1149}
1150
1151template <typename ProtoWithVarsAndCoeffs>
1152bool CpModelPresolver::CanonicalizeLinearExpressionInternal(
1153 const ConstraintProto& ct, ProtoWithVarsAndCoeffs* proto, int64* offset) {
1154 // First regroup the terms on the same variables and sum the fixed ones.
1155 //
1156 // TODO(user): Add a quick pass to skip most of the work below if the
1157 // constraint is already in canonical form?
1158 tmp_terms_.clear();
1159 int64 sum_of_fixed_terms = 0;
1160 bool remapped = false;
1161 const int old_size = proto->vars().size();
1162 DCHECK_EQ(old_size, proto->coeffs().size());
1163 for (int i = 0; i < old_size; ++i) {
1164 const int ref = proto->vars(i);
1165 const int var = PositiveRef(ref);
1166 const int64 coeff =
1167 RefIsPositive(ref) ? proto->coeffs(i) : -proto->coeffs(i);
1168 if (coeff == 0) continue;
1169
1170 if (context_->IsFixed(var)) {
1171 sum_of_fixed_terms += coeff * context_->MinOf(var);
1172 continue;
1173 }
1174
1175 // TODO(user): Avoid the quadratic loop for the corner case of many
1176 // enforcement literal (this should be pretty rare though).
1177 bool removed = false;
1178 for (const int enf : ct.enforcement_literal()) {
1179 if (var == PositiveRef(enf)) {
1180 if (RefIsPositive(enf)) {
1181 // If the constraint is enforced, we can assume the variable is at 1.
1182 sum_of_fixed_terms += coeff;
1183 } else {
1184 // We can assume the variable is at zero.
1185 }
1186 removed = true;
1187 break;
1188 }
1189 }
1190 if (removed) {
1191 context_->UpdateRuleStats("linear: enforcement literal in expression");
1192 continue;
1193 }
1194
1195 const AffineRelation::Relation r = context_->GetAffineRelation(var);
1196 if (r.representative != var) {
1197 remapped = true;
1198 sum_of_fixed_terms += coeff * r.offset;
1199 }
1200 tmp_terms_.push_back({r.representative, coeff * r.coeff});
1201 }
1202 proto->clear_vars();
1203 proto->clear_coeffs();
1204 std::sort(tmp_terms_.begin(), tmp_terms_.end());
1205 int current_var = 0;
1206 int64 current_coeff = 0;
1207 for (const auto entry : tmp_terms_) {
1208 CHECK(RefIsPositive(entry.first));
1209 if (entry.first == current_var) {
1210 current_coeff += entry.second;
1211 } else {
1212 if (current_coeff != 0) {
1213 proto->add_vars(current_var);
1214 proto->add_coeffs(current_coeff);
1215 }
1216 current_var = entry.first;
1217 current_coeff = entry.second;
1218 }
1219 }
1220 if (current_coeff != 0) {
1221 proto->add_vars(current_var);
1222 proto->add_coeffs(current_coeff);
1223 }
1224 if (remapped) {
1225 context_->UpdateRuleStats("linear: remapped using affine relations");
1226 }
1227 if (proto->vars().size() < old_size) {
1228 context_->UpdateRuleStats("linear: fixed or dup variables");
1229 }
1230 *offset = sum_of_fixed_terms;
1231 return remapped || proto->vars().size() < old_size;
1232}
1233
1234bool CpModelPresolver::CanonicalizeLinearExpression(
1235 const ConstraintProto& ct, LinearExpressionProto* exp) {
1236 int64 offset = 0;
1237 const bool result = CanonicalizeLinearExpressionInternal(ct, exp, &offset);
1238 exp->set_offset(exp->offset() + offset);
1239 return result;
1240}
1241
1242bool CpModelPresolver::CanonicalizeLinear(ConstraintProto* ct) {
1243 if (ct->constraint_case() != ConstraintProto::ConstraintCase::kLinear ||
1244 context_->ModelIsUnsat()) {
1245 return false;
1246 }
1247 int64 offset = 0;
1248 const bool result =
1249 CanonicalizeLinearExpressionInternal(*ct, ct->mutable_linear(), &offset);
1250 if (offset != 0) {
1252 ReadDomainFromProto(ct->linear()).AdditionWith(Domain(-offset)),
1253 ct->mutable_linear());
1254 }
1255 DivideLinearByGcd(ct);
1256 return result;
1257}
1258
1259bool CpModelPresolver::RemoveSingletonInLinear(ConstraintProto* ct) {
1260 if (ct->constraint_case() != ConstraintProto::ConstraintCase::kLinear ||
1261 context_->ModelIsUnsat()) {
1262 return false;
1263 }
1264
1265 std::set<int> index_to_erase;
1266 const int num_vars = ct->linear().vars().size();
1267 Domain rhs = ReadDomainFromProto(ct->linear());
1268
1269 // First pass. Process singleton column that are not in the objective. Note
1270 // that for postsolve, it is important that we process them in the same order
1271 // in which they will be removed.
1272 for (int i = 0; i < num_vars; ++i) {
1273 const int var = ct->linear().vars(i);
1274 const int64 coeff = ct->linear().coeffs(i);
1276 if (context_->VariableIsUniqueAndRemovable(var)) {
1277 bool exact;
1278 const auto term_domain =
1279 context_->DomainOf(var).MultiplicationBy(-coeff, &exact);
1280 if (!exact) continue;
1281
1282 // We do not do that if the domain of rhs becomes too complex.
1283 const Domain new_rhs = rhs.AdditionWith(term_domain);
1284 if (new_rhs.NumIntervals() > 100) continue;
1285
1286 // Note that we can't do that if we loose information in the
1287 // multiplication above because the new domain might not be as strict
1288 // as the initial constraint otherwise. TODO(user): because of the
1289 // addition, it might be possible to cover more cases though.
1290 context_->UpdateRuleStats("linear: singleton column");
1291 index_to_erase.insert(i);
1292 rhs = new_rhs;
1293 continue;
1294 }
1295 }
1296
1297 // If we didn't find any, look for the one appearing in the objective.
1298 if (index_to_erase.empty()) {
1299 // Note that we only do that if we have a non-reified equality.
1300 if (context_->params().presolve_substitution_level() <= 0) return false;
1301 if (!ct->enforcement_literal().empty()) return false;
1302
1303 // If it is possible to do so, note that we can transform constraint into
1304 // equalities in PropagateDomainsInLinear().
1305 if (rhs.Min() != rhs.Max()) return false;
1306
1307 for (int i = 0; i < num_vars; ++i) {
1308 const int var = ct->linear().vars(i);
1309 const int64 coeff = ct->linear().coeffs(i);
1311
1312 // If the variable appear only in the objective and we have an equality,
1313 // we can transfer the cost to the rest of the linear expression, and
1314 // remove that variable.
1315 //
1316 // Note that is similar to the substitution code in PresolveLinear() but
1317 // it doesn't require the variable to be implied free since we do not
1318 // remove the constraints afterwards, just the variable.
1319 if (!context_->VariableWithCostIsUniqueAndRemovable(var)) continue;
1320 DCHECK(context_->ObjectiveMap().contains(var));
1321
1322 // We only support substitution that does not require to multiply the
1323 // objective by some factor.
1324 //
1325 // TODO(user): If the objective is a single variable, we can actually
1326 // "absorb" any factor into the objective scaling.
1327 const int64 objective_coeff =
1328 gtl::FindOrDie(context_->ObjectiveMap(), var);
1329 CHECK_NE(coeff, 0);
1330 if (objective_coeff % coeff != 0) continue;
1331
1332 // We do not do that if the domain of rhs becomes too complex.
1333 bool exact;
1334 const auto term_domain =
1335 context_->DomainOf(var).MultiplicationBy(-coeff, &exact);
1336 if (!exact) continue;
1337 const Domain new_rhs = rhs.AdditionWith(term_domain);
1338 if (new_rhs.NumIntervals() > 100) continue;
1339
1340 // Special case: If the objective was a single variable, we can transfer
1341 // the domain of var to the objective, and just completely remove this
1342 // equality constraint like it is done in ExpandObjective().
1343 if (context_->ObjectiveMap().size() == 1) {
1344 if (!context_->IntersectDomainWith(
1346 objective_coeff))) {
1347 return true;
1348 }
1349
1350 // The intersection above might fix var, in which case, we just abort.
1351 if (context_->IsFixed(var)) continue;
1352
1353 // This makes sure the domain of var is propagated back to the
1354 // objective.
1355 if (!context_->CanonicalizeObjective()) {
1356 return context_->NotifyThatModelIsUnsat();
1357 }
1358
1359 // Normally, CanonicalizeObjective() shouldn't remove var because
1360 // we work on a linear constraint that has been canonicalized. We keep
1361 // the test here in case this ever happen so we are notified.
1362 if (!context_->ObjectiveMap().contains(var)) {
1363 LOG(WARNING) << "This was not supposed to happen and the presolve "
1364 "could be improved.";
1365 continue;
1366 }
1367 context_->UpdateRuleStats("linear: singleton column define objective.");
1368 context_->SubstituteVariableInObjective(var, coeff, *ct);
1369 context_->MarkVariableAsRemoved(var);
1370 *(context_->mapping_model->add_constraints()) = *ct;
1371 return RemoveConstraint(ct);
1372 }
1373
1374 // Update the objective and remove the variable from its equality
1375 // constraint by expanding its rhs. This might fail if the new linear
1376 // objective expression can lead to overflow.
1377 if (!context_->SubstituteVariableInObjective(var, coeff, *ct)) continue;
1378
1379 context_->UpdateRuleStats(
1380 "linear: singleton column in equality and in objective.");
1381 rhs = new_rhs;
1382 index_to_erase.insert(i);
1383 break;
1384 }
1385 }
1386 if (index_to_erase.empty()) return false;
1387
1388 // TODO(user): we could add the constraint to mapping_model only once
1389 // instead of adding a reduced version of it each time a new singleton
1390 // variable appear in the same constraint later. That would work but would
1391 // also force the postsolve to take search decisions...
1392 *(context_->mapping_model->add_constraints()) = *ct;
1393
1394 int new_size = 0;
1395 for (int i = 0; i < num_vars; ++i) {
1396 if (index_to_erase.count(i)) {
1397 context_->MarkVariableAsRemoved(ct->linear().vars(i));
1398 continue;
1399 }
1400 ct->mutable_linear()->set_coeffs(new_size, ct->linear().coeffs(i));
1401 ct->mutable_linear()->set_vars(new_size, ct->linear().vars(i));
1402 ++new_size;
1403 }
1404 ct->mutable_linear()->mutable_vars()->Truncate(new_size);
1405 ct->mutable_linear()->mutable_coeffs()->Truncate(new_size);
1406 FillDomainInProto(rhs, ct->mutable_linear());
1407 DivideLinearByGcd(ct);
1408 return true;
1409}
1410
1411bool CpModelPresolver::PresolveSmallLinear(ConstraintProto* ct) {
1412 if (ct->constraint_case() != ConstraintProto::ConstraintCase::kLinear ||
1413 context_->ModelIsUnsat()) {
1414 return false;
1415 }
1416
1417 // Empty constraint?
1418 if (ct->linear().vars().empty()) {
1419 context_->UpdateRuleStats("linear: empty");
1420 const Domain rhs = ReadDomainFromProto(ct->linear());
1421 if (rhs.Contains(0)) {
1422 return RemoveConstraint(ct);
1423 } else {
1424 return MarkConstraintAsFalse(ct);
1425 }
1426 }
1427
1428 // If the constraint is literal => x in domain and x = abs(abs_arg), we can
1429 // replace x by abs_arg and hopefully remove the variable x later.
1430 int abs_arg;
1431 if (ct->linear().vars_size() == 1 && ct->enforcement_literal_size() > 0 &&
1432 ct->linear().coeffs(0) == 1 &&
1433 context_->GetAbsRelation(ct->linear().vars(0), &abs_arg)) {
1434 // TODO(user): Deal with coeff = -1, here or during canonicalization.
1435 context_->UpdateRuleStats("linear: remove abs from abs(x) in domain");
1436 const Domain implied_abs_target_domain =
1437 ReadDomainFromProto(ct->linear())
1439 .IntersectionWith(context_->DomainOf(ct->linear().vars(0)));
1440
1441 if (implied_abs_target_domain.IsEmpty()) {
1442 return MarkConstraintAsFalse(ct);
1443 }
1444
1445 const Domain new_abs_var_domain =
1446 implied_abs_target_domain
1447 .UnionWith(implied_abs_target_domain.Negation())
1448 .IntersectionWith(context_->DomainOf(abs_arg));
1449
1450 if (new_abs_var_domain.IsEmpty()) {
1451 return MarkConstraintAsFalse(ct);
1452 }
1453
1454 ConstraintProto* new_ct = context_->working_model->add_constraints();
1455 new_ct->set_name(ct->name());
1456 for (const int literal : ct->enforcement_literal()) {
1457 new_ct->add_enforcement_literal(literal);
1458 }
1459 auto* arg = new_ct->mutable_linear();
1460 arg->add_vars(abs_arg);
1461 arg->add_coeffs(1);
1462 FillDomainInProto(new_abs_var_domain, new_ct->mutable_linear());
1464 return RemoveConstraint(ct);
1465 }
1466
1467 if (HasEnforcementLiteral(*ct)) {
1468 if (ct->enforcement_literal_size() != 1 || ct->linear().vars_size() != 1 ||
1469 (ct->linear().coeffs(0) != 1 && ct->linear().coeffs(0) == -1)) {
1470 return false;
1471 }
1472
1473 const int literal = ct->enforcement_literal(0);
1474 const LinearConstraintProto& linear = ct->linear();
1475 const int ref = linear.vars(0);
1476 const int var = PositiveRef(ref);
1477 const int64 coeff =
1478 RefIsPositive(ref) ? ct->linear().coeffs(0) : -ct->linear().coeffs(0);
1479
1480 if (linear.domain_size() == 2 && linear.domain(0) == linear.domain(1)) {
1481 const int64 value = RefIsPositive(ref) ? linear.domain(0) * coeff
1482 : -linear.domain(0) * coeff;
1484 // The domain is not actually modified, but we want to rescan the
1485 // constraints linked to this variable. See TODO below.
1486 context_->modified_domains.Set(var);
1487 }
1488 } else {
1489 const Domain complement = context_->DomainOf(ref).IntersectionWith(
1490 ReadDomainFromProto(linear).Complement());
1491 if (complement.Size() != 1) return false;
1492 const int64 value = RefIsPositive(ref) ? complement.Min() * coeff
1493 : -complement.Min() * coeff;
1495 // The domain is not actually modified, but we want to rescan the
1496 // constraints linked to this variable. See TODO below.
1497 context_->modified_domains.Set(var);
1498 }
1499 }
1500
1501 // TODO(user): if we have l1 <=> x == value && l2 => x == value, we
1502 // could rewrite the second constraint into l2 => l1.
1504 return false;
1505 }
1506
1507 // Size one constraint?
1508 if (ct->linear().vars().size() == 1) {
1509 const int64 coeff = RefIsPositive(ct->linear().vars(0))
1510 ? ct->linear().coeffs(0)
1511 : -ct->linear().coeffs(0);
1512 context_->UpdateRuleStats("linear: size one");
1513 const int var = PositiveRef(ct->linear().vars(0));
1514 const Domain rhs = ReadDomainFromProto(ct->linear());
1515 if (!context_->IntersectDomainWith(var,
1516 rhs.InverseMultiplicationBy(coeff))) {
1517 return true;
1518 }
1519 return RemoveConstraint(ct);
1520 }
1521
1522 // Detect affine relation.
1523 //
1524 // TODO(user): it might be better to first add only the affine relation with
1525 // a coefficient of magnitude 1, and later the one with larger coeffs.
1526 const LinearConstraintProto& arg = ct->linear();
1527 if (arg.vars_size() == 2) {
1528 const Domain rhs = ReadDomainFromProto(ct->linear());
1529 const int64 rhs_min = rhs.Min();
1530 const int64 rhs_max = rhs.Max();
1531 if (rhs_min == rhs_max) {
1532 const int v1 = arg.vars(0);
1533 const int v2 = arg.vars(1);
1534 const int64 coeff1 = arg.coeffs(0);
1535 const int64 coeff2 = arg.coeffs(1);
1536 bool added = false;
1537 if (coeff1 == 1) {
1538 added = context_->StoreAffineRelation(v1, v2, -coeff2, rhs_max);
1539 } else if (coeff2 == 1) {
1540 added = context_->StoreAffineRelation(v2, v1, -coeff1, rhs_max);
1541 } else if (coeff1 == -1) {
1542 added = context_->StoreAffineRelation(v1, v2, coeff2, -rhs_max);
1543 } else if (coeff2 == -1) {
1544 added = context_->StoreAffineRelation(v2, v1, coeff1, -rhs_max);
1545 }
1546 if (added) return RemoveConstraint(ct);
1547 }
1548 }
1549
1550 return false;
1551}
1552
1553namespace {
1554
1555// Return true if the given domain only restrict the values with an upper bound.
1556bool IsLeConstraint(const Domain& domain, const Domain& all_values) {
1557 return all_values.IntersectionWith(Domain(kint64min, domain.Max()))
1558 .IsIncludedIn(domain);
1559}
1560
1561// Same as IsLeConstraint() but in the other direction.
1562bool IsGeConstraint(const Domain& domain, const Domain& all_values) {
1563 return all_values.IntersectionWith(Domain(domain.Min(), kint64max))
1564 .IsIncludedIn(domain);
1565}
1566
1567// In the equation terms + coeff * var_domain \included rhs, returns true if can
1568// we always fix rhs to its min value for any value in terms. It is okay to
1569// not be as generic as possible here.
1570bool RhsCanBeFixedToMin(int64 coeff, const Domain& var_domain,
1571 const Domain& terms, const Domain& rhs) {
1572 if (var_domain.NumIntervals() != 1) return false;
1573 if (std::abs(coeff) != 1) return false;
1574
1575 // If for all values in terms, there is one value below rhs.Min(), then
1576 // because we add only one integer interval, if there is a feasible value, it
1577 // can be at rhs.Min().
1578 //
1579 // TODO(user): generalize to larger coeff magnitude if rhs is also a multiple
1580 // or if terms is a multiple.
1581 if (coeff == 1 && terms.Max() + var_domain.Min() <= rhs.Min()) {
1582 return true;
1583 }
1584 if (coeff == -1 && terms.Max() - var_domain.Max() <= rhs.Min()) {
1585 return true;
1586 }
1587 return false;
1588}
1589
1590bool RhsCanBeFixedToMax(int64 coeff, const Domain& var_domain,
1591 const Domain& terms, const Domain& rhs) {
1592 if (var_domain.NumIntervals() != 1) return false;
1593 if (std::abs(coeff) != 1) return false;
1594
1595 if (coeff == 1 && terms.Min() + var_domain.Max() >= rhs.Max()) {
1596 return true;
1597 }
1598 if (coeff == -1 && terms.Min() - var_domain.Min() >= rhs.Max()) {
1599 return true;
1600 }
1601 return false;
1602}
1603
1604// Remove from to_clear any entry not in current.
1605void TakeIntersectionWith(const absl::flat_hash_set<int>& current,
1606 absl::flat_hash_set<int>* to_clear) {
1607 std::vector<int> new_set;
1608 for (const int c : *to_clear) {
1609 if (current.contains(c)) new_set.push_back(c);
1610 }
1611 to_clear->clear();
1612 for (const int c : new_set) to_clear->insert(c);
1613}
1614
1615} // namespace
1616
1617bool CpModelPresolver::PropagateDomainsInLinear(int c, ConstraintProto* ct) {
1618 if (ct->constraint_case() != ConstraintProto::ConstraintCase::kLinear ||
1619 context_->ModelIsUnsat()) {
1620 return false;
1621 }
1622
1623 // Compute the implied rhs bounds from the variable ones.
1624 auto& term_domains = context_->tmp_term_domains;
1625 auto& left_domains = context_->tmp_left_domains;
1626 const int num_vars = ct->linear().vars_size();
1627 term_domains.resize(num_vars + 1);
1628 left_domains.resize(num_vars + 1);
1629 left_domains[0] = Domain(0);
1630 for (int i = 0; i < num_vars; ++i) {
1631 const int var = ct->linear().vars(i);
1632 const int64 coeff = ct->linear().coeffs(i);
1634 term_domains[i] = context_->DomainOf(var).MultiplicationBy(coeff);
1635 left_domains[i + 1] =
1636 left_domains[i].AdditionWith(term_domains[i]).RelaxIfTooComplex();
1637 }
1638 const Domain& implied_rhs = left_domains[num_vars];
1639
1640 // Abort if trivial.
1641 const Domain old_rhs = ReadDomainFromProto(ct->linear());
1642 if (implied_rhs.IsIncludedIn(old_rhs)) {
1643 context_->UpdateRuleStats("linear: always true");
1644 return RemoveConstraint(ct);
1645 }
1646
1647 // Incorporate the implied rhs information.
1648 Domain rhs = old_rhs.SimplifyUsingImpliedDomain(implied_rhs);
1649 if (rhs.IsEmpty()) {
1650 context_->UpdateRuleStats("linear: infeasible");
1651 return MarkConstraintAsFalse(ct);
1652 }
1653 if (rhs != old_rhs) {
1654 context_->UpdateRuleStats("linear: simplified rhs");
1655 }
1656 FillDomainInProto(rhs, ct->mutable_linear());
1657
1658 // Detect if it is always good for a term of this constraint to move towards
1659 // its lower (resp. upper) bound. This is the same as saying that this
1660 // constraint only bound in one direction.
1661 bool is_le_constraint = IsLeConstraint(rhs, implied_rhs);
1662 bool is_ge_constraint = IsGeConstraint(rhs, implied_rhs);
1663
1664 // Propagate the variable bounds.
1665 if (ct->enforcement_literal().size() > 1) return false;
1666
1667 bool new_bounds = false;
1668 bool recanonicalize = false;
1669 Domain negated_rhs = rhs.Negation();
1670 Domain right_domain(0);
1671 Domain new_domain;
1672 Domain implied_term_domain;
1673 term_domains[num_vars] = Domain(0);
1674 for (int i = num_vars - 1; i >= 0; --i) {
1675 const int var = ct->linear().vars(i);
1676 const int64 var_coeff = ct->linear().coeffs(i);
1677 right_domain =
1678 right_domain.AdditionWith(term_domains[i + 1]).RelaxIfTooComplex();
1679 implied_term_domain = left_domains[i].AdditionWith(right_domain);
1680 new_domain = implied_term_domain.AdditionWith(negated_rhs)
1681 .InverseMultiplicationBy(-var_coeff);
1682
1683 if (ct->enforcement_literal().empty()) {
1684 // Push the new domain.
1685 if (!context_->IntersectDomainWith(var, new_domain, &new_bounds)) {
1686 return true;
1687 }
1688 } else if (ct->enforcement_literal().size() == 1) {
1689 // We cannot push the new domain, but we can add some deduction.
1691 if (!context_->DomainOfVarIsIncludedIn(var, new_domain)) {
1692 context_->deductions.AddDeduction(ct->enforcement_literal(0), var,
1693 new_domain);
1694 }
1695 }
1696
1697 if (context_->IsFixed(var)) {
1698 // This will make sure we remove that fixed variable from the constraint.
1699 recanonicalize = true;
1700 continue;
1701 }
1702
1703 if (is_le_constraint || is_ge_constraint) {
1704 CHECK_NE(is_le_constraint, is_ge_constraint);
1705 if ((var_coeff > 0) == is_ge_constraint) {
1706 context_->var_to_lb_only_constraints[var].insert(c);
1707 } else {
1708 context_->var_to_ub_only_constraints[var].insert(c);
1709 }
1710
1711 // Simple dual fixing: If for any feasible solution, any solution with var
1712 // higher (resp. lower) is also valid, then we can fix that variable to
1713 // its bound if it also moves the objective in the good direction.
1714 //
1715 // A bit tricky. If a linear constraint was detected to not block a
1716 // variable in one direction, this shouldn't change later (expect in the
1717 // tightening code below, and we do take care of it). However variable can
1718 // appear in new constraints.
1719 if (!context_->keep_all_feasible_solutions) {
1720 const bool is_in_objective =
1721 context_->VarToConstraints(var).contains(-1);
1722 const int size =
1723 context_->VarToConstraints(var).size() - (is_in_objective ? 1 : 0);
1724 const int64 obj_coeff =
1725 is_in_objective ? gtl::FindOrDie(context_->ObjectiveMap(), var) : 0;
1726
1727 // We cannot fix anything if the domain of the objective is excluding
1728 // some objective values.
1729 if (obj_coeff != 0 && context_->ObjectiveDomainIsConstraining()) {
1730 continue;
1731 }
1732
1733 if (obj_coeff <= 0 &&
1734 context_->var_to_lb_only_constraints[var].size() >= size) {
1735 TakeIntersectionWith(context_->VarToConstraints(var),
1736 &(context_->var_to_lb_only_constraints[var]));
1737 if (context_->var_to_lb_only_constraints[var].size() >= size) {
1738 if (!context_->IntersectDomainWith(var,
1739 Domain(context_->MaxOf(var)))) {
1740 return false;
1741 }
1742 context_->UpdateRuleStats("linear: dual fixing");
1743 recanonicalize = true;
1744 continue;
1745 }
1746 }
1747 if (obj_coeff >= 0 &&
1748 context_->var_to_ub_only_constraints[var].size() >= size) {
1749 TakeIntersectionWith(context_->VarToConstraints(var),
1750 &(context_->var_to_ub_only_constraints[var]));
1751 if (context_->var_to_ub_only_constraints[var].size() >= size) {
1752 if (!context_->IntersectDomainWith(var,
1753 Domain(context_->MinOf(var)))) {
1754 return false;
1755 }
1756 context_->UpdateRuleStats("linear: dual fixing");
1757 recanonicalize = true;
1758 continue;
1759 }
1760 }
1761 }
1762 }
1763
1764 // The other transformations below require a non-reified constraint.
1765 if (!ct->enforcement_literal().empty()) continue;
1766
1767 // Given a variable that only appear in one constraint and in the
1768 // objective, for any feasible solution, it will be always better to move
1769 // this singleton variable as much as possible towards its good objective
1770 // direction. Sometimes, we can detect that we will always be able to do
1771 // this until the only constraint of this singleton variable is tight.
1772 //
1773 // When this happens, we can make the constraint an equality. Note that it
1774 // might not always be good to restrict constraint like this, but in this
1775 // case, the RemoveSingletonInLinear() code should be able to remove this
1776 // variable altogether.
1777 if (rhs.Min() != rhs.Max() &&
1779 const int64 obj_coeff = gtl::FindOrDie(context_->ObjectiveMap(), var);
1780 const bool same_sign = (var_coeff > 0) == (obj_coeff > 0);
1781 bool fixed = false;
1782 if (same_sign && RhsCanBeFixedToMin(var_coeff, context_->DomainOf(var),
1783 implied_term_domain, rhs)) {
1784 rhs = Domain(rhs.Min());
1785 fixed = true;
1786 }
1787 if (!same_sign && RhsCanBeFixedToMax(var_coeff, context_->DomainOf(var),
1788 implied_term_domain, rhs)) {
1789 rhs = Domain(rhs.Max());
1790 fixed = true;
1791 }
1792 if (fixed) {
1793 context_->UpdateRuleStats("linear: tightened into equality");
1794 FillDomainInProto(rhs, ct->mutable_linear());
1795 negated_rhs = rhs.Negation();
1796
1797 // Restart the loop.
1798 i = num_vars;
1799 right_domain = Domain(0);
1800
1801 // An equality is a >= (or <=) constraint iff all its term are fixed.
1802 // Since we restart the loop, we will detect that.
1803 is_le_constraint = false;
1804 is_ge_constraint = false;
1805 for (const int var : ct->linear().vars()) {
1806 context_->var_to_lb_only_constraints[var].erase(c);
1807 context_->var_to_ub_only_constraints[var].erase(c);
1808 }
1809 continue;
1810 }
1811 }
1812
1813 // Can we perform some substitution?
1814 //
1815 // TODO(user): there is no guarantee we will not miss some since we might
1816 // not reprocess a constraint once other have been deleted.
1817
1818 // Skip affine constraint. It is more efficient to substitute them lazily
1819 // when we process other constraints. Note that if we relax the fact that
1820 // we substitute only equalities, we can deal with inequality of size 2
1821 // here.
1822 if (ct->linear().vars().size() <= 2) continue;
1823
1824 // TODO(user): We actually do not need a strict equality when
1825 // keep_all_feasible_solutions is false, but that simplifies things as the
1826 // SubstituteVariable() function cannot fail this way.
1827 if (rhs.Min() != rhs.Max()) continue;
1828
1829 // Only consider "implied free" variables. Note that the coefficient of
1830 // magnitude 1 is important otherwise we can't easily remove the
1831 // constraint since the fact that the sum of the other terms must be a
1832 // multiple of coeff will not be enforced anymore.
1833 if (context_->DomainOf(var) != new_domain) continue;
1834 if (std::abs(var_coeff) != 1) continue;
1835 if (context_->params().presolve_substitution_level() <= 0) continue;
1836
1837 // NOTE: The mapping doesn't allow us to remove a variable if
1838 // 'keep_all_feasible_solutions' is true.
1839 if (context_->keep_all_feasible_solutions) continue;
1840
1841 bool is_in_objective = false;
1842 if (context_->VarToConstraints(var).contains(-1)) {
1843 is_in_objective = true;
1844 DCHECK(context_->ObjectiveMap().contains(var));
1845 }
1846
1847 // Only consider low degree columns.
1848 int col_size = context_->VarToConstraints(var).size();
1849 if (is_in_objective) col_size--;
1850 const int row_size = ct->linear().vars_size();
1851
1852 // This is actually an upper bound on the number of entries added since
1853 // some of them might already be present.
1854 const int num_entries_added = (row_size - 1) * (col_size - 1);
1855 const int num_entries_removed = col_size + row_size - 1;
1856
1857 if (num_entries_added > num_entries_removed) {
1858 continue;
1859 }
1860
1861 // Check pre-conditions on all the constraints in which this variable
1862 // appear. Basically they must all be linear.
1863 std::vector<int> others;
1864 bool abort = false;
1865 for (const int c : context_->VarToConstraints(var)) {
1866 if (c == kObjectiveConstraint) continue;
1867 if (c == kAffineRelationConstraint) {
1868 abort = true;
1869 break;
1870 }
1871 if (context_->working_model->mutable_constraints(c) == ct) continue;
1872 if (context_->working_model->constraints(c).constraint_case() !=
1873 ConstraintProto::ConstraintCase::kLinear) {
1874 abort = true;
1875 break;
1876 }
1877 for (const int ref :
1878 context_->working_model->constraints(c).enforcement_literal()) {
1879 if (PositiveRef(ref) == var) {
1880 abort = true;
1881 break;
1882 }
1883 }
1884 others.push_back(c);
1885 }
1886 if (abort) continue;
1887
1888 // Do the actual substitution.
1889 for (const int c : others) {
1890 // TODO(user): In some corner cases, this might create integer overflow
1891 // issues. The danger is limited since the range of the linear
1892 // expression used in the definition do not exceed the domain of the
1893 // variable we substitute.
1894 SubstituteVariable(var, var_coeff, *ct,
1895 context_->working_model->mutable_constraints(c));
1896
1897 // TODO(user): We should re-enqueue these constraints for presolve.
1898 context_->UpdateConstraintVariableUsage(c);
1899 }
1900
1901 // Substitute in objective.
1902 if (is_in_objective) {
1903 context_->SubstituteVariableInObjective(var, var_coeff, *ct);
1904 }
1905
1906 context_->UpdateRuleStats(
1907 absl::StrCat("linear: variable substitution ", others.size()));
1908
1909 // The variable now only appear in its definition and we can remove it
1910 // because it was implied free.
1911 //
1912 // Tricky: If the linear constraint contains other variables that are only
1913 // used here, then the postsolve needs more info. We do need to indicate
1914 // that whatever the value of those other variables, we will have a way to
1915 // assign var. We do that by putting it fist.
1916 CHECK_EQ(context_->VarToConstraints(var).size(), 1);
1917 context_->MarkVariableAsRemoved(var);
1918 const int ct_index = context_->mapping_model->constraints().size();
1919 *context_->mapping_model->add_constraints() = *ct;
1920 LinearConstraintProto* mapping_linear_ct =
1921 context_->mapping_model->mutable_constraints(ct_index)
1922 ->mutable_linear();
1923 std::swap(mapping_linear_ct->mutable_vars()->at(0),
1924 mapping_linear_ct->mutable_vars()->at(i));
1925 std::swap(mapping_linear_ct->mutable_coeffs()->at(0),
1926 mapping_linear_ct->mutable_coeffs()->at(i));
1927 return RemoveConstraint(ct);
1928 }
1929 if (new_bounds) {
1930 context_->UpdateRuleStats("linear: reduced variable domains");
1931 }
1932 if (recanonicalize) return CanonicalizeLinear(ct);
1933 return false;
1934}
1935
1936// Identify Boolean variable that makes the constraint always true when set to
1937// true or false. Moves such literal to the constraint enforcement literals
1938// list.
1939//
1940// We also generalize this to integer variable at one of their bound.
1941//
1942// This operation is similar to coefficient strengthening in the MIP world.
1943void CpModelPresolver::ExtractEnforcementLiteralFromLinearConstraint(
1944 int ct_index, ConstraintProto* ct) {
1945 if (ct->constraint_case() != ConstraintProto::ConstraintCase::kLinear ||
1946 context_->ModelIsUnsat()) {
1947 return;
1948 }
1949
1950 const LinearConstraintProto& arg = ct->linear();
1951 const int num_vars = arg.vars_size();
1952
1953 // No need to process size one constraints, they will be presolved separately.
1954 // We also do not want to split them in two.
1955 if (num_vars <= 1) return;
1956
1957 int64 min_sum = 0;
1958 int64 max_sum = 0;
1959 int64 max_coeff_magnitude = 0;
1960 for (int i = 0; i < num_vars; ++i) {
1961 const int ref = arg.vars(i);
1962 const int64 coeff = arg.coeffs(i);
1963 const int64 term_a = coeff * context_->MinOf(ref);
1964 const int64 term_b = coeff * context_->MaxOf(ref);
1965 max_coeff_magnitude = std::max(max_coeff_magnitude, std::abs(coeff));
1966 min_sum += std::min(term_a, term_b);
1967 max_sum += std::max(term_a, term_b);
1968 }
1969
1970 // We can only extract enforcement literals if the maximum coefficient
1971 // magnitude is large enough. Note that we handle complex domain.
1972 //
1973 // TODO(user): Depending on how we split below, the threshold are not the
1974 // same. This is maybe not too important, we just don't split as often as we
1975 // could, but it is still unclear if splitting is good.
1976 const auto& domain = ct->linear().domain();
1977 const int64 ub_threshold = domain[domain.size() - 2] - min_sum;
1978 const int64 lb_threshold = max_sum - domain[1];
1979 const Domain rhs_domain = ReadDomainFromProto(ct->linear());
1980 if (max_coeff_magnitude < std::max(ub_threshold, lb_threshold)) return;
1981
1982 // We need the constraint to be only bounded on one side in order to extract
1983 // enforcement literal.
1984 //
1985 // If it is boxed and we know that some coefficient are big enough (see test
1986 // above), then we split the constraint in two. That might not seems always
1987 // good, but for the CP propagation engine, we don't loose anything by doing
1988 // so, and for the LP we will regroup the constraints if they still have the
1989 // exact same coeff after the presolve.
1990 //
1991 // TODO(user): Creating two new constraints and removing the current one might
1992 // not be the most efficient, but it simplify the presolve code by not having
1993 // to do anything special to trigger a new presolving of these constraints.
1994 // Try to improve if this becomes a problem.
1995 //
1996 // TODO(user): At the end of the presolve we should probably remerge any
1997 // identical linear constraints. That also cover the corner cases where
1998 // constraints are just redundant...
1999 const bool lower_bounded = min_sum < rhs_domain.Min();
2000 const bool upper_bounded = max_sum > rhs_domain.Max();
2001 if (!lower_bounded && !upper_bounded) return;
2002 if (lower_bounded && upper_bounded) {
2003 context_->UpdateRuleStats("linear: split boxed constraint");
2004 ConstraintProto* new_ct1 = context_->working_model->add_constraints();
2005 *new_ct1 = *ct;
2006 if (!ct->name().empty()) {
2007 new_ct1->set_name(absl::StrCat(ct->name(), " (part 1)"));
2008 }
2009 FillDomainInProto(Domain(min_sum, rhs_domain.Max()),
2010 new_ct1->mutable_linear());
2011
2012 ConstraintProto* new_ct2 = context_->working_model->add_constraints();
2013 *new_ct2 = *ct;
2014 if (!ct->name().empty()) {
2015 new_ct2->set_name(absl::StrCat(ct->name(), " (part 2)"));
2016 }
2017 FillDomainInProto(rhs_domain.UnionWith(Domain(rhs_domain.Max(), max_sum)),
2018 new_ct2->mutable_linear());
2019
2021 return (void)RemoveConstraint(ct);
2022 }
2023
2024 // Any coefficient greater than this will cause the constraint to be trivially
2025 // satisfied when the variable move away from its bound. Note that as we
2026 // remove coefficient, the threshold do not change!
2027 const int64 threshold = lower_bounded ? ub_threshold : lb_threshold;
2028
2029 // Do we only extract Booleans?
2030 const bool only_booleans =
2031 !context_->params().presolve_extract_integer_enforcement();
2032
2033 // To avoid a quadratic loop, we will rewrite the linear expression at the
2034 // same time as we extract enforcement literals.
2035 int new_size = 0;
2036 int64 rhs_offset = 0;
2037 bool some_integer_encoding_were_extracted = false;
2038 LinearConstraintProto* mutable_arg = ct->mutable_linear();
2039 for (int i = 0; i < arg.vars_size(); ++i) {
2040 int ref = arg.vars(i);
2041 int64 coeff = arg.coeffs(i);
2042 if (coeff < 0) {
2043 ref = NegatedRef(ref);
2044 coeff = -coeff;
2045 }
2046
2047 const bool is_boolean = context_->CanBeUsedAsLiteral(ref);
2048 if (context_->IsFixed(ref) || coeff < threshold ||
2049 (only_booleans && !is_boolean)) {
2050 // We keep this term.
2051 mutable_arg->set_vars(new_size, mutable_arg->vars(i));
2052 mutable_arg->set_coeffs(new_size, mutable_arg->coeffs(i));
2053 ++new_size;
2054 continue;
2055 }
2056
2057 if (is_boolean) {
2058 context_->UpdateRuleStats("linear: extracted enforcement literal");
2059 } else {
2060 some_integer_encoding_were_extracted = true;
2061 context_->UpdateRuleStats(
2062 "linear: extracted integer enforcement literal");
2063 }
2064 if (lower_bounded) {
2065 ct->add_enforcement_literal(is_boolean
2066 ? NegatedRef(ref)
2067 : context_->GetOrCreateVarValueEncoding(
2068 ref, context_->MinOf(ref)));
2069 rhs_offset -= coeff * context_->MinOf(ref);
2070 } else {
2071 ct->add_enforcement_literal(is_boolean
2072 ? ref
2073 : context_->GetOrCreateVarValueEncoding(
2074 ref, context_->MaxOf(ref)));
2075 rhs_offset -= coeff * context_->MaxOf(ref);
2076 }
2077 }
2078 mutable_arg->mutable_vars()->Truncate(new_size);
2079 mutable_arg->mutable_coeffs()->Truncate(new_size);
2080 FillDomainInProto(rhs_domain.AdditionWith(Domain(rhs_offset)), mutable_arg);
2081 if (some_integer_encoding_were_extracted) {
2083 context_->UpdateConstraintVariableUsage(ct_index);
2084 }
2085}
2086
2087void CpModelPresolver::ExtractAtMostOneFromLinear(ConstraintProto* ct) {
2088 if (context_->ModelIsUnsat()) return;
2089 if (HasEnforcementLiteral(*ct)) return;
2090 const Domain rhs = ReadDomainFromProto(ct->linear());
2091
2092 const LinearConstraintProto& arg = ct->linear();
2093 const int num_vars = arg.vars_size();
2094 int64 min_sum = 0;
2095 int64 max_sum = 0;
2096 for (int i = 0; i < num_vars; ++i) {
2097 const int ref = arg.vars(i);
2098 const int64 coeff = arg.coeffs(i);
2099 const int64 term_a = coeff * context_->MinOf(ref);
2100 const int64 term_b = coeff * context_->MaxOf(ref);
2101 min_sum += std::min(term_a, term_b);
2102 max_sum += std::max(term_a, term_b);
2103 }
2104 for (const int type : {0, 1}) {
2105 std::vector<int> at_most_one;
2106 for (int i = 0; i < num_vars; ++i) {
2107 const int ref = arg.vars(i);
2108 const int64 coeff = arg.coeffs(i);
2109 if (context_->MinOf(ref) != 0) continue;
2110 if (context_->MaxOf(ref) != 1) continue;
2111
2112 if (type == 0) {
2113 // TODO(user): we could add one more Boolean with a lower coeff as long
2114 // as we have lower_coeff + min_of_other_coeff > rhs.Max().
2115 if (min_sum + 2 * std::abs(coeff) > rhs.Max()) {
2116 at_most_one.push_back(coeff > 0 ? ref : NegatedRef(ref));
2117 }
2118 } else {
2119 if (max_sum - 2 * std::abs(coeff) < rhs.Min()) {
2120 at_most_one.push_back(coeff > 0 ? NegatedRef(ref) : ref);
2121 }
2122 }
2123 }
2124 if (at_most_one.size() > 1) {
2125 if (type == 0) {
2126 context_->UpdateRuleStats("linear: extracted at most one (max).");
2127 } else {
2128 context_->UpdateRuleStats("linear: extracted at most one (min).");
2129 }
2130 ConstraintProto* new_ct = context_->working_model->add_constraints();
2131 new_ct->set_name(ct->name());
2132 for (const int ref : at_most_one) {
2133 new_ct->mutable_at_most_one()->add_literals(ref);
2134 }
2136 }
2137 }
2138}
2139
2140// Convert some linear constraint involving only Booleans to their Boolean
2141// form.
2142bool CpModelPresolver::PresolveLinearOnBooleans(ConstraintProto* ct) {
2143 if (ct->constraint_case() != ConstraintProto::ConstraintCase::kLinear ||
2144 context_->ModelIsUnsat()) {
2145 return false;
2146 }
2147
2148 const LinearConstraintProto& arg = ct->linear();
2149 const int num_vars = arg.vars_size();
2150 int64 min_coeff = kint64max;
2151 int64 max_coeff = 0;
2152 int64 min_sum = 0;
2153 int64 max_sum = 0;
2154 for (int i = 0; i < num_vars; ++i) {
2155 // We assume we already ran PresolveLinear().
2156 const int var = arg.vars(i);
2157 const int64 coeff = arg.coeffs(i);
2159 CHECK_NE(coeff, 0);
2160 if (context_->MinOf(var) != 0) return false;
2161 if (context_->MaxOf(var) != 1) return false;
2162
2163 if (coeff > 0) {
2164 max_sum += coeff;
2165 min_coeff = std::min(min_coeff, coeff);
2166 max_coeff = std::max(max_coeff, coeff);
2167 } else {
2168 // We replace the Boolean ref, by a ref to its negation (1 - x).
2169 min_sum += coeff;
2170 min_coeff = std::min(min_coeff, -coeff);
2171 max_coeff = std::max(max_coeff, -coeff);
2172 }
2173 }
2174 CHECK_LE(min_coeff, max_coeff);
2175
2176 // Detect trivially true/false constraints. Note that this is not necessarily
2177 // detected by PresolveLinear(). We do that here because we assume below
2178 // that this cannot happen.
2179 //
2180 // TODO(user): this could be generalized to constraint not containing only
2181 // Booleans.
2182 const Domain rhs_domain = ReadDomainFromProto(arg);
2183 if ((!rhs_domain.Contains(min_sum) &&
2184 min_sum + min_coeff > rhs_domain.Max()) ||
2185 (!rhs_domain.Contains(max_sum) &&
2186 max_sum - min_coeff < rhs_domain.Min())) {
2187 context_->UpdateRuleStats("linear: all booleans and trivially false");
2188 return MarkConstraintAsFalse(ct);
2189 }
2190 if (Domain(min_sum, max_sum).IsIncludedIn(rhs_domain)) {
2191 context_->UpdateRuleStats("linear: all booleans and trivially true");
2192 return RemoveConstraint(ct);
2193 }
2194
2195 // Detect clauses, reified ands, at_most_one.
2196 //
2197 // TODO(user): split a == 1 constraint or similar into a clause and an at
2198 // most one constraint?
2199 DCHECK(!rhs_domain.IsEmpty());
2200 if (min_sum + min_coeff > rhs_domain.Max()) {
2201 // All Boolean are false if the reified literal is true.
2202 context_->UpdateRuleStats("linear: negative reified and");
2203 const auto copy = arg;
2204 ct->mutable_bool_and()->clear_literals();
2205 for (int i = 0; i < num_vars; ++i) {
2206 ct->mutable_bool_and()->add_literals(
2207 copy.coeffs(i) > 0 ? NegatedRef(copy.vars(i)) : copy.vars(i));
2208 }
2209 return PresolveBoolAnd(ct);
2210 } else if (max_sum - min_coeff < rhs_domain.Min()) {
2211 // All Boolean are true if the reified literal is true.
2212 context_->UpdateRuleStats("linear: positive reified and");
2213 const auto copy = arg;
2214 ct->mutable_bool_and()->clear_literals();
2215 for (int i = 0; i < num_vars; ++i) {
2216 ct->mutable_bool_and()->add_literals(
2217 copy.coeffs(i) > 0 ? copy.vars(i) : NegatedRef(copy.vars(i)));
2218 }
2219 return PresolveBoolAnd(ct);
2220 } else if (min_sum + min_coeff >= rhs_domain.Min() &&
2221 rhs_domain.front().end >= max_sum) {
2222 // At least one Boolean is true.
2223 context_->UpdateRuleStats("linear: positive clause");
2224 const auto copy = arg;
2225 ct->mutable_bool_or()->clear_literals();
2226 for (int i = 0; i < num_vars; ++i) {
2227 ct->mutable_bool_or()->add_literals(
2228 copy.coeffs(i) > 0 ? copy.vars(i) : NegatedRef(copy.vars(i)));
2229 }
2230 return PresolveBoolOr(ct);
2231 } else if (max_sum - min_coeff <= rhs_domain.Max() &&
2232 rhs_domain.back().start <= min_sum) {
2233 // At least one Boolean is false.
2234 context_->UpdateRuleStats("linear: negative clause");
2235 const auto copy = arg;
2236 ct->mutable_bool_or()->clear_literals();
2237 for (int i = 0; i < num_vars; ++i) {
2238 ct->mutable_bool_or()->add_literals(
2239 copy.coeffs(i) > 0 ? NegatedRef(copy.vars(i)) : copy.vars(i));
2240 }
2241 return PresolveBoolOr(ct);
2242 } else if (!HasEnforcementLiteral(*ct) &&
2243 min_sum + max_coeff <= rhs_domain.Max() &&
2244 min_sum + 2 * min_coeff > rhs_domain.Max() &&
2245 rhs_domain.back().start <= min_sum) {
2246 // At most one Boolean is true.
2247 context_->UpdateRuleStats("linear: positive at most one");
2248 const auto copy = arg;
2249 ct->mutable_at_most_one()->clear_literals();
2250 for (int i = 0; i < num_vars; ++i) {
2251 ct->mutable_at_most_one()->add_literals(
2252 copy.coeffs(i) > 0 ? copy.vars(i) : NegatedRef(copy.vars(i)));
2253 }
2254 return true;
2255 } else if (!HasEnforcementLiteral(*ct) &&
2256 max_sum - max_coeff >= rhs_domain.Min() &&
2257 max_sum - 2 * min_coeff < rhs_domain.Min() &&
2258 rhs_domain.front().end >= max_sum) {
2259 // At most one Boolean is false.
2260 context_->UpdateRuleStats("linear: negative at most one");
2261 const auto copy = arg;
2262 ct->mutable_at_most_one()->clear_literals();
2263 for (int i = 0; i < num_vars; ++i) {
2264 ct->mutable_at_most_one()->add_literals(
2265 copy.coeffs(i) > 0 ? NegatedRef(copy.vars(i)) : copy.vars(i));
2266 }
2267 return true;
2268 } else if (!HasEnforcementLiteral(*ct) && rhs_domain.NumIntervals() == 1 &&
2269 min_sum < rhs_domain.Min() &&
2270 min_sum + min_coeff >= rhs_domain.Min() &&
2271 min_sum + 2 * min_coeff > rhs_domain.Max() &&
2272 min_sum + max_coeff <= rhs_domain.Max()) {
2273 context_->UpdateRuleStats("linear: positive equal one");
2274 ConstraintProto* exactly_one = context_->working_model->add_constraints();
2275 exactly_one->set_name(ct->name());
2276 for (int i = 0; i < num_vars; ++i) {
2277 exactly_one->mutable_exactly_one()->add_literals(
2278 arg.coeffs(i) > 0 ? arg.vars(i) : NegatedRef(arg.vars(i)));
2279 }
2281 return RemoveConstraint(ct);
2282 } else if (!HasEnforcementLiteral(*ct) && rhs_domain.NumIntervals() == 1 &&
2283 max_sum > rhs_domain.Max() &&
2284 max_sum - min_coeff <= rhs_domain.Max() &&
2285 max_sum - 2 * min_coeff < rhs_domain.Min() &&
2286 max_sum - max_coeff >= rhs_domain.Min()) {
2287 context_->UpdateRuleStats("linear: negative equal one");
2288 ConstraintProto* exactly_one = context_->working_model->add_constraints();
2289 exactly_one->set_name(ct->name());
2290 for (int i = 0; i < num_vars; ++i) {
2291 exactly_one->mutable_exactly_one()->add_literals(
2292 arg.coeffs(i) > 0 ? NegatedRef(arg.vars(i)) : arg.vars(i));
2293 }
2295 return RemoveConstraint(ct);
2296 }
2297
2298 // Expand small expression into clause.
2299 //
2300 // TODO(user): This is bad from a LP relaxation perspective. Do not do that
2301 // now? On another hand it is good for the SAT presolving.
2302 if (num_vars > 3) return false;
2303 context_->UpdateRuleStats("linear: small Boolean expression");
2304
2305 // Enumerate all possible value of the Booleans and add a clause if constraint
2306 // is false. TODO(user): the encoding could be made better in some cases.
2307 const int max_mask = (1 << arg.vars_size());
2308 for (int mask = 0; mask < max_mask; ++mask) {
2309 int64 value = 0;
2310 for (int i = 0; i < num_vars; ++i) {
2311 if ((mask >> i) & 1) value += arg.coeffs(i);
2312 }
2313 if (rhs_domain.Contains(value)) continue;
2314
2315 // Add a new clause to exclude this bad assignment.
2316 ConstraintProto* new_ct = context_->working_model->add_constraints();
2317 auto* new_arg = new_ct->mutable_bool_or();
2318 if (HasEnforcementLiteral(*ct)) {
2319 *new_ct->mutable_enforcement_literal() = ct->enforcement_literal();
2320 }
2321 for (int i = 0; i < num_vars; ++i) {
2322 new_arg->add_literals(((mask >> i) & 1) ? NegatedRef(arg.vars(i))
2323 : arg.vars(i));
2324 }
2325 }
2326
2328 return RemoveConstraint(ct);
2329}
2330
2331namespace {
2332
2333void AddLinearExpression(int64 coeff, const LinearExpressionProto& exp,
2334 LinearConstraintProto* linear_ct) {
2335 const int size = exp.vars().size();
2336 for (int i = 0; i < size; ++i) {
2337 linear_ct->add_vars(exp.vars(i));
2338 linear_ct->add_coeffs(coeff * exp.coeffs(i));
2339 }
2340 if (exp.offset() != 0) {
2342 .AdditionWith(Domain(-coeff * exp.offset())),
2343 linear_ct);
2344 }
2345}
2346
2347} // namespace
2348
2349bool CpModelPresolver::PresolveInterval(int c, ConstraintProto* ct) {
2350 if (context_->ModelIsUnsat()) return false;
2351
2352 if (ct->enforcement_literal().empty() && !ct->interval().has_start_view()) {
2353 bool changed = false;
2354 const int start = ct->interval().start();
2355 const int end = ct->interval().end();
2356 const int size = ct->interval().size();
2357 const Domain start_domain = context_->DomainOf(start);
2358 const Domain end_domain = context_->DomainOf(end);
2359 const Domain size_domain = context_->DomainOf(size);
2360 // Size can't be negative.
2361 if (!context_->IntersectDomainWith(size, Domain(0, context_->MaxOf(size)),
2362 &changed)) {
2363 return false;
2364 }
2365 if (!context_->IntersectDomainWith(
2366 end, start_domain.AdditionWith(size_domain), &changed)) {
2367 return false;
2368 }
2369 if (!context_->IntersectDomainWith(
2370 start, end_domain.AdditionWith(size_domain.Negation()), &changed)) {
2371 return false;
2372 }
2373 if (!context_->IntersectDomainWith(
2374 size, end_domain.AdditionWith(start_domain.Negation()), &changed)) {
2375 return false;
2376 }
2377 if (changed) {
2378 context_->UpdateRuleStats("interval: reduced domains");
2379 }
2380 }
2381
2382 if (context_->IntervalUsage(c) == 0) {
2383 if (!ct->interval().has_start_view()) {
2384 // Convert to linear.
2385 const int start = ct->interval().start();
2386 const int end = ct->interval().end();
2387 const int size = ct->interval().size();
2388 ConstraintProto* new_ct = context_->working_model->add_constraints();
2389 *(new_ct->mutable_enforcement_literal()) = ct->enforcement_literal();
2390 new_ct->mutable_linear()->add_domain(0);
2391 new_ct->mutable_linear()->add_domain(0);
2392 new_ct->mutable_linear()->add_vars(start);
2393 new_ct->mutable_linear()->add_coeffs(1);
2394 new_ct->mutable_linear()->add_vars(size);
2395 new_ct->mutable_linear()->add_coeffs(1);
2396 new_ct->mutable_linear()->add_vars(end);
2397 new_ct->mutable_linear()->add_coeffs(-1);
2399 }
2400 context_->UpdateRuleStats("interval: unused, converted to linear");
2401 return RemoveConstraint(ct);
2402 }
2403
2404 // TODO(user): Note that the conversion is not perfect for optional intervals.
2405 // because for a fixed size optional interval with a different start and end
2406 // variable, because of the optionality we will not be able to detect the
2407 // affine relation between start and end. So we will no remove a variable like
2408 // we do for non-optional fixed size intervals.
2409 if (context_->params().convert_intervals()) {
2410 bool changed = false;
2411 IntervalConstraintProto* interval = ct->mutable_interval();
2412 if (!ct->interval().has_start_view()) {
2413 changed = true;
2414
2415 // Fill the view fields.
2416 interval->mutable_start_view()->add_vars(interval->start());
2417 interval->mutable_start_view()->add_coeffs(1);
2418 interval->mutable_start_view()->set_offset(0);
2419 interval->mutable_size_view()->add_vars(interval->size());
2420 interval->mutable_size_view()->add_coeffs(1);
2421 interval->mutable_size_view()->set_offset(0);
2422 interval->mutable_end_view()->add_vars(interval->end());
2423 interval->mutable_end_view()->add_coeffs(1);
2424 interval->mutable_end_view()->set_offset(0);
2425
2426 // Create a new linear constraint to detect affine relation and propagate
2427 // the domain properly.
2428 //
2429 // Note(user): Aving an extra constraint is not super clean, but reusing
2430 // the code is a must.
2431 ConstraintProto* new_ct = context_->working_model->add_constraints();
2432 *(new_ct->mutable_enforcement_literal()) = ct->enforcement_literal();
2433 new_ct->mutable_linear()->add_domain(0);
2434 new_ct->mutable_linear()->add_domain(0);
2435 AddLinearExpression(1, interval->start_view(), new_ct->mutable_linear());
2436 AddLinearExpression(1, interval->size_view(), new_ct->mutable_linear());
2437 AddLinearExpression(-1, interval->end_view(), new_ct->mutable_linear());
2439
2440 // Set the old fields to their default. Not really needed.
2441 interval->set_start(0);
2442 interval->set_size(0);
2443 interval->set_end(0);
2444 }
2445
2446 changed |=
2447 CanonicalizeLinearExpression(*ct, interval->mutable_start_view());
2448 changed |= CanonicalizeLinearExpression(*ct, interval->mutable_size_view());
2449 changed |= CanonicalizeLinearExpression(*ct, interval->mutable_end_view());
2450 return changed;
2451 }
2452
2453 // If the interval is of fixed size, we can add the corresponsing affine
2454 // relation to our pool.
2455 //
2456 // TODO(user): This will currently add another linear relation to the proto
2457 // in addition to the interval at the end of the presolve though.
2458 if (/* DISABLES CODE */ (false) && ct->enforcement_literal().empty() &&
2459 context_->IsFixed(ct->interval().size())) {
2460 context_->StoreAffineRelation(ct->interval().end(), ct->interval().start(),
2461 1, context_->MinOf(ct->interval().size()));
2462 }
2463
2464 // This never change the constraint-variable graph.
2465 return false;
2466}
2467
2468bool CpModelPresolver::PresolveElement(ConstraintProto* ct) {
2469 if (context_->ModelIsUnsat()) return false;
2470
2471 const int index_ref = ct->element().index();
2472 const int target_ref = ct->element().target();
2473
2474 // TODO(user): think about this once we do have such constraint.
2475 if (HasEnforcementLiteral(*ct)) return false;
2476
2477 bool all_constants = true;
2478 absl::flat_hash_set<int64> constant_set;
2479 bool all_included_in_target_domain = true;
2480
2481 {
2482 bool reduced_index_domain = false;
2483 if (!context_->IntersectDomainWith(index_ref,
2484 Domain(0, ct->element().vars_size() - 1),
2485 &reduced_index_domain)) {
2486 return false;
2487 }
2488
2489 // Filter impossible index values if index == +/- target
2490 //
2491 // Note that this must be done before the unique_index/target rule.
2492 if (PositiveRef(target_ref) == PositiveRef(index_ref)) {
2493 std::vector<int64> possible_indices;
2494 const Domain& index_domain = context_->DomainOf(index_ref);
2495 for (const ClosedInterval& interval : index_domain) {
2496 for (int64 index_value = interval.start; index_value <= interval.end;
2497 ++index_value) {
2498 const int ref = ct->element().vars(index_value);
2499 const int64 target_value =
2500 target_ref == index_ref ? index_value : -index_value;
2501 if (context_->DomainContains(ref, target_value)) {
2502 possible_indices.push_back(target_value);
2503 }
2504 }
2505 }
2506 if (possible_indices.size() < index_domain.Size()) {
2507 if (!context_->IntersectDomainWith(
2508 index_ref, Domain::FromValues(possible_indices))) {
2509 return true;
2510 }
2511 context_->UpdateRuleStats(
2512 "element: reduced index domain when target equals index");
2513 }
2514 }
2515
2516 // Filter possible index values. Accumulate variable domains to build
2517 // a possible target domain.
2518 Domain infered_domain;
2519 const Domain& initial_index_domain = context_->DomainOf(index_ref);
2520 const Domain& target_domain = context_->DomainOf(target_ref);
2521 std::vector<int64> possible_indices;
2522 for (const ClosedInterval interval : initial_index_domain) {
2523 for (int value = interval.start; value <= interval.end; ++value) {
2524 CHECK_GE(value, 0);
2525 CHECK_LT(value, ct->element().vars_size());
2526 const int ref = ct->element().vars(value);
2527 const Domain& domain = context_->DomainOf(ref);
2528 if (domain.IntersectionWith(target_domain).IsEmpty()) continue;
2529 possible_indices.push_back(value);
2530 if (domain.IsFixed()) {
2531 constant_set.insert(domain.Min());
2532 } else {
2533 all_constants = false;
2534 }
2535 if (!domain.IsIncludedIn(target_domain)) {
2536 all_included_in_target_domain = false;
2537 }
2538 infered_domain = infered_domain.UnionWith(domain);
2539 }
2540 }
2541 if (possible_indices.size() < initial_index_domain.Size()) {
2542 if (!context_->IntersectDomainWith(
2543 index_ref, Domain::FromValues(possible_indices))) {
2544 return true;
2545 }
2546 context_->UpdateRuleStats("element: reduced index domain");
2547 }
2548 bool domain_modified = false;
2549 if (!context_->IntersectDomainWith(target_ref, infered_domain,
2550 &domain_modified)) {
2551 return true;
2552 }
2553 if (domain_modified) {
2554 context_->UpdateRuleStats("element: reduced target domain");
2555 }
2556 }
2557
2558 // If the index is fixed, this is a equality constraint.
2559 if (context_->IsFixed(index_ref)) {
2560 const int var = ct->element().vars(context_->MinOf(index_ref));
2561 if (var != target_ref) {
2562 LinearConstraintProto* const lin =
2563 context_->working_model->add_constraints()->mutable_linear();
2564 lin->add_vars(var);
2565 lin->add_coeffs(-1);
2566 lin->add_vars(target_ref);
2567 lin->add_coeffs(1);
2568 lin->add_domain(0);
2569 lin->add_domain(0);
2571 }
2572 context_->UpdateRuleStats("element: fixed index");
2573 return RemoveConstraint(ct);
2574 }
2575
2576 // If the accessible part of the array is made of a single constant value,
2577 // then we do not care about the index. And, because of the previous target
2578 // domain reduction, the target is also fixed.
2579 if (all_constants && constant_set.size() == 1) {
2580 CHECK(context_->IsFixed(target_ref));
2581 context_->UpdateRuleStats("element: one value array");
2582 return RemoveConstraint(ct);
2583 }
2584
2585 // Special case when the index is boolean, and the array does not contain
2586 // variables.
2587 if (context_->MinOf(index_ref) == 0 && context_->MaxOf(index_ref) == 1 &&
2588 all_constants) {
2589 const int64 v0 = context_->MinOf(ct->element().vars(0));
2590 const int64 v1 = context_->MinOf(ct->element().vars(1));
2591
2592 LinearConstraintProto* const lin =
2593 context_->working_model->add_constraints()->mutable_linear();
2594 lin->add_vars(target_ref);
2595 lin->add_coeffs(1);
2596 lin->add_vars(index_ref);
2597 lin->add_coeffs(v0 - v1);
2598 lin->add_domain(v0);
2599 lin->add_domain(v0);
2601 context_->UpdateRuleStats("element: linearize constant element of size 2");
2602 return RemoveConstraint(ct);
2603 }
2604
2605 // If the index has a canonical affine representative, rewrite the element.
2606 const AffineRelation::Relation r_index =
2607 context_->GetAffineRelation(index_ref);
2608 if (r_index.representative != index_ref) {
2609 // Checks the domains are synchronized.
2610 if (context_->DomainOf(r_index.representative).Size() >
2611 context_->DomainOf(index_ref).Size()) {
2612 // Postpone, we will come back later when domains are synchronized.
2613 return true;
2614 }
2615 const int r_ref = r_index.representative;
2616 const int64 r_min = context_->MinOf(r_ref);
2617 const int64 r_max = context_->MaxOf(r_ref);
2618 const int array_size = ct->element().vars_size();
2619 if (r_min != 0) {
2620 context_->UpdateRuleStats("TODO element: representative has bad domain");
2621 } else if (r_index.offset >= 0 && r_index.offset < array_size &&
2622 r_index.offset + r_max * r_index.coeff >= 0 &&
2623 r_index.offset + r_max * r_index.coeff < array_size) {
2624 // This will happen eventually when domains are synchronized.
2625 ElementConstraintProto* const element =
2626 context_->working_model->add_constraints()->mutable_element();
2627 for (int64 v = 0; v <= r_max; ++v) {
2628 const int64 scaled_index = v * r_index.coeff + r_index.offset;
2629 CHECK_GE(scaled_index, 0);
2630 CHECK_LT(scaled_index, array_size);
2631 element->add_vars(ct->element().vars(scaled_index));
2632 }
2633 element->set_index(r_ref);
2634 element->set_target(target_ref);
2635
2636 if (r_index.coeff == 1) {
2637 context_->UpdateRuleStats("element: shifed index ");
2638 } else {
2639 context_->UpdateRuleStats("element: scaled index");
2640 }
2642 return RemoveConstraint(ct);
2643 }
2644 }
2645
2646 const bool unique_index = context_->VariableIsUniqueAndRemovable(index_ref) ||
2647 context_->IsFixed(index_ref);
2648 if (all_constants && unique_index) {
2649 // This constraint is just here to reduce the domain of the target! We can
2650 // add it to the mapping_model to reconstruct the index value during
2651 // postsolve and get rid of it now.
2652 context_->UpdateRuleStats("element: trivial target domain reduction");
2653 context_->MarkVariableAsRemoved(index_ref);
2654 *(context_->mapping_model->add_constraints()) = *ct;
2655 return RemoveConstraint(ct);
2656 }
2657
2658 const bool unique_target =
2659 context_->VariableIsUniqueAndRemovable(target_ref) ||
2660 context_->IsFixed(target_ref);
2661 if (all_included_in_target_domain && unique_target) {
2662 context_->UpdateRuleStats("element: trivial index domain reduction");
2663 context_->MarkVariableAsRemoved(target_ref);
2664 *(context_->mapping_model->add_constraints()) = *ct;
2665 return RemoveConstraint(ct);
2666 }
2667
2668 if (unique_target && !context_->IsFixed(target_ref)) {
2669 context_->UpdateRuleStats("TODO element: target not used elsewhere");
2670 }
2671 if (unique_index) {
2672 context_->UpdateRuleStats("TODO element: index not used elsewhere");
2673 }
2674
2675 return false;
2676}
2677
2678bool CpModelPresolver::PresolveTable(ConstraintProto* ct) {
2679 if (context_->ModelIsUnsat()) return false;
2680 if (HasEnforcementLiteral(*ct)) return false;
2681 if (ct->table().vars().empty()) {
2682 context_->UpdateRuleStats("table: empty constraint");
2683 return RemoveConstraint(ct);
2684 }
2685
2686 // Filter the unreachable tuples.
2687 //
2688 // TODO(user): this is not supper efficient. Optimize if needed.
2689 const int num_vars = ct->table().vars_size();
2690 const int num_tuples = ct->table().values_size() / num_vars;
2691 std::vector<int64> tuple(num_vars);
2692 std::vector<std::vector<int64>> new_tuples;
2693 new_tuples.reserve(num_tuples);
2694 std::vector<absl::flat_hash_set<int64>> new_domains(num_vars);
2695 std::vector<AffineRelation::Relation> affine_relations;
2696
2697 absl::flat_hash_set<int> visited;
2698 for (const int ref : ct->table().vars()) {
2699 if (visited.contains(PositiveRef(ref))) {
2700 context_->UpdateRuleStats("TODO table: duplicate variables");
2701 } else {
2702 visited.insert(PositiveRef(ref));
2703 }
2704 }
2705
2706 bool modified_variables = false;
2707 for (int v = 0; v < num_vars; ++v) {
2708 const int ref = ct->table().vars(v);
2709 AffineRelation::Relation r = context_->GetAffineRelation(ref);
2710 affine_relations.push_back(r);
2711 if (r.representative != ref) {
2712 modified_variables = true;
2713 }
2714 }
2715
2716 for (int i = 0; i < num_tuples; ++i) {
2717 bool delete_row = false;
2718 std::string tmp;
2719 for (int j = 0; j < num_vars; ++j) {
2720 const int ref = ct->table().vars(j);
2721 int64 v = ct->table().values(i * num_vars + j);
2722 const AffineRelation::Relation& r = affine_relations[j];
2723 if (r.representative != ref) {
2724 const int64 inverse_value = (v - r.offset) / r.coeff;
2725 if (inverse_value * r.coeff + r.offset != v) {
2726 // Bad rounding.
2727 delete_row = true;
2728 break;
2729 }
2730 v = inverse_value;
2731 }
2732 tuple[j] = v;
2733 if (!context_->DomainContains(r.representative, v)) {
2734 delete_row = true;
2735 break;
2736 }
2737 }
2738 if (delete_row) continue;
2739 new_tuples.push_back(tuple);
2740 for (int j = 0; j < num_vars; ++j) {
2741 const int64 v = tuple[j];
2742 new_domains[j].insert(v);
2743 }
2744 }
2746
2747 // Update the list of tuples if needed.
2748 if (new_tuples.size() < num_tuples || modified_variables) {
2749 ct->mutable_table()->clear_values();
2750 for (const std::vector<int64>& t : new_tuples) {
2751 for (const int64 v : t) {
2752 ct->mutable_table()->add_values(v);
2753 }
2754 }
2755 if (new_tuples.size() < num_tuples) {
2756 context_->UpdateRuleStats("table: removed rows");
2757 }
2758 }
2759
2760 if (modified_variables) {
2761 for (int j = 0; j < num_vars; ++j) {
2762 const AffineRelation::Relation& r = affine_relations[j];
2763 if (r.representative != ct->table().vars(j)) {
2764 ct->mutable_table()->set_vars(j, r.representative);
2765 }
2766 }
2767 context_->UpdateRuleStats(
2768 "table: replace variable by canonical affine one");
2769 }
2770
2771 // Nothing more to do for negated tables.
2772 if (ct->table().negated()) return modified_variables;
2773
2774 // Filter the variable domains.
2775 bool changed = false;
2776 for (int j = 0; j < num_vars; ++j) {
2777 const int ref = ct->table().vars(j);
2778 if (!context_->IntersectDomainWith(
2779 PositiveRef(ref),
2780 Domain::FromValues(std::vector<int64>(new_domains[j].begin(),
2781 new_domains[j].end())),
2782 &changed)) {
2783 return true;
2784 }
2785 }
2786 if (changed) {
2787 context_->UpdateRuleStats("table: reduced variable domains");
2788 }
2789 if (num_vars == 1) {
2790 // Now that we properly update the domain, we can remove the constraint.
2791 context_->UpdateRuleStats("table: only one column!");
2792 return RemoveConstraint(ct);
2793 }
2794
2795 // Check that the table is not complete or just here to exclude a few tuples.
2796 double prod = 1.0;
2797 for (int j = 0; j < num_vars; ++j) prod *= new_domains[j].size();
2798 if (prod == new_tuples.size()) {
2799 context_->UpdateRuleStats("table: all tuples!");
2800 return RemoveConstraint(ct);
2801 }
2802
2803 // Convert to the negated table if we gain a lot of entries by doing so.
2804 // Note however that currently the negated table do not propagate as much as
2805 // it could.
2806 if (new_tuples.size() > 0.7 * prod) {
2807 // Enumerate all tuples.
2808 std::vector<std::vector<int64>> var_to_values(num_vars);
2809 for (int j = 0; j < num_vars; ++j) {
2810 var_to_values[j].assign(new_domains[j].begin(), new_domains[j].end());
2811 }
2812 std::vector<std::vector<int64>> all_tuples(prod);
2813 for (int i = 0; i < prod; ++i) {
2814 all_tuples[i].resize(num_vars);
2815 int index = i;
2816 for (int j = 0; j < num_vars; ++j) {
2817 all_tuples[i][j] = var_to_values[j][index % var_to_values[j].size()];
2818 index /= var_to_values[j].size();
2819 }
2820 }
2822
2823 // Compute the complement of new_tuples.
2824 std::vector<std::vector<int64>> diff(prod - new_tuples.size());
2825 std::set_difference(all_tuples.begin(), all_tuples.end(),
2826 new_tuples.begin(), new_tuples.end(), diff.begin());
2827
2828 // Negate the constraint.
2829 ct->mutable_table()->set_negated(!ct->table().negated());
2830 ct->mutable_table()->clear_values();
2831 for (const std::vector<int64>& t : diff) {
2832 for (const int64 v : t) ct->mutable_table()->add_values(v);
2833 }
2834 context_->UpdateRuleStats("table: negated");
2835 }
2836 return modified_variables;
2837}
2838
2839bool CpModelPresolver::PresolveAllDiff(ConstraintProto* ct) {
2840 if (context_->ModelIsUnsat()) return false;
2841 if (HasEnforcementLiteral(*ct)) return false;
2842
2843 AllDifferentConstraintProto& all_diff = *ct->mutable_all_diff();
2844
2845 bool constraint_has_changed = false;
2846 for (;;) {
2847 const int size = all_diff.vars_size();
2848 if (size == 0) {
2849 context_->UpdateRuleStats("all_diff: empty constraint");
2850 return RemoveConstraint(ct);
2851 }
2852 if (size == 1) {
2853 context_->UpdateRuleStats("all_diff: only one variable");
2854 return RemoveConstraint(ct);
2855 }
2856
2857 bool something_was_propagated = false;
2858 std::vector<int> new_variables;
2859 for (int i = 0; i < size; ++i) {
2860 if (!context_->IsFixed(all_diff.vars(i))) {
2861 new_variables.push_back(all_diff.vars(i));
2862 continue;
2863 }
2864
2865 const int64 value = context_->MinOf(all_diff.vars(i));
2866 bool propagated = false;
2867 for (int j = 0; j < size; ++j) {
2868 if (i == j) continue;
2869 if (context_->DomainContains(all_diff.vars(j), value)) {
2870 if (!context_->IntersectDomainWith(all_diff.vars(j),
2871 Domain(value).Complement())) {
2872 return true;
2873 }
2874 propagated = true;
2875 }
2876 }
2877 if (propagated) {
2878 context_->UpdateRuleStats("all_diff: propagated fixed variables");
2879 something_was_propagated = true;
2880 }
2881 }
2882
2883 std::sort(new_variables.begin(), new_variables.end());
2884 for (int i = 1; i < new_variables.size(); ++i) {
2885 if (new_variables[i] == new_variables[i - 1]) {
2886 return context_->NotifyThatModelIsUnsat(
2887 "Duplicate variable in all_diff");
2888 }
2889 }
2890
2891 if (new_variables.size() < all_diff.vars_size()) {
2892 all_diff.mutable_vars()->Clear();
2893 for (const int var : new_variables) {
2894 all_diff.add_vars(var);
2895 }
2896 context_->UpdateRuleStats("all_diff: removed fixed variables");
2897 something_was_propagated = true;
2898 constraint_has_changed = true;
2899 if (new_variables.size() <= 1) continue;
2900 }
2901
2902 // Propagate mandatory value if the all diff is actually a permutation.
2903 CHECK_GE(all_diff.vars_size(), 2);
2904 Domain domain = context_->DomainOf(all_diff.vars(0));
2905 for (int i = 1; i < all_diff.vars_size(); ++i) {
2906 domain = domain.UnionWith(context_->DomainOf(all_diff.vars(i)));
2907 }
2908 if (all_diff.vars_size() == domain.Size()) {
2909 absl::flat_hash_map<int64, std::vector<int>> value_to_refs;
2910 for (const int ref : all_diff.vars()) {
2911 for (const ClosedInterval& interval : context_->DomainOf(ref)) {
2912 for (int64 v = interval.start; v <= interval.end; ++v) {
2913 value_to_refs[v].push_back(ref);
2914 }
2915 }
2916 }
2917 bool propagated = false;
2918 for (const auto& it : value_to_refs) {
2919 if (it.second.size() == 1 &&
2920 context_->DomainOf(it.second.front()).Size() > 1) {
2921 const int ref = it.second.front();
2922 if (!context_->IntersectDomainWith(ref, Domain(it.first))) {
2923 return true;
2924 }
2925 propagated = true;
2926 }
2927 }
2928 if (propagated) {
2929 context_->UpdateRuleStats(
2930 "all_diff: propagated mandatory values in permutation");
2931 something_was_propagated = true;
2932 }
2933 }
2934 if (!something_was_propagated) break;
2935 }
2936
2937 return constraint_has_changed;
2938}
2939
2940namespace {
2941
2942// Returns the sorted list of literals for given bool_or or at_most_one
2943// constraint.
2944std::vector<int> GetLiteralsFromSetPPCConstraint(const ConstraintProto& ct) {
2945 std::vector<int> sorted_literals;
2946 if (ct.constraint_case() == ConstraintProto::kAtMostOne) {
2947 for (const int literal : ct.at_most_one().literals()) {
2948 sorted_literals.push_back(literal);
2949 }
2950 } else if (ct.constraint_case() == ConstraintProto::kBoolOr) {
2951 for (const int literal : ct.bool_or().literals()) {
2952 sorted_literals.push_back(literal);
2953 }
2954 } else if (ct.constraint_case() == ConstraintProto::kExactlyOne) {
2955 for (const int literal : ct.exactly_one().literals()) {
2956 sorted_literals.push_back(literal);
2957 }
2958 }
2959 std::sort(sorted_literals.begin(), sorted_literals.end());
2960 return sorted_literals;
2961}
2962
2963// Add the constraint (lhs => rhs) to the given proto. The hash map lhs ->
2964// bool_and constraint index is used to merge implications with the same lhs.
2965void AddImplication(int lhs, int rhs, CpModelProto* proto,
2966 absl::flat_hash_map<int, int>* ref_to_bool_and) {
2967 if (ref_to_bool_and->contains(lhs)) {
2968 const int ct_index = (*ref_to_bool_and)[lhs];
2969 proto->mutable_constraints(ct_index)->mutable_bool_and()->add_literals(rhs);
2970 } else if (ref_to_bool_and->contains(NegatedRef(rhs))) {
2971 const int ct_index = (*ref_to_bool_and)[NegatedRef(rhs)];
2972 proto->mutable_constraints(ct_index)->mutable_bool_and()->add_literals(
2973 NegatedRef(lhs));
2974 } else {
2975 (*ref_to_bool_and)[lhs] = proto->constraints_size();
2976 ConstraintProto* ct = proto->add_constraints();
2977 ct->add_enforcement_literal(lhs);
2978 ct->mutable_bool_and()->add_literals(rhs);
2979 }
2980}
2981
2982template <typename ClauseContainer>
2983void ExtractClauses(bool use_bool_and, const ClauseContainer& container,
2984 CpModelProto* proto) {
2985 // We regroup the "implication" into bool_and to have a more consise proto and
2986 // also for nicer information about the number of binary clauses.
2987 //
2988 // Important: however, we do not do that for the model used during presolving
2989 // since the order of the constraints might be important there depending on
2990 // how we perform the postsolve.
2991 absl::flat_hash_map<int, int> ref_to_bool_and;
2992 for (int i = 0; i < container.NumClauses(); ++i) {
2993 const std::vector<Literal>& clause = container.Clause(i);
2994 if (clause.empty()) continue;
2995
2996 // bool_and.
2997 if (use_bool_and && clause.size() == 2) {
2998 const int a = clause[0].IsPositive()
2999 ? clause[0].Variable().value()
3000 : NegatedRef(clause[0].Variable().value());
3001 const int b = clause[1].IsPositive()
3002 ? clause[1].Variable().value()
3003 : NegatedRef(clause[1].Variable().value());
3004 AddImplication(NegatedRef(a), b, proto, &ref_to_bool_and);
3005 continue;
3006 }
3007
3008 // bool_or.
3009 ConstraintProto* ct = proto->add_constraints();
3010 for (const Literal l : clause) {
3011 if (l.IsPositive()) {
3012 ct->mutable_bool_or()->add_literals(l.Variable().value());
3013 } else {
3014 ct->mutable_bool_or()->add_literals(NegatedRef(l.Variable().value()));
3015 }
3016 }
3017 }
3018}
3019
3020} // namespace
3021
3022int64 CpModelPresolver::StartMin(
3023 const IntervalConstraintProto& interval) const {
3024 if (interval.has_start_view()) return context_->MinOf(interval.start_view());
3025 return context_->MinOf(interval.start());
3026}
3027
3028int64 CpModelPresolver::EndMax(const IntervalConstraintProto& interval) const {
3029 if (interval.has_end_view()) return context_->MaxOf(interval.end_view());
3030 return context_->MaxOf(interval.end());
3031}
3032
3033int64 CpModelPresolver::SizeMin(const IntervalConstraintProto& interval) const {
3034 if (interval.has_size_view()) return context_->MinOf(interval.size_view());
3035 return context_->MinOf(interval.size());
3036}
3037
3038int64 CpModelPresolver::SizeMax(const IntervalConstraintProto& interval) const {
3039 if (interval.has_size_view()) return context_->MaxOf(interval.size_view());
3040 return context_->MaxOf(interval.size());
3041}
3042
3043bool CpModelPresolver::PresolveNoOverlap(ConstraintProto* ct) {
3044 if (context_->ModelIsUnsat()) return false;
3045 const NoOverlapConstraintProto& proto = ct->no_overlap();
3046 const int initial_num_intervals = proto.intervals_size();
3047
3048 // Filter absent intervals.
3049 int new_size = 0;
3050 for (int i = 0; i < proto.intervals_size(); ++i) {
3051 const int interval_index = proto.intervals(i);
3052 if (context_->working_model->constraints(interval_index)
3053 .constraint_case() ==
3054 ConstraintProto::ConstraintCase::CONSTRAINT_NOT_SET) {
3055 continue;
3056 }
3057 ct->mutable_no_overlap()->set_intervals(new_size++, interval_index);
3058 }
3059 ct->mutable_no_overlap()->mutable_intervals()->Truncate(new_size);
3060
3061 // Sort by start min.
3062 std::sort(
3063 ct->mutable_no_overlap()->mutable_intervals()->begin(),
3064 ct->mutable_no_overlap()->mutable_intervals()->end(),
3065 [this](int i1, int i2) {
3066 return StartMin(context_->working_model->constraints(i1).interval()) <
3067 StartMin(context_->working_model->constraints(i2).interval());
3068 });
3069
3070 // Remove intervals that cannot overlap any others.
3071 //
3072 // TODO(user): We might also want to split this constraints into many
3073 // independent no overlap constraints.
3074 int64 end_max_so_far = kint64min;
3075 new_size = 0;
3076 for (int i = 0; i < proto.intervals_size(); ++i) {
3077 const int interval_index = proto.intervals(i);
3078 const IntervalConstraintProto& interval =
3079 context_->working_model->constraints(interval_index).interval();
3080 const int64 end_max_of_previous_intervals = end_max_so_far;
3081 end_max_so_far = std::max(end_max_so_far, EndMax(interval));
3082 if (StartMin(interval) >= end_max_of_previous_intervals &&
3083 (i + 1 == proto.intervals_size() ||
3084 end_max_so_far <= StartMin(context_->working_model
3085 ->constraints(proto.intervals(i + 1))
3086 .interval()))) {
3087 context_->UpdateRuleStats("no_overlap: removed redundant intervals");
3088 continue;
3089 }
3090 ct->mutable_no_overlap()->set_intervals(new_size++, interval_index);
3091 }
3092 ct->mutable_no_overlap()->mutable_intervals()->Truncate(new_size);
3093
3094 if (proto.intervals_size() == 1) {
3095 context_->UpdateRuleStats("no_overlap: only one interval");
3096 return RemoveConstraint(ct);
3097 }
3098 if (proto.intervals().empty()) {
3099 context_->UpdateRuleStats("no_overlap: no intervals");
3100 return RemoveConstraint(ct);
3101 }
3102
3103 return new_size < initial_num_intervals;
3104}
3105
3106bool CpModelPresolver::PresolveCumulative(ConstraintProto* ct) {
3107 if (context_->ModelIsUnsat()) return false;
3108
3109 const CumulativeConstraintProto& proto = ct->cumulative();
3110
3111 // Filter absent intervals.
3112 int new_size = 0;
3113 bool changed = false;
3114 int num_zero_demand_removed = 0;
3115 int64 sum_of_max_demands = 0;
3116 for (int i = 0; i < proto.intervals_size(); ++i) {
3117 if (context_->working_model->constraints(proto.intervals(i))
3118 .constraint_case() ==
3119 ConstraintProto::ConstraintCase::CONSTRAINT_NOT_SET) {
3120 continue;
3121 }
3122
3123 const int demand_ref = proto.demands(i);
3124 const int64 demand_max = context_->MaxOf(demand_ref);
3125 if (demand_max == 0) {
3126 num_zero_demand_removed++;
3127 continue;
3128 }
3129 sum_of_max_demands += demand_max;
3130
3131 ct->mutable_cumulative()->set_intervals(new_size, proto.intervals(i));
3132 ct->mutable_cumulative()->set_demands(new_size, proto.demands(i));
3133 new_size++;
3134 }
3135
3136 const int64 capacity_min = context_->MinOf(proto.capacity());
3137 const int64 capacity_max = context_->MaxOf(proto.capacity());
3138
3139 if (sum_of_max_demands <= capacity_min) {
3140 context_->UpdateRuleStats("cumulative: capacity exceeds sum of demands");
3141 return RemoveConstraint(ct);
3142 }
3143
3144 if (new_size < proto.intervals_size()) {
3145 changed = true;
3146 ct->mutable_cumulative()->mutable_intervals()->Truncate(new_size);
3147 ct->mutable_cumulative()->mutable_demands()->Truncate(new_size);
3148 }
3149
3150 if (num_zero_demand_removed > 0) {
3151 context_->UpdateRuleStats("cumulative: removed intervals with no demands");
3152 }
3153
3154 if (new_size == 0) {
3155 context_->UpdateRuleStats("cumulative: no intervals");
3156 return RemoveConstraint(ct);
3157 }
3158
3159 if (HasEnforcementLiteral(*ct)) return changed;
3160
3161 const int num_intervals = proto.intervals_size();
3162 bool with_start_view = false;
3163 std::vector<int> start_refs(num_intervals, -1);
3164
3165 int num_duration_one = 0;
3166 int num_greater_half_capacity = 0;
3167
3168 bool has_optional_interval = false;
3169 for (int i = 0; i < num_intervals; ++i) {
3170 // TODO(user): adapt in the presence of optional intervals.
3171 const ConstraintProto& ct =
3172 context_->working_model->constraints(proto.intervals(i));
3173 if (!ct.enforcement_literal().empty()) has_optional_interval = true;
3174 const IntervalConstraintProto& interval = ct.interval();
3175 if (interval.has_start_view()) with_start_view = true;
3176 start_refs[i] = interval.start();
3177 const int demand_ref = proto.demands(i);
3178 if (SizeMin(interval) == 1 && SizeMax(interval) == 1) {
3179 num_duration_one++;
3180 }
3181 if (SizeMin(interval) == 0) {
3182 // The behavior for zero-duration interval is currently not the same in
3183 // the no-overlap and the cumulative constraint.
3184 return changed;
3185 }
3186
3187 const int64 demand_min = context_->MinOf(demand_ref);
3188 const int64 demand_max = context_->MaxOf(demand_ref);
3189 if (demand_min > capacity_max / 2) {
3190 num_greater_half_capacity++;
3191 }
3192 if (demand_min > capacity_max) {
3193 context_->UpdateRuleStats("cumulative: demand_min exceeds capacity max");
3194 if (ct.enforcement_literal().empty()) {
3195 return context_->NotifyThatModelIsUnsat();
3196 } else {
3197 CHECK_EQ(ct.enforcement_literal().size(), 1);
3198 if (!context_->SetLiteralToFalse(ct.enforcement_literal(0))) {
3199 return true;
3200 }
3201 }
3202 return changed;
3203 } else if (demand_max > capacity_max) {
3204 if (ct.enforcement_literal().empty()) {
3205 context_->UpdateRuleStats(
3206 "cumulative: demand_max exceeds capacity max.");
3207 if (!context_->IntersectDomainWith(demand_ref,
3208 Domain(kint64min, capacity_max))) {
3209 return true;
3210 }
3211 } else {
3212 // TODO(user): we abort because we cannot convert this to a no_overlap
3213 // for instance.
3214 context_->UpdateRuleStats(
3215 "cumulative: demand_max of optional interval exceeds capacity.");
3216 return changed;
3217 }
3218 }
3219 }
3220
3221 if (num_greater_half_capacity == num_intervals) {
3222 if (num_duration_one == num_intervals && !has_optional_interval &&
3223 !with_start_view) {
3224 context_->UpdateRuleStats("cumulative: convert to all_different");
3225 ConstraintProto* new_ct = context_->working_model->add_constraints();
3226 auto* arg = new_ct->mutable_all_diff();
3227 for (const int var : start_refs) arg->add_vars(var);
3229 return RemoveConstraint(ct);
3230 } else {
3231 context_->UpdateRuleStats("cumulative: convert to no_overlap");
3232 ConstraintProto* new_ct = context_->working_model->add_constraints();
3233 auto* arg = new_ct->mutable_no_overlap();
3234 for (const int interval : proto.intervals()) {
3235 arg->add_intervals(interval);
3236 }
3238 return RemoveConstraint(ct);
3239 }
3240 }
3241
3242 return changed;
3243}
3244
3245bool CpModelPresolver::PresolveRoutes(ConstraintProto* ct) {
3246 if (context_->ModelIsUnsat()) return false;
3247 if (HasEnforcementLiteral(*ct)) return false;
3248 RoutesConstraintProto& proto = *ct->mutable_routes();
3249
3250 int new_size = 0;
3251 const int num_arcs = proto.literals_size();
3252 for (int i = 0; i < num_arcs; ++i) {
3253 const int ref = proto.literals(i);
3254 const int tail = proto.tails(i);
3255 const int head = proto.heads(i);
3256 if (context_->LiteralIsFalse(ref)) {
3257 context_->UpdateRuleStats("routes: removed false arcs");
3258 continue;
3259 }
3260 proto.set_literals(new_size, ref);
3261 proto.set_tails(new_size, tail);
3262 proto.set_heads(new_size, head);
3263 ++new_size;
3264 }
3265 if (new_size < num_arcs) {
3266 proto.mutable_literals()->Truncate(new_size);
3267 proto.mutable_tails()->Truncate(new_size);
3268 proto.mutable_heads()->Truncate(new_size);
3269 return true;
3270 }
3271 return false;
3272}
3273
3274bool CpModelPresolver::PresolveCircuit(ConstraintProto* ct) {
3275 if (context_->ModelIsUnsat()) return false;
3276 if (HasEnforcementLiteral(*ct)) return false;
3277 CircuitConstraintProto& proto = *ct->mutable_circuit();
3278
3279 // The indexing might not be dense, so fix that first.
3280 ReindexArcs(ct->mutable_circuit()->mutable_tails(),
3281 ct->mutable_circuit()->mutable_heads());
3282
3283 // Convert the flat structure to a graph, note that we includes all the arcs
3284 // here (even if they are at false).
3285 std::vector<std::vector<int>> incoming_arcs;
3286 std::vector<std::vector<int>> outgoing_arcs;
3287 int num_nodes = 0;
3288 const int num_arcs = proto.literals_size();
3289 for (int i = 0; i < num_arcs; ++i) {
3290 const int ref = proto.literals(i);
3291 const int tail = proto.tails(i);
3292 const int head = proto.heads(i);
3293 num_nodes = std::max(num_nodes, std::max(tail, head) + 1);
3294 if (std::max(tail, head) >= incoming_arcs.size()) {
3295 incoming_arcs.resize(std::max(tail, head) + 1);
3296 outgoing_arcs.resize(std::max(tail, head) + 1);
3297 }
3298 incoming_arcs[head].push_back(ref);
3299 outgoing_arcs[tail].push_back(ref);
3300 }
3301
3302 // Note that it is important to reach the fixed point here:
3303 // One arc at true, then all other arc at false. This is because we rely
3304 // on this in case the circuit is fully specified below.
3305 //
3306 // TODO(user): Use a better complexity if needed.
3307 bool loop_again = true;
3308 int num_fixed_at_true = 0;
3309 while (loop_again) {
3310 loop_again = false;
3311 for (const auto* node_to_refs : {&incoming_arcs, &outgoing_arcs}) {
3312 for (const std::vector<int>& refs : *node_to_refs) {
3313 if (refs.size() == 1) {
3314 if (!context_->LiteralIsTrue(refs.front())) {
3315 ++num_fixed_at_true;
3316 if (!context_->SetLiteralToTrue(refs.front())) return true;
3317 }
3318 continue;
3319 }
3320
3321 // At most one true, so if there is one, mark all the other to false.
3322 int num_true = 0;
3323 int true_ref;
3324 for (const int ref : refs) {
3325 if (context_->LiteralIsTrue(ref)) {
3326 ++num_true;
3327 true_ref = ref;
3328 break;
3329 }
3330 }
3331 if (num_true > 1) {
3332 return context_->NotifyThatModelIsUnsat();
3333 }
3334 if (num_true == 1) {
3335 for (const int ref : refs) {
3336 if (ref != true_ref) {
3337 if (!context_->IsFixed(ref)) {
3338 context_->UpdateRuleStats("circuit: set literal to false.");
3339 loop_again = true;
3340 }
3341 if (!context_->SetLiteralToFalse(ref)) return true;
3342 }
3343 }
3344 }
3345 }
3346 }
3347 }
3348 if (num_fixed_at_true > 0) {
3349 context_->UpdateRuleStats("circuit: fixed singleton arcs.");
3350 }
3351
3352 // Remove false arcs.
3353 int new_size = 0;
3354 int num_true = 0;
3355 int circuit_start = -1;
3356 std::vector<int> next(num_nodes, -1);
3357 std::vector<int> new_in_degree(num_nodes, 0);
3358 std::vector<int> new_out_degree(num_nodes, 0);
3359 for (int i = 0; i < num_arcs; ++i) {
3360 const int ref = proto.literals(i);
3361 if (context_->LiteralIsFalse(ref)) continue;
3362 if (context_->LiteralIsTrue(ref)) {
3363 if (next[proto.tails(i)] != -1) {
3364 return context_->NotifyThatModelIsUnsat();
3365 }
3366 next[proto.tails(i)] = proto.heads(i);
3367 if (proto.tails(i) != proto.heads(i)) {
3368 circuit_start = proto.tails(i);
3369 }
3370 ++num_true;
3371 }
3372 ++new_out_degree[proto.tails(i)];
3373 ++new_in_degree[proto.heads(i)];
3374 proto.set_tails(new_size, proto.tails(i));
3375 proto.set_heads(new_size, proto.heads(i));
3376 proto.set_literals(new_size, proto.literals(i));
3377 ++new_size;
3378 }
3379
3380 // Detect infeasibility due to a node having no more incoming or outgoing arc.
3381 // This is a bit tricky because for now the meaning of the constraint says
3382 // that all nodes that appear in at least one of the arcs must be in the
3383 // circuit or have a self-arc. So if any such node ends up with an incoming or
3384 // outgoing degree of zero once we remove false arcs then the constraint is
3385 // infeasible!
3386 for (int i = 0; i < num_nodes; ++i) {
3387 // Skip initially ignored node.
3388 if (incoming_arcs[i].empty() && outgoing_arcs[i].empty()) continue;
3389
3390 if (new_in_degree[i] == 0 || new_out_degree[i] == 0) {
3391 return context_->NotifyThatModelIsUnsat();
3392 }
3393 }
3394
3395 // Test if a subcircuit is already present.
3396 if (circuit_start != -1) {
3397 std::vector<bool> visited(num_nodes, false);
3398 int current = circuit_start;
3399 while (current != -1 && !visited[current]) {
3400 visited[current] = true;
3401 current = next[current];
3402 }
3403 if (current == circuit_start) {
3404 // We have a sub-circuit! mark all other arc false except self-loop not in
3405 // circuit.
3406 for (int i = 0; i < num_arcs; ++i) {
3407 if (visited[proto.tails(i)]) continue;
3408 if (proto.tails(i) == proto.heads(i)) {
3409 if (!context_->SetLiteralToTrue(proto.literals(i))) return true;
3410 } else {
3411 if (!context_->SetLiteralToFalse(proto.literals(i))) return true;
3412 }
3413 }
3414 context_->UpdateRuleStats("circuit: fully specified.");
3415 return RemoveConstraint(ct);
3416 }
3417 } else {
3418 // All self loop?
3419 if (num_true == new_size) {
3420 context_->UpdateRuleStats("circuit: empty circuit.");
3421 return RemoveConstraint(ct);
3422 }
3423 }
3424
3425 // Look for in/out-degree of two, this will imply that one of the indicator
3426 // Boolean is equal to the negation of the other.
3427 for (int i = 0; i < num_nodes; ++i) {
3428 for (const std::vector<int>* arc_literals :
3429 {&incoming_arcs[i], &outgoing_arcs[i]}) {
3430 std::vector<int> literals;
3431 for (const int ref : *arc_literals) {
3432 if (context_->LiteralIsFalse(ref)) continue;
3433 if (context_->LiteralIsTrue(ref)) {
3434 literals.clear();
3435 break;
3436 }
3437 literals.push_back(ref);
3438 }
3439 if (literals.size() == 2 && literals[0] != NegatedRef(literals[1])) {
3440 context_->UpdateRuleStats("circuit: degree 2");
3441 context_->StoreBooleanEqualityRelation(literals[0],
3442 NegatedRef(literals[1]));
3443 }
3444 }
3445 }
3446
3447 // Truncate the circuit and return.
3448 if (new_size < num_arcs) {
3449 proto.mutable_tails()->Truncate(new_size);
3450 proto.mutable_heads()->Truncate(new_size);
3451 proto.mutable_literals()->Truncate(new_size);
3452 context_->UpdateRuleStats("circuit: removed false arcs.");
3453 return true;
3454 }
3455 return false;
3456}
3457
3458bool CpModelPresolver::PresolveAutomaton(ConstraintProto* ct) {
3459 if (context_->ModelIsUnsat()) return false;
3460 if (HasEnforcementLiteral(*ct)) return false;
3461 AutomatonConstraintProto& proto = *ct->mutable_automaton();
3462 if (proto.vars_size() == 0 || proto.transition_label_size() == 0) {
3463 return false;
3464 }
3465
3466 bool all_affine = true;
3467 std::vector<AffineRelation::Relation> affine_relations;
3468 for (int v = 0; v < proto.vars_size(); ++v) {
3469 const int var = ct->automaton().vars(v);
3470 AffineRelation::Relation r = context_->GetAffineRelation(PositiveRef(var));
3471 affine_relations.push_back(r);
3472 if (r.representative == var) {
3473 all_affine = false;
3474 break;
3475 }
3476 if (v > 0 && (r.coeff != affine_relations[v - 1].coeff ||
3477 r.offset != affine_relations[v - 1].offset)) {
3478 all_affine = false;
3479 break;
3480 }
3481 }
3482
3483 if (all_affine) { // Unscale labels.
3484 for (int v = 0; v < proto.vars_size(); ++v) {
3485 proto.set_vars(v, affine_relations[v].representative);
3486 }
3487 const AffineRelation::Relation rep = affine_relations.front();
3488 int new_size = 0;
3489 for (int t = 0; t < proto.transition_tail_size(); ++t) {
3490 const int64 label = proto.transition_label(t);
3491 int64 inverse_label = (label - rep.offset) / rep.coeff;
3492 if (inverse_label * rep.coeff + rep.offset == label) {
3493 if (new_size != t) {
3494 proto.set_transition_tail(new_size, proto.transition_tail(t));
3495 proto.set_transition_head(new_size, proto.transition_head(t));
3496 }
3497 proto.set_transition_label(new_size, inverse_label);
3498 new_size++;
3499 }
3500 }
3501 if (new_size < proto.transition_tail_size()) {
3502 proto.mutable_transition_tail()->Truncate(new_size);
3503 proto.mutable_transition_label()->Truncate(new_size);
3504 proto.mutable_transition_head()->Truncate(new_size);
3505 context_->UpdateRuleStats("automaton: remove invalid transitions");
3506 }
3507 context_->UpdateRuleStats("automaton: unscale all affine labels");
3508 return true;
3509 }
3510
3511 Domain hull = context_->DomainOf(proto.vars(0));
3512 for (int v = 1; v < proto.vars_size(); ++v) {
3513 hull = hull.UnionWith(context_->DomainOf(proto.vars(v)));
3514 }
3515
3516 int new_size = 0;
3517 for (int t = 0; t < proto.transition_tail_size(); ++t) {
3518 const int64 label = proto.transition_label(t);
3519 if (hull.Contains(label)) {
3520 if (new_size != t) {
3521 proto.set_transition_tail(new_size, proto.transition_tail(t));
3522 proto.set_transition_label(new_size, label);
3523 proto.set_transition_head(new_size, proto.transition_head(t));
3524 }
3525 new_size++;
3526 }
3527 }
3528 if (new_size < proto.transition_tail_size()) {
3529 proto.mutable_transition_tail()->Truncate(new_size);
3530 proto.mutable_transition_label()->Truncate(new_size);
3531 proto.mutable_transition_head()->Truncate(new_size);
3532 context_->UpdateRuleStats("automaton: remove invalid transitions");
3533 return false;
3534 }
3535
3536 const int n = proto.vars_size();
3537 const std::vector<int> vars = {proto.vars().begin(), proto.vars().end()};
3538
3539 // Compute the set of reachable state at each time point.
3540 std::vector<std::set<int64>> reachable_states(n + 1);
3541 reachable_states[0].insert(proto.starting_state());
3542 reachable_states[n] = {proto.final_states().begin(),
3543 proto.final_states().end()};
3544
3545 // Forward.
3546 //
3547 // TODO(user): filter using the domain of vars[time] that may not contain
3548 // all the possible transitions.
3549 for (int time = 0; time + 1 < n; ++time) {
3550 for (int t = 0; t < proto.transition_tail_size(); ++t) {
3551 const int64 tail = proto.transition_tail(t);
3552 const int64 label = proto.transition_label(t);
3553 const int64 head = proto.transition_head(t);
3554 if (!gtl::ContainsKey(reachable_states[time], tail)) continue;
3555 if (!context_->DomainContains(vars[time], label)) continue;
3556 reachable_states[time + 1].insert(head);
3557 }
3558 }
3559
3560 std::vector<std::set<int64>> reached_values(n);
3561
3562 // Backward.
3563 for (int time = n - 1; time >= 0; --time) {
3564 std::set<int64> new_set;
3565 for (int t = 0; t < proto.transition_tail_size(); ++t) {
3566 const int64 tail = proto.transition_tail(t);
3567 const int64 label = proto.transition_label(t);
3568 const int64 head = proto.transition_head(t);
3569
3570 if (!gtl::ContainsKey(reachable_states[time], tail)) continue;
3571 if (!context_->DomainContains(vars[time], label)) continue;
3572 if (!gtl::ContainsKey(reachable_states[time + 1], head)) continue;
3573 new_set.insert(tail);
3574 reached_values[time].insert(label);
3575 }
3576 reachable_states[time].swap(new_set);
3577 }
3578
3579 bool removed_values = false;
3580 for (int time = 0; time < n; ++time) {
3581 if (!context_->IntersectDomainWith(
3582 vars[time],
3584 {reached_values[time].begin(), reached_values[time].end()}),
3585 &removed_values)) {
3586 return false;
3587 }
3588 }
3589 if (removed_values) {
3590 context_->UpdateRuleStats("automaton: reduced variable domains");
3591 }
3592 return false;
3593}
3594
3595bool CpModelPresolver::PresolveReservoir(ConstraintProto* ct) {
3596 if (context_->ModelIsUnsat()) return false;
3597 if (HasEnforcementLiteral(*ct)) return false;
3598
3599 bool changed = false;
3600
3601 ReservoirConstraintProto& mutable_reservoir = *ct->mutable_reservoir();
3602 if (mutable_reservoir.actives().empty()) {
3603 const int true_literal = context_->GetOrCreateConstantVar(1);
3604 for (int i = 0; i < mutable_reservoir.times_size(); ++i) {
3605 mutable_reservoir.add_actives(true_literal);
3606 }
3607 changed = true;
3608 }
3609
3610 const auto& demand_is_null = [&](int i) {
3611 return mutable_reservoir.demands(i) == 0 ||
3612 context_->LiteralIsFalse(mutable_reservoir.actives(i));
3613 };
3614
3615 // Remove zero demands, and inactive events.
3616 int num_zeros = 0;
3617 for (int i = 0; i < mutable_reservoir.demands_size(); ++i) {
3618 if (demand_is_null(i)) num_zeros++;
3619 }
3620
3621 if (num_zeros > 0) { // Remove null events
3622 changed = true;
3623 int new_size = 0;
3624 for (int i = 0; i < mutable_reservoir.demands_size(); ++i) {
3625 if (demand_is_null(i)) continue;
3626 mutable_reservoir.set_demands(new_size, mutable_reservoir.demands(i));
3627 mutable_reservoir.set_times(new_size, mutable_reservoir.times(i));
3628 mutable_reservoir.set_actives(new_size, mutable_reservoir.actives(i));
3629 new_size++;
3630 }
3631
3632 mutable_reservoir.mutable_demands()->Truncate(new_size);
3633 mutable_reservoir.mutable_times()->Truncate(new_size);
3634 mutable_reservoir.mutable_actives()->Truncate(new_size);
3635
3636 context_->UpdateRuleStats(
3637 "reservoir: remove zero demands or inactive events.");
3638 }
3639
3640 const int num_events = mutable_reservoir.demands_size();
3641 int64 gcd = mutable_reservoir.demands().empty()
3642 ? 0
3643 : std::abs(mutable_reservoir.demands(0));
3644 int num_positives = 0;
3645 int num_negatives = 0;
3646 int64 sum_of_demands = 0;
3647 int64 max_sum_of_positive_demands = 0;
3648 int64 min_sum_of_negative_demands = 0;
3649 for (int i = 0; i < num_events; ++i) {
3650 const int64 demand = mutable_reservoir.demands(i);
3651 sum_of_demands += demand;
3652 gcd = MathUtil::GCD64(gcd, std::abs(demand));
3653 if (demand > 0) {
3654 num_positives++;
3655 max_sum_of_positive_demands += demand;
3656 } else {
3657 DCHECK_LT(demand, 0);
3658 num_negatives++;
3659 min_sum_of_negative_demands += demand;
3660 }
3661 }
3662
3663 if (min_sum_of_negative_demands >= mutable_reservoir.min_level() &&
3664 max_sum_of_positive_demands <= mutable_reservoir.max_level()) {
3665 context_->UpdateRuleStats("reservoir: always feasible");
3666 return RemoveConstraint(ct);
3667 }
3668
3669 if (min_sum_of_negative_demands > mutable_reservoir.max_level() ||
3670 max_sum_of_positive_demands < mutable_reservoir.min_level()) {
3671 context_->UpdateRuleStats("reservoir: trivially infeasible");
3672 return context_->NotifyThatModelIsUnsat();
3673 }
3674
3675 if (min_sum_of_negative_demands > mutable_reservoir.min_level()) {
3676 mutable_reservoir.set_min_level(min_sum_of_negative_demands);
3677 context_->UpdateRuleStats(
3678 "reservoir: increase min_level to reachable value");
3679 }
3680
3681 if (max_sum_of_positive_demands < mutable_reservoir.max_level()) {
3682 mutable_reservoir.set_max_level(max_sum_of_positive_demands);
3683 context_->UpdateRuleStats("reservoir: reduce max_level to reachable value");
3684 }
3685
3686 if (mutable_reservoir.min_level() <= 0 &&
3687 mutable_reservoir.max_level() >= 0 &&
3688 (num_positives == 0 || num_negatives == 0)) {
3689 // If all demands have the same sign, and if the initial state is
3690 // always feasible, we do not care about the order, just the sum.
3691 auto* const sum =
3692 context_->working_model->add_constraints()->mutable_linear();
3693 int64 fixed_contrib = 0;
3694 for (int i = 0; i < mutable_reservoir.demands_size(); ++i) {
3695 const int64 demand = mutable_reservoir.demands(i);
3696 DCHECK_NE(demand, 0);
3697
3698 const int active = mutable_reservoir.actives(i);
3699 if (RefIsPositive(active)) {
3700 sum->add_vars(active);
3701 sum->add_coeffs(demand);
3702 } else {
3703 sum->add_vars(PositiveRef(active));
3704 sum->add_coeffs(-demand);
3705 fixed_contrib += demand;
3706 }
3707 }
3708 sum->add_domain(mutable_reservoir.min_level() - fixed_contrib);
3709 sum->add_domain(mutable_reservoir.max_level() - fixed_contrib);
3710 context_->UpdateRuleStats("reservoir: converted to linear");
3711 return RemoveConstraint(ct);
3712 }
3713
3714 if (gcd > 1) {
3715 for (int i = 0; i < mutable_reservoir.demands_size(); ++i) {
3716 mutable_reservoir.set_demands(i, mutable_reservoir.demands(i) / gcd);
3717 }
3718
3719 // Adjust min and max levels.
3720 // max level is always rounded down.
3721 // min level is always rounded up.
3722 const Domain reduced_domain =
3723 Domain({mutable_reservoir.min_level(), mutable_reservoir.max_level()})
3724 .InverseMultiplicationBy(gcd);
3725 mutable_reservoir.set_min_level(reduced_domain.Min());
3726 mutable_reservoir.set_max_level(reduced_domain.Max());
3727 context_->UpdateRuleStats("reservoir: simplify demands and levels by gcd.");
3728 }
3729
3730 if (num_positives == 1 && num_negatives > 0) {
3731 context_->UpdateRuleStats(
3732 "TODO reservoir: one producer, multiple consumers.");
3733 }
3734
3735 absl::flat_hash_set<std::pair<int, int>> time_active_set;
3736 for (int i = 0; i < mutable_reservoir.demands_size(); ++i) {
3737 const std::pair<int, int> key = std::make_pair(
3738 mutable_reservoir.times(i), mutable_reservoir.actives(i));
3739 if (time_active_set.contains(key)) {
3740 context_->UpdateRuleStats("TODO reservoir: merge synchronized events.");
3741 break;
3742 } else {
3743 time_active_set.insert(key);
3744 }
3745 }
3746
3747 return changed;
3748}
3749
3750// TODO(user): It is probably more efficient to keep all the bool_and in a
3751// global place during all the presolve, and just output them at the end
3752// rather than modifying more than once the proto.
3753void CpModelPresolver::ExtractBoolAnd() {
3754 absl::flat_hash_map<int, int> ref_to_bool_and;
3755 const int num_constraints = context_->working_model->constraints_size();
3756 std::vector<int> to_remove;
3757 for (int c = 0; c < num_constraints; ++c) {
3758 const ConstraintProto& ct = context_->working_model->constraints(c);
3759 if (HasEnforcementLiteral(ct)) continue;
3760
3761 if (ct.constraint_case() == ConstraintProto::ConstraintCase::kBoolOr &&
3762 ct.bool_or().literals().size() == 2) {
3763 AddImplication(NegatedRef(ct.bool_or().literals(0)),
3764 ct.bool_or().literals(1), context_->working_model,
3765 &ref_to_bool_and);
3766 to_remove.push_back(c);
3767 continue;
3768 }
3769
3770 if (ct.constraint_case() == ConstraintProto::ConstraintCase::kAtMostOne &&
3771 ct.at_most_one().literals().size() == 2) {
3772 AddImplication(ct.at_most_one().literals(0),
3773 NegatedRef(ct.at_most_one().literals(1)),
3774 context_->working_model, &ref_to_bool_and);
3775 to_remove.push_back(c);
3776 continue;
3777 }
3778 }
3779
3781 for (const int c : to_remove) {
3782 ConstraintProto* ct = context_->working_model->mutable_constraints(c);
3783 CHECK(RemoveConstraint(ct));
3784 context_->UpdateConstraintVariableUsage(c);
3785 }
3786}
3787
3788void CpModelPresolver::Probe() {
3789 if (context_->ModelIsUnsat()) return;
3790
3791 // Update the domain in the current CpModelProto.
3792 for (int i = 0; i < context_->working_model->variables_size(); ++i) {
3793 FillDomainInProto(context_->DomainOf(i),
3794 context_->working_model->mutable_variables(i));
3795 }
3796 const CpModelProto& model_proto = *(context_->working_model);
3797
3798 // Load the constraints in a local model.
3799 //
3800 // TODO(user): The model we load does not contain affine relations! But
3801 // ideally we should be able to remove all of them once we allow more complex
3802 // constraints to contains linear expression.
3803 //
3804 // TODO(user): remove code duplication with cp_model_solver. Here we also do
3805 // not run the heuristic to decide which variable to fully encode.
3806 //
3807 // TODO(user): Maybe do not load slow to propagate constraints? for instance
3808 // we do not use any linear relaxation here.
3809 Model model;
3810
3811 // Adapt some of the parameters during this probing phase.
3812 auto* local_param = model.GetOrCreate<SatParameters>();
3813 *local_param = context_->params();
3814 local_param->set_use_implied_bounds(false);
3815
3816 model.GetOrCreate<TimeLimit>()->MergeWithGlobalTimeLimit(
3817 context_->time_limit());
3818 model.Register<ModelRandomGenerator>(context_->random());
3819 auto* encoder = model.GetOrCreate<IntegerEncoder>();
3820 encoder->DisableImplicationBetweenLiteral();
3821 auto* mapping = model.GetOrCreate<CpModelMapping>();
3822 mapping->CreateVariables(model_proto, false, &model);
3823 mapping->DetectOptionalVariables(model_proto, &model);
3824 mapping->ExtractEncoding(model_proto, &model);
3825 auto* sat_solver = model.GetOrCreate<SatSolver>();
3826 for (const ConstraintProto& ct : model_proto.constraints()) {
3827 if (mapping->ConstraintIsAlreadyLoaded(&ct)) continue;
3829 if (sat_solver->IsModelUnsat()) {
3830 return (void)context_->NotifyThatModelIsUnsat();
3831 }
3832 }
3833 encoder->AddAllImplicationsBetweenAssociatedLiterals();
3834 if (!sat_solver->Propagate()) {
3835 return (void)context_->NotifyThatModelIsUnsat();
3836 }
3837
3838 // Probe.
3839 //
3840 // TODO(user): Compute the transitive reduction instead of just the
3841 // equivalences, and use the newly learned binary clauses?
3842 auto* implication_graph = model.GetOrCreate<BinaryImplicationGraph>();
3843 auto* prober = model.GetOrCreate<Prober>();
3844 prober->ProbeBooleanVariables(/*deterministic_time_limit=*/1.0);
3846 model.GetOrCreate<TimeLimit>()->GetElapsedDeterministicTime());
3847 if (sat_solver->IsModelUnsat() || !implication_graph->DetectEquivalences()) {
3848 return (void)context_->NotifyThatModelIsUnsat();
3849 }
3850
3851 // Update the presolve context with fixed Boolean variables.
3852 CHECK_EQ(sat_solver->CurrentDecisionLevel(), 0);
3853 for (int i = 0; i < sat_solver->LiteralTrail().Index(); ++i) {
3854 const Literal l = sat_solver->LiteralTrail()[i];
3855 const int var = mapping->GetProtoVariableFromBooleanVariable(l.Variable());
3856 if (var >= 0) {
3857 const int ref = l.IsPositive() ? var : NegatedRef(var);
3858 if (!context_->SetLiteralToTrue(ref)) return;
3859 }
3860 }
3861
3862 const int num_variables = context_->working_model->variables().size();
3863 auto* integer_trail = model.GetOrCreate<IntegerTrail>();
3864 for (int var = 0; var < num_variables; ++var) {
3865 // Restrict IntegerVariable domain.
3866 // Note that Boolean are already dealt with above.
3867 if (!mapping->IsBoolean(var)) {
3868 if (!context_->IntersectDomainWith(
3869 var,
3870 integer_trail->InitialVariableDomain(mapping->Integer(var)))) {
3871 return;
3872 }
3873 continue;
3874 }
3875
3876 // Add Boolean equivalence relations.
3877 const Literal l = mapping->Literal(var);
3878 const Literal r = implication_graph->RepresentativeOf(l);
3879 if (r != l) {
3880 const int r_var =
3881 mapping->GetProtoVariableFromBooleanVariable(r.Variable());
3882 CHECK_GE(r_var, 0);
3884 var, r.IsPositive() ? r_var : NegatedRef(r_var));
3885 }
3886 }
3887}
3888
3889// TODO(user): What to do with the at_most_one/exactly_one constraints?
3890// currently we do not take them into account here.
3891void CpModelPresolver::PresolvePureSatPart() {
3892 // TODO(user,user): Reenable some SAT presolve with
3893 // keep_all_feasible_solutions set to true.
3894 if (context_->ModelIsUnsat() || context_->keep_all_feasible_solutions) return;
3895
3896 const int num_variables = context_->working_model->variables_size();
3897 SatPostsolver sat_postsolver(num_variables);
3898 SatPresolver sat_presolver(&sat_postsolver);
3899 sat_presolver.SetNumVariables(num_variables);
3900 sat_presolver.SetTimeLimit(context_->time_limit());
3901
3902 SatParameters params = context_->params();
3903
3904 // The "full solver" postsolve does not support changing the value of a
3905 // variable from the solution of the presolved problem, and we do need this
3906 // for blocked clause. It should be possible to allow for this by adding extra
3907 // variable to the mapping model at presolve and some linking constraints, but
3908 // this is messy.
3909 if (params.cp_model_postsolve_with_full_solver()) {
3910 params.set_presolve_blocked_clause(false);
3911 }
3912
3913 // TODO(user): BVA takes times and do not seems to help on the minizinc
3914 // benchmarks. That said, it was useful on pure sat problems, so we may
3915 // want to enable it.
3916 params.set_presolve_use_bva(false);
3917 sat_presolver.SetParameters(params);
3918
3919 // Converts a cp_model literal ref to a sat::Literal used by SatPresolver.
3920 absl::flat_hash_set<int> used_variables;
3921 auto convert = [&used_variables](int ref) {
3922 used_variables.insert(PositiveRef(ref));
3923 if (RefIsPositive(ref)) return Literal(BooleanVariable(ref), true);
3924 return Literal(BooleanVariable(NegatedRef(ref)), false);
3925 };
3926
3927 // We need all Boolean constraints to be presolved before loading them below.
3928 // Otherwise duplicate literals might result in a wrong outcome.
3929 //
3930 // TODO(user): Be a bit more efficient, and enforce this invariant before we
3931 // reach this point?
3932 for (int c = 0; c < context_->working_model->constraints_size(); ++c) {
3933 const ConstraintProto& ct = context_->working_model->constraints(c);
3934 if (ct.constraint_case() == ConstraintProto::ConstraintCase::kBoolOr ||
3935 ct.constraint_case() == ConstraintProto::ConstraintCase::kBoolAnd) {
3936 if (PresolveOneConstraint(c)) {
3937 context_->UpdateConstraintVariableUsage(c);
3938 }
3939 if (context_->ModelIsUnsat()) return;
3940 }
3941 }
3942
3943 // Load all Clauses into the presolver and remove them from the current model.
3944 //
3945 // TODO(user): The removing and adding back of the same clause when nothing
3946 // happens in the presolve "seems" bad. That said, complexity wise, it is
3947 // a lot faster that what happens in the presolve though.
3948 //
3949 // TODO(user): Add the "small" at most one constraints to the SAT presolver by
3950 // expanding them to implications? that could remove a lot of clauses. Do that
3951 // when we are sure we don't load duplicates at_most_one/implications in the
3952 // solver. Ideally, the pure sat presolve could be improved to handle at most
3953 // one, and we could merge this with what the ProcessSetPPC() is doing.
3954 std::vector<Literal> clause;
3955 int num_removed_constraints = 0;
3956 for (int i = 0; i < context_->working_model->constraints_size(); ++i) {
3957 const ConstraintProto& ct = context_->working_model->constraints(i);
3958
3959 if (ct.constraint_case() == ConstraintProto::ConstraintCase::kBoolOr) {
3960 ++num_removed_constraints;
3961 clause.clear();
3962 for (const int ref : ct.bool_or().literals()) {
3963 clause.push_back(convert(ref));
3964 }
3965 for (const int ref : ct.enforcement_literal()) {
3966 clause.push_back(convert(ref).Negated());
3967 }
3968 sat_presolver.AddClause(clause);
3969
3970 context_->working_model->mutable_constraints(i)->Clear();
3971 context_->UpdateConstraintVariableUsage(i);
3972 continue;
3973 }
3974
3975 if (ct.constraint_case() == ConstraintProto::ConstraintCase::kBoolAnd) {
3976 ++num_removed_constraints;
3977 std::vector<Literal> clause;
3978 for (const int ref : ct.enforcement_literal()) {
3979 clause.push_back(convert(ref).Negated());
3980 }
3981 clause.push_back(Literal(kNoLiteralIndex)); // will be replaced below.
3982 for (const int ref : ct.bool_and().literals()) {
3983 clause.back() = convert(ref);
3984 sat_presolver.AddClause(clause);
3985 }
3986
3987 context_->working_model->mutable_constraints(i)->Clear();
3988 context_->UpdateConstraintVariableUsage(i);
3989 continue;
3990 }
3991 }
3992
3993 // Abort early if there was no Boolean constraints.
3994 if (num_removed_constraints == 0) return;
3995
3996 // Mark the variables appearing elsewhere or in the objective as non-removable
3997 // by the sat presolver.
3998 //
3999 // TODO(user): do not remove variable that appear in the decision heuristic?
4000 // TODO(user): We could go further for variable with only one polarity by
4001 // removing variable from the objective if they can be set to their "low"
4002 // objective value, and also removing enforcement literal that can be set to
4003 // false and don't appear elsewhere.
4004 std::vector<bool> can_be_removed(num_variables, false);
4005 for (int i = 0; i < num_variables; ++i) {
4006 if (context_->VarToConstraints(i).empty()) {
4007 can_be_removed[i] = true;
4008 }
4009
4010 // Because we might not have reached the presove "fixed point" above, some
4011 // variable in the added clauses might be fixed. We need to indicate this to
4012 // the SAT presolver.
4013 if (used_variables.contains(i) && context_->IsFixed(i)) {
4014 if (context_->LiteralIsTrue(i)) {
4015 sat_presolver.AddClause({convert(i)});
4016 } else {
4017 sat_presolver.AddClause({convert(NegatedRef(i))});
4018 }
4019 }
4020 }
4021
4022 // Run the presolve for a small number of passes.
4023 // TODO(user): Add probing like we do in the pure sat solver presolve loop?
4024 // TODO(user): Add a time limit, this can be slow on big SAT problem.
4025 const int num_passes = params.presolve_use_bva() ? 4 : 1;
4026 for (int i = 0; i < num_passes; ++i) {
4027 const int old_num_clause = sat_postsolver.NumClauses();
4028 if (!sat_presolver.Presolve(can_be_removed, context_->log_info())) {
4029 VLOG(1) << "UNSAT during SAT presolve.";
4030 return (void)context_->NotifyThatModelIsUnsat();
4031 }
4032 if (old_num_clause == sat_postsolver.NumClauses()) break;
4033 }
4034
4035 // Add any new variables to our internal structure.
4036 const int new_num_variables = sat_presolver.NumVariables();
4037 if (new_num_variables > context_->working_model->variables_size()) {
4038 VLOG(1) << "New variables added by the SAT presolver.";
4039 for (int i = context_->working_model->variables_size();
4040 i < new_num_variables; ++i) {
4041 IntegerVariableProto* var_proto =
4042 context_->working_model->add_variables();
4043 var_proto->add_domain(0);
4044 var_proto->add_domain(1);
4045 }
4046 context_->InitializeNewDomains();
4047 }
4048
4049 // Add the presolver clauses back into the model.
4050 ExtractClauses(/*use_bool_and=*/true, sat_presolver, context_->working_model);
4051
4052 // Update the constraints <-> variables graph.
4054
4055 // Add the sat_postsolver clauses to mapping_model.
4056 //
4057 // TODO(user): Mark removed variable as removed to detect any potential bugs.
4058 ExtractClauses(/*use_bool_and=*/false, sat_postsolver,
4059 context_->mapping_model);
4060}
4061
4062// TODO(user): The idea behind this was that it is better to have an objective
4063// as spreaded as possible. However on some problems this have the opposite
4064// effect. Like on a triangular matrix where each expansion reduced the size
4065// of the objective by one. Investigate and fix?
4066void CpModelPresolver::ExpandObjective() {
4067 if (context_->ModelIsUnsat()) return;
4068
4069 // The objective is already loaded in the constext, but we re-canonicalize
4070 // it with the latest information.
4071 if (!context_->CanonicalizeObjective()) {
4072 (void)context_->NotifyThatModelIsUnsat();
4073 return;
4074 }
4075
4076 if (context_->time_limit()->LimitReached()) {
4077 context_->WriteObjectiveToProto();
4078 return;
4079 }
4080
4081 // If the objective is a single variable, then we can usually remove this
4082 // variable if it is only used in one linear equality constraint and we do
4083 // just one expansion. This is because the domain of the variable will be
4084 // transferred to our objective_domain.
4085 int unique_expanded_constraint = -1;
4086 const bool objective_was_a_single_variable =
4087 context_->ObjectiveMap().size() == 1;
4088
4089 // To avoid a bad complexity, we need to compute the number of relevant
4090 // constraints for each variables.
4091 const int num_variables = context_->working_model->variables_size();
4092 const int num_constraints = context_->working_model->constraints_size();
4093 absl::flat_hash_set<int> relevant_constraints;
4094 std::vector<int> var_to_num_relevant_constraints(num_variables, 0);
4095 for (int ct_index = 0; ct_index < num_constraints; ++ct_index) {
4096 const ConstraintProto& ct = context_->working_model->constraints(ct_index);
4097 // Skip everything that is not a linear equality constraint.
4098 if (!ct.enforcement_literal().empty() ||
4099 ct.constraint_case() != ConstraintProto::ConstraintCase::kLinear ||
4100 ct.linear().domain().size() != 2 ||
4101 ct.linear().domain(0) != ct.linear().domain(1)) {
4102 continue;
4103 }
4104
4105 relevant_constraints.insert(ct_index);
4106 const int num_terms = ct.linear().vars_size();
4107 for (int i = 0; i < num_terms; ++i) {
4108 var_to_num_relevant_constraints[PositiveRef(ct.linear().vars(i))]++;
4109 }
4110 }
4111
4112 std::set<int> var_to_process;
4113 for (const auto entry : context_->ObjectiveMap()) {
4114 const int var = entry.first;
4116 if (var_to_num_relevant_constraints[var] != 0) {
4117 var_to_process.insert(var);
4118 }
4119 }
4120
4121 // We currently never expand a variable more than once.
4122 int num_expansions = 0;
4123 absl::flat_hash_set<int> processed_vars;
4124 std::vector<int> new_vars_in_objective;
4125 while (!relevant_constraints.empty()) {
4126 // Find a not yet expanded var.
4127 int objective_var = -1;
4128 while (!var_to_process.empty()) {
4129 const int var = *var_to_process.begin();
4130 CHECK(!processed_vars.contains(var));
4131 if (var_to_num_relevant_constraints[var] == 0) {
4132 processed_vars.insert(var);
4133 var_to_process.erase(var);
4134 continue;
4135 }
4136 if (!context_->ObjectiveMap().contains(var)) {
4137 // We do not mark it as processed in case it re-appear later.
4138 var_to_process.erase(var);
4139 continue;
4140 }
4141 objective_var = var;
4142 break;
4143 }
4144
4145 if (objective_var == -1) break;
4146 CHECK(RefIsPositive(objective_var));
4147 processed_vars.insert(objective_var);
4148 var_to_process.erase(objective_var);
4149
4150 int expanded_linear_index = -1;
4151 int64 objective_coeff_in_expanded_constraint;
4152 int64 size_of_expanded_constraint = 0;
4153 const auto& non_deterministic_list =
4154 context_->VarToConstraints(objective_var);
4155 std::vector<int> constraints_with_objective(non_deterministic_list.begin(),
4156 non_deterministic_list.end());
4157 std::sort(constraints_with_objective.begin(),
4158 constraints_with_objective.end());
4159 for (const int ct_index : constraints_with_objective) {
4160 if (relevant_constraints.count(ct_index) == 0) continue;
4161 const ConstraintProto& ct =
4162 context_->working_model->constraints(ct_index);
4163
4164 // This constraint is relevant now, but it will never be later because
4165 // it will contain the objective_var which is already processed!
4166 relevant_constraints.erase(ct_index);
4167 const int num_terms = ct.linear().vars_size();
4168 for (int i = 0; i < num_terms; ++i) {
4169 var_to_num_relevant_constraints[PositiveRef(ct.linear().vars(i))]--;
4170 }
4171
4172 // Find the coefficient of objective_var in this constraint, and perform
4173 // various checks.
4174 //
4175 // TODO(user): This can crash the program if for some reason the linear
4176 // constraint was not canonicalized and contains the objective variable
4177 // twice. Currently this can only happen if the time limit was reached
4178 // before all constraints where processed, but because we abort at the
4179 // beginning of the function when this is the case we should be safe.
4180 // However, it might be more robust to just handle this case properly.
4181 bool is_present = false;
4182 int64 objective_coeff;
4183 for (int i = 0; i < num_terms; ++i) {
4184 const int ref = ct.linear().vars(i);
4185 const int64 coeff = ct.linear().coeffs(i);
4186 if (PositiveRef(ref) == objective_var) {
4187 CHECK(!is_present) << "Duplicate variables not supported.";
4188 is_present = true;
4189 objective_coeff = (ref == objective_var) ? coeff : -coeff;
4190 } else {
4191 // This is not possible since we only consider relevant constraints.
4192 CHECK(!processed_vars.contains(PositiveRef(ref)));
4193 }
4194 }
4195 CHECK(is_present);
4196
4197 // We use the longest equality we can find.
4198 //
4199 // TODO(user): Deal with objective_coeff with a magnitude greater than
4200 // 1? This will only be possible if we change the objective coeff type
4201 // to double.
4202 if (std::abs(objective_coeff) == 1 &&
4203 num_terms > size_of_expanded_constraint) {
4204 expanded_linear_index = ct_index;
4205 size_of_expanded_constraint = num_terms;
4206 objective_coeff_in_expanded_constraint = objective_coeff;
4207 }
4208 }
4209
4210 if (expanded_linear_index != -1) {
4211 context_->UpdateRuleStats("objective: expanded objective constraint.");
4212
4213 // Update the objective map. Note that the division is possible because
4214 // currently we only expand with coeff with a magnitude of 1.
4215 CHECK_EQ(std::abs(objective_coeff_in_expanded_constraint), 1);
4216 const ConstraintProto& ct =
4217 context_->working_model->constraints(expanded_linear_index);
4219 objective_var, objective_coeff_in_expanded_constraint, ct,
4220 &new_vars_in_objective);
4221
4222 // Add not yet processed new variables.
4223 for (const int var : new_vars_in_objective) {
4224 if (!processed_vars.contains(var)) var_to_process.insert(var);
4225 }
4226
4227 // If the objective variable wasn't used in other constraints and it can
4228 // be reconstructed whatever the value of the other variables, we can
4229 // remove the constraint.
4230 //
4231 // TODO(user): It should be possible to refactor the code so this is
4232 // automatically done by the linear constraint singleton presolve rule.
4233 if (context_->VarToConstraints(objective_var).size() == 1 &&
4234 !context_->keep_all_feasible_solutions) {
4235 // Compute implied domain on objective_var.
4236 Domain implied_domain = ReadDomainFromProto(ct.linear());
4237 for (int i = 0; i < size_of_expanded_constraint; ++i) {
4238 const int ref = ct.linear().vars(i);
4239 if (PositiveRef(ref) == objective_var) continue;
4240 implied_domain =
4241 implied_domain
4242 .AdditionWith(context_->DomainOf(ref).MultiplicationBy(
4243 -ct.linear().coeffs(i)))
4244 .RelaxIfTooComplex();
4245 }
4246 implied_domain = implied_domain.InverseMultiplicationBy(
4247 objective_coeff_in_expanded_constraint);
4248
4249 // Remove the constraint if the implied domain is included in the
4250 // domain of the objective_var term.
4251 if (implied_domain.IsIncludedIn(context_->DomainOf(objective_var))) {
4252 context_->UpdateRuleStats("objective: removed objective constraint.");
4253 *(context_->mapping_model->add_constraints()) = ct;
4254 context_->working_model->mutable_constraints(expanded_linear_index)
4255 ->Clear();
4256 context_->UpdateConstraintVariableUsage(expanded_linear_index);
4257 } else {
4258 unique_expanded_constraint = expanded_linear_index;
4259 }
4260 }
4261
4262 ++num_expansions;
4263 }
4264 }
4265
4266 // Special case: If we just did one expansion of a single variable, then we
4267 // can remove the expanded constraints if the objective wasn't used elsewhere.
4268 if (num_expansions == 1 && objective_was_a_single_variable &&
4269 unique_expanded_constraint != -1) {
4270 context_->UpdateRuleStats(
4271 "objective: removed unique objective constraint.");
4272 ConstraintProto* mutable_ct = context_->working_model->mutable_constraints(
4273 unique_expanded_constraint);
4274 *(context_->mapping_model->add_constraints()) = *mutable_ct;
4275 mutable_ct->Clear();
4276 context_->UpdateConstraintVariableUsage(unique_expanded_constraint);
4277 }
4278
4279 // We re-do a canonicalization with the final linear expression.
4280 if (!context_->CanonicalizeObjective()) {
4281 (void)context_->NotifyThatModelIsUnsat();
4282 return;
4283 }
4284 context_->WriteObjectiveToProto();
4285}
4286
4287void CpModelPresolver::MergeNoOverlapConstraints() {
4288 if (context_->ModelIsUnsat()) return;
4289
4290 const int num_constraints = context_->working_model->constraints_size();
4291 int old_num_no_overlaps = 0;
4292 int old_num_intervals = 0;
4293
4294 // Extract the no-overlap constraints.
4295 std::vector<int> disjunctive_index;
4296 std::vector<std::vector<Literal>> cliques;
4297 for (int c = 0; c < num_constraints; ++c) {
4298 const ConstraintProto& ct = context_->working_model->constraints(c);
4299 if (ct.constraint_case() != ConstraintProto::ConstraintCase::kNoOverlap) {
4300 continue;
4301 }
4302 std::vector<Literal> clique;
4303 for (const int i : ct.no_overlap().intervals()) {
4304 clique.push_back(Literal(BooleanVariable(i), true));
4305 }
4306 cliques.push_back(clique);
4307 disjunctive_index.push_back(c);
4308
4309 old_num_no_overlaps++;
4310 old_num_intervals += clique.size();
4311 }
4312 if (old_num_no_overlaps == 0) return;
4313
4314 // We reuse the max-clique code from sat.
4315 Model local_model;
4316 local_model.GetOrCreate<Trail>()->Resize(num_constraints);
4317 auto* graph = local_model.GetOrCreate<BinaryImplicationGraph>();
4318 graph->Resize(num_constraints);
4319 for (const std::vector<Literal>& clique : cliques) {
4320 // All variables at false is always a valid solution of the local model,
4321 // so this should never return UNSAT.
4322 CHECK(graph->AddAtMostOne(clique));
4323 }
4324 CHECK(graph->DetectEquivalences());
4325 graph->TransformIntoMaxCliques(
4326 &cliques, context_->params().merge_no_overlap_work_limit());
4327
4328 // Replace each no-overlap with an extended version, or remove if empty.
4329 int new_num_no_overlaps = 0;
4330 int new_num_intervals = 0;
4331 for (int i = 0; i < cliques.size(); ++i) {
4332 const int ct_index = disjunctive_index[i];
4333 ConstraintProto* ct =
4334 context_->working_model->mutable_constraints(ct_index);
4335 ct->Clear();
4336 if (cliques[i].empty()) continue;
4337 for (const Literal l : cliques[i]) {
4338 CHECK(l.IsPositive());
4339 ct->mutable_no_overlap()->add_intervals(l.Variable().value());
4340 }
4341 new_num_no_overlaps++;
4342 new_num_intervals += cliques[i].size();
4343 }
4344 if (old_num_intervals != new_num_intervals ||
4345 old_num_no_overlaps != new_num_no_overlaps) {
4346 VLOG(1) << absl::StrCat("Merged ", old_num_no_overlaps, " no-overlaps (",
4347 old_num_intervals, " intervals) into ",
4348 new_num_no_overlaps, " no-overlaps (",
4349 new_num_intervals, " intervals).");
4350 context_->UpdateRuleStats("no_overlap: merged constraints");
4351 }
4352}
4353
4354// TODO(user): Should we take into account the exactly_one constraints? note
4355// that such constraint cannot be extended. If if a literal implies two literals
4356// at one inside an exactly one constraint then it must be false. Similarly if
4357// it implies all literals at zero inside the exactly one.
4358void CpModelPresolver::TransformIntoMaxCliques() {
4359 if (context_->ModelIsUnsat()) return;
4360
4361 auto convert = [](int ref) {
4362 if (RefIsPositive(ref)) return Literal(BooleanVariable(ref), true);
4363 return Literal(BooleanVariable(NegatedRef(ref)), false);
4364 };
4365 const int num_constraints = context_->working_model->constraints_size();
4366
4367 // Extract the bool_and and at_most_one constraints.
4368 std::vector<std::vector<Literal>> cliques;
4369
4370 for (int c = 0; c < num_constraints; ++c) {
4371 ConstraintProto* ct = context_->working_model->mutable_constraints(c);
4372 if (ct->constraint_case() == ConstraintProto::ConstraintCase::kAtMostOne) {
4373 std::vector<Literal> clique;
4374 for (const int ref : ct->at_most_one().literals()) {
4375 clique.push_back(convert(ref));
4376 }
4377 cliques.push_back(clique);
4378 if (RemoveConstraint(ct)) {
4379 context_->UpdateConstraintVariableUsage(c);
4380 }
4381 } else if (ct->constraint_case() ==
4382 ConstraintProto::ConstraintCase::kBoolAnd) {
4383 if (ct->enforcement_literal().size() != 1) continue;
4384 const Literal enforcement = convert(ct->enforcement_literal(0));
4385 for (const int ref : ct->bool_and().literals()) {
4386 if (ref == ct->enforcement_literal(0)) continue;
4387 cliques.push_back({enforcement, convert(ref).Negated()});
4388 }
4389 if (RemoveConstraint(ct)) {
4390 context_->UpdateConstraintVariableUsage(c);
4391 }
4392 }
4393 }
4394
4395 const int num_old_cliques = cliques.size();
4396
4397 // We reuse the max-clique code from sat.
4398 Model local_model;
4399 const int num_variables = context_->working_model->variables().size();
4400 local_model.GetOrCreate<Trail>()->Resize(num_variables);
4401 auto* graph = local_model.GetOrCreate<BinaryImplicationGraph>();
4402 graph->Resize(num_variables);
4403 for (const std::vector<Literal>& clique : cliques) {
4404 if (!graph->AddAtMostOne(clique)) {
4405 return (void)context_->NotifyThatModelIsUnsat();
4406 }
4407 }
4408 if (!graph->DetectEquivalences()) {
4409 return (void)context_->NotifyThatModelIsUnsat();
4410 }
4411 graph->TransformIntoMaxCliques(
4412 &cliques, context_->params().merge_at_most_one_work_limit());
4413
4414 // Add the Boolean variable equivalence detected by DetectEquivalences().
4415 // Those are needed because TransformIntoMaxCliques() will replace all
4416 // variable by its representative.
4417 for (int var = 0; var < num_variables; ++var) {
4418 const Literal l = Literal(BooleanVariable(var), true);
4419 if (graph->RepresentativeOf(l) != l) {
4420 const Literal r = graph->RepresentativeOf(l);
4422 var, r.IsPositive() ? r.Variable().value()
4423 : NegatedRef(r.Variable().value()));
4424 }
4425 }
4426
4427 int num_new_cliques = 0;
4428 for (const std::vector<Literal>& clique : cliques) {
4429 if (clique.empty()) continue;
4430 num_new_cliques++;
4431 ConstraintProto* ct = context_->working_model->add_constraints();
4432 for (const Literal literal : clique) {
4433 if (literal.IsPositive()) {
4434 ct->mutable_at_most_one()->add_literals(literal.Variable().value());
4435 } else {
4436 ct->mutable_at_most_one()->add_literals(
4437 NegatedRef(literal.Variable().value()));
4438 }
4439 }
4440
4441 // Make sure we do not have duplicate variable reference.
4442 PresolveAtMostOne(ct);
4443 }
4445 if (num_new_cliques != num_old_cliques) {
4446 context_->UpdateRuleStats("at_most_one: transformed into max clique.");
4447 }
4448
4449 if (context_->log_info()) {
4450 LOG(INFO) << "Merged " << num_old_cliques << " into " << num_new_cliques
4451 << " cliques.";
4452 }
4453}
4454
4456 if (context_->ModelIsUnsat()) return false;
4457 ConstraintProto* ct = context_->working_model->mutable_constraints(c);
4458
4459 // Generic presolve to exploit variable/literal equivalence.
4460 if (ExploitEquivalenceRelations(c, ct)) {
4461 context_->UpdateConstraintVariableUsage(c);
4462 }
4463
4464 // Generic presolve for reified constraint.
4465 if (PresolveEnforcementLiteral(ct)) {
4466 context_->UpdateConstraintVariableUsage(c);
4467 }
4468
4469 // Call the presolve function for this constraint if any.
4470 switch (ct->constraint_case()) {
4471 case ConstraintProto::ConstraintCase::kBoolOr:
4472 return PresolveBoolOr(ct);
4473 case ConstraintProto::ConstraintCase::kBoolAnd:
4474 return PresolveBoolAnd(ct);
4475 case ConstraintProto::ConstraintCase::kAtMostOne:
4476 return PresolveAtMostOne(ct);
4477 case ConstraintProto::ConstraintCase::kExactlyOne:
4478 return PresolveExactlyOne(ct);
4479 case ConstraintProto::ConstraintCase::kBoolXor:
4480 return PresolveBoolXor(ct);
4481 case ConstraintProto::ConstraintCase::kIntMax:
4482 if (ct->int_max().vars_size() == 2 &&
4483 NegatedRef(ct->int_max().vars(0)) == ct->int_max().vars(1)) {
4484 return PresolveIntAbs(ct);
4485 } else {
4486 return PresolveIntMax(ct);
4487 }
4488 case ConstraintProto::ConstraintCase::kIntMin:
4489 return PresolveIntMin(ct);
4490 case ConstraintProto::ConstraintCase::kLinMax:
4491 return PresolveLinMax(ct);
4492 case ConstraintProto::ConstraintCase::kLinMin:
4493 return PresolveLinMin(ct);
4494 case ConstraintProto::ConstraintCase::kIntProd:
4495 return PresolveIntProd(ct);
4496 case ConstraintProto::ConstraintCase::kIntDiv:
4497 return PresolveIntDiv(ct);
4498 case ConstraintProto::ConstraintCase::kLinear: {
4499 if (CanonicalizeLinear(ct)) {
4500 context_->UpdateConstraintVariableUsage(c);
4501 }
4502 if (PresolveSmallLinear(ct)) {
4503 context_->UpdateConstraintVariableUsage(c);
4504 }
4505 if (PropagateDomainsInLinear(c, ct)) {
4506 context_->UpdateConstraintVariableUsage(c);
4507 }
4508 if (PresolveSmallLinear(ct)) {
4509 context_->UpdateConstraintVariableUsage(c);
4510 }
4511 // We first propagate the domains before calling this presolve rule.
4512 if (RemoveSingletonInLinear(ct)) {
4513 context_->UpdateConstraintVariableUsage(c);
4514
4515 // There is no need to re-do a propagation here, but the constraint
4516 // size might have been reduced.
4517 if (PresolveSmallLinear(ct)) {
4518 context_->UpdateConstraintVariableUsage(c);
4519 }
4520 }
4521 if (PresolveLinearOnBooleans(ct)) {
4522 context_->UpdateConstraintVariableUsage(c);
4523 }
4524 if (ct->constraint_case() == ConstraintProto::ConstraintCase::kLinear) {
4525 const int old_num_enforcement_literals = ct->enforcement_literal_size();
4526 ExtractEnforcementLiteralFromLinearConstraint(c, ct);
4527 if (ct->constraint_case() ==
4528 ConstraintProto::ConstraintCase::CONSTRAINT_NOT_SET) {
4529 context_->UpdateConstraintVariableUsage(c);
4530 return true;
4531 }
4532 if (ct->enforcement_literal_size() > old_num_enforcement_literals &&
4533 PresolveSmallLinear(ct)) {
4534 context_->UpdateConstraintVariableUsage(c);
4535 }
4536 PresolveLinearEqualityModuloTwo(ct);
4537 }
4538 return false;
4539 }
4540 case ConstraintProto::ConstraintCase::kInterval:
4541 return PresolveInterval(c, ct);
4542 case ConstraintProto::ConstraintCase::kElement:
4543 return PresolveElement(ct);
4544 case ConstraintProto::ConstraintCase::kTable:
4545 return PresolveTable(ct);
4546 case ConstraintProto::ConstraintCase::kAllDiff:
4547 return PresolveAllDiff(ct);
4548 case ConstraintProto::ConstraintCase::kNoOverlap:
4549 return PresolveNoOverlap(ct);
4550 case ConstraintProto::ConstraintCase::kCumulative:
4551 return PresolveCumulative(ct);
4552 case ConstraintProto::ConstraintCase::kCircuit:
4553 return PresolveCircuit(ct);
4554 case ConstraintProto::ConstraintCase::kRoutes:
4555 return PresolveRoutes(ct);
4556 case ConstraintProto::ConstraintCase::kAutomaton:
4557 return PresolveAutomaton(ct);
4558 case ConstraintProto::ConstraintCase::kReservoir:
4559 return PresolveReservoir(ct);
4560 default:
4561 return false;
4562 }
4563}
4564
4565// We deal with all the all 3 x 3 possible types of inclusion.
4566//
4567// Returns false iff the model is UNSAT.
4568bool CpModelPresolver::ProcessSetPPCSubset(
4569 int c1, int c2, const std::vector<int>& c2_minus_c1,
4570 const std::vector<int>& original_constraint_index,
4571 std::vector<bool>* removed) {
4572 if (context_->ModelIsUnsat()) return false;
4573
4574 CHECK(!(*removed)[c1]);
4575 CHECK(!(*removed)[c2]);
4576
4577 ConstraintProto* ct1 = context_->working_model->mutable_constraints(
4578 original_constraint_index[c1]);
4579 ConstraintProto* ct2 = context_->working_model->mutable_constraints(
4580 original_constraint_index[c2]);
4581
4582 if ((ct1->constraint_case() == ConstraintProto::kBoolOr ||
4583 ct1->constraint_case() == ConstraintProto::kExactlyOne) &&
4584 (ct2->constraint_case() == ConstraintProto::kAtMostOne ||
4585 ct2->constraint_case() == ConstraintProto::kExactlyOne)) {
4586 context_->UpdateRuleStats("setppc: bool_or in at_most_one.");
4587
4588 // Fix extras in c2 to 0, note that these will be removed from the
4589 // constraint later.
4590 for (const int literal : c2_minus_c1) {
4591 if (!context_->SetLiteralToFalse(literal)) return false;
4592 context_->UpdateRuleStats("setppc: fixed variables");
4593 }
4594
4595 // Change c2 to exactly_one if not already.
4596 if (ct2->constraint_case() != ConstraintProto::kExactlyOne) {
4597 ConstraintProto copy = *ct2;
4598 (*ct2->mutable_exactly_one()->mutable_literals()) =
4599 copy.at_most_one().literals();
4600 }
4601
4602 // Remove c1.
4603 (*removed)[c1] = true;
4604 ct1->Clear();
4605 context_->UpdateConstraintVariableUsage(original_constraint_index[c1]);
4606 return true;
4607 }
4608
4609 if ((ct1->constraint_case() == ConstraintProto::kBoolOr ||
4610 ct1->constraint_case() == ConstraintProto::kExactlyOne) &&
4611 ct2->constraint_case() == ConstraintProto::kBoolOr) {
4612 context_->UpdateRuleStats("setppc: removed dominated constraints");
4613
4614 (*removed)[c2] = true;
4615 ct2->Clear();
4616 context_->UpdateConstraintVariableUsage(original_constraint_index[c2]);
4617 return true;
4618 }
4619
4620 if (ct1->constraint_case() == ConstraintProto::kAtMostOne &&
4621 (ct2->constraint_case() == ConstraintProto::kAtMostOne ||
4622 ct2->constraint_case() == ConstraintProto::kExactlyOne)) {
4623 context_->UpdateRuleStats("setppc: removed dominated constraints");
4624 (*removed)[c1] = true;
4625 ct1->Clear();
4626 context_->UpdateConstraintVariableUsage(original_constraint_index[c1]);
4627 return true;
4628 }
4629
4630 // We can't deduce anything in the last remaining case:
4631 // ct1->constraint_case() == ConstraintProto::kAtMostOne &&
4632 // ct2->constraint_case() == ConstraintProto::kBoolOr
4633
4634 return true;
4635}
4636
4637// TODO(user,user): TransformIntoMaxCliques() convert the bool_and to
4638// at_most_one, but maybe also duplicating them into bool_or would allow this
4639// function to do more presolving.
4640bool CpModelPresolver::ProcessSetPPC() {
4641 const int num_constraints = context_->working_model->constraints_size();
4642
4643 // Signatures of all the constraints. In the signature the bit i is 1 if it
4644 // contains a literal l such that l%64 = i.
4645 std::vector<uint64> signatures;
4646
4647 // Graph of constraints to literals. constraint_literals[c] contains all the
4648 // literals in constraint indexed by 'c' in sorted order.
4649 std::vector<std::vector<int>> constraint_literals;
4650
4651 // Graph of literals to constraints. literals_to_constraints[l] contains the
4652 // vector of constraint indices in which literal 'l' or 'neg(l)' appears.
4653 std::vector<std::vector<int>> literals_to_constraints;
4654
4655 // vector of booleans indicating if the constraint was already removed.
4656 std::vector<bool> removed;
4657
4658 // The containers above use the local indices for setppc constraints. We store
4659 // the original constraint indices corresponding to those local indices here.
4660 std::vector<int> original_constraint_index;
4661
4662 // Fill the initial constraint <-> literal graph, compute signatures and
4663 // initialize other containers defined above.
4664 int num_setppc_constraints = 0;
4665 for (int c = 0; c < num_constraints; ++c) {
4666 ConstraintProto* ct = context_->working_model->mutable_constraints(c);
4667 if (ct->constraint_case() == ConstraintProto::ConstraintCase::kBoolOr ||
4668 ct->constraint_case() == ConstraintProto::ConstraintCase::kAtMostOne ||
4669 ct->constraint_case() == ConstraintProto::ConstraintCase::kExactlyOne) {
4670 // Because TransformIntoMaxCliques() can detect literal equivalence
4671 // relation, we make sure the constraints are presolved before being
4672 // inspected.
4673 if (PresolveOneConstraint(c)) {
4674 context_->UpdateConstraintVariableUsage(c);
4675 }
4676 if (context_->ModelIsUnsat()) return false;
4677 }
4678 if (ct->constraint_case() == ConstraintProto::ConstraintCase::kBoolOr ||
4679 ct->constraint_case() == ConstraintProto::ConstraintCase::kAtMostOne ||
4680 ct->constraint_case() == ConstraintProto::ConstraintCase::kExactlyOne) {
4681 constraint_literals.push_back(GetLiteralsFromSetPPCConstraint(*ct));
4682
4683 uint64 signature = 0;
4684 for (const int literal : constraint_literals.back()) {
4685 const int positive_literal = PositiveRef(literal);
4686 signature |= (int64{1} << (positive_literal % 64));
4687 DCHECK_GE(positive_literal, 0);
4688 if (positive_literal >= literals_to_constraints.size()) {
4689 literals_to_constraints.resize(positive_literal + 1);
4690 }
4691 literals_to_constraints[positive_literal].push_back(
4692 num_setppc_constraints);
4693 }
4694 signatures.push_back(signature);
4695 removed.push_back(false);
4696 original_constraint_index.push_back(c);
4697 num_setppc_constraints++;
4698 }
4699 }
4700 VLOG(1) << "#setppc constraints: " << num_setppc_constraints;
4701
4702 // Set of constraint pairs which are already compared.
4703 absl::flat_hash_set<std::pair<int, int>> compared_constraints;
4704 for (const std::vector<int>& literal_to_constraints :
4705 literals_to_constraints) {
4706 for (int index1 = 0; index1 < literal_to_constraints.size(); ++index1) {
4707 if (context_->time_limit()->LimitReached()) return true;
4708
4709 const int c1 = literal_to_constraints[index1];
4710 if (removed[c1]) continue;
4711 const std::vector<int>& c1_literals = constraint_literals[c1];
4712
4713 for (int index2 = index1 + 1; index2 < literal_to_constraints.size();
4714 ++index2) {
4715 const int c2 = literal_to_constraints[index2];
4716 if (removed[c2]) continue;
4717 if (removed[c1]) break;
4718
4719 // TODO(user): This should not happen. Investigate.
4720 if (c1 == c2) continue;
4721
4722 CHECK_LT(c1, c2);
4723 if (gtl::ContainsKey(compared_constraints,
4724 std::pair<int, int>(c1, c2))) {
4725 continue;
4726 }
4727 compared_constraints.insert({c1, c2});
4728
4729 // Hard limit on number of comparisons to avoid spending too much time
4730 // here.
4731 if (compared_constraints.size() >= 50000) return true;
4732
4733 const bool smaller = (signatures[c1] & ~signatures[c2]) == 0;
4734 const bool larger = (signatures[c2] & ~signatures[c1]) == 0;
4735 if (!(smaller || larger)) continue;
4736
4737 // Check if literals in c1 is subset of literals in c2 or vice versa.
4738 const std::vector<int>& c2_literals = constraint_literals[c2];
4739
4740 // TODO(user): Try avoiding computation of set differences if
4741 // possible.
4742 std::vector<int> c1_minus_c2;
4743 gtl::STLSetDifference(c1_literals, c2_literals, &c1_minus_c2);
4744 std::vector<int> c2_minus_c1;
4745 gtl::STLSetDifference(c2_literals, c1_literals, &c2_minus_c1);
4746
4747 if (c1_minus_c2.empty()) { // c1 included in c2.
4748 if (!ProcessSetPPCSubset(c1, c2, c2_minus_c1,
4749 original_constraint_index, &removed)) {
4750 return false;
4751 }
4752 } else if (c2_minus_c1.empty()) { // c2 included in c1.
4753 if (!ProcessSetPPCSubset(c2, c1, c1_minus_c2,
4754 original_constraint_index, &removed)) {
4755 return false;
4756 }
4757 }
4758 }
4759 }
4760 }
4761
4762 return true;
4763}
4764
4765void CpModelPresolver::TryToSimplifyDomain(int var) {
4768 if (context_->ModelIsUnsat()) return;
4769 if (context_->IsFixed(var)) return;
4770 if (context_->VariableIsNotUsedAnymore(var)) return;
4771
4772 const AffineRelation::Relation r = context_->GetAffineRelation(var);
4773 if (r.representative != var) return;
4774
4775 if (context_->VariableIsOnlyUsedInEncoding(var)) {
4776 // TODO(user): We can remove such variable and its constraints by:
4777 // - Adding proper constraints on the enforcement literals. Simple case is
4778 // exactly one (or at most one) in case of a full (or partial) encoding.
4779 // - Moving all the old constraints to the mapping model.
4780 // - Updating the search heuristics/hint properly.
4781 context_->UpdateRuleStats("TODO variables: only used in encoding.");
4782 }
4783
4784 // Only process discrete domain.
4785 const Domain& domain = context_->DomainOf(var);
4786
4787 // Special case for non-Boolean domain of size 2.
4788 if (domain.Size() == 2 && (domain.Min() != 0 || domain.Max() != 1)) {
4790 return;
4791 }
4792
4793 if (domain.NumIntervals() != domain.Size()) return;
4794
4795 const int64 var_min = domain.Min();
4796 int64 gcd = domain[1].start - var_min;
4797 for (int index = 2; index < domain.NumIntervals(); ++index) {
4798 const ClosedInterval& i = domain[index];
4799 CHECK_EQ(i.start, i.end);
4800 const int64 shifted_value = i.start - var_min;
4801 CHECK_GE(shifted_value, 0);
4802
4803 gcd = MathUtil::GCD64(gcd, shifted_value);
4804 if (gcd == 1) break;
4805 }
4806 if (gcd == 1) return;
4807
4808 int new_var_index;
4809 {
4810 std::vector<int64> scaled_values;
4811 for (int index = 0; index < domain.NumIntervals(); ++index) {
4812 const ClosedInterval& i = domain[index];
4813 CHECK_EQ(i.start, i.end);
4814 const int64 shifted_value = i.start - var_min;
4815 scaled_values.push_back(shifted_value / gcd);
4816 }
4817 new_var_index = context_->NewIntVar(Domain::FromValues(scaled_values));
4818 }
4819 if (context_->ModelIsUnsat()) return;
4820
4821 CHECK(context_->StoreAffineRelation(var, new_var_index, gcd, var_min));
4822 context_->UpdateRuleStats("variables: canonicalize affine domain");
4824}
4825
4826// Adds all affine relations to our model for the variables that are still used.
4827void CpModelPresolver::EncodeAllAffineRelations() {
4828 int64 num_added = 0;
4829 for (int var = 0; var < context_->working_model->variables_size(); ++var) {
4830 if (context_->IsFixed(var)) continue;
4831
4832 const AffineRelation::Relation r = context_->GetAffineRelation(var);
4833 if (r.representative == var) continue;
4834
4835 if (!context_->keep_all_feasible_solutions) {
4836 // TODO(user): It seems some affine relation are still removable at this
4837 // stage even though they should be removed inside PresolveToFixPoint().
4838 // Investigate. For now, we just remove such relations.
4839 if (context_->VariableIsNotUsedAnymore(var)) continue;
4840 if (!PresolveAffineRelationIfAny(var)) break;
4841 if (context_->VariableIsNotUsedAnymore(var)) continue;
4842 if (context_->IsFixed(var)) continue;
4843 }
4844
4845 ++num_added;
4846 ConstraintProto* ct = context_->working_model->add_constraints();
4847 auto* arg = ct->mutable_linear();
4848 arg->add_vars(var);
4849 arg->add_coeffs(1);
4850 arg->add_vars(r.representative);
4851 arg->add_coeffs(-r.coeff);
4852 arg->add_domain(r.offset);
4853 arg->add_domain(r.offset);
4855 }
4856
4857 // Now that we encoded all remaining affine relation with constraints, we
4858 // remove the special marker to have a proper constraint variable graph.
4860
4861 if (context_->log_info() && num_added > 0) {
4862 LOG(INFO) << num_added << " affine relations still in the model.";
4863 }
4864}
4865
4866// Presolve a variable in relation with its representative.
4867bool CpModelPresolver::PresolveAffineRelationIfAny(int var) {
4868 if (context_->VariableIsNotUsedAnymore(var)) return true;
4869
4870 const AffineRelation::Relation r = context_->GetAffineRelation(var);
4871 if (r.representative == var) return true;
4872
4873 // Propagate domains.
4874 if (!context_->PropagateAffineRelation(var)) return false;
4875
4876 // Once an affine relation is detected, the variables should be added to
4877 // the kAffineRelationConstraint. The only way to be unmarked is if the
4878 // variable do not appear in any other constraint and is not a representative,
4879 // in which case it should never be added back.
4880 if (context_->IsFixed(var)) return true;
4882 CHECK(!context_->VariableIsNotUsedAnymore(r.representative));
4883
4884 // If var is no longer used, remove. Note that we can always do that since we
4885 // propagated the domain above and so we can find a feasible value for a for
4886 // any value of the representative.
4887 if (context_->VariableIsUniqueAndRemovable(var)) {
4888 // Add relation with current representative to the mapping model.
4889 ConstraintProto* ct = context_->mapping_model->add_constraints();
4890 auto* arg = ct->mutable_linear();
4891 arg->add_vars(var);
4892 arg->add_coeffs(1);
4893 arg->add_vars(r.representative);
4894 arg->add_coeffs(-r.coeff);
4895 arg->add_domain(r.offset);
4896 arg->add_domain(r.offset);
4898 }
4899 return true;
4900}
4901
4902void CpModelPresolver::PresolveToFixPoint() {
4903 if (context_->ModelIsUnsat()) return;
4904
4905 // Limit on number of operations.
4906 const int64 max_num_operations =
4907 context_->params().cp_model_max_num_presolve_operations() > 0
4908 ? context_->params().cp_model_max_num_presolve_operations()
4909 : kint64max;
4910
4911 // This is used for constraint having unique variables in them (i.e. not
4912 // appearing anywhere else) to not call the presolve more than once for this
4913 // reason.
4914 absl::flat_hash_set<std::pair<int, int>> var_constraint_pair_already_called;
4915
4916 TimeLimit* time_limit = context_->time_limit();
4917
4918 // The queue of "active" constraints, initialized to the non-empty ones.
4919 std::vector<bool> in_queue(context_->working_model->constraints_size(),
4920 false);
4921 std::deque<int> queue;
4922 for (int c = 0; c < in_queue.size(); ++c) {
4923 if (context_->working_model->constraints(c).constraint_case() !=
4924 ConstraintProto::ConstraintCase::CONSTRAINT_NOT_SET) {
4925 in_queue[c] = true;
4926 queue.push_back(c);
4927 }
4928 }
4929
4930 // When thinking about how the presolve works, it seems like a good idea to
4931 // process the "simple" constraints first in order to be more efficient.
4932 // In September 2019, experiment on the flatzinc problems shows no changes in
4933 // the results. We should actually count the number of rules triggered.
4934 if (context_->params().permute_presolve_constraint_order()) {
4935 std::shuffle(queue.begin(), queue.end(), *context_->random());
4936 } else {
4937 std::sort(queue.begin(), queue.end(), [this](int a, int b) {
4938 const int score_a = context_->ConstraintToVars(a).size();
4939 const int score_b = context_->ConstraintToVars(b).size();
4940 return score_a < score_b || (score_a == score_b && a < b);
4941 });
4942 }
4943
4944 while (!queue.empty() && !context_->ModelIsUnsat()) {
4945 if (time_limit->LimitReached()) break;
4946 if (context_->num_presolve_operations > max_num_operations) break;
4947 while (!queue.empty() && !context_->ModelIsUnsat()) {
4948 if (time_limit->LimitReached()) break;
4949 if (context_->num_presolve_operations > max_num_operations) break;
4950 const int c = queue.front();
4951 in_queue[c] = false;
4952 queue.pop_front();
4953
4954 const int old_num_constraint =
4955 context_->working_model->constraints_size();
4956 const bool changed = PresolveOneConstraint(c);
4957 if (context_->ModelIsUnsat() && context_->log_info()) {
4958 LOG(INFO) << "Unsat after presolving constraint #" << c
4959 << " (warning, dump might be inconsistent): "
4960 << context_->working_model->constraints(c).ShortDebugString();
4961 }
4962
4963 // Add to the queue any newly created constraints.
4964 const int new_num_constraints =
4965 context_->working_model->constraints_size();
4966 if (new_num_constraints > old_num_constraint) {
4968 in_queue.resize(new_num_constraints, true);
4969 for (int c = old_num_constraint; c < new_num_constraints; ++c) {
4970 queue.push_back(c);
4971 }
4972 }
4973
4974 // TODO(user): Is seems safer to simply remove the changed Boolean.
4975 // We loose a bit of performance, but the code is simpler.
4976 if (changed) {
4977 context_->UpdateConstraintVariableUsage(c);
4978 }
4979 }
4980
4981 // We also make sure all affine relations are propagated and any not
4982 // yet canonicalized domain is.
4983 //
4984 // TODO(user): maybe we can avoid iterating over all variables, but then
4985 // we already do that below.
4986 const int current_num_variables = context_->working_model->variables_size();
4987 for (int v = 0; v < current_num_variables; ++v) {
4988 if (context_->ModelIsUnsat()) return;
4989 if (!PresolveAffineRelationIfAny(v)) return;
4990
4991 // Try to canonicalize the domain, note that we should have detected all
4992 // affine relations before, so we don't recreate "canononical" variables
4993 // if they already exist in the model.
4994 TryToSimplifyDomain(v);
4996 }
4997
4998 // Re-add to the queue the constraints that touch a variable that changed.
4999 //
5000 // TODO(user): Avoid reprocessing the constraints that changed the variables
5001 // with the use of timestamp.
5002 if (context_->ModelIsUnsat()) return;
5003 in_queue.resize(context_->working_model->constraints_size(), false);
5004 for (const int v : context_->modified_domains.PositionsSetAtLeastOnce()) {
5005 if (context_->VariableIsNotUsedAnymore(v)) continue;
5006 if (context_->IsFixed(v)) context_->ExploitFixedDomain(v);
5007 for (const int c : context_->VarToConstraints(v)) {
5008 if (c >= 0 && !in_queue[c]) {
5009 in_queue[c] = true;
5010 queue.push_back(c);
5011 }
5012 }
5013 }
5014
5015 // Re-add to the queue constraints that have unique variables. Note that to
5016 // not enter an infinite loop, we call each (var, constraint) pair at most
5017 // once.
5018 const int num_vars = context_->working_model->variables_size();
5019 for (int v = 0; v < num_vars; ++v) {
5020 const auto& constraints = context_->VarToConstraints(v);
5021 if (constraints.size() != 1) continue;
5022 const int c = *constraints.begin();
5023 if (c < 0) continue;
5024
5025 // Note that to avoid bad complexity in problem like a TSP with just one
5026 // big constraint. we mark all the singleton variables of a constraint
5027 // even if this constraint is already in the queue.
5028 if (gtl::ContainsKey(var_constraint_pair_already_called,
5029 std::pair<int, int>(v, c))) {
5030 continue;
5031 }
5032 var_constraint_pair_already_called.insert({v, c});
5033
5034 if (!in_queue[c]) {
5035 in_queue[c] = true;
5036 queue.push_back(c);
5037 }
5038 }
5039
5040 // Make sure the order is deterministic! because var_to_constraints[]
5041 // order changes from one run to the next.
5042 std::sort(queue.begin(), queue.end());
5043 context_->modified_domains.SparseClearAll();
5044 }
5045
5046 if (context_->ModelIsUnsat()) return;
5047
5048 // Detect & exploit dominance between variables, or variables that can move
5049 // freely in one direction. Or variables that are just blocked by one
5050 // constraint in one direction.
5051 //
5052 // TODO(user): We can support assumptions but we need to not cut them out of
5053 // the feasible region.
5054 if (!context_->keep_all_feasible_solutions &&
5055 context_->working_model->assumptions().empty()) {
5056 VarDomination var_dom;
5057 DualBoundStrengthening dual_bound_strengthening;
5058 DetectDominanceRelations(*context_, &var_dom, &dual_bound_strengthening);
5059 if (!dual_bound_strengthening.Strengthen(context_)) return;
5060
5061 // TODO(user): The Strengthen() function above might make some inequality
5062 // tight. Currently, because we only do that for implication, this will not
5063 // change who dominate who, but in general we should process the new
5064 // constraint direction before calling this.
5065 if (!ExploitDominanceRelations(var_dom, context_)) return;
5066 }
5067
5068 // Second "pass" for transformation better done after all of the above and
5069 // that do not need a fix-point loop.
5070 //
5071 // TODO(user): Also add deductions achieved during probing!
5072 //
5073 // TODO(user): ideally we should "wake-up" any constraint that contains an
5074 // absent interval in the main propagation loop above. But we currently don't
5075 // maintain such list.
5076 const int num_constraints = context_->working_model->constraints_size();
5077 for (int c = 0; c < num_constraints; ++c) {
5078 ConstraintProto* ct = context_->working_model->mutable_constraints(c);
5079 switch (ct->constraint_case()) {
5080 case ConstraintProto::ConstraintCase::kNoOverlap:
5081 // Filter out absent intervals.
5082 if (PresolveNoOverlap(ct)) {
5083 context_->UpdateConstraintVariableUsage(c);
5084 }
5085 break;
5086 case ConstraintProto::ConstraintCase::kNoOverlap2D:
5087 // TODO(user): Implement if we ever support optional intervals in
5088 // this constraint. Currently we do not.
5089 break;
5090 case ConstraintProto::ConstraintCase::kCumulative:
5091 // Filter out absent intervals.
5092 if (PresolveCumulative(ct)) {
5093 context_->UpdateConstraintVariableUsage(c);
5094 }
5095 break;
5096 case ConstraintProto::ConstraintCase::kBoolOr: {
5097 // Try to infer domain reductions from clauses and the saved "implies in
5098 // domain" relations.
5099 for (const auto& pair :
5100 context_->deductions.ProcessClause(ct->bool_or().literals())) {
5101 bool modified = false;
5102 if (!context_->IntersectDomainWith(pair.first, pair.second,
5103 &modified)) {
5104 return;
5105 }
5106 if (modified) {
5107 context_->UpdateRuleStats("deductions: reduced variable domain");
5108 }
5109 }
5110 break;
5111 }
5112 default:
5113 break;
5114 }
5115 }
5116
5118}
5119
5121 LOG(INFO) << "- " << context->NumAffineRelations()
5122 << " affine relations were detected.";
5123 LOG(INFO) << "- " << context->NumEquivRelations()
5124 << " variable equivalence relations were detected.";
5125 std::map<std::string, int> sorted_rules(context->stats_by_rule_name.begin(),
5126 context->stats_by_rule_name.end());
5127 for (const auto& entry : sorted_rules) {
5128 if (entry.second == 1) {
5129 LOG(INFO) << "- rule '" << entry.first << "' was applied 1 time.";
5130 } else {
5131 LOG(INFO) << "- rule '" << entry.first << "' was applied " << entry.second
5132 << " times.";
5133 }
5134 }
5135}
5136
5137// =============================================================================
5138// Public API.
5139// =============================================================================
5140
5142 std::vector<int>* postsolve_mapping) {
5143 CpModelPresolver presolver(context, postsolve_mapping);
5144 return presolver.Presolve();
5145}
5146
5148 std::vector<int>* postsolve_mapping)
5149 : postsolve_mapping_(postsolve_mapping), context_(context) {
5150 // TODO(user): move in the context.
5151 context_->keep_all_feasible_solutions =
5152 context_->params().keep_all_feasible_solutions_in_presolve() ||
5153 context_->params().enumerate_all_solutions() ||
5154 context_->params().fill_tightened_domains_in_response() ||
5155 !context_->params().cp_model_presolve();
5156
5157 // We copy the search strategy to the mapping_model.
5158 for (const auto& decision_strategy :
5159 context_->working_model->search_strategy()) {
5160 *(context_->mapping_model->add_search_strategy()) = decision_strategy;
5161 }
5162
5163 // Initialize the initial context.working_model domains.
5164 context_->InitializeNewDomains();
5165
5166 // Initialize the objective.
5167 context_->ReadObjectiveFromProto();
5168 if (!context_->CanonicalizeObjective()) {
5169 (void)context_->NotifyThatModelIsUnsat();
5170 }
5171
5172 // Note that we delay the call to UpdateNewConstraintsVariableUsage() for
5173 // efficiency during LNS presolve.
5174}
5175
5176// The presolve works as follow:
5177//
5178// First stage:
5179// We will process all active constraints until a fix point is reached. During
5180// this stage:
5181// - Variable will never be deleted, but their domain will be reduced.
5182// - Constraint will never be deleted (they will be marked as empty if needed).
5183// - New variables and new constraints can be added after the existing ones.
5184// - Constraints are added only when needed to the mapping_problem if they are
5185// needed during the postsolve.
5186//
5187// Second stage:
5188// - All the variables domain will be copied to the mapping_model.
5189// - Everything will be remapped so that only the variables appearing in some
5190// constraints will be kept and their index will be in [0, num_new_variables).
5192 context_->enable_stats = context_->log_info();
5193
5194 // If presolve is false, just run expansion.
5195 if (!context_->params().cp_model_presolve()) {
5197 ExpandCpModel(context_);
5198 if (context_->log_info()) LogInfoFromContext(context_);
5199 return true;
5200 }
5201
5202 // Before initializing the constraint <-> variable graph (which is costly), we
5203 // run a bunch of simple presolve rules. Note that these function should NOT
5204 // use the graph, or the part that uses it should properly check for
5205 // context_->ConstraintVariableGraphIsUpToDate() before doing anything that
5206 // depends on the graph.
5207 for (int c = 0; c < context_->working_model->constraints_size(); ++c) {
5208 ConstraintProto* ct = context_->working_model->mutable_constraints(c);
5209 PresolveEnforcementLiteral(ct);
5210 switch (ct->constraint_case()) {
5211 case ConstraintProto::ConstraintCase::kBoolOr:
5212 PresolveBoolOr(ct);
5213 break;
5214 case ConstraintProto::ConstraintCase::kBoolAnd:
5215 PresolveBoolAnd(ct);
5216 break;
5217 case ConstraintProto::ConstraintCase::kAtMostOne:
5218 PresolveAtMostOne(ct);
5219 break;
5220 case ConstraintProto::ConstraintCase::kExactlyOne:
5221 PresolveExactlyOne(ct);
5222 break;
5223 case ConstraintProto::ConstraintCase::kLinear:
5224 CanonicalizeLinear(ct);
5225 break;
5226 default:
5227 break;
5228 }
5229 if (context_->ModelIsUnsat()) break;
5230 }
5234
5235 // Main propagation loop.
5236 for (int iter = 0; iter < context_->params().max_presolve_iterations();
5237 ++iter) {
5238 context_->UpdateRuleStats("presolve: iteration");
5239 // Save some quantities to decide if we abort early the iteration loop.
5240 const int64 old_num_presolve_op = context_->num_presolve_operations;
5241 int old_num_non_empty_constraints = 0;
5242 for (int c = 0; c < context_->working_model->constraints_size(); ++c) {
5243 const auto type =
5244 context_->working_model->constraints(c).constraint_case();
5245 if (type == ConstraintProto::ConstraintCase::CONSTRAINT_NOT_SET) continue;
5246 old_num_non_empty_constraints++;
5247 }
5248
5249 // TODO(user): The presolve transformations we do after this is called might
5250 // result in even more presolve if we were to call this again! improve the
5251 // code. See for instance plusexample_6_sat.fzn were represolving the
5252 // presolved problem reduces it even more.
5253 PresolveToFixPoint();
5254
5255 // Call expansion.
5256 ExpandCpModel(context_);
5258
5259 // TODO(user): do that and the pure-SAT part below more than once.
5260 if (context_->params().cp_model_probing_level() > 0) {
5261 if (!context_->time_limit()->LimitReached()) {
5262 Probe();
5263 PresolveToFixPoint();
5264 }
5265 }
5266
5267 // Runs SAT specific presolve on the pure-SAT part of the problem.
5268 // Note that because this can only remove/fix variable not used in the other
5269 // part of the problem, there is no need to redo more presolve afterwards.
5270 if (context_->params().cp_model_use_sat_presolve()) {
5271 if (!context_->time_limit()->LimitReached()) {
5272 PresolvePureSatPart();
5273 }
5274 }
5275
5276 // Extract redundant at most one constraint form the linear ones.
5277 //
5278 // TODO(user): more generally if we do some probing, the same relation will
5279 // be detected (and more). Also add an option to turn this off?
5280 //
5281 // TODO(user): instead of extracting at most one, extract pairwise conflicts
5282 // and add them to bool_and clauses? this is some sort of small scale
5283 // probing, but good for sat presolve and clique later?
5284 if (!context_->ModelIsUnsat() && iter == 0) {
5285 const int old_size = context_->working_model->constraints_size();
5286 for (int c = 0; c < old_size; ++c) {
5287 ConstraintProto* ct = context_->working_model->mutable_constraints(c);
5288 if (ct->constraint_case() != ConstraintProto::ConstraintCase::kLinear) {
5289 continue;
5290 }
5291 ExtractAtMostOneFromLinear(ct);
5292 }
5294 }
5295
5296 if (iter == 0) TransformIntoMaxCliques();
5297
5298 // TODO(user): Decide where is the best place for this. Fow now we do it
5299 // after max clique to get all the bool_and converted to at most ones.
5300 if (context_->params().symmetry_level() > 0 && !context_->ModelIsUnsat() &&
5301 !context_->time_limit()->LimitReached() &&
5302 !context_->keep_all_feasible_solutions) {
5304 }
5305
5306 // Process set packing, partitioning and covering constraint.
5307 if (!context_->time_limit()->LimitReached()) {
5308 ProcessSetPPC();
5309 ExtractBoolAnd();
5310
5311 // Call the main presolve to remove the fixed variables and do more
5312 // deductions.
5313 PresolveToFixPoint();
5314 }
5315
5316 // Exit the loop if the reduction is not so large.
5317 if (context_->num_presolve_operations - old_num_presolve_op <
5318 0.8 * (context_->working_model->variables_size() +
5319 old_num_non_empty_constraints)) {
5320 break;
5321 }
5322 }
5323
5324 // Regroup no-overlaps into max-cliques.
5325 if (!context_->ModelIsUnsat()) {
5326 MergeNoOverlapConstraints();
5327 }
5328
5329 // Tries to spread the objective amongst many variables.
5330 if (context_->working_model->has_objective() && !context_->ModelIsUnsat()) {
5331 ExpandObjective();
5332 }
5333
5334 // Adds all needed affine relation to context_->working_model.
5335 if (!context_->ModelIsUnsat()) {
5336 EncodeAllAffineRelations();
5337 }
5338
5339 // Remove duplicate constraints.
5340 //
5341 // TODO(user): We might want to do that earlier so that our count of variable
5342 // usage is not biased by duplicate constraints.
5343 const std::vector<int> duplicates =
5345 if (!duplicates.empty()) {
5346 for (const int c : duplicates) {
5347 const int type =
5348 context_->working_model->constraints(c).constraint_case();
5349 if (type == ConstraintProto::ConstraintCase::kInterval) {
5350 // TODO(user): we could delete duplicate identical interval, but we need
5351 // to make sure reference to them are updated.
5352 continue;
5353 }
5354
5355 context_->working_model->mutable_constraints(c)->Clear();
5356 context_->UpdateConstraintVariableUsage(c);
5357 context_->UpdateRuleStats("removed duplicate constraints");
5358 }
5359 }
5360
5361 if (context_->ModelIsUnsat()) {
5362 if (context_->log_info()) LogInfoFromContext(context_);
5363
5364 // Set presolved_model to the simplest UNSAT problem (empty clause).
5365 context_->working_model->Clear();
5366 context_->working_model->add_constraints()->mutable_bool_or();
5367 return true;
5368 }
5369
5370 // The strategy variable indices will be remapped in ApplyVariableMapping()
5371 // but first we use the representative of the affine relations for the
5372 // variables that are not present anymore.
5373 //
5374 // Note that we properly take into account the sign of the coefficient which
5375 // will result in the same domain reduction strategy. Moreover, if the
5376 // variable order is not CHOOSE_FIRST, then we also encode the associated
5377 // affine transformation in order to preserve the order.
5378 absl::flat_hash_set<int> used_variables;
5379 for (DecisionStrategyProto& strategy :
5380 *context_->working_model->mutable_search_strategy()) {
5381 DecisionStrategyProto copy = strategy;
5382 strategy.clear_variables();
5383 for (const int ref : copy.variables()) {
5384 const int var = PositiveRef(ref);
5385
5386 // Remove fixed variables.
5387 if (context_->IsFixed(var)) continue;
5388
5389 // There is not point having a variable appear twice, so we only keep
5390 // the first occurrence in the first strategy in which it occurs.
5391 if (gtl::ContainsKey(used_variables, var)) continue;
5392 used_variables.insert(var);
5393
5394 if (context_->VarToConstraints(var).empty()) {
5395 const AffineRelation::Relation r = context_->GetAffineRelation(var);
5396 if (!context_->VarToConstraints(r.representative).empty()) {
5397 const int rep = (r.coeff > 0) == RefIsPositive(ref)
5398 ? r.representative
5400 strategy.add_variables(rep);
5401 if (strategy.variable_selection_strategy() !=
5402 DecisionStrategyProto::CHOOSE_FIRST) {
5403 DecisionStrategyProto::AffineTransformation* t =
5404 strategy.add_transformations();
5405 t->set_var(rep);
5406 t->set_offset(r.offset);
5407 t->set_positive_coeff(std::abs(r.coeff));
5408 }
5409 } else {
5410 // TODO(user): this variable was removed entirely by the presolve (no
5411 // equivalent variable present). We simply ignore it entirely which
5412 // might result in a different search...
5413 }
5414 } else {
5415 strategy.add_variables(ref);
5416 }
5417 }
5418 }
5419
5420 // Sync the domains.
5421 for (int i = 0; i < context_->working_model->variables_size(); ++i) {
5422 FillDomainInProto(context_->DomainOf(i),
5423 context_->working_model->mutable_variables(i));
5424 DCHECK_GT(context_->working_model->variables(i).domain_size(), 0);
5425 }
5426
5427 // Set the variables of the mapping_model.
5428 context_->mapping_model->mutable_variables()->CopyFrom(
5429 context_->working_model->variables());
5430
5431 // Remove all the unused variables from the presolved model.
5432 postsolve_mapping_->clear();
5433 std::vector<int> mapping(context_->working_model->variables_size(), -1);
5434 for (int i = 0; i < context_->working_model->variables_size(); ++i) {
5435 if (context_->VariableIsNotUsedAnymore(i) &&
5436 !context_->keep_all_feasible_solutions) {
5437 continue;
5438 }
5439 mapping[i] = postsolve_mapping_->size();
5440 postsolve_mapping_->push_back(i);
5441 }
5442
5443 if (context_->params().permute_variable_randomly()) {
5444 std::shuffle(postsolve_mapping_->begin(), postsolve_mapping_->end(),
5445 *context_->random());
5446 for (int i = 0; i < postsolve_mapping_->size(); ++i) {
5447 mapping[(*postsolve_mapping_)[i]] = i;
5448 }
5449 }
5450
5452 ApplyVariableMapping(mapping, *context_);
5453
5454 // Compact all non-empty constraint at the beginning.
5456
5457 // Hack to display the number of deductions stored.
5458 if (context_->deductions.NumDeductions() > 0) {
5459 context_->UpdateRuleStats(absl::StrCat(
5460 "deductions: ", context_->deductions.NumDeductions(), " stored"));
5461 }
5462
5463 // Stats and checks.
5464 if (context_->log_info()) LogInfoFromContext(context_);
5465
5466 // One possible error that is difficult to avoid here: because of our
5467 // objective expansion, we might detect a possible overflow...
5468 //
5469 // TODO(user): We could abort the expansion when this happen.
5470 {
5471 const std::string error = ValidateCpModel(*context_->working_model);
5472 if (!error.empty()) {
5473 if (context_->log_info()) {
5474 LOG(INFO) << "Error while validating postsolved model: " << error;
5475 }
5476 return false;
5477 }
5478 }
5479 {
5480 const std::string error = ValidateCpModel(*context_->mapping_model);
5481 if (!error.empty()) {
5482 if (context_->log_info()) {
5483 LOG(INFO) << "Error while validating mapping_model model: " << error;
5484 }
5485 return false;
5486 }
5487 }
5488 return true;
5489}
5490
5491void ApplyVariableMapping(const std::vector<int>& mapping,
5492 const PresolveContext& context) {
5493 CpModelProto* proto = context.working_model;
5494
5495 // Remap all the variable/literal references in the constraints and the
5496 // enforcement literals in the variables.
5497 auto mapping_function = [&mapping](int* ref) {
5498 const int image = mapping[PositiveRef(*ref)];
5499 CHECK_GE(image, 0);
5500 *ref = RefIsPositive(*ref) ? image : NegatedRef(image);
5501 };
5502 for (ConstraintProto& ct_ref : *proto->mutable_constraints()) {
5503 ApplyToAllVariableIndices(mapping_function, &ct_ref);
5504 ApplyToAllLiteralIndices(mapping_function, &ct_ref);
5505 }
5506
5507 // Remap the objective variables.
5508 if (proto->has_objective()) {
5509 for (int& mutable_ref : *proto->mutable_objective()->mutable_vars()) {
5510 mapping_function(&mutable_ref);
5511 }
5512 }
5513
5514 // Remap the assumptions.
5515 for (int& mutable_ref : *proto->mutable_assumptions()) {
5516 mapping_function(&mutable_ref);
5517 }
5518
5519 // Remap the search decision heuristic.
5520 // Note that we delete any heuristic related to a removed variable.
5521 for (DecisionStrategyProto& strategy : *proto->mutable_search_strategy()) {
5522 DecisionStrategyProto copy = strategy;
5523 strategy.clear_variables();
5524 for (const int ref : copy.variables()) {
5525 const int image = mapping[PositiveRef(ref)];
5526 if (image >= 0) {
5527 strategy.add_variables(RefIsPositive(ref) ? image : NegatedRef(image));
5528 }
5529 }
5530 strategy.clear_transformations();
5531 for (const auto& transform : copy.transformations()) {
5532 const int ref = transform.var();
5533 const int image = mapping[PositiveRef(ref)];
5534 if (image >= 0) {
5535 auto* new_transform = strategy.add_transformations();
5536 *new_transform = transform;
5537 new_transform->set_var(RefIsPositive(ref) ? image : NegatedRef(image));
5538 }
5539 }
5540 }
5541
5542 // Remap the solution hint.
5543 if (proto->has_solution_hint()) {
5544 auto* mutable_hint = proto->mutable_solution_hint();
5545 int new_size = 0;
5546 for (int i = 0; i < mutable_hint->vars_size(); ++i) {
5547 const int old_ref = mutable_hint->vars(i);
5548 const int64 old_value = mutable_hint->values(i);
5549
5550 // Note that if (old_value - r.offset) is not divisible by r.coeff, then
5551 // the hint is clearly infeasible, but we still set it to a "close" value.
5552 const AffineRelation::Relation r = context.GetAffineRelation(old_ref);
5553 const int var = r.representative;
5554 const int64 value = (old_value - r.offset) / r.coeff;
5555
5556 const int image = mapping[var];
5557 if (image >= 0) {
5558 mutable_hint->set_vars(new_size, image);
5559 mutable_hint->set_values(new_size, value);
5560 ++new_size;
5561 }
5562 }
5563 if (new_size > 0) {
5564 mutable_hint->mutable_vars()->Truncate(new_size);
5565 mutable_hint->mutable_values()->Truncate(new_size);
5566 } else {
5567 proto->clear_solution_hint();
5568 }
5569 }
5570
5571 // Move the variable definitions.
5572 std::vector<IntegerVariableProto> new_variables;
5573 for (int i = 0; i < mapping.size(); ++i) {
5574 const int image = mapping[i];
5575 if (image < 0) continue;
5576 if (image >= new_variables.size()) {
5577 new_variables.resize(image + 1, IntegerVariableProto());
5578 }
5579 new_variables[image].Swap(proto->mutable_variables(i));
5580 }
5581 proto->clear_variables();
5582 for (IntegerVariableProto& proto_ref : new_variables) {
5583 proto->add_variables()->Swap(&proto_ref);
5584 }
5585
5586 // Check that all variables are used.
5587 for (const IntegerVariableProto& v : proto->variables()) {
5588 CHECK_GT(v.domain_size(), 0);
5589 }
5590}
5591
5592std::vector<int> FindDuplicateConstraints(const CpModelProto& model_proto) {
5593 std::vector<int> result;
5594
5595 // We use a map hash: serialized_constraint_proto hash -> constraint index.
5596 ConstraintProto copy;
5597 absl::flat_hash_map<int64, int> equiv_constraints;
5598
5599 std::string s;
5600 const int num_constraints = model_proto.constraints().size();
5601 for (int c = 0; c < num_constraints; ++c) {
5602 if (model_proto.constraints(c).constraint_case() ==
5603 ConstraintProto::ConstraintCase::CONSTRAINT_NOT_SET) {
5604 continue;
5605 }
5606
5607 // We ignore names when comparing constraints.
5608 //
5609 // TODO(user): This is not particularly efficient.
5610 copy = model_proto.constraints(c);
5611 copy.clear_name();
5612 s = copy.SerializeAsString();
5613
5614 const int64 hash = std::hash<std::string>()(s);
5615 const auto insert = equiv_constraints.insert({hash, c});
5616 if (!insert.second) {
5617 // Already present!
5618 const int other_c_with_same_hash = insert.first->second;
5619 copy = model_proto.constraints(other_c_with_same_hash);
5620 copy.clear_name();
5621 if (s == copy.SerializeAsString()) {
5622 result.push_back(c);
5623 }
5624 }
5625 }
5626
5627 return result;
5628}
5629
5630} // namespace sat
5631} // namespace operations_research
int64 min
Definition: alldiff_cst.cc:138
int64 max
Definition: alldiff_cst.cc:139
#define CHECK(condition)
Definition: base/logging.h:495
#define DCHECK_NE(val1, val2)
Definition: base/logging.h:886
#define CHECK_LT(val1, val2)
Definition: base/logging.h:700
#define CHECK_EQ(val1, val2)
Definition: base/logging.h:697
#define CHECK_GE(val1, val2)
Definition: base/logging.h:701
#define CHECK_GT(val1, val2)
Definition: base/logging.h:702
#define DCHECK_GE(val1, val2)
Definition: base/logging.h:889
#define CHECK_NE(val1, val2)
Definition: base/logging.h:698
#define DCHECK_GT(val1, val2)
Definition: base/logging.h:890
#define DCHECK_LT(val1, val2)
Definition: base/logging.h:888
#define LOG(severity)
Definition: base/logging.h:420
#define DCHECK(condition)
Definition: base/logging.h:884
#define CHECK_LE(val1, val2)
Definition: base/logging.h:699
#define DCHECK_EQ(val1, val2)
Definition: base/logging.h:885
#define VLOG(verboselevel)
Definition: base/logging.h:978
bool IsIncludedIn(const Domain &domain) const
Returns true iff D is included in the given domain.
Domain MultiplicationBy(int64 coeff, bool *exact=nullptr) const
Returns {x ∈ Int64, ∃ e ∈ D, x = e * coeff}.
static Domain FromValues(std::vector< int64 > values)
Creates a domain from the union of an unsorted list of integer values.
int64 Size() const
Returns the number of elements in the domain.
Domain InverseMultiplicationBy(const int64 coeff) const
Returns {x ∈ Int64, ∃ e ∈ D, x * coeff = e}.
Domain AdditionWith(const Domain &domain) const
Returns {x ∈ Int64, ∃ a ∈ D, ∃ b ∈ domain, x = a + b}.
ClosedInterval front() const
Domain UnionWith(const Domain &domain) const
Returns the union of D and domain.
Domain IntersectionWith(const Domain &domain) const
Returns the intersection of D and domain.
Domain RelaxIfTooComplex() const
If NumIntervals() is too large, this return a superset of the domain.
Domain DivisionBy(int64 coeff) const
Returns {x ∈ Int64, ∃ e ∈ D, x = e / coeff}.
static int64 GCD64(int64 x, int64 y)
Definition: mathutil.h:107
void Set(IntegerType index)
Definition: bitset.h:803
const std::vector< IntegerType > & PositionsSetAtLeastOnce() const
Definition: bitset.h:813
bool LimitReached()
Returns true when the external limit is true, or the deterministic time is over the deterministic lim...
Definition: time_limit.h:532
void AdvanceDeterministicTime(double deterministic_duration)
Advances the deterministic time.
Definition: time_limit.h:226
CpModelPresolver(PresolveContext *context, std::vector< int > *postsolve_mapping)
std::vector< std::pair< int, Domain > > ProcessClause(absl::Span< const int > clause)
void AddDeduction(int literal_ref, int var, Domain domain)
bool StoreAbsRelation(int target_ref, int ref)
const absl::flat_hash_set< int > & VarToConstraints(int var) const
ABSL_MUST_USE_RESULT bool IntersectDomainWith(int ref, const Domain &domain, bool *domain_modified=nullptr)
std::vector< absl::flat_hash_set< int > > var_to_lb_only_constraints
int GetOrCreateVarValueEncoding(int ref, int64 value)
bool DomainContains(int ref, int64 value) const
bool StoreLiteralImpliesVarNEqValue(int literal, int var, int64 value)
bool DomainOfVarIsIncludedIn(int var, const Domain &domain)
bool VariableWithCostIsUniqueAndRemovable(int ref) const
const absl::flat_hash_map< int, int64 > & ObjectiveMap() const
ABSL_MUST_USE_RESULT bool SetLiteralToTrue(int lit)
bool StoreLiteralImpliesVarEqValue(int literal, int var, int64 value)
std::vector< absl::flat_hash_set< int > > var_to_ub_only_constraints
ABSL_MUST_USE_RESULT bool NotifyThatModelIsUnsat(const std::string &message="")
const std::vector< int > & ConstraintToVars(int c) const
void StoreBooleanEqualityRelation(int ref_a, int ref_b)
bool SubstituteVariableInObjective(int var_in_equality, int64 coeff_in_equality, const ConstraintProto &equality, std::vector< int > *new_vars_in_objective=nullptr)
void UpdateRuleStats(const std::string &name, int num_times=1)
ABSL_MUST_USE_RESULT bool CanonicalizeObjective()
AffineRelation::Relation GetAffineRelation(int ref) const
ABSL_MUST_USE_RESULT bool SetLiteralToFalse(int lit)
const SatParameters & params() const
absl::flat_hash_set< int > tmp_literal_set
bool GetAbsRelation(int target_ref, int *ref)
bool StoreAffineRelation(int ref_x, int ref_y, int64 coeff, int64 offset)
Block * next
CpModelProto proto
CpModelProto const * model_proto
SharedTimeLimit * time_limit
const std::string name
const Constraint * ct
int64 value
IntVar * var
Definition: expr_array.cc:1858
GRBmodel * model
GurobiMPCallbackContext * context
static const int64 kint64max
int64_t int64
static const int32 kint32min
uint64_t uint64
static const int64 kint64min
const int WARNING
Definition: log_severity.h:31
const int INFO
Definition: log_severity.h:31
int64 hash
Definition: matrix_utils.cc:60
const Collection::value_type::second_type & FindOrDie(const Collection &collection, const typename Collection::value_type::first_type &key)
Definition: map_util.h:176
bool ContainsKey(const Collection &collection, const Key &key)
Definition: map_util.h:170
void STLSetDifference(const In1 &a, const In2 &b, Out *out, Compare compare)
Definition: stl_util.h:595
void STLSortAndRemoveDuplicates(T *v, const LessFunc &less_func)
Definition: stl_util.h:58
int ReindexArcs(IntContainer *tails, IntContainer *heads)
Definition: circuit.h:168
bool HasEnforcementLiteral(const ConstraintProto &ct)
bool DetectAndExploitSymmetriesInPresolve(PresolveContext *context)
bool LoadConstraint(const ConstraintProto &ct, Model *m)
void SetToNegatedLinearExpression(const LinearExpressionProto &input_expr, LinearExpressionProto *output_negated_expr)
Domain ReadDomainFromProto(const ProtoWithDomain &proto)
const LiteralIndex kNoLiteralIndex(-1)
void ExpandCpModel(PresolveContext *context)
constexpr int kAffineRelationConstraint
bool PresolveCpModel(PresolveContext *context, std::vector< int > *postsolve_mapping)
void ApplyToAllLiteralIndices(const std::function< void(int *)> &f, ConstraintProto *ct)
void SubstituteVariable(int var, int64 var_coeff_in_definition, const ConstraintProto &definition, ConstraintProto *ct)
void DetectDominanceRelations(const PresolveContext &context, VarDomination *var_domination, DualBoundStrengthening *dual_bound_strengthening)
void ApplyToAllIntervalIndices(const std::function< void(int *)> &f, ConstraintProto *ct)
std::vector< int > FindDuplicateConstraints(const CpModelProto &model_proto)
void ApplyToAllVariableIndices(const std::function< void(int *)> &f, ConstraintProto *ct)
void LogInfoFromContext(const PresolveContext *context)
constexpr int kObjectiveConstraint
bool ExploitDominanceRelations(const VarDomination &var_domination, PresolveContext *context)
bool RefIsPositive(int ref)
void ApplyVariableMapping(const std::vector< int > &mapping, const PresolveContext &context)
std::string ValidateCpModel(const CpModelProto &model)
void FillDomainInProto(const Domain &domain, ProtoWithDomain *proto)
The vehicle routing library lets one model and solve generic vehicle routing problems ranging from th...
Literal literal
Definition: optimization.cc:84
int index
Definition: pack.cc:508
if(!yyg->yy_init)
Definition: parser.yy.cc:965
ColIndex representative
int64 time
Definition: resource.cc:1683
int64 demand
Definition: resource.cc:123
IntervalVar * interval
Definition: resource.cc:98
int64 tail
int64 head