[CP-SAT] more work on clause sharing; optimize linear code; better pr…

…ocessing of hints; fill_tightened_domains_in_response now works after a normal search; do not use inprocessing in hint search

ortools/sat/BUILD.bazel

@@ -600,6 +600,7 @@ cc_library(
         "//ortools/base",
         "//ortools/base:types",
         "//ortools/port:proto_utils",
+        "//ortools/util:logging",
         "//ortools/util:sorted_interval_list",
         "@com_google_absl//absl/log:check",
     ],
@@ -1509,6 +1510,7 @@ cc_library(
         "@com_google_absl//absl/container:btree",
         "@com_google_absl//absl/container:flat_hash_map",
         "@com_google_absl//absl/container:flat_hash_set",
+        "@com_google_absl//absl/log",
         "@com_google_absl//absl/log:check",
         "@com_google_absl//absl/meta:type_traits",
         "@com_google_absl//absl/numeric:int128",
@@ -2072,6 +2074,7 @@ cc_library(
         "@com_google_absl//absl/strings",
         "@com_google_absl//absl/synchronization",
         "@com_google_absl//absl/time",
+        "@com_google_absl//absl/types:span",
     ],
 )
 
@@ -2146,6 +2149,7 @@ cc_binary(
         "//ortools/base:path",
         "//ortools/util:file_util",
         "//ortools/util:filelineiter",
+        "//ortools/util:sorted_interval_list",
         "@com_google_absl//absl/flags:flag",
         "@com_google_absl//absl/log",
         "@com_google_absl//absl/log:check",

ortools/sat/cp_model.proto

@@ -698,7 +698,7 @@ enum CpSolverStatus {
   // An optimal feasible solution has been found.
   //
   // More generally, this status represent a success. So we also return OPTIMAL
-  // if we find a solution for a pure feasiblity problem or if a gap limit has
+  // if we find a solution for a pure feasibility problem or if a gap limit has
   // been specified and we return a solution within this limit. In the case
   // where we need to return all the feasible solution, this status will only be
   // returned if we enumerated all of them; If we stopped before, we will return
@@ -752,10 +752,14 @@ message CpSolverResponse {
   // is only filled with the info derived during a normal search and we do not
   // have any dedicated algorithm to improve it.
   //
-  // If the problem is a feasibility problem, then these bounds will be valid
-  // for any feasible solution. If the problem is an optimization problem, then
-  // these bounds will only be valid for any OPTIMAL solutions, it can exclude
-  // sub-optimal feasible ones.
+  // Warning: if you didn't set keep_all_feasible_solutions_in_presolve, then
+  // these domains might exclude valid feasible solution. Otherwise for a
+  // feasibility problem, all feasible solution should be there.
+  //
+  // Warning: For an optimization problem, these will correspond to valid bounds
+  // for the problem of finding an improving solution to the best one found so
+  // far. It might be better to solve a feasibility version if one just want to
+  // explore the feasible region.
   repeated IntegerVariableProto tightened_variables = 21;
 
   // A subset of the model "assumptions" field. This will only be filled if the

ortools/sat/cp_model_expand.cc

@@ -466,7 +466,7 @@ void ExpandLinMax(ConstraintProto* ct, PresolveContext* context) {
   const int num_exprs = ct->lin_max().exprs().size();
   if (num_exprs < 2) return;
 
-  // We have a special treatment for Abs, Earlyness, Tardiness, and all
+  // We have a special treatment for Abs, Earliness, Tardiness, and all
   // affine_max where there is only one variable present in all the expressions.
   if (ExpressionsContainsOnlyOneVar(ct->lin_max().exprs())) return;
 
@@ -476,7 +476,7 @@ void ExpandLinMax(ConstraintProto* ct, PresolveContext* context) {
   // - target >= ai
   for (const LinearExpressionProto& expr : ct->lin_max().exprs()) {
     ConstraintProto* new_ct = context->working_model->add_constraints();
-    LinearConstraintProto* lin = ct->mutable_linear();
+    LinearConstraintProto* lin = new_ct->mutable_linear();
     lin->add_domain(0);
     lin->add_domain(std::numeric_limits<int64_t>::max());
     AddLinearExpressionToLinearConstraint(ct->lin_max().target(), 1, lin);

ortools/sat/cp_model_postsolve.cc

@@ -23,6 +23,7 @@
 #include "ortools/port/proto_utils.h"
 #include "ortools/sat/cp_model.pb.h"
 #include "ortools/sat/cp_model_utils.h"
+#include "ortools/util/logging.h"
 #include "ortools/util/sorted_interval_list.h"
 
 namespace operations_research {
@@ -406,5 +407,127 @@ void PostsolveResponse(const int64_t num_variables_in_original_model,
   }
 }
 
+void FillTightenedDomainInResponse(const CpModelProto& original_model,
+                                   const CpModelProto& mapping_proto,
+                                   const std::vector<int>& postsolve_mapping,
+                                   const std::vector<Domain>& search_domains,
+                                   CpSolverResponse* response,
+                                   SolverLogger* logger) {
+  // The [0, num_vars) part will contain the tightened domains.
+  const int num_original_vars = original_model.variables().size();
+  const int num_expanded_vars = mapping_proto.variables().size();
+  CHECK_LE(num_original_vars, num_expanded_vars);
+  std::vector<Domain> domains(num_expanded_vars);
+
+  // Start with the domain from the mapping proto. Note that by construction
+  // this should be tighter than the original variable domains.
+  for (int i = 0; i < num_expanded_vars; ++i) {
+    domains[i] = ReadDomainFromProto(mapping_proto.variables(i));
+    if (i < num_original_vars) {
+      CHECK(domains[i].IsIncludedIn(
+          ReadDomainFromProto(original_model.variables(i))));
+    }
+  }
+
+  // The first test is for the corner case of presolve closing the problem,
+  // in which case there is no more info to process.
+  int num_common_vars = 0;
+  int num_affine_reductions = 0;
+  if (!search_domains.empty()) {
+    if (postsolve_mapping.empty()) {
+      // Currently no mapping should mean all variables are in common. This
+      // happen when presolve is disabled, but we might still have more
+      // variables due to expansion for instance.
+      //
+      // There is also the corner case of presolve closing the problem,
+      CHECK_GE(search_domains.size(), num_original_vars);
+      num_common_vars = num_original_vars;
+      for (int i = 0; i < num_original_vars; ++i) {
+        domains[i] = domains[i].IntersectionWith(search_domains[i]);
+      }
+    } else {
+      // This is the normal presolve case.
+      // Intersect the domain of the variables in common.
+      CHECK_EQ(postsolve_mapping.size(), search_domains.size());
+      for (int search_i = 0; search_i < postsolve_mapping.size(); ++search_i) {
+        const int i_in_mapping_model = postsolve_mapping[search_i];
+        if (i_in_mapping_model < num_original_vars) {
+          ++num_common_vars;
+        }
+        domains[i_in_mapping_model] =
+            domains[i_in_mapping_model].IntersectionWith(
+                search_domains[search_i]);
+      }
+
+      // Look for affine relation, and do more intersection.
+      for (const ConstraintProto& ct : mapping_proto.constraints()) {
+        if (ct.constraint_case() != ConstraintProto::kLinear) continue;
+        const LinearConstraintProto& lin = ct.linear();
+        if (lin.vars().size() != 2) continue;
+        if (lin.domain().size() != 2) continue;
+        if (lin.domain(0) != lin.domain(1)) continue;
+        int v1 = lin.vars(0);
+        int v2 = lin.vars(1);
+        int c1 = lin.coeffs(0);
+        int c2 = lin.coeffs(1);
+        if (v2 < num_original_vars && v1 >= num_original_vars) {
+          std::swap(v1, v2);
+          std::swap(c1, c2);
+        }
+        if (v1 < num_original_vars && v2 >= num_original_vars) {
+          // We can reduce the domain of v1 by using the affine relation
+          // and the domain of v2.
+          // We have c1 * v2 + c2 * v2 = offset;
+          const int64_t offset = lin.domain(0);
+          const Domain restriction =
+              Domain(offset)
+                  .AdditionWith(domains[v2].ContinuousMultiplicationBy(-c2))
+                  .InverseMultiplicationBy(c1);
+          if (!domains[v1].IsIncludedIn(restriction)) {
+            ++num_affine_reductions;
+            domains[v1] = domains[v1].IntersectionWith(restriction);
+          }
+        }
+      }
+    }
+  }
+
+  // Copy the names and replace domains.
+  *response->mutable_tightened_variables() = original_model.variables();
+  int num_tigher_domains = 0;
+  int num_empty = 0;
+  int num_fixed = 0;
+  for (int i = 0; i < num_original_vars; ++i) {
+    FillDomainInProto(domains[i], response->mutable_tightened_variables(i));
+    if (domains[i].IsEmpty()) {
+      ++num_empty;
+      continue;
+    }
+
+    if (domains[i].IsFixed()) num_fixed++;
+    const Domain original = ReadDomainFromProto(original_model.variables(i));
+    if (domains[i] != original) {
+      DCHECK(domains[i].IsIncludedIn(original));
+      ++num_tigher_domains;
+    }
+  }
+
+  // Some stats.
+  if (num_empty > 0) {
+    SOLVER_LOG(logger, num_empty,
+               " tightened domains are empty. This should not happen except if "
+               "we proven infeasibility or optimality.");
+  }
+  SOLVER_LOG(logger, "Filled tightened domains in the response.");
+  SOLVER_LOG(logger, "[TighteningInfo] num_tighter:", num_tigher_domains,
+             " num_fixed:", num_fixed,
+             " num_affine_reductions:", num_affine_reductions);
+  SOLVER_LOG(logger,
+             "[TighteningInfo] original_num_variables:", num_original_vars,
+             " during_presolve:", num_expanded_vars,
+             " after:", search_domains.size(), " in_common:", num_common_vars);
+  SOLVER_LOG(logger, "");
+}
+
 }  // namespace sat
 }  // namespace operations_research

ortools/sat/cp_model_postsolve.h

@@ -19,6 +19,8 @@
 
 #include "ortools/base/types.h"
 #include "ortools/sat/cp_model.pb.h"
+#include "ortools/util/logging.h"
+#include "ortools/util/sorted_interval_list.h"
 
 namespace operations_research {
 namespace sat {
@@ -46,6 +48,16 @@ void PostsolveResponse(int64_t num_variables_in_original_model,
                        const std::vector<int>& postsolve_mapping,
                        std::vector<int64_t>* solution);
 
+// Try to postsolve with a "best-effort" the reduced domain from the presolved
+// model to the user given model. See the documentation of the CpSolverResponse
+// tightened_variables field for more information on the caveats.
+void FillTightenedDomainInResponse(const CpModelProto& original_model,
+                                   const CpModelProto& mapping_proto,
+                                   const std::vector<int>& postsolve_mapping,
+                                   const std::vector<Domain>& search_domains,
+                                   CpSolverResponse* response,
+                                   SolverLogger* logger);
+
 }  // namespace sat
 }  // namespace operations_research
 

ortools/sat/cp_model_presolve.cc

@@ -7284,8 +7284,13 @@ bool CpModelPresolver::PresolvePureSatPart() {
   // for blocked clause. It should be possible to allow for this by adding extra
   // variable to the mapping model at presolve and some linking constraints, but
   // this is messy.
+  //
+  // We also disable this if the user asked for tightened domain as this might
+  // fix variable to a potentially infeasible value, and just correct them later
+  // during postsolve of a particular solution.
   SatParameters params = context_->params();
-  if (params.debug_postsolve_with_full_solver()) {
+  if (params.debug_postsolve_with_full_solver() ||
+      params.fill_tightened_domains_in_response()) {
     params.set_presolve_blocked_clause(false);
   }
 
@@ -12352,6 +12357,19 @@ CpSolverStatus CpModelPresolver::InfeasibleStatus() {
   return CpSolverStatus::INFEASIBLE;
 }
 
+void CpModelPresolver::InitializeMappingModelVariables() {
+  // Sync the domains.
+  for (int i = 0; i < context_->working_model->variables_size(); ++i) {
+    FillDomainInProto(context_->DomainOf(i),
+                      context_->working_model->mutable_variables(i));
+    DCHECK_GT(context_->working_model->variables(i).domain_size(), 0);
+  }
+
+  // Set the variables of the mapping_model.
+  context_->mapping_model->mutable_variables()->CopyFrom(
+      context_->working_model->variables());
+}
+
 // The presolve works as follow:
 //
 // First stage:
@@ -12431,6 +12449,13 @@ CpSolverStatus CpModelPresolver::Presolve() {
 
     // We need to append all the variable equivalence that are still used!
     EncodeAllAffineRelations();
+
+    // Make sure we also have an initialized mapping model as we use this for
+    // filling the tightened variables. Even without presolve, we do some
+    // trivial presolving during the initial copy of the model, and expansion
+    // might do more.
+    InitializeMappingModelVariables();
+
     if (logger_->LoggingIsEnabled()) context_->LogInfo();
     return CpSolverStatus::UNKNOWN;
   }
@@ -12660,16 +12685,8 @@ CpSolverStatus CpModelPresolver::Presolve() {
     google::protobuf::util::Truncate(strategy.mutable_exprs(), new_size);
   }
 
-  // Sync the domains.
-  for (int i = 0; i < context_->working_model->variables_size(); ++i) {
-    FillDomainInProto(context_->DomainOf(i),
-                      context_->working_model->mutable_variables(i));
-    DCHECK_GT(context_->working_model->variables(i).domain_size(), 0);
-  }
-
-  // Set the variables of the mapping_model.
-  context_->mapping_model->mutable_variables()->CopyFrom(
-      context_->working_model->variables());
+  // Sync the domains and initialize the mapping model variables.
+  InitializeMappingModelVariables();
 
   // Remove all the unused variables from the presolved model.
   postsolve_mapping_->clear();

ortools/sat/cp_model_presolve.h

@@ -91,6 +91,11 @@ class CpModelPresolver {
   // A simple helper that logs the rules applied so far and return INFEASIBLE.
   CpSolverStatus InfeasibleStatus();
 
+  // At the end of presolve, the mapping model is initialized to contains all
+  // the variable from the original model + the one created during presolve
+  // expand. It also contains the tightened domains.
+  void InitializeMappingModelVariables();
+
   // Runs the inner loop of the presolver.
   bool ProcessChangedVariables(std::vector<bool>* in_queue,
                                std::deque<int>* queue);

ortools/sat/cp_model_search.cc

@@ -687,9 +687,13 @@ absl::flat_hash_map<std::string, SatParameters> GetNamedParameters(
     SatParameters new_params = base_params;
     new_params.set_stop_after_first_solution(false);
     new_params.set_cp_model_presolve(true);
+
+    // We disable costly presolve/inprocessing.
+    new_params.set_use_sat_inprocessing(false);
     new_params.set_cp_model_probing_level(0);
     new_params.set_symmetry_level(0);
     new_params.set_find_big_linear_overlap(false);
+
     new_params.set_log_search_progress(false);
     new_params.set_debug_crash_on_bad_hint(false);  // Can happen in lns.
     new_params.set_solution_pool_size(1);  // Keep the best solution found.