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"
35 #include "ortools/base/logging.h"
36 #include "ortools/base/map_util.h"
37 #include "ortools/base/mathutil.h"
38 #include "ortools/base/stl_util.h"
40 #include "ortools/sat/circuit.h"
41 #include "ortools/sat/cp_model.pb.h"
49 #include "ortools/sat/probing.h"
50 #include "ortools/sat/sat_base.h"
51 #include "ortools/sat/sat_parameters.pb.h"
54 
55 namespace operations_research {
56 namespace sat {
57 
58 bool 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 
93 bool 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 
138 bool 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 
186 bool 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.
191  if (HasEnforcementLiteral(*ct)) {
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 
262 ABSL_MUST_USE_RESULT bool CpModelPresolver::MarkConstraintAsFalse(
263  ConstraintProto* ct) {
264  if (HasEnforcementLiteral(*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 
278 bool 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 
345 bool 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 
457 bool 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 
479 bool 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 
509 bool 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 
661 bool 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 
674 bool 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 
724 bool 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 
801 bool 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 
812 bool 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 
976 bool 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 
1041 bool 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 
1094 void 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 
1118 void 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 
1151 template <typename ProtoWithVarsAndCoeffs>
1152 bool 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 
1234 bool 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 
1242 bool 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 
1259 bool 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 
1411 bool 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;
1483  if (context_->StoreLiteralImpliesVarEqValue(literal, var, value)) {
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;
1494  if (context_->StoreLiteralImpliesVarNEqValue(literal, var, value)) {
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 
1553 namespace {
1554 
1555 // Return true if the given domain only restrict the values with an upper bound.
1556 bool 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.
1562 bool 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.
1570 bool 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 
1590 bool 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.
1605 void 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 
1617 bool 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.
1943 void 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 
2087 void 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.
2142 bool 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 
2331 namespace {
2332 
2333 void 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 
2349 bool 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 
2468 bool 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 
2678 bool 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  }
2745  gtl::STLSortAndRemoveDuplicates(&new_tuples);
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  }
2821  gtl::STLSortAndRemoveDuplicates(&all_tuples);
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 
2839 bool 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 
2940 namespace {
2941 
2942 // Returns the sorted list of literals for given bool_or or at_most_one
2943 // constraint.
2944 std::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.
2965 void 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 
2982 template <typename ClauseContainer>
2983 void 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 
3022 int64 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 
3028 int64 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 
3033 int64 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 
3038 int64 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 
3043 bool 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 
3106 bool 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 
3245 bool 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 
3274 bool 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 
3458 bool 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 
3595 bool 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.
3753 void 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 
3788 void 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);
3845  context_->time_limit()->AdvanceDeterministicTime(
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);
3883  context_->StoreBooleanEqualityRelation(
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.
3891 void 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?
4066 void 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 
4287 void 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.
4358 void 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);
4421  context_->StoreBooleanEqualityRelation(
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.
4568 bool 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.
4640 bool 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 
4765 void 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)) {
4789  context_->CanonicalizeDomainOfSizeTwo(var);
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.
4827 void 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.
4867 bool 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;
4881  CHECK(context_->VarToConstraints(var).contains(kAffineRelationConstraint));
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 
4902 void 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 
5491 void 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 
5592 std::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
const std::vector< IntegerType > & PositionsSetAtLeastOnce() const
Definition: bitset.h:813
void Set(IntegerType index)
Definition: bitset.h:803
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)
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
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
const std::vector< int > & ConstraintToVars(int c) const
ABSL_MUST_USE_RESULT bool NotifyThatModelIsUnsat(const std::string &message="")
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)
const absl::flat_hash_map< int, int64 > & ObjectiveMap() const
void UpdateRuleStats(const std::string &name, int num_times=1)
ABSL_MUST_USE_RESULT bool CanonicalizeObjective()
const SatParameters & params() const
AffineRelation::Relation GetAffineRelation(int ref) const
const absl::flat_hash_set< int > & VarToConstraints(int var) const
ABSL_MUST_USE_RESULT bool SetLiteralToFalse(int lit)
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
void STLSortAndRemoveDuplicates(T *v, const LessFunc &less_func)
Definition: stl_util.h:58
void STLSetDifference(const In1 &a, const In2 &b, Out *out, Compare compare)
Definition: stl_util.h:595
bool ContainsKey(const Collection &collection, const Key &key)
Definition: map_util.h:170
const Collection::value_type::second_type & FindOrDie(const Collection &collection, const typename Collection::value_type::first_type &key)
Definition: map_util.h:176
bool DetectAndExploitSymmetriesInPresolve(PresolveContext *context)
bool LoadConstraint(const ConstraintProto &ct, Model *m)
bool RefIsPositive(int ref)
void SetToNegatedLinearExpression(const LinearExpressionProto &input_expr, LinearExpressionProto *output_negated_expr)
const LiteralIndex kNoLiteralIndex(-1)
bool HasEnforcementLiteral(const ConstraintProto &ct)
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)
int ReindexArcs(IntContainer *tails, IntContainer *heads)
Definition: circuit.h:168
void FillDomainInProto(const Domain &domain, ProtoWithDomain *proto)
std::vector< int > FindDuplicateConstraints(const CpModelProto &model_proto)
void LogInfoFromContext(const PresolveContext *context)
Domain ReadDomainFromProto(const ProtoWithDomain &proto)
void ApplyToAllVariableIndices(const std::function< void(int *)> &f, ConstraintProto *ct)
constexpr int kObjectiveConstraint
bool ExploitDominanceRelations(const VarDomination &var_domination, PresolveContext *context)
void ApplyVariableMapping(const std::vector< int > &mapping, const PresolveContext &context)
std::string ValidateCpModel(const CpModelProto &model)
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