Google OR-Tools v9.12
a fast and portable software suite for combinatorial optimization
Loading...
Searching...
No Matches
integer_expr.cc
Go to the documentation of this file.
1// Copyright 2010-2025 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 <algorithm>
17#include <cstdint>
18#include <cstdlib>
19#include <cstring>
20#include <limits>
21#include <utility>
22#include <vector>
23
24#include "absl/log/check.h"
25#include "absl/numeric/int128.h"
26#include "absl/types/span.h"
29#include "ortools/sat/integer.h"
32#include "ortools/sat/model.h"
34#include "ortools/sat/util.h"
37
38namespace operations_research {
39namespace sat {
40
41template <bool use_int128>
43 absl::Span<const Literal> enforcement_literals,
44 absl::Span<const IntegerVariable> vars,
45 absl::Span<const IntegerValue> coeffs, IntegerValue upper, Model* model)
46 : upper_bound_(upper),
47 shared_(
48 model->GetOrCreate<LinearConstraintPropagator<use_int128>::Shared>()),
49 size_(vars.size()),
50 vars_(new IntegerVariable[size_]),
51 coeffs_(new IntegerValue[size_]),
52 max_variations_(new IntegerValue[size_]) {
53 // TODO(user): deal with this corner case.
54 CHECK(!vars.empty());
55
56 // Copy data.
57 memcpy(vars_.get(), vars.data(), size_ * sizeof(IntegerVariable));
58 memcpy(coeffs_.get(), coeffs.data(), size_ * sizeof(IntegerValue));
59
60 // Handle negative coefficients.
61 for (int i = 0; i < size_; ++i) {
62 if (coeffs_[i] < 0) {
63 vars_[i] = NegationOf(vars_[i]);
64 coeffs_[i] = -coeffs_[i];
65 }
66 }
67
68 // Literal reason will only be used with the negation of enforcement_literals.
69 // It will stay constant. We also do not store enforcement_literals, but
70 // retrieve them from there.
71 literal_reason_.reserve(enforcement_literals.size());
72 for (const Literal literal : enforcement_literals) {
73 literal_reason_.push_back(literal.Negated());
74 }
76 // Initialize the reversible numbers.
77 rev_num_fixed_vars_ = 0;
78 rev_lb_fixed_vars_ = IntegerValue(0);
79}
80
81// TODO(user): Avoid duplication with other constructor.
82template <bool use_int128>
84 LinearConstraint ct, Model* model)
85 : upper_bound_(ct.ub),
86 shared_(
87 model->GetOrCreate<LinearConstraintPropagator<use_int128>::Shared>()),
88 size_(ct.num_terms),
89 vars_(std::move(ct.vars)),
90 coeffs_(std::move(ct.coeffs)),
91 max_variations_(new IntegerValue[size_]) {
92 // TODO(user): deal with this corner case.
93 CHECK_GT(size_, 0);
94
95 // Handle negative coefficients.
96 for (int i = 0; i < size_; ++i) {
97 if (coeffs_[i] < 0) {
98 vars_[i] = NegationOf(vars_[i]);
99 coeffs_[i] = -coeffs_[i];
100 }
102
103 // Initialize the reversible numbers.
104 rev_num_fixed_vars_ = 0;
105 rev_lb_fixed_vars_ = IntegerValue(0);
106}
107
108template <bool use_int128>
109void LinearConstraintPropagator<use_int128>::FillIntegerReason() {
110 shared_->integer_reason.clear();
111 shared_->reason_coeffs.clear();
112 for (int i = 0; i < size_; ++i) {
113 const IntegerVariable var = vars_[i];
114 if (!shared_->integer_trail->VariableLowerBoundIsFromLevelZero(var)) {
115 shared_->integer_reason.push_back(
116 shared_->integer_trail->LowerBoundAsLiteral(var));
117 shared_->reason_coeffs.push_back(coeffs_[i]);
118 }
119 }
120}
121
122namespace {
123IntegerValue CappedCast(absl::int128 input, IntegerValue cap) {
124 if (input >= absl::int128(cap.value())) {
125 return cap;
126 }
127 return IntegerValue(static_cast<int64_t>(input));
128}
129
130} // namespace
131
132// NOTE(user): This is only used with int128, so we code only a single version.
133template <bool use_int128>
134std::pair<IntegerValue, IntegerValue>
136 IntegerLiteral integer_literal, IntegerVariable target_var) const {
137 // The code below is wrong if integer_literal and target_var are the same.
138 // In this case we return the trivial bounds.
139 if (PositiveVariable(integer_literal.var) == PositiveVariable(target_var)) {
140 if (integer_literal.var == target_var) {
141 return {kMinIntegerValue, integer_literal.bound};
142 } else {
143 return {integer_literal.Negated().bound, kMinIntegerValue};
144 }
145 }
146
147 // Recall that all our coefficient are positive.
148 bool literal_var_present = false;
149 bool literal_var_present_positively = false;
150 IntegerValue var_coeff;
151
152 bool target_var_present_negatively = false;
153 absl::int128 target_coeff;
154
155 // Warning: It is important to do the computation like the propagation is
156 // doing it to be sure we don't have overflow, since this is what we check
157 // when creating constraints.
158 absl::int128 lb_128 = 0;
159 for (int i = 0; i < size_; ++i) {
160 const IntegerVariable var = vars_[i];
161 const IntegerValue coeff = coeffs_[i];
162 if (var == NegationOf(target_var)) {
163 target_coeff = absl::int128(coeff.value());
164 target_var_present_negatively = true;
165 }
166
167 const IntegerValue lb = shared_->integer_trail->LowerBound(var);
168 lb_128 += absl::int128(coeff.value()) * absl::int128(lb.value());
169 if (PositiveVariable(var) == PositiveVariable(integer_literal.var)) {
170 var_coeff = coeff;
171 literal_var_present = true;
172 literal_var_present_positively = (var == integer_literal.var);
173 }
174 }
175
176 if (!literal_var_present || !target_var_present_negatively) {
178 }
179
180 // The upper bound on NegationOf(target_var) are lb(-target) + slack / coeff.
181 // So the lower bound on target_var is ub - slack / coeff.
182 const absl::int128 slack128 = absl::int128(upper_bound_.value()) - lb_128;
183 const IntegerValue target_lb = shared_->integer_trail->LowerBound(target_var);
184 const IntegerValue target_ub = shared_->integer_trail->UpperBound(target_var);
185 if (slack128 <= 0) {
186 // TODO(user): If there is a conflict (negative slack) we can be more
187 // precise.
188 return {target_ub, target_ub};
189 }
190
191 const IntegerValue target_diff = target_ub - target_lb;
192 const IntegerValue delta = CappedCast(slack128 / target_coeff, target_diff);
193
194 // A literal means var >= bound.
195 if (literal_var_present_positively) {
196 // We have var_coeff * var in the expression, the literal is var >= bound.
197 // When it is false, it is not relevant as implied_lb used var >= lb.
198 // When it is true, the diff is bound - lb.
199 const IntegerValue diff =
200 std::max(IntegerValue(0),
201 integer_literal.bound -
202 shared_->integer_trail->LowerBound(integer_literal.var));
203 const absl::int128 tighter_slack =
204 std::max(absl::int128(0), slack128 - absl::int128(var_coeff.value()) *
205 absl::int128(diff.value()));
206 const IntegerValue tighter_delta =
207 CappedCast(tighter_slack / target_coeff, target_diff);
208 return {target_ub - delta, target_ub - tighter_delta};
209 } else {
210 // We have var_coeff * -var in the expression, the literal is var >= bound.
211 // When it is true, it is not relevant as implied_lb used -var >= -ub.
212 // And when it is false it means var < bound, so -var >= -bound + 1
213 const IntegerValue diff =
214 std::max(IntegerValue(0),
215 shared_->integer_trail->UpperBound(integer_literal.var) -
216 integer_literal.bound + 1);
217 const absl::int128 tighter_slack =
218 std::max(absl::int128(0), slack128 - absl::int128(var_coeff.value()) *
219 absl::int128(diff.value()));
220 const IntegerValue tighter_delta =
221 CappedCast(tighter_slack / target_coeff, target_diff);
222 return {target_ub - tighter_delta, target_ub - delta};
223 }
224}
225
226template <bool use_int128>
228 int /*id*/, IntegerValue propagation_slack, IntegerVariable var_to_explain,
229 int trail_index, std::vector<Literal>* literals_reason,
230 std::vector<int>* trail_indices_reason) {
231 *literals_reason = literal_reason_;
232 trail_indices_reason->clear();
233 shared_->reason_coeffs.clear();
234 for (int i = 0; i < size_; ++i) {
235 const IntegerVariable var = vars_[i];
236 if (PositiveVariable(var) == PositiveVariable(var_to_explain)) {
237 continue;
238 }
239 const int index =
240 shared_->integer_trail->FindTrailIndexOfVarBefore(var, trail_index);
241 if (index >= 0) {
242 trail_indices_reason->push_back(index);
243 if (propagation_slack > 0) {
244 shared_->reason_coeffs.push_back(coeffs_[i]);
245 }
246 }
247 }
248 if (propagation_slack > 0) {
249 shared_->integer_trail->RelaxLinearReason(
250 propagation_slack, shared_->reason_coeffs, trail_indices_reason);
251 }
252}
253
254template <bool use_int128>
256 // Reified case: If any of the enforcement_literals are false, we ignore the
257 // constraint.
258 int num_unassigned_enforcement_literal = 0;
259 LiteralIndex unique_unnasigned_literal = kNoLiteralIndex;
260 for (const Literal negated_enforcement : literal_reason_) {
261 const Literal literal = negated_enforcement.Negated();
262 if (shared_->assignment.LiteralIsFalse(literal)) return true;
263 if (!shared_->assignment.LiteralIsTrue(literal)) {
264 ++num_unassigned_enforcement_literal;
265 unique_unnasigned_literal = literal.Index();
266 }
267 }
268
269 // Unfortunately, we can't propagate anything if we have more than one
270 // unassigned enforcement literal.
271 if (num_unassigned_enforcement_literal > 1) return true;
272
273 int num_fixed_vars = rev_num_fixed_vars_;
274 IntegerValue lb_fixed_vars = rev_lb_fixed_vars_;
275
276 // Compute the new lower bound and update the reversible structures.
277 absl::int128 lb_128 = 0;
278 IntegerValue lb_unfixed_vars = IntegerValue(0);
279 for (int i = num_fixed_vars; i < size_; ++i) {
280 const IntegerVariable var = vars_[i];
281 const IntegerValue coeff = coeffs_[i];
282 const IntegerValue lb = shared_->integer_trail->LowerBound(var);
283 const IntegerValue ub = shared_->integer_trail->UpperBound(var);
284 if (use_int128) {
285 lb_128 += absl::int128(lb.value()) * absl::int128(coeff.value());
286 continue;
287 }
288
289 if (lb != ub) {
290 max_variations_[i] = (ub - lb) * coeff;
291 lb_unfixed_vars += lb * coeff;
292 } else {
293 // Update the set of fixed variables.
294 std::swap(vars_[i], vars_[num_fixed_vars]);
295 std::swap(coeffs_[i], coeffs_[num_fixed_vars]);
296 std::swap(max_variations_[i], max_variations_[num_fixed_vars]);
297 num_fixed_vars++;
298 lb_fixed_vars += lb * coeff;
299 }
300 }
301 shared_->time_limit->AdvanceDeterministicTime(
302 static_cast<double>(size_ - num_fixed_vars) * 5e-9);
303
304 // Save the current sum of fixed variables.
305 if (is_registered_ && num_fixed_vars != rev_num_fixed_vars_) {
306 CHECK(!use_int128);
307 shared_->rev_integer_value_repository->SaveState(&rev_lb_fixed_vars_);
308 shared_->rev_int_repository->SaveState(&rev_num_fixed_vars_);
309 rev_num_fixed_vars_ = num_fixed_vars;
310 rev_lb_fixed_vars_ = lb_fixed_vars;
311 }
312
313 // If use_int128 is true, the slack or propagation slack can be larger than
314 // this. To detect overflow with capped arithmetic, it is important the slack
315 // used in our algo never exceed this value.
316 const absl::int128 max_slack = std::numeric_limits<int64_t>::max() - 1;
317
318 // Conflict?
319 IntegerValue slack;
320 absl::int128 slack128;
321 if (use_int128) {
322 slack128 = absl::int128(upper_bound_.value()) - lb_128;
323 if (slack128 < 0) {
324 // It is fine if we don't relax as much as possible.
325 // Note that RelaxLinearReason() is overflow safe.
326 slack = static_cast<int64_t>(std::max(-max_slack, slack128));
327 }
328 } else {
329 slack = upper_bound_ - (lb_fixed_vars + lb_unfixed_vars);
330 }
331 if (slack < 0) {
332 FillIntegerReason();
333 shared_->integer_trail->RelaxLinearReason(
334 -slack - 1, shared_->reason_coeffs, &shared_->integer_reason);
335
336 if (num_unassigned_enforcement_literal == 1) {
337 // Propagate the only non-true literal to false.
338 const Literal to_propagate = Literal(unique_unnasigned_literal).Negated();
339 std::vector<Literal> tmp = literal_reason_;
340 tmp.erase(std::find(tmp.begin(), tmp.end(), to_propagate));
341 shared_->integer_trail->EnqueueLiteral(to_propagate, tmp,
342 shared_->integer_reason);
343 return true;
344 }
345 return shared_->integer_trail->ReportConflict(literal_reason_,
346 shared_->integer_reason);
347 }
348
349 // We can only propagate more if all the enforcement literals are true.
350 if (num_unassigned_enforcement_literal > 0) return true;
351
352 // The lower bound of all the variables except one can be used to update the
353 // upper bound of the last one.
354 for (int i = num_fixed_vars; i < size_; ++i) {
355 if (!use_int128 && max_variations_[i] <= slack) continue;
356
357 // TODO(user): If the new ub fall into an hole of the variable, we can
358 // actually relax the reason more by computing a better slack.
359 const IntegerVariable var = vars_[i];
360 const IntegerValue coeff = coeffs_[i];
361 const IntegerValue lb = shared_->integer_trail->LowerBound(var);
362
363 IntegerValue new_ub;
364 IntegerValue propagation_slack;
365 if (use_int128) {
366 const absl::int128 coeff128 = absl::int128(coeff.value());
367 const absl::int128 div128 = slack128 / coeff128;
368 const IntegerValue ub = shared_->integer_trail->UpperBound(var);
369 if (absl::int128(lb.value()) + div128 >= absl::int128(ub.value())) {
370 continue;
371 }
372 new_ub = lb + IntegerValue(static_cast<int64_t>(div128));
373 propagation_slack = static_cast<int64_t>(
374 std::min(max_slack, (div128 + 1) * coeff128 - slack128 - 1));
375 } else {
376 const IntegerValue div = slack / coeff;
377 new_ub = lb + div;
378 propagation_slack = (div + 1) * coeff - slack - 1;
379 }
380 if (!shared_->integer_trail->EnqueueWithLazyReason(
381 IntegerLiteral::LowerOrEqual(var, new_ub), 0, propagation_slack,
382 this)) {
383 // TODO(user): this is never supposed to happen since if we didn't have a
384 // conflict above, we should be able to reduce the upper bound. It might
385 // indicate an issue with our Boolean <-> integer encoding.
386 return false;
387 }
388 }
389
390 return true;
391}
392
393template <bool use_int128>
395 // TODO(user): Deal with enforcements. It is just a bit of code to read the
396 // value of the literals at level zero.
397 if (!literal_reason_.empty()) return true;
398
399 // Compute the new lower bound and update the reversible structures.
400 absl::int128 lb_128 = 0;
401 IntegerValue min_activity = IntegerValue(0);
402 for (int i = 0; i < size_; ++i) {
403 const IntegerVariable var = vars_[i];
404 const IntegerValue coeff = coeffs_[i];
405 const IntegerValue lb = shared_->integer_trail->LevelZeroLowerBound(var);
406 if (use_int128) {
407 lb_128 += absl::int128(lb.value()) * absl::int128(coeff.value());
408 } else {
409 const IntegerValue ub = shared_->integer_trail->LevelZeroUpperBound(var);
410 max_variations_[i] = (ub - lb) * coeff;
411 min_activity += lb * coeff;
412 }
413 }
414 shared_->time_limit->AdvanceDeterministicTime(
415 static_cast<double>(size_ * 1e-9));
416
417 // Conflict?
418 IntegerValue slack;
419 absl::int128 slack128;
420 if (use_int128) {
421 slack128 = absl::int128(upper_bound_.value()) - lb_128;
422 if (slack128 < 0) {
423 return shared_->integer_trail->ReportConflict({}, {});
424 }
425 } else {
426 slack = upper_bound_ - min_activity;
427 if (slack < 0) {
428 return shared_->integer_trail->ReportConflict({}, {});
429 }
430 }
431
432 // The lower bound of all the variables except one can be used to update the
433 // upper bound of the last one.
434 for (int i = 0; i < size_; ++i) {
435 if (!use_int128 && max_variations_[i] <= slack) continue;
436
437 const IntegerVariable var = vars_[i];
438 const IntegerValue coeff = coeffs_[i];
439 const IntegerValue lb = shared_->integer_trail->LevelZeroLowerBound(var);
440
441 IntegerValue new_ub;
442 if (use_int128) {
443 const IntegerValue ub = shared_->integer_trail->LevelZeroUpperBound(var);
444 if (absl::int128((ub - lb).value()) * absl::int128(coeff.value()) <=
445 slack128) {
446 continue;
447 }
448 const absl::int128 div128 = slack128 / absl::int128(coeff.value());
449 new_ub = lb + IntegerValue(static_cast<int64_t>(div128));
450 } else {
451 const IntegerValue div = slack / coeff;
452 new_ub = lb + div;
453 }
454 if (!shared_->integer_trail->Enqueue(
455 IntegerLiteral::LowerOrEqual(var, new_ub), {}, {})) {
456 return false;
457 }
458 }
459
460 return true;
461}
462
463template <bool use_int128>
465 GenericLiteralWatcher* watcher) {
466 is_registered_ = true;
467 const int id = watcher->Register(this);
468 for (int i = 0; i < size_; ++i) {
469 watcher->WatchLowerBound(vars_[i], id);
470 }
471 for (const Literal negated_enforcement : literal_reason_) {
472 // We only watch the true direction.
473 //
474 // TODO(user): if there is more than one, maybe we should watch more to
475 // propagate a "conflict" as soon as only one is unassigned?
476 watcher->WatchLiteral(negated_enforcement.Negated(), id);
477 }
478}
479
480// Explicit declaration.
483
485 const std::vector<IntegerVariable>& vars,
486 const std::vector<IntegerValue>& coeffs,
487 Model* model)
488 : target_(target),
489 vars_(vars),
490 coeffs_(coeffs),
491 trail_(model->GetOrCreate<Trail>()),
492 integer_trail_(model->GetOrCreate<IntegerTrail>()) {
493 auto* watcher = model->GetOrCreate<GenericLiteralWatcher>();
494 const int id = watcher->Register(this);
495 watcher->SetPropagatorPriority(id, 2);
496 watcher->WatchIntegerVariable(target, id);
497 for (const IntegerVariable& var : vars_) {
498 watcher->WatchIntegerVariable(var, id);
499 }
500}
501
502// TODO(user): We could go even further than just the GCD, and do more
503// arithmetic to tighten the target bounds. See for instance a problem like
504// ej.mps.gz that we don't solve easily, but has just 3 variables! the goal is
505// to minimize X, given 31013 X - 41014 Y - 51015 Z = -31013 (all >=0, Y and Z
506// bounded with high values). I know some MIP solvers have a basic linear
507// diophantine equation support.
509 // TODO(user): Once the GCD is not 1, we could at any level make sure the
510 // objective is of the correct form. For now, this only happen in a few
511 // miplib problem that we close quickly, so I didn't add the extra code yet.
512 if (trail_->CurrentDecisionLevel() != 0) return true;
513
514 int64_t gcd = 0;
515 IntegerValue sum(0);
516 for (int i = 0; i < vars_.size(); ++i) {
517 if (integer_trail_->IsFixed(vars_[i])) {
518 sum += coeffs_[i] * integer_trail_->LowerBound(vars_[i]);
519 continue;
520 }
521 gcd = MathUtil::GCD64(gcd, std::abs(coeffs_[i].value()));
522 if (gcd == 1) break;
523 }
524 if (gcd == 0) return true; // All fixed.
525
526 if (gcd > gcd_) {
527 VLOG(1) << "Objective gcd: " << gcd;
528 }
529 CHECK_GE(gcd, gcd_);
530 gcd_ = IntegerValue(gcd);
531
532 const IntegerValue lb = integer_trail_->LowerBound(target_);
533 const IntegerValue lb_remainder = PositiveRemainder(lb - sum, gcd_);
534 if (lb_remainder != 0) {
535 if (!integer_trail_->Enqueue(
536 IntegerLiteral::GreaterOrEqual(target_, lb + gcd_ - lb_remainder),
537 {}, {}))
538 return false;
539 }
540
541 const IntegerValue ub = integer_trail_->UpperBound(target_);
542 const IntegerValue ub_remainder =
543 PositiveRemainder(ub - sum, IntegerValue(gcd));
544 if (ub_remainder != 0) {
545 if (!integer_trail_->Enqueue(
546 IntegerLiteral::LowerOrEqual(target_, ub - ub_remainder), {}, {}))
547 return false;
548 }
549
550 return true;
551}
552
553MinPropagator::MinPropagator(std::vector<AffineExpression> vars,
554 AffineExpression min_var,
555 IntegerTrail* integer_trail)
556 : vars_(std::move(vars)),
557 min_var_(min_var),
558 integer_trail_(integer_trail) {}
559
561 if (vars_.empty()) return true;
562
563 // Count the number of interval that are possible candidate for the min.
564 // Only the intervals for which lb > current_min_ub cannot.
565 const IntegerLiteral min_ub_literal =
566 integer_trail_->UpperBoundAsLiteral(min_var_);
567 const IntegerValue current_min_ub = integer_trail_->UpperBound(min_var_);
568 int num_intervals_that_can_be_min = 0;
569 int last_possible_min_interval = 0;
570
571 IntegerValue min = kMaxIntegerValue;
572 for (int i = 0; i < vars_.size(); ++i) {
573 const IntegerValue lb = integer_trail_->LowerBound(vars_[i]);
574 min = std::min(min, lb);
575 if (lb <= current_min_ub) {
576 ++num_intervals_that_can_be_min;
577 last_possible_min_interval = i;
578 }
579 }
580
581 // Propagation a)
582 if (min > integer_trail_->LowerBound(min_var_)) {
583 integer_reason_.clear();
584 for (const AffineExpression& var : vars_) {
585 integer_reason_.push_back(var.GreaterOrEqual(min));
586 }
587 if (!integer_trail_->SafeEnqueue(min_var_.GreaterOrEqual(min),
588 integer_reason_)) {
589 return false;
590 }
591 }
592
593 // Propagation b)
594 if (num_intervals_that_can_be_min == 1) {
595 const IntegerValue ub_of_only_candidate =
596 integer_trail_->UpperBound(vars_[last_possible_min_interval]);
597 if (current_min_ub < ub_of_only_candidate) {
598 integer_reason_.clear();
599
600 // The reason is that all the other interval start after current_min_ub.
601 // And that min_ub has its current value.
602 integer_reason_.push_back(min_ub_literal);
603 for (const AffineExpression& var : vars_) {
604 if (var == vars_[last_possible_min_interval]) continue;
605 integer_reason_.push_back(var.GreaterOrEqual(current_min_ub + 1));
606 }
607 if (!integer_trail_->SafeEnqueue(
608 vars_[last_possible_min_interval].LowerOrEqual(current_min_ub),
609 integer_reason_)) {
610 return false;
611 }
612 }
613 }
614
615 // Conflict.
616 //
617 // TODO(user): Not sure this code is useful since this will be detected
618 // by the fact that the [lb, ub] of the min is empty. It depends on the
619 // propagation order though, but probably the precedences propagator would
620 // propagate before this one. So change this to a CHECK?
621 if (num_intervals_that_can_be_min == 0) {
622 integer_reason_.clear();
623
624 // Almost the same as propagation b).
625 integer_reason_.push_back(min_ub_literal);
626 for (const AffineExpression& var : vars_) {
627 IntegerLiteral lit = var.GreaterOrEqual(current_min_ub + 1);
628 if (lit != IntegerLiteral::TrueLiteral()) {
629 integer_reason_.push_back(lit);
630 }
631 }
632 return integer_trail_->ReportConflict(integer_reason_);
633 }
634
635 return true;
636}
637
639 const int id = watcher->Register(this);
640 for (const AffineExpression& var : vars_) {
641 watcher->WatchLowerBound(var, id);
642 }
643 watcher->WatchUpperBound(min_var_, id);
644}
645
646LinMinPropagator::LinMinPropagator(std::vector<LinearExpression> exprs,
647 IntegerVariable min_var, Model* model)
648 : exprs_(std::move(exprs)),
649 min_var_(min_var),
650 model_(model),
651 integer_trail_(model_->GetOrCreate<IntegerTrail>()) {}
652
653void LinMinPropagator::Explain(int id, IntegerValue propagation_slack,
654 IntegerVariable var_to_explain, int trail_index,
655 std::vector<Literal>* literals_reason,
656 std::vector<int>* trail_indices_reason) {
657 const auto& vars = exprs_[id].vars;
658 const auto& coeffs = exprs_[id].coeffs;
659 literals_reason->clear();
660 trail_indices_reason->clear();
661 std::vector<IntegerValue> reason_coeffs;
662 const int size = vars.size();
663 for (int i = 0; i < size; ++i) {
664 const IntegerVariable var = vars[i];
665 if (PositiveVariable(var) == PositiveVariable(var_to_explain)) {
666 continue;
667 }
668 const int index =
669 integer_trail_->FindTrailIndexOfVarBefore(var, trail_index);
670 if (index >= 0) {
671 trail_indices_reason->push_back(index);
672 if (propagation_slack > 0) {
673 reason_coeffs.push_back(coeffs[i]);
674 }
675 }
676 }
677 if (propagation_slack > 0) {
678 integer_trail_->RelaxLinearReason(propagation_slack, reason_coeffs,
679 trail_indices_reason);
680 }
681 // Now add the old integer_reason that triggered this propagation.
682 for (IntegerLiteral reason_lit : integer_reason_for_unique_candidate_) {
683 const int index =
684 integer_trail_->FindTrailIndexOfVarBefore(reason_lit.var, trail_index);
685 if (index >= 0) {
686 trail_indices_reason->push_back(index);
687 }
688 }
689}
690
691bool LinMinPropagator::PropagateLinearUpperBound(
692 int id, absl::Span<const IntegerVariable> vars,
693 absl::Span<const IntegerValue> coeffs, const IntegerValue upper_bound) {
694 IntegerValue sum_lb = IntegerValue(0);
695 const int num_vars = vars.size();
696 max_variations_.resize(num_vars);
697 for (int i = 0; i < num_vars; ++i) {
698 const IntegerVariable var = vars[i];
699 const IntegerValue coeff = coeffs[i];
700 // The coefficients are assumed to be positive for this to work properly.
701 DCHECK_GE(coeff, 0);
702 const IntegerValue lb = integer_trail_->LowerBound(var);
703 const IntegerValue ub = integer_trail_->UpperBound(var);
704 max_variations_[i] = (ub - lb) * coeff;
705 sum_lb += lb * coeff;
706 }
707
708 model_->GetOrCreate<TimeLimit>()->AdvanceDeterministicTime(
709 static_cast<double>(num_vars) * 1e-9);
710
711 const IntegerValue slack = upper_bound - sum_lb;
712 if (slack < 0) {
713 // Conflict.
714 local_reason_.clear();
715 reason_coeffs_.clear();
716 for (int i = 0; i < num_vars; ++i) {
717 const IntegerVariable var = vars[i];
718 if (!integer_trail_->VariableLowerBoundIsFromLevelZero(var)) {
719 local_reason_.push_back(integer_trail_->LowerBoundAsLiteral(var));
720 reason_coeffs_.push_back(coeffs[i]);
721 }
722 }
723 integer_trail_->RelaxLinearReason(-slack - 1, reason_coeffs_,
724 &local_reason_);
725 local_reason_.insert(local_reason_.end(),
726 integer_reason_for_unique_candidate_.begin(),
727 integer_reason_for_unique_candidate_.end());
728 return integer_trail_->ReportConflict({}, local_reason_);
729 }
730
731 // The lower bound of all the variables except one can be used to update the
732 // upper bound of the last one.
733 for (int i = 0; i < num_vars; ++i) {
734 if (max_variations_[i] <= slack) continue;
735
736 const IntegerVariable var = vars[i];
737 const IntegerValue coeff = coeffs[i];
738 const IntegerValue div = slack / coeff;
739 const IntegerValue new_ub = integer_trail_->LowerBound(var) + div;
740 const IntegerValue propagation_slack = (div + 1) * coeff - slack - 1;
741 if (!integer_trail_->EnqueueWithLazyReason(
742 IntegerLiteral::LowerOrEqual(var, new_ub), id, propagation_slack,
743 this)) {
744 return false;
745 }
746 }
747 return true;
748}
749
751 if (exprs_.empty()) return true;
752
753 // Count the number of interval that are possible candidate for the min.
754 // Only the intervals for which lb > current_min_ub cannot.
755 const IntegerValue current_min_ub = integer_trail_->UpperBound(min_var_);
756 int num_intervals_that_can_be_min = 0;
757 int last_possible_min_interval = 0;
758
759 expr_lbs_.clear();
760 IntegerValue min_of_linear_expression_lb = kMaxIntegerValue;
761 for (int i = 0; i < exprs_.size(); ++i) {
762 const IntegerValue lb = exprs_[i].Min(*integer_trail_);
763 expr_lbs_.push_back(lb);
764 min_of_linear_expression_lb = std::min(min_of_linear_expression_lb, lb);
765 if (lb <= current_min_ub) {
766 ++num_intervals_that_can_be_min;
767 last_possible_min_interval = i;
768 }
769 }
770
771 // Propagation a) lb(min) >= lb(MIN(exprs)) = MIN(lb(exprs));
772
773 // Conflict will be detected by the fact that the [lb, ub] of the min is
774 // empty. In case of conflict, we just need the reason for pushing UB + 1.
775 if (min_of_linear_expression_lb > current_min_ub) {
776 min_of_linear_expression_lb = current_min_ub + 1;
777 }
778 if (min_of_linear_expression_lb > integer_trail_->LowerBound(min_var_)) {
779 local_reason_.clear();
780 for (int i = 0; i < exprs_.size(); ++i) {
781 const IntegerValue slack = expr_lbs_[i] - min_of_linear_expression_lb;
782 integer_trail_->AppendRelaxedLinearReason(slack, exprs_[i].coeffs,
783 exprs_[i].vars, &local_reason_);
784 }
785 if (!integer_trail_->Enqueue(IntegerLiteral::GreaterOrEqual(
786 min_var_, min_of_linear_expression_lb),
787 {}, local_reason_)) {
788 return false;
789 }
790 }
791
792 // Propagation b) ub(min) >= ub(MIN(exprs)) and we can't propagate anything
793 // here unless there is just one possible expression 'e' that can be the min:
794 // for all u != e, lb(u) > ub(min);
795 // In this case, ub(min) >= ub(e).
796 if (num_intervals_that_can_be_min == 1) {
797 const IntegerValue ub_of_only_candidate =
798 exprs_[last_possible_min_interval].Max(*integer_trail_);
799 if (current_min_ub < ub_of_only_candidate) {
800 // For this propagation, we only need to fill the integer reason once at
801 // the lowest level. At higher levels this reason still remains valid.
802 if (rev_unique_candidate_ == 0) {
803 integer_reason_for_unique_candidate_.clear();
804
805 // The reason is that all the other interval start after current_min_ub.
806 // And that min_ub has its current value.
807 integer_reason_for_unique_candidate_.push_back(
808 integer_trail_->UpperBoundAsLiteral(min_var_));
809 for (int i = 0; i < exprs_.size(); ++i) {
810 if (i == last_possible_min_interval) continue;
811 const IntegerValue slack = expr_lbs_[i] - (current_min_ub + 1);
812 integer_trail_->AppendRelaxedLinearReason(
813 slack, exprs_[i].coeffs, exprs_[i].vars,
814 &integer_reason_for_unique_candidate_);
815 }
816 rev_unique_candidate_ = 1;
817 }
818
819 return PropagateLinearUpperBound(
820 last_possible_min_interval, exprs_[last_possible_min_interval].vars,
821 exprs_[last_possible_min_interval].coeffs,
822 current_min_ub - exprs_[last_possible_min_interval].offset);
823 }
824 }
825
826 return true;
827}
828
830 const int id = watcher->Register(this);
831 bool has_var_also_in_exprs = false;
832 for (const LinearExpression& expr : exprs_) {
833 for (int i = 0; i < expr.vars.size(); ++i) {
834 const IntegerVariable& var = expr.vars[i];
835 const IntegerValue coeff = expr.coeffs[i];
836 if (coeff > 0) {
837 watcher->WatchLowerBound(var, id);
838 } else {
839 watcher->WatchUpperBound(var, id);
840 }
841 has_var_also_in_exprs |= (var == min_var_);
842 }
843 }
844 watcher->WatchUpperBound(min_var_, id);
845 watcher->RegisterReversibleInt(id, &rev_unique_candidate_);
846 if (has_var_also_in_exprs) {
848 }
849}
850
853 IntegerTrail* integer_trail)
854 : a_(a), b_(b), p_(p), integer_trail_(integer_trail) {}
855
856// We want all affine expression to be either non-negative or across zero.
857bool ProductPropagator::CanonicalizeCases() {
858 if (integer_trail_->UpperBound(a_) <= 0) {
859 a_ = a_.Negated();
860 p_ = p_.Negated();
861 }
862 if (integer_trail_->UpperBound(b_) <= 0) {
863 b_ = b_.Negated();
864 p_ = p_.Negated();
865 }
866
867 // If both a and b positive, p must be too.
868 if (integer_trail_->LowerBound(a_) >= 0 &&
869 integer_trail_->LowerBound(b_) >= 0) {
870 return integer_trail_->SafeEnqueue(
871 p_.GreaterOrEqual(0), {a_.GreaterOrEqual(0), b_.GreaterOrEqual(0)});
872 }
873
874 // Otherwise, make sure p is non-negative or across zero.
875 if (integer_trail_->UpperBound(p_) <= 0) {
876 if (integer_trail_->LowerBound(a_) < 0) {
877 DCHECK_GT(integer_trail_->UpperBound(a_), 0);
878 a_ = a_.Negated();
879 p_ = p_.Negated();
880 } else {
881 DCHECK_LT(integer_trail_->LowerBound(b_), 0);
882 DCHECK_GT(integer_trail_->UpperBound(b_), 0);
883 b_ = b_.Negated();
884 p_ = p_.Negated();
885 }
886 }
887
888 return true;
889}
890
891// Note that this propagation is exact, except on the domain of p as this
892// involves more complex arithmetic.
893//
894// TODO(user): We could tighten the bounds on p by removing extreme value that
895// do not contains divisor in the domains of a or b. There is an algo in O(
896// smallest domain size between a or b).
897bool ProductPropagator::PropagateWhenAllNonNegative() {
898 {
899 const IntegerValue max_a = integer_trail_->UpperBound(a_);
900 const IntegerValue max_b = integer_trail_->UpperBound(b_);
901 const IntegerValue new_max = CapProdI(max_a, max_b);
902 if (new_max < integer_trail_->UpperBound(p_)) {
903 if (!integer_trail_->SafeEnqueue(
904 p_.LowerOrEqual(new_max),
905 {integer_trail_->UpperBoundAsLiteral(a_),
906 integer_trail_->UpperBoundAsLiteral(b_), a_.GreaterOrEqual(0),
907 b_.GreaterOrEqual(0)})) {
908 return false;
909 }
910 }
911 }
912
913 {
914 const IntegerValue min_a = integer_trail_->LowerBound(a_);
915 const IntegerValue min_b = integer_trail_->LowerBound(b_);
916 const IntegerValue new_min = CapProdI(min_a, min_b);
917
918 // The conflict test is needed because when new_min is large, we could
919 // have an overflow in p_.GreaterOrEqual(new_min);
920 if (new_min > integer_trail_->UpperBound(p_)) {
921 return integer_trail_->ReportConflict(
922 {integer_trail_->UpperBoundAsLiteral(p_),
923 integer_trail_->LowerBoundAsLiteral(a_),
924 integer_trail_->LowerBoundAsLiteral(b_)});
925 }
926 if (new_min > integer_trail_->LowerBound(p_)) {
927 if (!integer_trail_->SafeEnqueue(
928 p_.GreaterOrEqual(new_min),
929 {integer_trail_->LowerBoundAsLiteral(a_),
930 integer_trail_->LowerBoundAsLiteral(b_)})) {
931 return false;
932 }
933 }
934 }
935
936 for (int i = 0; i < 2; ++i) {
937 const AffineExpression a = i == 0 ? a_ : b_;
938 const AffineExpression b = i == 0 ? b_ : a_;
939 const IntegerValue max_a = integer_trail_->UpperBound(a);
940 const IntegerValue min_b = integer_trail_->LowerBound(b);
941 const IntegerValue min_p = integer_trail_->LowerBound(p_);
942 const IntegerValue max_p = integer_trail_->UpperBound(p_);
943 const IntegerValue prod = CapProdI(max_a, min_b);
944 if (prod > max_p) {
945 if (!integer_trail_->SafeEnqueue(a.LowerOrEqual(FloorRatio(max_p, min_b)),
946 {integer_trail_->LowerBoundAsLiteral(b),
947 integer_trail_->UpperBoundAsLiteral(p_),
948 p_.GreaterOrEqual(0)})) {
949 return false;
950 }
951 } else if (prod < min_p && max_a != 0) {
952 if (!integer_trail_->SafeEnqueue(
953 b.GreaterOrEqual(CeilRatio(min_p, max_a)),
954 {integer_trail_->UpperBoundAsLiteral(a),
955 integer_trail_->LowerBoundAsLiteral(p_), a.GreaterOrEqual(0)})) {
956 return false;
957 }
958 }
959 }
960
961 return true;
962}
963
964// This assumes p > 0, p = a * X, and X can take any value.
965// We can propagate max of a by computing a bound on the min b when positive.
966// The expression b is just used to detect when there is no solution given the
967// upper bound of b.
968bool ProductPropagator::PropagateMaxOnPositiveProduct(AffineExpression a,
970 IntegerValue min_p,
971 IntegerValue max_p) {
972 const IntegerValue max_a = integer_trail_->UpperBound(a);
973 if (max_a <= 0) return true;
974 DCHECK_GT(min_p, 0);
975
976 if (max_a >= min_p) {
977 if (max_p < max_a) {
978 if (!integer_trail_->SafeEnqueue(
979 a.LowerOrEqual(max_p),
980 {p_.LowerOrEqual(max_p), p_.GreaterOrEqual(1)})) {
981 return false;
982 }
983 }
984 return true;
985 }
986
987 const IntegerValue min_pos_b = CeilRatio(min_p, max_a);
988 if (min_pos_b > integer_trail_->UpperBound(b)) {
989 if (!integer_trail_->SafeEnqueue(
990 b.LowerOrEqual(0), {integer_trail_->LowerBoundAsLiteral(p_),
991 integer_trail_->UpperBoundAsLiteral(a),
992 integer_trail_->UpperBoundAsLiteral(b)})) {
993 return false;
994 }
995 return true;
996 }
997
998 const IntegerValue new_max_a = FloorRatio(max_p, min_pos_b);
999 if (new_max_a < integer_trail_->UpperBound(a)) {
1000 if (!integer_trail_->SafeEnqueue(
1001 a.LowerOrEqual(new_max_a),
1002 {integer_trail_->LowerBoundAsLiteral(p_),
1003 integer_trail_->UpperBoundAsLiteral(a),
1004 integer_trail_->UpperBoundAsLiteral(p_)})) {
1005 return false;
1006 }
1007 }
1008 return true;
1009}
1010
1012 if (!CanonicalizeCases()) return false;
1013
1014 // In the most common case, we use better reasons even though the code
1015 // below would propagate the same.
1016 const int64_t min_a = integer_trail_->LowerBound(a_).value();
1017 const int64_t min_b = integer_trail_->LowerBound(b_).value();
1018 if (min_a >= 0 && min_b >= 0) {
1019 // This was done by CanonicalizeCases().
1020 DCHECK_GE(integer_trail_->LowerBound(p_), 0);
1021 return PropagateWhenAllNonNegative();
1022 }
1023
1024 // Lets propagate on p_ first, the max/min is given by one of: max_a * max_b,
1025 // max_a * min_b, min_a * max_b, min_a * min_b. This is true, because any
1026 // product x * y, depending on the sign, is dominated by one of these.
1027 //
1028 // TODO(user): In the reasons, including all 4 bounds is always correct, but
1029 // we might be able to relax some of them.
1030 const IntegerValue max_a = integer_trail_->UpperBound(a_);
1031 const IntegerValue max_b = integer_trail_->UpperBound(b_);
1032 const IntegerValue p1 = CapProdI(max_a, max_b);
1033 const IntegerValue p2 = CapProdI(max_a, min_b);
1034 const IntegerValue p3 = CapProdI(min_a, max_b);
1035 const IntegerValue p4 = CapProdI(min_a, min_b);
1036 const IntegerValue new_max_p = std::max({p1, p2, p3, p4});
1037 if (new_max_p < integer_trail_->UpperBound(p_)) {
1038 if (!integer_trail_->SafeEnqueue(
1039 p_.LowerOrEqual(new_max_p),
1040 {integer_trail_->LowerBoundAsLiteral(a_),
1041 integer_trail_->LowerBoundAsLiteral(b_),
1042 integer_trail_->UpperBoundAsLiteral(a_),
1043 integer_trail_->UpperBoundAsLiteral(b_)})) {
1044 return false;
1045 }
1046 }
1047 const IntegerValue new_min_p = std::min({p1, p2, p3, p4});
1048 if (new_min_p > integer_trail_->LowerBound(p_)) {
1049 if (!integer_trail_->SafeEnqueue(
1050 p_.GreaterOrEqual(new_min_p),
1051 {integer_trail_->LowerBoundAsLiteral(a_),
1052 integer_trail_->LowerBoundAsLiteral(b_),
1053 integer_trail_->UpperBoundAsLiteral(a_),
1054 integer_trail_->UpperBoundAsLiteral(b_)})) {
1055 return false;
1056 }
1057 }
1058
1059 // Lets propagate on a and b.
1060 const IntegerValue min_p = integer_trail_->LowerBound(p_);
1061 const IntegerValue max_p = integer_trail_->UpperBound(p_);
1062
1063 // We need a bit more propagation to avoid bad cases below.
1064 const bool zero_is_possible = min_p <= 0;
1065 if (!zero_is_possible) {
1066 if (integer_trail_->LowerBound(a_) == 0) {
1067 if (!integer_trail_->SafeEnqueue(
1068 a_.GreaterOrEqual(1),
1069 {p_.GreaterOrEqual(1), a_.GreaterOrEqual(0)})) {
1070 return false;
1071 }
1072 }
1073 if (integer_trail_->LowerBound(b_) == 0) {
1074 if (!integer_trail_->SafeEnqueue(
1075 b_.GreaterOrEqual(1),
1076 {p_.GreaterOrEqual(1), b_.GreaterOrEqual(0)})) {
1077 return false;
1078 }
1079 }
1080 if (integer_trail_->LowerBound(a_) >= 0 &&
1081 integer_trail_->LowerBound(b_) <= 0) {
1082 return integer_trail_->SafeEnqueue(
1083 b_.GreaterOrEqual(1), {a_.GreaterOrEqual(0), p_.GreaterOrEqual(1)});
1084 }
1085 if (integer_trail_->LowerBound(b_) >= 0 &&
1086 integer_trail_->LowerBound(a_) <= 0) {
1087 return integer_trail_->SafeEnqueue(
1088 a_.GreaterOrEqual(1), {b_.GreaterOrEqual(0), p_.GreaterOrEqual(1)});
1089 }
1090 }
1091
1092 for (int i = 0; i < 2; ++i) {
1093 // p = a * b, what is the min/max of a?
1094 const AffineExpression a = i == 0 ? a_ : b_;
1095 const AffineExpression b = i == 0 ? b_ : a_;
1096 const IntegerValue max_b = integer_trail_->UpperBound(b);
1097 const IntegerValue min_b = integer_trail_->LowerBound(b);
1098
1099 // If the domain of b contain zero, we can't propagate anything on a.
1100 // Because of CanonicalizeCases(), we just deal with min_b > 0 here.
1101 if (zero_is_possible && min_b <= 0) continue;
1102
1103 // Here both a and b are across zero, but zero is not possible.
1104 if (min_b < 0 && max_b > 0) {
1105 CHECK_GT(min_p, 0); // Because zero is not possible.
1106
1107 // If a is not across zero, we will deal with this on the next
1108 // Propagate() call.
1109 if (!PropagateMaxOnPositiveProduct(a, b, min_p, max_p)) {
1110 return false;
1111 }
1112 if (!PropagateMaxOnPositiveProduct(a.Negated(), b.Negated(), min_p,
1113 max_p)) {
1114 return false;
1115 }
1116 continue;
1117 }
1118
1119 // This shouldn't happen here.
1120 // If it does, we should reach the fixed point on the next iteration.
1121 if (min_b <= 0) continue;
1122 if (min_p >= 0) {
1123 return integer_trail_->SafeEnqueue(
1124 a.GreaterOrEqual(0), {p_.GreaterOrEqual(0), b.GreaterOrEqual(1)});
1125 }
1126 if (max_p <= 0) {
1127 return integer_trail_->SafeEnqueue(
1128 a.LowerOrEqual(0), {p_.LowerOrEqual(0), b.GreaterOrEqual(1)});
1129 }
1130
1131 // So min_b > 0 and p is across zero: min_p < 0 and max_p > 0.
1132 const IntegerValue new_max_a = FloorRatio(max_p, min_b);
1133 if (new_max_a < integer_trail_->UpperBound(a)) {
1134 if (!integer_trail_->SafeEnqueue(
1135 a.LowerOrEqual(new_max_a),
1136 {integer_trail_->UpperBoundAsLiteral(p_),
1137 integer_trail_->LowerBoundAsLiteral(b)})) {
1138 return false;
1139 }
1140 }
1141 const IntegerValue new_min_a = CeilRatio(min_p, min_b);
1142 if (new_min_a > integer_trail_->LowerBound(a)) {
1143 if (!integer_trail_->SafeEnqueue(
1144 a.GreaterOrEqual(new_min_a),
1145 {integer_trail_->LowerBoundAsLiteral(p_),
1146 integer_trail_->LowerBoundAsLiteral(b)})) {
1147 return false;
1148 }
1149 }
1150 }
1151
1152 return true;
1153}
1154
1156 const int id = watcher->Register(this);
1157 watcher->WatchAffineExpression(a_, id);
1158 watcher->WatchAffineExpression(b_, id);
1159 watcher->WatchAffineExpression(p_, id);
1161}
1162
1164 IntegerTrail* integer_trail)
1165 : x_(x), s_(s), integer_trail_(integer_trail) {
1166 CHECK_GE(integer_trail->LevelZeroLowerBound(x), 0);
1167}
1168
1169// Propagation from x to s: s in [min_x * min_x, max_x * max_x].
1170// Propagation from s to x: x in [ceil(sqrt(min_s)), floor(sqrt(max_s))].
1172 const IntegerValue min_x = integer_trail_->LowerBound(x_);
1173 const IntegerValue min_s = integer_trail_->LowerBound(s_);
1174 const IntegerValue min_x_square = CapProdI(min_x, min_x);
1175 if (min_x_square > min_s) {
1176 if (!integer_trail_->SafeEnqueue(s_.GreaterOrEqual(min_x_square),
1177 {x_.GreaterOrEqual(min_x)})) {
1178 return false;
1179 }
1180 } else if (min_x_square < min_s) {
1181 const IntegerValue new_min(CeilSquareRoot(min_s.value()));
1182 if (!integer_trail_->SafeEnqueue(
1183 x_.GreaterOrEqual(new_min),
1184 {s_.GreaterOrEqual((new_min - 1) * (new_min - 1) + 1)})) {
1185 return false;
1186 }
1187 }
1188
1189 const IntegerValue max_x = integer_trail_->UpperBound(x_);
1190 const IntegerValue max_s = integer_trail_->UpperBound(s_);
1191 const IntegerValue max_x_square = CapProdI(max_x, max_x);
1192 if (max_x_square < max_s) {
1193 if (!integer_trail_->SafeEnqueue(s_.LowerOrEqual(max_x_square),
1194 {x_.LowerOrEqual(max_x)})) {
1195 return false;
1196 }
1197 } else if (max_x_square > max_s) {
1198 const IntegerValue new_max(FloorSquareRoot(max_s.value()));
1199 if (!integer_trail_->SafeEnqueue(
1200 x_.LowerOrEqual(new_max),
1201 {s_.LowerOrEqual(CapProdI(new_max + 1, new_max + 1) - 1)})) {
1202 return false;
1203 }
1204 }
1205
1206 return true;
1207}
1208
1210 const int id = watcher->Register(this);
1211 watcher->WatchAffineExpression(x_, id);
1212 watcher->WatchAffineExpression(s_, id);
1214}
1215
1217 AffineExpression denom,
1218 AffineExpression div,
1219 IntegerTrail* integer_trail)
1220 : num_(num),
1221 denom_(denom),
1222 div_(div),
1223 negated_denom_(denom.Negated()),
1224 negated_num_(num.Negated()),
1225 negated_div_(div.Negated()),
1226 integer_trail_(integer_trail) {}
1227
1228// TODO(user): We can propagate more, especially in the case where denom
1229// spans across 0.
1230// TODO(user): We can propagate a bit more if min_div = 0:
1231// (min_num > -min_denom).
1233 if (integer_trail_->LowerBound(denom_) < 0 &&
1234 integer_trail_->UpperBound(denom_) > 0) {
1235 return true;
1236 }
1237
1238 AffineExpression num = num_;
1239 AffineExpression negated_num = negated_num_;
1240 AffineExpression denom = denom_;
1241 AffineExpression negated_denom = negated_denom_;
1242
1243 if (integer_trail_->UpperBound(denom) < 0) {
1244 std::swap(num, negated_num);
1245 std::swap(denom, negated_denom);
1246 }
1247
1248 if (!PropagateSigns(num, denom, div_)) return false;
1249
1250 if (integer_trail_->UpperBound(num) >= 0 &&
1251 integer_trail_->UpperBound(div_) >= 0 &&
1252 !PropagateUpperBounds(num, denom, div_)) {
1253 return false;
1254 }
1255
1256 if (integer_trail_->UpperBound(negated_num) >= 0 &&
1257 integer_trail_->UpperBound(negated_div_) >= 0 &&
1258 !PropagateUpperBounds(negated_num, denom, negated_div_)) {
1259 return false;
1260 }
1261
1262 if (integer_trail_->LowerBound(num) >= 0 &&
1263 integer_trail_->LowerBound(div_) >= 0) {
1264 return PropagatePositiveDomains(num, denom, div_);
1265 }
1266
1267 if (integer_trail_->LowerBound(negated_num) >= 0 &&
1268 integer_trail_->LowerBound(negated_div_) >= 0) {
1269 return PropagatePositiveDomains(negated_num, denom, negated_div_);
1270 }
1271
1272 return true;
1273}
1274
1275bool DivisionPropagator::PropagateSigns(AffineExpression num,
1276 AffineExpression denom,
1277 AffineExpression div) {
1278 const IntegerValue min_num = integer_trail_->LowerBound(num);
1279 const IntegerValue max_num = integer_trail_->UpperBound(num);
1280 const IntegerValue min_div = integer_trail_->LowerBound(div);
1281 const IntegerValue max_div = integer_trail_->UpperBound(div);
1282
1283 // If num >= 0, as denom > 0, then div must be >= 0.
1284 if (min_num >= 0 && min_div < 0) {
1285 if (!integer_trail_->SafeEnqueue(
1286 div.GreaterOrEqual(0),
1287 {num.GreaterOrEqual(0), denom.GreaterOrEqual(1)})) {
1288 return false;
1289 }
1290 }
1291
1292 // If div > 0, as denom > 0, then num must be > 0.
1293 if (min_num <= 0 && min_div > 0) {
1294 if (!integer_trail_->SafeEnqueue(
1295 num.GreaterOrEqual(1),
1296 {div.GreaterOrEqual(1), denom.GreaterOrEqual(1)})) {
1297 return false;
1298 }
1299 }
1300
1301 // If num <= 0, as denom > 0, then div must be <= 0.
1302 if (max_num <= 0 && max_div > 0) {
1303 if (!integer_trail_->SafeEnqueue(
1304 div.LowerOrEqual(0),
1305 {num.LowerOrEqual(0), denom.GreaterOrEqual(1)})) {
1306 return false;
1307 }
1308 }
1309
1310 // If div < 0, as denom > 0, then num must be < 0.
1311 if (max_num >= 0 && max_div < 0) {
1312 if (!integer_trail_->SafeEnqueue(
1313 num.LowerOrEqual(-1),
1314 {div.LowerOrEqual(-1), denom.GreaterOrEqual(1)})) {
1315 return false;
1316 }
1317 }
1318
1319 return true;
1320}
1321
1322bool DivisionPropagator::PropagateUpperBounds(AffineExpression num,
1323 AffineExpression denom,
1324 AffineExpression div) {
1325 const IntegerValue max_num = integer_trail_->UpperBound(num);
1326 const IntegerValue min_denom = integer_trail_->LowerBound(denom);
1327 const IntegerValue max_denom = integer_trail_->UpperBound(denom);
1328 const IntegerValue max_div = integer_trail_->UpperBound(div);
1329
1330 const IntegerValue new_max_div = max_num / min_denom;
1331 if (max_div > new_max_div) {
1332 if (!integer_trail_->SafeEnqueue(
1333 div.LowerOrEqual(new_max_div),
1334 {num.LowerOrEqual(max_num), denom.GreaterOrEqual(min_denom)})) {
1335 return false;
1336 }
1337 }
1338
1339 // We start from num / denom <= max_div.
1340 // num < (max_div + 1) * denom
1341 // num + 1 <= (max_div + 1) * max_denom.
1342 const IntegerValue new_max_num =
1343 CapAddI(CapProdI(max_div + 1, max_denom), -1);
1344 if (max_num > new_max_num) {
1345 if (!integer_trail_->SafeEnqueue(
1346 num.LowerOrEqual(new_max_num),
1347 {denom.LowerOrEqual(max_denom), denom.GreaterOrEqual(1),
1348 div.LowerOrEqual(max_div)})) {
1349 return false;
1350 }
1351 }
1352
1353 return true;
1354}
1355
1356bool DivisionPropagator::PropagatePositiveDomains(AffineExpression num,
1357 AffineExpression denom,
1358 AffineExpression div) {
1359 const IntegerValue min_num = integer_trail_->LowerBound(num);
1360 const IntegerValue max_num = integer_trail_->UpperBound(num);
1361 const IntegerValue min_denom = integer_trail_->LowerBound(denom);
1362 const IntegerValue max_denom = integer_trail_->UpperBound(denom);
1363 const IntegerValue min_div = integer_trail_->LowerBound(div);
1364 const IntegerValue max_div = integer_trail_->UpperBound(div);
1365
1366 const IntegerValue new_min_div = min_num / max_denom;
1367 if (min_div < new_min_div) {
1368 if (!integer_trail_->SafeEnqueue(
1369 div.GreaterOrEqual(new_min_div),
1370 {num.GreaterOrEqual(min_num), denom.LowerOrEqual(max_denom),
1371 denom.GreaterOrEqual(1)})) {
1372 return false;
1373 }
1374 }
1375
1376 // We start from num / denom >= min_div.
1377 // num >= min_div * denom.
1378 // num >= min_div * min_denom.
1379 const IntegerValue new_min_num = CapProdI(min_denom, min_div);
1380 if (min_num < new_min_num) {
1381 if (!integer_trail_->SafeEnqueue(
1382 num.GreaterOrEqual(new_min_num),
1383 {denom.GreaterOrEqual(min_denom), div.GreaterOrEqual(min_div)})) {
1384 return false;
1385 }
1386 }
1387
1388 // We start with num / denom >= min_div.
1389 // So num >= min_div * denom
1390 // If min_div == 0 we can't deduce anything.
1391 // Otherwise, denom <= num / min_div and denom <= max_num / min_div.
1392 if (min_div > 0) {
1393 const IntegerValue new_max_denom = max_num / min_div;
1394 if (max_denom > new_max_denom) {
1395 if (!integer_trail_->SafeEnqueue(
1396 denom.LowerOrEqual(new_max_denom),
1397 {num.LowerOrEqual(max_num), num.GreaterOrEqual(0),
1398 div.GreaterOrEqual(min_div), denom.GreaterOrEqual(1)})) {
1399 return false;
1400 }
1401 }
1402 }
1403
1404 // denom >= CeilRatio(num + 1, max_div + 1)
1405 // >= CeilRatio(min_num + 1, max_div + 1).
1406 const IntegerValue new_min_denom = CeilRatio(min_num + 1, max_div + 1);
1407 if (min_denom < new_min_denom) {
1408 if (!integer_trail_->SafeEnqueue(
1409 denom.GreaterOrEqual(new_min_denom),
1410 {num.GreaterOrEqual(min_num), div.LowerOrEqual(max_div),
1411 div.GreaterOrEqual(0), denom.GreaterOrEqual(1)})) {
1412 return false;
1413 }
1414 }
1415
1416 return true;
1417}
1418
1420 const int id = watcher->Register(this);
1421 watcher->WatchAffineExpression(num_, id);
1422 watcher->WatchAffineExpression(denom_, id);
1423 watcher->WatchAffineExpression(div_, id);
1425}
1426
1428 IntegerValue b,
1430 IntegerTrail* integer_trail)
1431 : a_(a), b_(b), c_(c), integer_trail_(integer_trail) {
1432 CHECK_GT(b_, 0);
1433}
1434
1436 const IntegerValue min_a = integer_trail_->LowerBound(a_);
1437 const IntegerValue max_a = integer_trail_->UpperBound(a_);
1438 IntegerValue min_c = integer_trail_->LowerBound(c_);
1439 IntegerValue max_c = integer_trail_->UpperBound(c_);
1440
1441 if (max_a / b_ < max_c) {
1442 max_c = max_a / b_;
1443 if (!integer_trail_->SafeEnqueue(
1444 c_.LowerOrEqual(max_c),
1445 {integer_trail_->UpperBoundAsLiteral(a_)})) {
1446 return false;
1447 }
1448 } else if (max_a / b_ > max_c) {
1449 const IntegerValue new_max_a =
1450 max_c >= 0 ? max_c * b_ + b_ - 1 : CapProdI(max_c, b_);
1451 CHECK_LT(new_max_a, max_a);
1452 if (!integer_trail_->SafeEnqueue(
1453 a_.LowerOrEqual(new_max_a),
1454 {integer_trail_->UpperBoundAsLiteral(c_)})) {
1455 return false;
1456 }
1457 }
1458
1459 if (min_a / b_ > min_c) {
1460 min_c = min_a / b_;
1461 if (!integer_trail_->SafeEnqueue(
1462 c_.GreaterOrEqual(min_c),
1463 {integer_trail_->LowerBoundAsLiteral(a_)})) {
1464 return false;
1465 }
1466 } else if (min_a / b_ < min_c) {
1467 const IntegerValue new_min_a =
1468 min_c > 0 ? CapProdI(min_c, b_) : min_c * b_ - b_ + 1;
1469 CHECK_GT(new_min_a, min_a);
1470 if (!integer_trail_->SafeEnqueue(
1471 a_.GreaterOrEqual(new_min_a),
1472 {integer_trail_->LowerBoundAsLiteral(c_)})) {
1473 return false;
1474 }
1475 }
1476
1477 return true;
1478}
1479
1481 const int id = watcher->Register(this);
1482 watcher->WatchAffineExpression(a_, id);
1483 watcher->WatchAffineExpression(c_, id);
1484}
1485
1487 IntegerValue mod,
1488 AffineExpression target,
1489 IntegerTrail* integer_trail)
1490 : expr_(expr), mod_(mod), target_(target), integer_trail_(integer_trail) {
1491 CHECK_GT(mod_, 0);
1492}
1493
1495 if (!PropagateSignsAndTargetRange()) return false;
1496 if (!PropagateOuterBounds()) return false;
1497
1498 if (integer_trail_->LowerBound(expr_) >= 0) {
1499 if (!PropagateBoundsWhenExprIsPositive(expr_, target_)) return false;
1500 } else if (integer_trail_->UpperBound(expr_) <= 0) {
1501 if (!PropagateBoundsWhenExprIsPositive(expr_.Negated(),
1502 target_.Negated())) {
1503 return false;
1504 }
1505 }
1506
1507 return true;
1508}
1509
1510bool FixedModuloPropagator::PropagateSignsAndTargetRange() {
1511 // Initial domain reduction on the target.
1512 if (integer_trail_->UpperBound(target_) >= mod_) {
1513 if (!integer_trail_->SafeEnqueue(target_.LowerOrEqual(mod_ - 1), {})) {
1514 return false;
1515 }
1516 }
1517
1518 if (integer_trail_->LowerBound(target_) <= -mod_) {
1519 if (!integer_trail_->SafeEnqueue(target_.GreaterOrEqual(1 - mod_), {})) {
1520 return false;
1521 }
1522 }
1523
1524 // The sign of target_ is fixed by the sign of expr_.
1525 if (integer_trail_->LowerBound(expr_) >= 0 &&
1526 integer_trail_->LowerBound(target_) < 0) {
1527 if (!integer_trail_->SafeEnqueue(target_.GreaterOrEqual(0),
1528 {expr_.GreaterOrEqual(0)})) {
1529 return false;
1530 }
1531 }
1532
1533 if (integer_trail_->UpperBound(expr_) <= 0 &&
1534 integer_trail_->UpperBound(target_) > 0) {
1535 if (!integer_trail_->SafeEnqueue(target_.LowerOrEqual(0),
1536 {expr_.LowerOrEqual(0)})) {
1537 return false;
1538 }
1539 }
1540
1541 return true;
1542}
1543
1544bool FixedModuloPropagator::PropagateOuterBounds() {
1545 const IntegerValue min_expr = integer_trail_->LowerBound(expr_);
1546 const IntegerValue max_expr = integer_trail_->UpperBound(expr_);
1547 const IntegerValue min_target = integer_trail_->LowerBound(target_);
1548 const IntegerValue max_target = integer_trail_->UpperBound(target_);
1549
1550 if (max_expr % mod_ > max_target) {
1551 if (!integer_trail_->SafeEnqueue(
1552 expr_.LowerOrEqual((max_expr / mod_) * mod_ + max_target),
1553 {integer_trail_->UpperBoundAsLiteral(target_),
1554 integer_trail_->UpperBoundAsLiteral(expr_)})) {
1555 return false;
1556 }
1557 }
1558
1559 if (min_expr % mod_ < min_target) {
1560 if (!integer_trail_->SafeEnqueue(
1561 expr_.GreaterOrEqual((min_expr / mod_) * mod_ + min_target),
1562 {integer_trail_->LowerBoundAsLiteral(expr_),
1563 integer_trail_->LowerBoundAsLiteral(target_)})) {
1564 return false;
1565 }
1566 }
1567
1568 if (min_expr / mod_ == max_expr / mod_) {
1569 if (min_target < min_expr % mod_) {
1570 if (!integer_trail_->SafeEnqueue(
1571 target_.GreaterOrEqual(min_expr - (min_expr / mod_) * mod_),
1572 {integer_trail_->LowerBoundAsLiteral(target_),
1573 integer_trail_->UpperBoundAsLiteral(target_),
1574 integer_trail_->LowerBoundAsLiteral(expr_),
1575 integer_trail_->UpperBoundAsLiteral(expr_)})) {
1576 return false;
1577 }
1578 }
1579
1580 if (max_target > max_expr % mod_) {
1581 if (!integer_trail_->SafeEnqueue(
1582 target_.LowerOrEqual(max_expr - (max_expr / mod_) * mod_),
1583 {integer_trail_->LowerBoundAsLiteral(target_),
1584 integer_trail_->UpperBoundAsLiteral(target_),
1585 integer_trail_->LowerBoundAsLiteral(expr_),
1586 integer_trail_->UpperBoundAsLiteral(expr_)})) {
1587 return false;
1588 }
1589 }
1590 } else if (min_expr / mod_ == 0 && min_target < 0) {
1591 // expr == target when expr <= 0.
1592 if (min_target < min_expr) {
1593 if (!integer_trail_->SafeEnqueue(
1594 target_.GreaterOrEqual(min_expr),
1595 {integer_trail_->LowerBoundAsLiteral(target_),
1596 integer_trail_->LowerBoundAsLiteral(expr_)})) {
1597 return false;
1598 }
1599 }
1600 } else if (max_expr / mod_ == 0 && max_target > 0) {
1601 // expr == target when expr >= 0.
1602 if (max_target > max_expr) {
1603 if (!integer_trail_->SafeEnqueue(
1604 target_.LowerOrEqual(max_expr),
1605 {integer_trail_->UpperBoundAsLiteral(target_),
1606 integer_trail_->UpperBoundAsLiteral(expr_)})) {
1607 return false;
1608 }
1609 }
1610 }
1611
1612 return true;
1613}
1614
1615bool FixedModuloPropagator::PropagateBoundsWhenExprIsPositive(
1616 AffineExpression expr, AffineExpression target) {
1617 const IntegerValue min_target = integer_trail_->LowerBound(target);
1618 DCHECK_GE(min_target, 0);
1619 const IntegerValue max_target = integer_trail_->UpperBound(target);
1620
1621 // The propagation rules below will not be triggered if the domain of target
1622 // covers [0..mod_ - 1].
1623 if (min_target == 0 && max_target == mod_ - 1) return true;
1624
1625 const IntegerValue min_expr = integer_trail_->LowerBound(expr);
1626 const IntegerValue max_expr = integer_trail_->UpperBound(expr);
1627
1628 if (max_expr % mod_ < min_target) {
1629 DCHECK_GE(max_expr, 0);
1630 if (!integer_trail_->SafeEnqueue(
1631 expr.LowerOrEqual((max_expr / mod_ - 1) * mod_ + max_target),
1632 {integer_trail_->UpperBoundAsLiteral(expr),
1633 integer_trail_->LowerBoundAsLiteral(target),
1634 integer_trail_->UpperBoundAsLiteral(target)})) {
1635 return false;
1636 }
1637 }
1638
1639 if (min_expr % mod_ > max_target) {
1640 DCHECK_GE(min_expr, 0);
1641 if (!integer_trail_->SafeEnqueue(
1642 expr.GreaterOrEqual((min_expr / mod_ + 1) * mod_ + min_target),
1643 {integer_trail_->LowerBoundAsLiteral(target),
1644 integer_trail_->UpperBoundAsLiteral(target),
1645 integer_trail_->LowerBoundAsLiteral(expr)})) {
1646 return false;
1647 }
1648 }
1649
1650 return true;
1651}
1652
1654 const int id = watcher->Register(this);
1655 watcher->WatchAffineExpression(expr_, id);
1656 watcher->WatchAffineExpression(target_, id);
1658}
1659
1660} // namespace sat
1661} // namespace operations_research
static int64_t GCD64(int64_t x, int64_t y)
Definition mathutil.h:108
void RegisterWith(GenericLiteralWatcher *watcher)
DivisionPropagator(AffineExpression num, AffineExpression denom, AffineExpression div, IntegerTrail *integer_trail)
void RegisterWith(GenericLiteralWatcher *watcher)
FixedDivisionPropagator(AffineExpression a, IntegerValue b, AffineExpression c, IntegerTrail *integer_trail)
FixedModuloPropagator(AffineExpression expr, IntegerValue mod, AffineExpression target, IntegerTrail *integer_trail)
void RegisterWith(GenericLiteralWatcher *watcher)
void WatchAffineExpression(AffineExpression e, int id)
Definition integer.h:1154
void WatchLiteral(Literal l, int id, int watch_index=-1)
Definition integer.h:1428
void WatchUpperBound(IntegerVariable var, int id, int watch_index=-1)
Definition integer.h:1454
void WatchLowerBound(IntegerVariable var, int id, int watch_index=-1)
Definition integer.h:1436
int Register(PropagatorInterface *propagator)
Registers a propagator and returns its unique ids.
Definition integer.cc:2327
IntegerValue LowerBound(IntegerVariable i) const
Returns the current lower/upper bound of the given integer variable.
Definition integer.h:1301
bool VariableLowerBoundIsFromLevelZero(IntegerVariable var) const
Returns true if the variable lower bound is still the one from level zero.
Definition integer.h:722
ABSL_MUST_USE_RESULT bool SafeEnqueue(IntegerLiteral i_lit, absl::Span< const IntegerLiteral > integer_reason)
Definition integer.cc:1253
IntegerValue UpperBound(IntegerVariable i) const
Definition integer.h:1305
IntegerLiteral LowerBoundAsLiteral(IntegerVariable i) const
Definition integer.h:1344
IntegerValue LevelZeroLowerBound(IntegerVariable var) const
Returns globally valid lower/upper bound on the given integer variable.
Definition integer.h:1396
LevelZeroEquality(IntegerVariable target, const std::vector< IntegerVariable > &vars, const std::vector< IntegerValue > &coeffs, Model *model)
void RegisterWith(GenericLiteralWatcher *watcher)
void Explain(int id, IntegerValue propagation_slack, IntegerVariable var_to_explain, int trail_index, std::vector< Literal > *literals_reason, std::vector< int > *trail_indices_reason) final
For LazyReasonInterface.
LinMinPropagator(std::vector< LinearExpression > exprs, IntegerVariable min_var, Model *model)
std::pair< IntegerValue, IntegerValue > ConditionalLb(IntegerLiteral integer_literal, IntegerVariable target_var) const
NOTE(user): This is only used with int128, so we code only a single version.
void RegisterWith(GenericLiteralWatcher *watcher)
LinearConstraintPropagator(absl::Span< const Literal > enforcement_literals, absl::Span< const IntegerVariable > vars, absl::Span< const IntegerValue > coeffs, IntegerValue upper_bound, Model *model)
void Explain(int id, IntegerValue propagation_slack, IntegerVariable var_to_explain, int trail_index, std::vector< Literal > *literals_reason, std::vector< int > *trail_indices_reason) final
For LazyReasonInterface.
LiteralIndex Index() const
Definition sat_base.h:91
void RegisterWith(GenericLiteralWatcher *watcher)
MinPropagator(std::vector< AffineExpression > vars, AffineExpression min_var, IntegerTrail *integer_trail)
ProductPropagator(AffineExpression a, AffineExpression b, AffineExpression p, IntegerTrail *integer_trail)
void RegisterWith(GenericLiteralWatcher *watcher)
SquarePropagator(AffineExpression x, AffineExpression s, IntegerTrail *integer_trail)
void RegisterWith(GenericLiteralWatcher *watcher)
IntegerValue FloorRatio(IntegerValue dividend, IntegerValue positive_divisor)
constexpr IntegerValue kMaxIntegerValue(std::numeric_limits< IntegerValue::ValueType >::max() - 1)
IntegerValue CeilRatio(IntegerValue dividend, IntegerValue positive_divisor)
const LiteralIndex kNoLiteralIndex(-1)
int64_t FloorSquareRoot(int64_t a)
The argument must be non-negative.
Definition util.cc:300
int64_t CeilSquareRoot(int64_t a)
Definition util.cc:309
std::vector< IntegerVariable > NegationOf(absl::Span< const IntegerVariable > vars)
Returns the vector of the negated variables.
Definition integer.cc:52
constexpr IntegerValue kMinIntegerValue(-kMaxIntegerValue.value())
IntegerVariable PositiveVariable(IntegerVariable i)
IntegerValue PositiveRemainder(IntegerValue dividend, IntegerValue positive_divisor)
IntegerValue CapAddI(IntegerValue a, IntegerValue b)
std::function< int64_t(const Model &)> UpperBound(IntegerVariable v)
Definition integer.h:1545
IntegerValue CapProdI(IntegerValue a, IntegerValue b)
Overflows and saturated arithmetic.
In SWIG mode, we don't want anything besides these top-level includes.
STL namespace.
static int input(yyscan_t yyscanner)
IntegerLiteral GreaterOrEqual(IntegerValue bound) const
var * coeff + constant >= bound.
IntegerLiteral LowerOrEqual(IntegerValue bound) const
var * coeff + constant <= bound.
static IntegerLiteral GreaterOrEqual(IntegerVariable i, IntegerValue bound)
IntegerLiteral Negated() const
The negation of x >= bound is x <= bound - 1.
static IntegerLiteral LowerOrEqual(IntegerVariable i, IntegerValue bound)