[CP-SAT] fixes and speed improvements

ortools/sat/BUILD.bazel

@@ -969,6 +969,7 @@ cc_library(
         ":implied_bounds",
         ":integer",
         ":intervals",
+        ":linear_constraint_manager",
         ":linear_programming_constraint",
         ":model",
         ":probing",

ortools/sat/clause.cc

@@ -498,6 +498,7 @@ void ClauseManager::DeleteRemovedClauses() {
 void BinaryImplicationGraph::Resize(int num_variables) {
   SCOPED_TIME_STAT(&stats_);
   implications_.resize(num_variables << 1);
+  implies_something_.resize(num_variables << 1);
   might_have_dups_.resize(num_variables << 1);
   is_redundant_.resize(implications_.size());
   is_removed_.resize(implications_.size(), false);
@@ -550,6 +551,8 @@ bool BinaryImplicationGraph::AddBinaryClause(Literal a, Literal b) {
   estimated_sizes_[b.NegatedIndex()]++;
   implications_[a.NegatedIndex()].push_back(b);
   implications_[b.NegatedIndex()].push_back(a);
+  implies_something_.Set(a.NegatedIndex());
+  implies_something_.Set(b.NegatedIndex());
   NotifyPossibleDuplicate(a);
   NotifyPossibleDuplicate(b);
   is_dag_ = false;
@@ -724,6 +727,7 @@ bool BinaryImplicationGraph::CleanUpAndAddAtMostOnes(int base_index) {
         for (const Literal b : at_most_one) {
           if (a == b) continue;
           implications_[a].push_back(b.Negated());
+          implies_something_.Set(a);
           NotifyPossibleDuplicate(a);
         }
       }
@@ -741,6 +745,7 @@ bool BinaryImplicationGraph::CleanUpAndAddAtMostOnes(int base_index) {
       }
       DCHECK(!is_redundant_[l]);
       at_most_ones_[l].push_back(local_start);
+      implies_something_.Set(l);
     }
 
     // Add sentinel.
@@ -2056,6 +2061,8 @@ void BinaryImplicationGraph::MarkDescendants(Literal root) {
   is_marked_.Set(root);
   for (int j = 0; j < stack_size; ++j) {
     const Literal current = stack[j];
+    if (!implies_something_[current]) continue;
+
     for (const Literal l : implications_[current]) {
       if (!is_marked[l] && !is_redundant[l]) {
         is_marked_.SetUnsafe(l);

ortools/sat/clause.h

@@ -862,6 +862,14 @@ class BinaryImplicationGraph : public SatPropagator {
   const int at_most_one_max_expansion_size_;
   int at_most_one_iterator_ = 0;
 
+  // Invariant: implies_something_[l] should be true iff implications_[l] or
+  // at_most_ones_[l] might be non-empty.
+  //
+  // For problems with a large number of variables and sparse implications_ or
+  // at_most_ones_ entries, checking this is way faster during
+  // MarkDescendants(). See for instance proteindesign122trx11p8.pb.gz.
+  Bitset64<LiteralIndex> implies_something_;
+
   // Used by GenerateAtMostOnesWithLargeWeight().
   std::vector<std::vector<Literal>> tmp_cuts_;
 

ortools/sat/cp_model_presolve.cc

@@ -9806,7 +9806,7 @@ void CpModelPresolver::FindBigVerticalLinearOverlap(
       // Skip bad constraint.
       if (common_part.size() < 2) continue;
 
-      // Update common part.
+      // Update coeff_map.
       block.push_back({c, x_coeff});
       coeff_map.clear();
       for (const auto [var, coeff] : common_part) {
@@ -9816,8 +9816,8 @@ void CpModelPresolver::FindBigVerticalLinearOverlap(
 
     // We have a candidate.
     const int64_t saved_nz =
-        ComputeNonZeroReduction(block.size(), common_part.size());
-    if (saved_nz < 100) continue;
+        ComputeNonZeroReduction(block.size(), coeff_map.size());
+    if (saved_nz < 30) continue;
 
     // Fix multiples, currently this contain the coeff of x for each constraint.
     const int64_t base_x = coeff_map.at(x);
@@ -9826,7 +9826,7 @@ void CpModelPresolver::FindBigVerticalLinearOverlap(
       multipier /= base_x;
     }
 
-    // Introduce new_var = common_part and perform the substitution.
+    // Introduce new_var = coeff_map and perform the substitution.
     if (!RemoveCommonPart(coeff_map, block, helper)) continue;
     ++num_blocks;
     nz_reduction += saved_nz;
@@ -12306,8 +12306,12 @@ CpSolverStatus CpModelPresolver::Presolve() {
       ActivityBoundHelper activity_amo_helper;
       activity_amo_helper.AddAllAtMostOnes(*context_->working_model);
       FindBigAtMostOneAndLinearOverlap(&activity_amo_helper);
-      FindBigHorizontalLinearOverlap(&activity_amo_helper);
+
+      // Heuristic: vertical introduce smaller defining constraint and appear in
+      // many constraints, so might be more constrained. We might also still
+      // make horizontal rectangle with the variable introduced.
       FindBigVerticalLinearOverlap(&activity_amo_helper);
+      FindBigHorizontalLinearOverlap(&activity_amo_helper);
     }
     if (context_->ModelIsUnsat()) return InfeasibleStatus();
 

ortools/sat/cp_model_solver.cc

@@ -4013,7 +4013,7 @@ CpSolverResponse SolveCpModel(const CpModelProto& model_proto, Model* model) {
     SOLVER_LOG(
         logger,
         "Warning: solving with assumptions was requested in a non-fully "
-        "supported setting.\nWe will assume these assumptions true while "
+        "supported setting.\nWe will assumes these assumptions true while "
         "solving, but if the model is infeasible, you will not get a useful "
         "'sufficient_assumptions_for_infeasibility' field in the response, it "
         "will include all assumptions.");

ortools/sat/cp_model_symmetries.cc

@@ -323,6 +323,40 @@ std::unique_ptr<Graph> GenerateGraphForSymmetryDetection(
         }
         break;
       }
+      case ConstraintProto::kLinMax: {
+        const LinearExpressionProto& target_expr =
+            constraint.lin_max().target();
+
+        std::vector<int64_t> target_color = color;
+        target_color.push_back(target_expr.offset());
+        const int target_node = new_node(target_color);
+
+        for (int j = 0; j < target_expr.vars_size(); ++j) {
+          const int var = target_expr.vars(j);
+          DCHECK(RefIsPositive(var));
+          const int64_t coeff = target_expr.coeffs(j);
+          graph->AddArc(get_coefficient_node(var, coeff), target_node);
+        }
+
+        for (int i = 0; i < constraint.lin_max().exprs_size(); ++i) {
+          const LinearExpressionProto& expr = constraint.lin_max().exprs(i);
+
+          std::vector<int64_t> local_color = color;
+          local_color.push_back(expr.offset());
+          const int local_node = new_node(local_color);
+
+          for (int j = 0; j < expr.vars().size(); ++j) {
+            const int var = expr.vars(j);
+            DCHECK(RefIsPositive(var));
+            const int64_t coeff = expr.coeffs(j);
+            graph->AddArc(get_coefficient_node(var, coeff), local_node);
+          }
+
+          graph->AddArc(local_node, target_node);
+        }
+
+        break;
+      }
       case ConstraintProto::kInterval: {
         // We create 3 constraint nodes (for start, size and end) including the
         // offset. We connect these to their terms like for a linear constraint.

ortools/sat/integer_search.cc

@@ -34,6 +34,7 @@
 #include "ortools/sat/implied_bounds.h"
 #include "ortools/sat/integer.h"
 #include "ortools/sat/intervals.h"
+#include "ortools/sat/linear_constraint_manager.h"
 #include "ortools/sat/linear_programming_constraint.h"
 #include "ortools/sat/model.h"
 #include "ortools/sat/probing.h"
@@ -177,6 +178,28 @@ std::function<BooleanOrIntegerLiteral()> FirstUnassignedVarAtItsMinHeuristic(
   };
 }
 
+std::function<BooleanOrIntegerLiteral()> MostFractionalHeuristic(Model* model) {
+  auto* lp_values = model->GetOrCreate<ModelLpValues>();
+  auto* integer_trail = model->GetOrCreate<IntegerTrail>();
+  return [lp_values, integer_trail, model]() {
+    double best_fractionality = 0.0;
+    BooleanOrIntegerLiteral decision;
+    for (IntegerVariable var(0); var < lp_values->size(); var += 2) {
+      if (integer_trail->IsFixed(var)) continue;
+      const double lp_value = (*lp_values)[var];
+      const double fractionality = std::abs(lp_value - std::round(lp_value));
+      if (fractionality > best_fractionality) {
+        best_fractionality = fractionality;
+
+        // This choose <= value if possible.
+        decision = BooleanOrIntegerLiteral(SplitAroundGivenValue(
+            var, IntegerValue(std::floor(lp_value)), model));
+      }
+    }
+    return decision;
+  };
+}
+
 std::function<BooleanOrIntegerLiteral()>
 UnassignedVarWithLowestMinAtItsMinHeuristic(
     const std::vector<IntegerVariable>& vars, Model* model) {

ortools/sat/integer_search.h

@@ -171,6 +171,9 @@ IntegerLiteral SplitDomainUsingBestSolutionValue(IntegerVariable var,
 std::function<BooleanOrIntegerLiteral()> FirstUnassignedVarAtItsMinHeuristic(
     const std::vector<IntegerVariable>& vars, Model* model);
 
+// Choose the variable with most fractional LP value.
+std::function<BooleanOrIntegerLiteral()> MostFractionalHeuristic(Model* model);
+
 // Decision heuristic for SolveIntegerProblemWithLazyEncoding(). Like
 // FirstUnassignedVarAtItsMinHeuristic() but the function will return the
 // literal corresponding to the fact that the currently non-assigned variable

ortools/sat/lb_tree_search.cc

@@ -80,9 +80,10 @@ LbTreeSearch::LbTreeSearch(Model* model)
   // We use the normal SAT search but we will bump the variable activity
   // slightly differently. In addition to the conflicts, we also bump it each
   // time the objective lower bound increase in a sub-node.
-  search_heuristic_ =
-      SequentialSearch({SatSolverHeuristic(model),
-                        model->GetOrCreate<SearchHeuristics>()->fixed_search});
+  search_heuristic_ = SequentialSearch(
+      {SatSolverHeuristic(model), MostFractionalHeuristic(model),
+       IntegerValueSelectionHeuristic(
+           model->GetOrCreate<SearchHeuristics>()->fixed_search, model)});
 }
 
 void LbTreeSearch::UpdateParentObjective(int level) {