[CP-SAT] better management of cuts; better overflow detection; do not…

… run feasibility_jump if search is globally finished; fix #3842

ortools/sat/BUILD.bazel

@@ -1229,6 +1229,7 @@ cc_library(
         ":linear_constraint",
         ":model",
         "//ortools/base",
+        "//ortools/base:cleanup",
         "//ortools/base:mathutil",
         "//ortools/graph",
         "//ortools/graph:max_flow",

ortools/sat/cuts.cc

@@ -274,15 +274,14 @@ void CutDataBuilder::AddOrMergeTerm(const CutTerm& term, IntegerValue t,
     //
     // If we cannot merge the term, we will keep them separate. The produced cut
     // will be less strong, but can still be used.
-    const int64_t new_coeff =
-        CapAdd(cut->terms[entry_index].coeff.value(), term.coeff.value());
-    const int64_t overflow_check = CapProd(t.value(), new_coeff);
-    if (AtMinOrMaxInt64(new_coeff) || AtMinOrMaxInt64(overflow_check)) {
+    const IntegerValue new_coeff =
+        CapAddI(cut->terms[entry_index].coeff, term.coeff);
+    if (AtMinOrMaxInt64I(new_coeff) || ProdOverflow(t, new_coeff)) {
       // If we cannot merge the term, we keep them separate.
       cut->terms.push_back(term);
     } else {
       ++num_merges_;
-      cut->terms[entry_index].coeff = IntegerValue(new_coeff);
+      cut->terms[entry_index].coeff = new_coeff;
     }
   }
 }
@@ -328,22 +327,6 @@ namespace {
 // is needed to avoid numerical issues and adding cuts with minor effect.
 const double kMinCutViolation = 1e-4;
 
-IntegerValue CapProdI(IntegerValue a, IntegerValue b) {
-  return IntegerValue(CapProd(a.value(), b.value()));
-}
-
-IntegerValue CapSubI(IntegerValue a, IntegerValue b) {
-  return IntegerValue(CapSub(a.value(), b.value()));
-}
-
-IntegerValue CapAddI(IntegerValue a, IntegerValue b) {
-  return IntegerValue(CapAdd(a.value(), b.value()));
-}
-
-bool ProdOverflow(IntegerValue t, IntegerValue value) {
-  return AtMinOrMaxInt64(CapProd(t.value(), value.value()));
-}
-
 IntegerValue PositiveRemainder(absl::int128 dividend,
                                IntegerValue positive_divisor) {
   DCHECK_GT(positive_divisor, 0);
@@ -933,7 +916,7 @@ bool IntegerRoundingCutHelper::ComputeCut(
   }
 
   // Use implied bounds to "lift" Booleans into the cut.
-  // This should lead to stronger cuts even if the norms migth be worse.
+  // This should lead to stronger cuts even if the norms might be worse.
   num_ib_used_ = 0;
   if (ib_processor != nullptr) {
     const auto [num_lb, num_ub] = ib_processor->PostprocessWithImpliedBound(
@@ -2643,8 +2626,7 @@ bool BuildMaxAffineUpConstraint(
   // -delta_y * var + delta_x * target <= delta_x * y_at_min - delta_y * x_min
   //
   // Checks the rhs for overflows.
-  if (AtMinOrMaxInt64(CapProd(delta_x.value(), y_at_min.value())) ||
-      AtMinOrMaxInt64(CapProd(delta_y.value(), x_min.value()))) {
+  if (ProdOverflow(delta_x, y_at_min) || ProdOverflow(delta_y, x_min)) {
     return false;
   }
 

ortools/sat/feasibility_jump.cc

@@ -348,6 +348,9 @@ std::function<void()> FeasibilityJumpSolver::GenerateTask(int64_t /*task_id*/) {
     if (shared_bounds_ != nullptr) {
       shared_bounds_->UpdateDomains(&var_domains_);
       for (int var = 0; var < var_domains_.size(); ++var) {
+        // We abort if the problem is trivially UNSAT. This might happen while
+        // we are cleaning up all workers at the end of a search.
+        if (var_domains_[var].IsEmpty()) return;
         var_has_two_values_[var] = var_domains_[var].HasTwoValues();
       }
     }

ortools/sat/integer.h

@@ -105,6 +105,28 @@ inline IntegerValue FloorRatio(IntegerValue dividend,
   return result - adjust;
 }
 
+// Overflows and saturated arithmetic.
+
+inline IntegerValue CapProdI(IntegerValue a, IntegerValue b) {
+  return IntegerValue(CapProd(a.value(), b.value()));
+}
+
+inline IntegerValue CapSubI(IntegerValue a, IntegerValue b) {
+  return IntegerValue(CapSub(a.value(), b.value()));
+}
+
+inline IntegerValue CapAddI(IntegerValue a, IntegerValue b) {
+  return IntegerValue(CapAdd(a.value(), b.value()));
+}
+
+inline bool ProdOverflow(IntegerValue t, IntegerValue value) {
+  return AtMinOrMaxInt64(CapProd(t.value(), value.value()));
+}
+
+inline bool AtMinOrMaxInt64I(IntegerValue t) {
+  return AtMinOrMaxInt64(t.value());
+}
+
 // Returns dividend - FloorRatio(dividend, divisor) * divisor;
 //
 // This function is around the same speed than the computation above, but it
@@ -118,17 +140,21 @@ inline IntegerValue PositiveRemainder(IntegerValue dividend,
   return m < 0 ? m + positive_divisor : m;
 }
 
+inline bool AddTo(IntegerValue a, IntegerValue* result) {
+  if (AtMinOrMaxInt64I(a)) return false;
+  const IntegerValue add = CapAddI(a, *result);
+  if (AtMinOrMaxInt64I(add)) return false;
+  *result = add;
+  return true;
+}
+
 // Computes result += a * b, and return false iff there is an overflow.
 inline bool AddProductTo(IntegerValue a, IntegerValue b, IntegerValue* result) {
-  const int64_t prod = CapProd(a.value(), b.value());
-  if (prod == std::numeric_limits<int64_t>::min() ||
-      prod == std::numeric_limits<int64_t>::max())
-    return false;
-  const int64_t add = CapAdd(prod, result->value());
-  if (add == std::numeric_limits<int64_t>::min() ||
-      add == std::numeric_limits<int64_t>::max())
-    return false;
-  *result = IntegerValue(add);
+  const IntegerValue prod = CapProdI(a, b);
+  if (AtMinOrMaxInt64I(prod)) return false;
+  const IntegerValue add = CapAddI(prod, *result);
+  if (AtMinOrMaxInt64I(add)) return false;
+  *result = add;
   return true;
 }
 

ortools/sat/integer_expr.cc

@@ -851,7 +851,7 @@ bool ProductPropagator::PropagateWhenAllNonNegative() {
   {
     const IntegerValue max_a = integer_trail_->UpperBound(a_);
     const IntegerValue max_b = integer_trail_->UpperBound(b_);
-    const IntegerValue new_max(CapProd(max_a.value(), max_b.value()));
+    const IntegerValue new_max = CapProdI(max_a, max_b);
     if (new_max < integer_trail_->UpperBound(p_)) {
       if (!integer_trail_->SafeEnqueue(
               p_.LowerOrEqual(new_max),
@@ -866,7 +866,7 @@ bool ProductPropagator::PropagateWhenAllNonNegative() {
   {
     const IntegerValue min_a = integer_trail_->LowerBound(a_);
     const IntegerValue min_b = integer_trail_->LowerBound(b_);
-    const IntegerValue new_min(CapProd(min_a.value(), min_b.value()));
+    const IntegerValue new_min = CapProdI(min_a, min_b);
 
     // The conflict test is needed because when new_min is large, we could
     // have an overflow in p_.GreaterOrEqual(new_min);
@@ -893,7 +893,7 @@ bool ProductPropagator::PropagateWhenAllNonNegative() {
     const IntegerValue min_b = integer_trail_->LowerBound(b);
     const IntegerValue min_p = integer_trail_->LowerBound(p_);
     const IntegerValue max_p = integer_trail_->UpperBound(p_);
-    const IntegerValue prod(CapProd(max_a.value(), min_b.value()));
+    const IntegerValue prod = CapProdI(max_a, min_b);
     if (prod > max_p) {
       if (!integer_trail_->SafeEnqueue(a.LowerOrEqual(FloorRatio(max_p, min_b)),
                                        {integer_trail_->LowerBoundAsLiteral(b),
@@ -980,12 +980,12 @@ bool ProductPropagator::Propagate() {
   //
   // TODO(user): In the reasons, including all 4 bounds is always correct, but
   // we might be able to relax some of them.
-  const int64_t max_a = integer_trail_->UpperBound(a_).value();
-  const int64_t max_b = integer_trail_->UpperBound(b_).value();
-  const IntegerValue p1(CapProd(max_a, max_b));
-  const IntegerValue p2(CapProd(max_a, min_b));
-  const IntegerValue p3(CapProd(min_a, max_b));
-  const IntegerValue p4(CapProd(min_a, min_b));
+  const IntegerValue max_a = integer_trail_->UpperBound(a_);
+  const IntegerValue max_b = integer_trail_->UpperBound(b_);
+  const IntegerValue p1 = CapProdI(max_a, max_b);
+  const IntegerValue p2 = CapProdI(max_a, min_b);
+  const IntegerValue p3 = CapProdI(min_a, max_b);
+  const IntegerValue p4 = CapProdI(min_a, min_b);
   const IntegerValue new_max_p = std::max({p1, p2, p3, p4});
   if (new_max_p < integer_trail_->UpperBound(p_)) {
     if (!integer_trail_->SafeEnqueue(
@@ -1124,7 +1124,7 @@ SquarePropagator::SquarePropagator(AffineExpression x, AffineExpression s,
 bool SquarePropagator::Propagate() {
   const IntegerValue min_x = integer_trail_->LowerBound(x_);
   const IntegerValue min_s = integer_trail_->LowerBound(s_);
-  const IntegerValue min_x_square(CapProd(min_x.value(), min_x.value()));
+  const IntegerValue min_x_square = CapProdI(min_x, min_x);
   if (min_x_square > min_s) {
     if (!integer_trail_->SafeEnqueue(s_.GreaterOrEqual(min_x_square),
                                      {x_.GreaterOrEqual(min_x)})) {
@@ -1141,7 +1141,7 @@ bool SquarePropagator::Propagate() {
 
   const IntegerValue max_x = integer_trail_->UpperBound(x_);
   const IntegerValue max_s = integer_trail_->UpperBound(s_);
-  const IntegerValue max_x_square(CapProd(max_x.value(), max_x.value()));
+  const IntegerValue max_x_square = CapProdI(max_x, max_x);
   if (max_x_square < max_s) {
     if (!integer_trail_->SafeEnqueue(s_.LowerOrEqual(max_x_square),
                                      {x_.LowerOrEqual(max_x)})) {
@@ -1151,9 +1151,7 @@ bool SquarePropagator::Propagate() {
     const IntegerValue new_max(FloorSquareRoot(max_s.value()));
     if (!integer_trail_->SafeEnqueue(
             x_.LowerOrEqual(new_max),
-            {s_.LowerOrEqual(IntegerValue(CapProd(new_max.value() + 1,
-                                                  new_max.value() + 1)) -
-                             1)})) {
+            {s_.LowerOrEqual(CapProdI(new_max + 1, new_max + 1) - 1)})) {
       return false;
     }
   }
@@ -1273,7 +1271,7 @@ bool DivisionPropagator::PropagateUpperBounds(AffineExpression num,
   // num < (max_div + 1) * denom
   // num + 1 <= (max_div + 1) * max_denom.
   const IntegerValue new_max_num =
-      IntegerValue(CapAdd(CapProd(max_div.value() + 1, max_denom.value()), -1));
+      CapAddI(CapProdI(max_div + 1, max_denom), -1);
   if (max_num > new_max_num) {
     if (!integer_trail_->SafeEnqueue(
             num.LowerOrEqual(new_max_num),
@@ -1309,8 +1307,7 @@ bool DivisionPropagator::PropagatePositiveDomains(AffineExpression num,
   // We start from num / denom >= min_div.
   // num >= min_div * denom.
   // num >= min_div * min_denom.
-  const IntegerValue new_min_num =
-      IntegerValue(CapProd(min_denom.value(), min_div.value()));
+  const IntegerValue new_min_num = CapProdI(min_denom, min_div);
   if (min_num < new_min_num) {
     if (!integer_trail_->SafeEnqueue(
             num.GreaterOrEqual(new_min_num),
@@ -1382,8 +1379,7 @@ bool FixedDivisionPropagator::Propagate() {
     }
   } else if (max_a / b_ > max_c) {
     const IntegerValue new_max_a =
-        max_c >= 0 ? max_c * b_ + b_ - 1
-                   : IntegerValue(CapProd(max_c.value(), b_.value()));
+        max_c >= 0 ? max_c * b_ + b_ - 1 : CapProdI(max_c, b_);
     CHECK_LT(new_max_a, max_a);
     if (!integer_trail_->SafeEnqueue(
             a_.LowerOrEqual(new_max_a),
@@ -1401,8 +1397,7 @@ bool FixedDivisionPropagator::Propagate() {
     }
   } else if (min_a / b_ < min_c) {
     const IntegerValue new_min_a =
-        min_c > 0 ? IntegerValue(CapProd(min_c.value(), b_.value()))
-                  : min_c * b_ - b_ + 1;
+        min_c > 0 ? CapProdI(min_c, b_) : min_c * b_ - b_ + 1;
     CHECK_GT(new_min_a, min_a);
     if (!integer_trail_->SafeEnqueue(
             a_.GreaterOrEqual(new_min_a),

ortools/sat/linear_constraint_manager.cc

@@ -46,6 +46,28 @@
 namespace operations_research {
 namespace sat {
 
+bool PossibleOverflow(const IntegerTrail& integer_trail,
+                      const LinearConstraint& constraint) {
+  IntegerValue min_activity(0);
+  IntegerValue max_activity(0);
+  const int size = constraint.vars.size();
+  for (int i = 0; i < size; ++i) {
+    const IntegerVariable var = constraint.vars[i];
+    const IntegerValue coeff = constraint.coeffs[i];
+    CHECK_NE(coeff, 0);
+    const IntegerValue lb = integer_trail.LevelZeroLowerBound(var);
+    const IntegerValue ub = integer_trail.LevelZeroUpperBound(var);
+    if (coeff > 0) {
+      if (!AddProductTo(lb, coeff, &min_activity)) return true;
+      if (!AddProductTo(ub, coeff, &max_activity)) return true;
+    } else {
+      if (!AddProductTo(ub, coeff, &min_activity)) return true;
+      if (!AddProductTo(lb, coeff, &max_activity)) return true;
+    }
+  }
+  return AtMinOrMaxInt64(CapSub(max_activity.value(), min_activity.value()));
+}
+
 namespace {
 
 const LinearConstraintManager::ConstraintIndex kInvalidConstraintIndex(-1);
@@ -290,6 +312,12 @@ bool LinearConstraintManager::AddCut(const LinearConstraint& ct,
     return false;
   }
 
+  // TODO(user): We could prevent overflow by dividing more. Note that mainly
+  // happen with super large variable domain since we usually restrict the size
+  // of the generated coefficients in our cuts. So it shouldn't be that
+  // important.
+  if (PossibleOverflow(integer_trail_, ct)) return false;
+
   bool added = false;
   const ConstraintIndex ct_index = Add(ct, &added);
 

ortools/sat/linear_constraint_manager.h

@@ -41,6 +41,11 @@
 namespace operations_research {
 namespace sat {
 
+// Tests for possible overflow in the simplification of the given linear
+// constraint.
+bool PossibleOverflow(const IntegerTrail& integer_trail,
+                      const LinearConstraint& constraint);
+
 // Stores for each IntegerVariable its temporary LP solution.
 //
 // This is shared between all LinearProgrammingConstraint because in the corner