1 // Copyright 2010-2021 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 14 #include "ortools/sat/cp_model_solver.h" 15 16 #include <algorithm> 17 #include <atomic> 18 #include <cmath> 19 #include <cstdint> 20 #include <functional> 21 #include <limits> 22 #include <map> 23 #include <memory> 24 #include <random> 25 #include <set> 26 #include <string> 27 #include <thread> 28 #include <utility> 29 #include <vector> 30 31 #if !defined(__PORTABLE_PLATFORM__) 32 #include "absl/synchronization/notification.h" 33 #include "google/protobuf/text_format.h" 34 #include "ortools/base/file.h" 35 #include "ortools/util/sigint.h" 36 #endif // __PORTABLE_PLATFORM__ 37 38 #include "absl/base/attributes.h" 39 #include "absl/container/flat_hash_set.h" 40 #include "absl/memory/memory.h" 41 #include "absl/status/status.h" 42 #include "absl/strings/str_cat.h" 43 #include "absl/strings/str_format.h" 44 #include "absl/strings/str_join.h" 45 #include "absl/strings/str_split.h" 46 #include "absl/synchronization/mutex.h" 47 #include "ortools/base/cleanup.h" 48 #include "ortools/base/commandlineflags.h" 49 #include "ortools/base/int_type.h" 50 #include "ortools/base/integral_types.h" 51 #include "ortools/base/logging.h" 52 #include "ortools/base/map_util.h" 53 #include "ortools/base/threadpool.h" 54 #include "ortools/base/timer.h" 55 #include "ortools/base/version.h" 56 #include "ortools/base/vlog_is_on.h" 57 #include "ortools/graph/connected_components.h" 58 #include "ortools/port/proto_utils.h" 59 #include "ortools/sat/circuit.h" 60 #include "ortools/sat/clause.h" 61 #include "ortools/sat/cp_model.pb.h" 62 #include "ortools/sat/cp_model_checker.h" 63 #include "ortools/sat/cp_model_lns.h" 64 #include "ortools/sat/cp_model_loader.h" 65 #include "ortools/sat/cp_model_postsolve.h" 66 #include "ortools/sat/cp_model_presolve.h" 67 #include "ortools/sat/cp_model_search.h" 68 #include "ortools/sat/cp_model_symmetries.h" 69 #include "ortools/sat/cp_model_utils.h" 70 #include "ortools/sat/cuts.h" 71 #include "ortools/sat/drat_checker.h" 72 #include "ortools/sat/drat_proof_handler.h" 73 #include "ortools/sat/feasibility_pump.h" 74 #include "ortools/sat/integer.h" 75 #include "ortools/sat/integer_expr.h" 76 #include "ortools/sat/integer_search.h" 77 #include "ortools/sat/intervals.h" 78 #include "ortools/sat/lb_tree_search.h" 79 #include "ortools/sat/linear_programming_constraint.h" 80 #include "ortools/sat/linear_relaxation.h" 81 #include "ortools/sat/lp_utils.h" 82 #include "ortools/sat/optimization.h" 83 #include "ortools/sat/parameters_validation.h" 84 #include "ortools/sat/precedences.h" 85 #include "ortools/sat/probing.h" 86 #include "ortools/sat/rins.h" 87 #include "ortools/sat/sat_base.h" 88 #include "ortools/sat/sat_inprocessing.h" 89 #include "ortools/sat/sat_parameters.pb.h" 90 #include "ortools/sat/sat_solver.h" 91 #include "ortools/sat/simplification.h" 92 #include "ortools/sat/subsolver.h" 93 #include "ortools/sat/synchronization.h" 94 #include "ortools/util/logging.h" 95 #include "ortools/util/sorted_interval_list.h" 96 #include "ortools/util/time_limit.h" 97 98 #if defined(_MSC_VER) 99 ABSL_FLAG(std::string, cp_model_dump_prefix, ".\\", 100 "Prefix filename for all dumped files"); 101 #else 102 ABSL_FLAG(std::string, cp_model_dump_prefix, "/tmp/", 103 "Prefix filename for all dumped files"); 104 #endif 105 ABSL_FLAG(bool, cp_model_dump_models, false, 106 "DEBUG ONLY. When set to true, SolveCpModel() will dump its model " 107 "protos (original model, presolved model, mapping model) in text " 108 "format to 'FLAGS_cp_model_dump_prefix'{model|presolved_model|" 109 "mapping_model}.pb.txt."); 110 111 ABSL_FLAG(bool, cp_model_dump_lns, false, 112 "DEBUG ONLY. When set to true, solve will dump all " 113 "lns models proto in text format to " 114 "'FLAGS_cp_model_dump_prefix'lns_xxx.pb.txt."); 115 116 ABSL_FLAG( 117 bool, cp_model_dump_problematic_lns, false, 118 "DEBUG ONLY. Similar to --cp_model_dump_lns, but only dump fragment for " 119 "which we got an issue while validating the postsolved solution. This " 120 "allows to debug presolve issues without dumping all the models."); 121 122 ABSL_FLAG(bool, cp_model_dump_response, false, 123 "DEBUG ONLY. If true, the final response of each solve will be " 124 "dumped to 'FLAGS_cp_model_dump_prefix'response.pb.txt"); 125 126 ABSL_FLAG(std::string, cp_model_params, "", 127 "This is interpreted as a text SatParameters proto. The " 128 "specified fields will override the normal ones for all solves."); 129 130 ABSL_FLAG(std::string, drat_output, "", 131 "If non-empty, a proof in DRAT format will be written to this file. " 132 "This will only be used for pure-SAT problems."); 133 134 ABSL_FLAG(bool, drat_check, false, 135 "If true, a proof in DRAT format will be stored in memory and " 136 "checked if the problem is UNSAT. This will only be used for " 137 "pure-SAT problems."); 138 139 ABSL_FLAG(double, max_drat_time_in_seconds, 140 std::numeric_limits<double>::infinity(), 141 "Maximum time in seconds to check the DRAT proof. This will only " 142 "be used is the drat_check flag is enabled."); 143 144 ABSL_FLAG(bool, cp_model_check_intermediate_solutions, false, 145 "When true, all intermediate solutions found by the solver will be " 146 "checked. This can be expensive, therefore it is off by default."); 147 148 ABSL_FLAG(std::string, contention_profile, "", 149 "If non-empty, dump a contention pprof proto to the specified " 150 "destination at the end of the solve."); 151 152 namespace operations_research { 153 namespace sat { 154 155 namespace { 156 157 // Makes the string fit in one line by cutting it in the middle if necessary. 158 std::string Summarize(const std::string& input) { 159 if (input.size() < 105) return input; 160 const int half = 50; 161 return absl::StrCat(input.substr(0, half), " ... ", 162 input.substr(input.size() - half, half)); 163 } 164 165 } // namespace. 166 167 // ============================================================================= 168 // Public API. 169 // ============================================================================= 170 171 std::string CpModelStats(const CpModelProto& model_proto) { 172 std::map<std::string, int> num_constraints_by_name; 173 std::map<std::string, int> num_reif_constraints_by_name; 174 std::map<std::string, int> name_to_num_literals; 175 std::map<std::string, int> name_to_num_terms; 176 177 int no_overlap_2d_num_rectangles = 0; 178 int no_overlap_2d_num_optional_rectangles = 0; 179 int no_overlap_2d_num_linear_areas = 0; 180 int no_overlap_2d_num_quadratic_areas = 0; 181 182 int cumulative_num_intervals = 0; 183 int cumulative_num_optional_intervals = 0; 184 int cumulative_num_variable_sizes = 0; 185 int cumulative_num_variable_demands = 0; 186 187 int no_overlap_num_intervals = 0; 188 int no_overlap_num_optional_intervals = 0; 189 int no_overlap_num_variable_sizes = 0; 190 191 for (const ConstraintProto& ct : model_proto.constraints()) { 192 std::string name = ConstraintCaseName(ct.constraint_case()); 193 194 // We split the linear constraints into 3 buckets has it gives more insight 195 // on the type of problem we are facing. 196 if (ct.constraint_case() == ConstraintProto::ConstraintCase::kLinear) { 197 if (ct.linear().vars_size() == 1) name += "1"; 198 if (ct.linear().vars_size() == 2) name += "2"; 199 if (ct.linear().vars_size() == 3) name += "3"; 200 if (ct.linear().vars_size() > 3) name += "N"; 201 } 202 203 num_constraints_by_name[name]++; 204 if (!ct.enforcement_literal().empty()) { 205 num_reif_constraints_by_name[name]++; 206 } 207 208 auto variable_is_fixed = [&model_proto](int ref) { 209 const IntegerVariableProto& proto = 210 model_proto.variables(PositiveRef(ref)); 211 return proto.domain_size() == 2 && proto.domain(0) == proto.domain(1); 212 }; 213 214 auto expression_is_fixed = 215 [&variable_is_fixed](const LinearExpressionProto& expr) { 216 for (const int ref : expr.vars()) { 217 if (!variable_is_fixed(ref)) { 218 return false; 219 } 220 } 221 return true; 222 }; 223 224 auto interval_has_fixed_size = [&model_proto, &expression_is_fixed](int c) { 225 return expression_is_fixed(model_proto.constraints(c).interval().size()); 226 }; 227 228 auto constraint_is_optional = [&model_proto](int i) { 229 return !model_proto.constraints(i).enforcement_literal().empty(); 230 }; 231 232 // For pure Boolean constraints, we also display the total number of literal 233 // involved as this gives a good idea of the problem size. 234 if (ct.constraint_case() == ConstraintProto::ConstraintCase::kBoolOr) { 235 name_to_num_literals[name] += ct.bool_or().literals().size(); 236 } else if (ct.constraint_case() == 237 ConstraintProto::ConstraintCase::kBoolAnd) { 238 name_to_num_literals[name] += ct.bool_and().literals().size(); 239 } else if (ct.constraint_case() == 240 ConstraintProto::ConstraintCase::kAtMostOne) { 241 name_to_num_literals[name] += ct.at_most_one().literals().size(); 242 } else if (ct.constraint_case() == 243 ConstraintProto::ConstraintCase::kExactlyOne) { 244 name_to_num_literals[name] += ct.exactly_one().literals().size(); 245 } else if (ct.constraint_case() == 246 ConstraintProto::ConstraintCase::kNoOverlap2D) { 247 const int num_boxes = ct.no_overlap_2d().x_intervals_size(); 248 no_overlap_2d_num_rectangles += num_boxes; 249 for (int i = 0; i < num_boxes; ++i) { 250 const int x_interval = ct.no_overlap_2d().x_intervals(i); 251 const int y_interval = ct.no_overlap_2d().y_intervals(i); 252 if (constraint_is_optional(x_interval) || 253 constraint_is_optional(y_interval)) { 254 no_overlap_2d_num_optional_rectangles++; 255 } 256 const int num_fixed = interval_has_fixed_size(x_interval) + 257 interval_has_fixed_size(y_interval); 258 if (num_fixed == 0) { 259 no_overlap_2d_num_quadratic_areas++; 260 } else if (num_fixed == 1) { 261 no_overlap_2d_num_linear_areas++; 262 } 263 } 264 } else if (ct.constraint_case() == 265 ConstraintProto::ConstraintCase::kNoOverlap) { 266 const int num_intervals = ct.no_overlap().intervals_size(); 267 no_overlap_num_intervals += num_intervals; 268 for (int i = 0; i < num_intervals; ++i) { 269 const int interval = ct.no_overlap().intervals(i); 270 if (constraint_is_optional(interval)) { 271 no_overlap_num_optional_intervals++; 272 } 273 if (!interval_has_fixed_size(interval)) { 274 no_overlap_num_variable_sizes++; 275 } 276 } 277 } else if (ct.constraint_case() == 278 ConstraintProto::ConstraintCase::kCumulative) { 279 const int num_intervals = ct.cumulative().intervals_size(); 280 cumulative_num_intervals += num_intervals; 281 for (int i = 0; i < num_intervals; ++i) { 282 const int interval = ct.cumulative().intervals(i); 283 if (constraint_is_optional(interval)) { 284 cumulative_num_optional_intervals++; 285 } 286 if (!interval_has_fixed_size(interval)) { 287 cumulative_num_variable_sizes++; 288 } 289 if (!expression_is_fixed(ct.cumulative().demands(i))) { 290 cumulative_num_variable_demands++; 291 } 292 } 293 } 294 295 if (ct.constraint_case() == ConstraintProto::ConstraintCase::kLinear && 296 ct.linear().vars_size() > 3) { 297 name_to_num_terms[name] += ct.linear().vars_size(); 298 } 299 } 300 301 int num_constants = 0; 302 std::set<int64_t> constant_values; 303 std::map<Domain, int> num_vars_per_domains; 304 for (const IntegerVariableProto& var : model_proto.variables()) { 305 if (var.domain_size() == 2 && var.domain(0) == var.domain(1)) { 306 ++num_constants; 307 constant_values.insert(var.domain(0)); 308 } else { 309 num_vars_per_domains[ReadDomainFromProto(var)]++; 310 } 311 } 312 313 std::string result; 314 if (model_proto.has_objective() || 315 model_proto.has_floating_point_objective()) { 316 absl::StrAppend(&result, "optimization model '", model_proto.name(), 317 "':\n"); 318 } else { 319 absl::StrAppend(&result, "satisfaction model '", model_proto.name(), 320 "':\n"); 321 } 322 323 for (const DecisionStrategyProto& strategy : model_proto.search_strategy()) { 324 absl::StrAppend( 325 &result, "Search strategy: on ", strategy.variables_size(), 326 " variables, ", 327 ProtoEnumToString<DecisionStrategyProto::VariableSelectionStrategy>( 328 strategy.variable_selection_strategy()), 329 ", ", 330 ProtoEnumToString<DecisionStrategyProto::DomainReductionStrategy>( 331 strategy.domain_reduction_strategy()), 332 "\n"); 333 } 334 335 const std::string objective_string = 336 model_proto.has_objective() 337 ? absl::StrCat(" (", model_proto.objective().vars_size(), 338 " in objective)") 339 : (model_proto.has_floating_point_objective() 340 ? absl::StrCat( 341 " (", model_proto.floating_point_objective().vars_size(), 342 " in floating point objective)") 343 : ""); 344 absl::StrAppend(&result, "#Variables: ", model_proto.variables_size(), 345 objective_string, "\n"); 346 if (num_vars_per_domains.size() < 100) { 347 for (const auto& entry : num_vars_per_domains) { 348 const std::string temp = absl::StrCat(" - ", entry.second, " in ", 349 entry.first.ToString(), "\n"); 350 absl::StrAppend(&result, Summarize(temp)); 351 } 352 } else { 353 int64_t max_complexity = 0; 354 int64_t min = std::numeric_limits<int64_t>::max(); 355 int64_t max = std::numeric_limits<int64_t>::min(); 356 for (const auto& entry : num_vars_per_domains) { 357 min = std::min(min, entry.first.Min()); 358 max = std::max(max, entry.first.Max()); 359 max_complexity = std::max( 360 max_complexity, static_cast<int64_t>(entry.first.NumIntervals())); 361 } 362 absl::StrAppend(&result, " - ", num_vars_per_domains.size(), 363 " different domains in [", min, ",", max, 364 "] with a largest complexity of ", max_complexity, ".\n"); 365 } 366 367 if (num_constants > 0) { 368 const std::string temp = 369 absl::StrCat(" - ", num_constants, " constants in {", 370 absl::StrJoin(constant_values, ","), "} \n"); 371 absl::StrAppend(&result, Summarize(temp)); 372 } 373 374 std::vector<std::string> constraints; 375 constraints.reserve(num_constraints_by_name.size()); 376 for (const auto& entry : num_constraints_by_name) { 377 const std::string& name = entry.first; 378 constraints.push_back(absl::StrCat("#", name, ": ", entry.second)); 379 if (gtl::ContainsKey(num_reif_constraints_by_name, name)) { 380 absl::StrAppend(&constraints.back(), 381 " (#enforced: ", num_reif_constraints_by_name[name], ")"); 382 } 383 if (gtl::ContainsKey(name_to_num_literals, name)) { 384 absl::StrAppend(&constraints.back(), 385 " (#literals: ", name_to_num_literals[name], ")"); 386 } 387 if (gtl::ContainsKey(name_to_num_terms, name)) { 388 absl::StrAppend(&constraints.back(), 389 " (#terms: ", name_to_num_terms[name], ")"); 390 } 391 if (name == "kNoOverlap2D") { 392 absl::StrAppend(&constraints.back(), 393 " (#rectangles: ", no_overlap_2d_num_rectangles); 394 if (no_overlap_2d_num_optional_rectangles > 0) { 395 absl::StrAppend(&constraints.back(), 396 ", #optional: ", no_overlap_2d_num_optional_rectangles); 397 } 398 if (no_overlap_2d_num_linear_areas > 0) { 399 absl::StrAppend(&constraints.back(), 400 ", #linear_areas: ", no_overlap_2d_num_linear_areas); 401 } 402 if (no_overlap_2d_num_quadratic_areas > 0) { 403 absl::StrAppend(&constraints.back(), ", #quadratic_areas: ", 404 no_overlap_2d_num_quadratic_areas); 405 } 406 absl::StrAppend(&constraints.back(), ")"); 407 } else if (name == "kCumulative") { 408 absl::StrAppend(&constraints.back(), 409 " (#intervals: ", cumulative_num_intervals); 410 if (cumulative_num_optional_intervals > 0) { 411 absl::StrAppend(&constraints.back(), 412 ", #optional: ", cumulative_num_optional_intervals); 413 } 414 if (cumulative_num_variable_sizes > 0) { 415 absl::StrAppend(&constraints.back(), 416 ", #variable_sizes: ", cumulative_num_variable_sizes); 417 } 418 if (cumulative_num_variable_demands > 0) { 419 absl::StrAppend(&constraints.back(), ", #variable_demands: ", 420 cumulative_num_variable_demands); 421 } 422 absl::StrAppend(&constraints.back(), ")"); 423 } else if (name == "kNoOverlap") { 424 absl::StrAppend(&constraints.back(), 425 " (#intervals: ", no_overlap_num_intervals); 426 if (no_overlap_num_optional_intervals > 0) { 427 absl::StrAppend(&constraints.back(), 428 ", #optional: ", no_overlap_num_optional_intervals); 429 } 430 if (no_overlap_num_variable_sizes > 0) { 431 absl::StrAppend(&constraints.back(), 432 ", #variable_sizes: ", no_overlap_num_variable_sizes); 433 } 434 absl::StrAppend(&constraints.back(), ")"); 435 } 436 } 437 std::sort(constraints.begin(), constraints.end()); 438 absl::StrAppend(&result, absl::StrJoin(constraints, "\n")); 439 440 return result; 441 } 442 443 std::string CpSolverResponseStats(const CpSolverResponse& response, 444 bool has_objective) { 445 std::string result; 446 absl::StrAppend(&result, "CpSolverResponse summary:"); 447 absl::StrAppend(&result, "\nstatus: ", 448 ProtoEnumToString<CpSolverStatus>(response.status())); 449 450 if (has_objective && response.status() != CpSolverStatus::INFEASIBLE) { 451 absl::StrAppendFormat(&result, "\nobjective: %.16g", 452 response.objective_value()); 453 absl::StrAppendFormat(&result, "\nbest_bound: %.16g", 454 response.best_objective_bound()); 455 } else { 456 absl::StrAppend(&result, "\nobjective: NA"); 457 absl::StrAppend(&result, "\nbest_bound: NA"); 458 } 459 460 absl::StrAppend(&result, "\nbooleans: ", response.num_booleans()); 461 absl::StrAppend(&result, "\nconflicts: ", response.num_conflicts()); 462 absl::StrAppend(&result, "\nbranches: ", response.num_branches()); 463 464 // TODO(user): This is probably better named "binary_propagation", but we just 465 // output "propagations" to be consistent with sat/analyze.sh. 466 absl::StrAppend(&result, 467 "\npropagations: ", response.num_binary_propagations()); 468 absl::StrAppend( 469 &result, "\ninteger_propagations: ", response.num_integer_propagations()); 470 471 absl::StrAppend(&result, "\nrestarts: ", response.num_restarts()); 472 absl::StrAppend(&result, "\nlp_iterations: ", response.num_lp_iterations()); 473 absl::StrAppend(&result, "\nwalltime: ", response.wall_time()); 474 absl::StrAppend(&result, "\nusertime: ", response.user_time()); 475 absl::StrAppend(&result, 476 "\ndeterministic_time: ", response.deterministic_time()); 477 absl::StrAppend(&result, "\ngap_integral: ", response.gap_integral()); 478 absl::StrAppend(&result, "\n"); 479 return result; 480 } 481 482 namespace { 483 484 #if !defined(__PORTABLE_PLATFORM__) 485 #endif // __PORTABLE_PLATFORM__ 486 487 void FillSolutionInResponse(const CpModelProto& model_proto, const Model& model, 488 CpSolverResponse* response) { 489 response->clear_solution(); 490 491 auto* mapping = model.Get<CpModelMapping>(); 492 auto* trail = model.Get<Trail>(); 493 494 std::vector<int64_t> solution; 495 for (int i = 0; i < model_proto.variables_size(); ++i) { 496 if (mapping->IsInteger(i)) { 497 const IntegerVariable var = mapping->Integer(i); 498 499 // For ignored or not fully instanciated variable, we just use the 500 // lower bound. 501 solution.push_back(model.Get(LowerBound(var))); 502 } else { 503 DCHECK(mapping->IsBoolean(i)); 504 const Literal literal = mapping->Literal(i); 505 if (trail->Assignment().LiteralIsAssigned(literal)) { 506 solution.push_back(model.Get(Value(literal))); 507 } else { 508 // Just use the lower bound if the variable is not fully instantiated. 509 solution.push_back(0); 510 } 511 } 512 } 513 514 if (DEBUG_MODE || 515 absl::GetFlag(FLAGS_cp_model_check_intermediate_solutions)) { 516 // TODO(user): Checks against initial model. 517 CHECK(SolutionIsFeasible(model_proto, solution)); 518 } 519 for (const int64_t value : solution) { 520 response->add_solution(value); 521 } 522 } 523 524 namespace { 525 526 IntegerVariable GetOrCreateVariableWithTightBound( 527 const std::vector<std::pair<IntegerVariable, int64_t>>& terms, 528 Model* model) { 529 if (terms.empty()) return model->Add(ConstantIntegerVariable(0)); 530 if (terms.size() == 1 && terms.front().second == 1) { 531 return terms.front().first; 532 } 533 if (terms.size() == 1 && terms.front().second == -1) { 534 return NegationOf(terms.front().first); 535 } 536 537 int64_t sum_min = 0; 538 int64_t sum_max = 0; 539 for (const std::pair<IntegerVariable, int64_t> var_coeff : terms) { 540 const int64_t min_domain = model->Get(LowerBound(var_coeff.first)); 541 const int64_t max_domain = model->Get(UpperBound(var_coeff.first)); 542 const int64_t coeff = var_coeff.second; 543 const int64_t prod1 = min_domain * coeff; 544 const int64_t prod2 = max_domain * coeff; 545 sum_min += std::min(prod1, prod2); 546 sum_max += std::max(prod1, prod2); 547 } 548 return model->Add(NewIntegerVariable(sum_min, sum_max)); 549 } 550 551 IntegerVariable GetOrCreateVariableGreaterOrEqualToSumOf( 552 const std::vector<std::pair<IntegerVariable, int64_t>>& terms, 553 Model* model) { 554 if (terms.empty()) return model->Add(ConstantIntegerVariable(0)); 555 if (terms.size() == 1 && terms.front().second == 1) { 556 return terms.front().first; 557 } 558 if (terms.size() == 1 && terms.front().second == -1) { 559 return NegationOf(terms.front().first); 560 } 561 562 // Create a new variable and link it with the linear terms. 563 const IntegerVariable new_var = 564 GetOrCreateVariableWithTightBound(terms, model); 565 std::vector<IntegerVariable> vars; 566 std::vector<int64_t> coeffs; 567 for (const auto& term : terms) { 568 vars.push_back(term.first); 569 coeffs.push_back(term.second); 570 } 571 vars.push_back(new_var); 572 coeffs.push_back(-1); 573 model->Add(WeightedSumLowerOrEqual(vars, coeffs, 0)); 574 return new_var; 575 } 576 577 } // namespace 578 579 // Adds one LinearProgrammingConstraint per connected component of the model. 580 IntegerVariable AddLPConstraints(const CpModelProto& model_proto, Model* m) { 581 const LinearRelaxation relaxation = ComputeLinearRelaxation(model_proto, m); 582 583 // The bipartite graph of LP constraints might be disconnected: 584 // make a partition of the variables into connected components. 585 // Constraint nodes are indexed by [0..num_lp_constraints), 586 // variable nodes by [num_lp_constraints..num_lp_constraints+num_variables). 587 // 588 // TODO(user): look into biconnected components. 589 const int num_lp_constraints = relaxation.linear_constraints.size(); 590 const int num_lp_cut_generators = relaxation.cut_generators.size(); 591 const int num_integer_variables = 592 m->GetOrCreate<IntegerTrail>()->NumIntegerVariables().value(); 593 DenseConnectedComponentsFinder components; 594 components.SetNumberOfNodes(num_lp_constraints + num_lp_cut_generators + 595 num_integer_variables); 596 auto get_constraint_index = [](int ct_index) { return ct_index; }; 597 auto get_cut_generator_index = [num_lp_constraints](int cut_index) { 598 return num_lp_constraints + cut_index; 599 }; 600 auto get_var_index = [num_lp_constraints, 601 num_lp_cut_generators](IntegerVariable var) { 602 return num_lp_constraints + num_lp_cut_generators + var.value(); 603 }; 604 for (int i = 0; i < num_lp_constraints; i++) { 605 for (const IntegerVariable var : relaxation.linear_constraints[i].vars) { 606 components.AddEdge(get_constraint_index(i), get_var_index(var)); 607 } 608 } 609 for (int i = 0; i < num_lp_cut_generators; ++i) { 610 for (const IntegerVariable var : relaxation.cut_generators[i].vars) { 611 components.AddEdge(get_cut_generator_index(i), get_var_index(var)); 612 } 613 } 614 615 const int num_components = components.GetNumberOfComponents(); 616 std::vector<int> component_sizes(num_components, 0); 617 const std::vector<int> index_to_component = components.GetComponentIds(); 618 for (int i = 0; i < num_lp_constraints; i++) { 619 ++component_sizes[index_to_component[get_constraint_index(i)]]; 620 } 621 for (int i = 0; i < num_lp_cut_generators; i++) { 622 ++component_sizes[index_to_component[get_cut_generator_index(i)]]; 623 } 624 625 // Make sure any constraint that touch the objective is not discarded even 626 // if it is the only one in its component. This is important to propagate 627 // as much as possible the objective bound by using any bounds the LP give 628 // us on one of its components. This is critical on the zephyrus problems for 629 // instance. 630 auto* mapping = m->GetOrCreate<CpModelMapping>(); 631 for (int i = 0; i < model_proto.objective().coeffs_size(); ++i) { 632 const IntegerVariable var = 633 mapping->Integer(model_proto.objective().vars(i)); 634 ++component_sizes[index_to_component[get_var_index(var)]]; 635 } 636 637 // Dispatch every constraint to its LinearProgrammingConstraint. 638 std::vector<LinearProgrammingConstraint*> lp_constraints(num_components, 639 nullptr); 640 std::vector<std::vector<LinearConstraint>> component_to_constraints( 641 num_components); 642 for (int i = 0; i < num_lp_constraints; i++) { 643 const int c = index_to_component[get_constraint_index(i)]; 644 if (component_sizes[c] <= 1) continue; 645 component_to_constraints[c].push_back(relaxation.linear_constraints[i]); 646 if (lp_constraints[c] == nullptr) { 647 lp_constraints[c] = m->Create<LinearProgrammingConstraint>(); 648 } 649 // Load the constraint. 650 lp_constraints[c]->AddLinearConstraint(relaxation.linear_constraints[i]); 651 } 652 653 // Dispatch every cut generator to its LinearProgrammingConstraint. 654 for (int i = 0; i < num_lp_cut_generators; i++) { 655 const int c = index_to_component[get_cut_generator_index(i)]; 656 if (lp_constraints[c] == nullptr) { 657 lp_constraints[c] = m->Create<LinearProgrammingConstraint>(); 658 } 659 lp_constraints[c]->AddCutGenerator(std::move(relaxation.cut_generators[i])); 660 } 661 662 // Add the objective. 663 std::vector<std::vector<std::pair<IntegerVariable, int64_t>>> 664 component_to_cp_terms(num_components); 665 std::vector<std::pair<IntegerVariable, int64_t>> top_level_cp_terms; 666 int num_components_containing_objective = 0; 667 if (model_proto.has_objective()) { 668 // First pass: set objective coefficients on the lp constraints, and store 669 // the cp terms in one vector per component. 670 for (int i = 0; i < model_proto.objective().coeffs_size(); ++i) { 671 const IntegerVariable var = 672 mapping->Integer(model_proto.objective().vars(i)); 673 const int64_t coeff = model_proto.objective().coeffs(i); 674 const int c = index_to_component[get_var_index(var)]; 675 if (lp_constraints[c] != nullptr) { 676 lp_constraints[c]->SetObjectiveCoefficient(var, IntegerValue(coeff)); 677 component_to_cp_terms[c].push_back(std::make_pair(var, coeff)); 678 } else { 679 // Component is too small. We still need to store the objective term. 680 top_level_cp_terms.push_back(std::make_pair(var, coeff)); 681 } 682 } 683 // Second pass: Build the cp sub-objectives per component. 684 for (int c = 0; c < num_components; ++c) { 685 if (component_to_cp_terms[c].empty()) continue; 686 const IntegerVariable sub_obj_var = 687 GetOrCreateVariableGreaterOrEqualToSumOf(component_to_cp_terms[c], m); 688 top_level_cp_terms.push_back(std::make_pair(sub_obj_var, 1)); 689 lp_constraints[c]->SetMainObjectiveVariable(sub_obj_var); 690 num_components_containing_objective++; 691 } 692 } 693 694 const IntegerVariable main_objective_var = 695 model_proto.has_objective() 696 ? GetOrCreateVariableGreaterOrEqualToSumOf(top_level_cp_terms, m) 697 : kNoIntegerVariable; 698 699 // Register LP constraints. Note that this needs to be done after all the 700 // constraints have been added. 701 for (LinearProgrammingConstraint* lp_constraint : lp_constraints) { 702 if (lp_constraint == nullptr) continue; 703 lp_constraint->RegisterWith(m); 704 VLOG(3) << "LP constraint: " << lp_constraint->DimensionString() << "."; 705 } 706 707 VLOG(3) << top_level_cp_terms.size() 708 << " terms in the main objective linear equation (" 709 << num_components_containing_objective << " from LP constraints)."; 710 return main_objective_var; 711 } 712 713 } // namespace 714 715 // Used by NewFeasibleSolutionObserver to register observers. 716 struct SolutionObservers { 717 explicit SolutionObservers(Model* model) {} 718 std::vector<std::function<void(const CpSolverResponse& response)>> observers; 719 }; 720 721 std::function<void(Model*)> NewFeasibleSolutionObserver( 722 const std::function<void(const CpSolverResponse& response)>& observer) { 723 return [=](Model* model) { 724 model->GetOrCreate<SolutionObservers>()->observers.push_back(observer); 725 }; 726 } 727 728 #if !defined(__PORTABLE_PLATFORM__) 729 // TODO(user): Support it on android. 730 std::function<SatParameters(Model*)> NewSatParameters( 731 const std::string& params) { 732 sat::SatParameters parameters; 733 if (!params.empty()) { 734 CHECK(google::protobuf::TextFormat::ParseFromString(params, ¶meters)) 735 << params; 736 } 737 return NewSatParameters(parameters); 738 } 739 #endif // __PORTABLE_PLATFORM__ 740 741 std::function<SatParameters(Model*)> NewSatParameters( 742 const sat::SatParameters& parameters) { 743 return [=](Model* model) { 744 // Tricky: It is important to initialize the model parameters before any 745 // of the solver object are created, so that by default they use the given 746 // parameters. 747 // 748 // TODO(user): A notable exception to this is the TimeLimit which is 749 // currently not initializing itself from the SatParameters in the model. It 750 // will also starts counting from the time of its creation. It will be good 751 // to find a solution that is less error prone. 752 *model->GetOrCreate<SatParameters>() = parameters; 753 return parameters; 754 }; 755 } 756 757 namespace { 758 759 // Registers a callback that will export variables bounds fixed at level 0 of 760 // the search. This should not be registered to a LNS search. 761 void RegisterVariableBoundsLevelZeroExport( 762 const CpModelProto& model_proto, SharedBoundsManager* shared_bounds_manager, 763 Model* model) { 764 CHECK(shared_bounds_manager != nullptr); 765 int saved_trail_index = 0; 766 auto broadcast_level_zero_bounds = 767 [&model_proto, saved_trail_index, model, shared_bounds_manager]( 768 const std::vector<IntegerVariable>& modified_vars) mutable { 769 CpModelMapping* const mapping = model->GetOrCreate<CpModelMapping>(); 770 771 std::vector<int> model_variables; 772 std::vector<int64_t> new_lower_bounds; 773 std::vector<int64_t> new_upper_bounds; 774 absl::flat_hash_set<int> visited_variables; 775 776 // Inspect the modified IntegerVariables. 777 auto* integer_trail = model->Get<IntegerTrail>(); 778 for (const IntegerVariable& var : modified_vars) { 779 const IntegerVariable positive_var = PositiveVariable(var); 780 const int model_var = 781 mapping->GetProtoVariableFromIntegerVariable(positive_var); 782 if (model_var == -1 || visited_variables.contains(model_var)) { 783 // TODO(user): I don't think we should see the same model_var twice 784 // here so maybe we don't need the visited_variables.contains() 785 // part. 786 continue; 787 } 788 789 visited_variables.insert(model_var); 790 const int64_t new_lb = 791 integer_trail->LevelZeroLowerBound(positive_var).value(); 792 const int64_t new_ub = 793 integer_trail->LevelZeroUpperBound(positive_var).value(); 794 // TODO(user): We could imagine an API based on atomic<int64_t> 795 // that could preemptively check if this new bounds are improving. 796 model_variables.push_back(model_var); 797 new_lower_bounds.push_back(new_lb); 798 new_upper_bounds.push_back(new_ub); 799 } 800 801 // Inspect the newly modified Booleans. 802 auto* trail = model->Get<Trail>(); 803 for (; saved_trail_index < trail->Index(); ++saved_trail_index) { 804 const Literal fixed_literal = (*trail)[saved_trail_index]; 805 const int model_var = mapping->GetProtoVariableFromBooleanVariable( 806 fixed_literal.Variable()); 807 if (model_var == -1 || visited_variables.contains(model_var)) { 808 // If the variable is already visited, it should mean that this 809 // Boolean also has an IntegerVariable view, and we should already 810 // have set its bound correctly. 811 continue; 812 } 813 814 visited_variables.insert(model_var); 815 model_variables.push_back(model_var); 816 if (fixed_literal.IsPositive()) { 817 new_lower_bounds.push_back(1); 818 new_upper_bounds.push_back(1); 819 } else { 820 new_lower_bounds.push_back(0); 821 new_upper_bounds.push_back(0); 822 } 823 } 824 825 if (!model_variables.empty()) { 826 shared_bounds_manager->ReportPotentialNewBounds( 827 model_proto, model->Name(), model_variables, new_lower_bounds, 828 new_upper_bounds); 829 830 // If we are not in interleave_search we synchronize right away. 831 if (!model->Get<SatParameters>()->interleave_search()) { 832 shared_bounds_manager->Synchronize(); 833 } 834 } 835 }; 836 837 // The callback will just be called on NEWLY modified var. So initially, 838 // we do want to read all variables. 839 // 840 // TODO(user): Find a better way? It seems nicer to register this before 841 // any variable is modified. But then we don't want to call it each time 842 // we reach level zero during probing. It should be better to only call 843 // it when a new variable has been fixed. 844 const IntegerVariable num_vars = 845 model->GetOrCreate<IntegerTrail>()->NumIntegerVariables(); 846 std::vector<IntegerVariable> all_variables; 847 all_variables.reserve(num_vars.value()); 848 for (IntegerVariable var(0); var < num_vars; ++var) { 849 all_variables.push_back(var); 850 } 851 broadcast_level_zero_bounds(all_variables); 852 853 model->GetOrCreate<GenericLiteralWatcher>() 854 ->RegisterLevelZeroModifiedVariablesCallback(broadcast_level_zero_bounds); 855 } 856 857 // Registers a callback to import new variables bounds stored in the 858 // shared_bounds_manager. These bounds are imported at level 0 of the search 859 // in the linear scan minimize function. 860 void RegisterVariableBoundsLevelZeroImport( 861 const CpModelProto& model_proto, SharedBoundsManager* shared_bounds_manager, 862 Model* model) { 863 CHECK(shared_bounds_manager != nullptr); 864 auto* integer_trail = model->GetOrCreate<IntegerTrail>(); 865 CpModelMapping* const mapping = model->GetOrCreate<CpModelMapping>(); 866 const int id = shared_bounds_manager->RegisterNewId(); 867 868 const auto& import_level_zero_bounds = [&model_proto, shared_bounds_manager, 869 model, integer_trail, id, mapping]() { 870 std::vector<int> model_variables; 871 std::vector<int64_t> new_lower_bounds; 872 std::vector<int64_t> new_upper_bounds; 873 shared_bounds_manager->GetChangedBounds( 874 id, &model_variables, &new_lower_bounds, &new_upper_bounds); 875 bool new_bounds_have_been_imported = false; 876 for (int i = 0; i < model_variables.size(); ++i) { 877 const int model_var = model_variables[i]; 878 // This can happen if a boolean variables is forced to have an 879 // integer view in one thread, and not in another thread. 880 if (!mapping->IsInteger(model_var)) continue; 881 const IntegerVariable var = mapping->Integer(model_var); 882 const IntegerValue new_lb(new_lower_bounds[i]); 883 const IntegerValue new_ub(new_upper_bounds[i]); 884 const IntegerValue old_lb = integer_trail->LowerBound(var); 885 const IntegerValue old_ub = integer_trail->UpperBound(var); 886 const bool changed_lb = new_lb > old_lb; 887 const bool changed_ub = new_ub < old_ub; 888 if (!changed_lb && !changed_ub) continue; 889 890 new_bounds_have_been_imported = true; 891 if (VLOG_IS_ON(3)) { 892 const IntegerVariableProto& var_proto = 893 model_proto.variables(model_var); 894 const std::string& var_name = 895 var_proto.name().empty() 896 ? absl::StrCat("anonymous_var(", model_var, ")") 897 : var_proto.name(); 898 LOG(INFO) << " '" << model->Name() << "' imports new bounds for " 899 << var_name << ": from [" << old_lb << ", " << old_ub 900 << "] to [" << new_lb << ", " << new_ub << "]"; 901 } 902 903 if (changed_lb && 904 !integer_trail->Enqueue(IntegerLiteral::GreaterOrEqual(var, new_lb), 905 {}, {})) { 906 return false; 907 } 908 if (changed_ub && 909 !integer_trail->Enqueue(IntegerLiteral::LowerOrEqual(var, new_ub), {}, 910 {})) { 911 return false; 912 } 913 } 914 if (new_bounds_have_been_imported && 915 !model->GetOrCreate<SatSolver>()->FinishPropagation()) { 916 return false; 917 } 918 return true; 919 }; 920 model->GetOrCreate<LevelZeroCallbackHelper>()->callbacks.push_back( 921 import_level_zero_bounds); 922 } 923 924 // Registers a callback that will report improving objective best bound. 925 // It will be called each time new objective bound are propagated at level zero. 926 void RegisterObjectiveBestBoundExport( 927 IntegerVariable objective_var, 928 SharedResponseManager* shared_response_manager, Model* model) { 929 auto* integer_trail = model->Get<IntegerTrail>(); 930 const auto broadcast_objective_lower_bound = 931 [objective_var, integer_trail, shared_response_manager, 932 model](const std::vector<IntegerVariable>& unused) { 933 shared_response_manager->UpdateInnerObjectiveBounds( 934 model->Name(), integer_trail->LevelZeroLowerBound(objective_var), 935 integer_trail->LevelZeroUpperBound(objective_var)); 936 // If we are not in interleave_search we synchronize right away. 937 if (!model->Get<SatParameters>()->interleave_search()) { 938 shared_response_manager->Synchronize(); 939 } 940 }; 941 model->GetOrCreate<GenericLiteralWatcher>() 942 ->RegisterLevelZeroModifiedVariablesCallback( 943 broadcast_objective_lower_bound); 944 } 945 946 // Registers a callback to import new objective bounds. It will be called each 947 // time the search main loop is back to level zero. Note that it the presence of 948 // assumptions, this will not happen until the set of assumptions is changed. 949 void RegisterObjectiveBoundsImport( 950 SharedResponseManager* shared_response_manager, Model* model) { 951 auto* solver = model->GetOrCreate<SatSolver>(); 952 auto* integer_trail = model->GetOrCreate<IntegerTrail>(); 953 auto* objective = model->GetOrCreate<ObjectiveDefinition>(); 954 const std::string name = model->Name(); 955 const auto import_objective_bounds = [name, solver, integer_trail, objective, 956 shared_response_manager]() { 957 if (solver->AssumptionLevel() != 0) return true; 958 bool propagate = false; 959 960 const IntegerValue external_lb = 961 shared_response_manager->SynchronizedInnerObjectiveLowerBound(); 962 const IntegerValue current_lb = 963 integer_trail->LowerBound(objective->objective_var); 964 if (external_lb > current_lb) { 965 if (!integer_trail->Enqueue(IntegerLiteral::GreaterOrEqual( 966 objective->objective_var, external_lb), 967 {}, {})) { 968 return false; 969 } 970 propagate = true; 971 } 972 973 const IntegerValue external_ub = 974 shared_response_manager->SynchronizedInnerObjectiveUpperBound(); 975 const IntegerValue current_ub = 976 integer_trail->UpperBound(objective->objective_var); 977 if (external_ub < current_ub) { 978 if (!integer_trail->Enqueue(IntegerLiteral::LowerOrEqual( 979 objective->objective_var, external_ub), 980 {}, {})) { 981 return false; 982 } 983 propagate = true; 984 } 985 986 if (!propagate) return true; 987 988 VLOG(3) << "'" << name << "' imports objective bounds: external [" 989 << objective->ScaleIntegerObjective(external_lb) << ", " 990 << objective->ScaleIntegerObjective(external_ub) << "], current [" 991 << objective->ScaleIntegerObjective(current_lb) << ", " 992 << objective->ScaleIntegerObjective(current_ub) << "]"; 993 994 return solver->FinishPropagation(); 995 }; 996 997 model->GetOrCreate<LevelZeroCallbackHelper>()->callbacks.push_back( 998 import_objective_bounds); 999 } 1000 1001 void LoadBaseModel(const CpModelProto& model_proto, Model* model) { 1002 auto* shared_response_manager = model->GetOrCreate<SharedResponseManager>(); 1003 CHECK(shared_response_manager != nullptr); 1004 auto* sat_solver = model->GetOrCreate<SatSolver>(); 1005 1006 // Simple function for the few places where we do "return unsat()". 1007 const auto unsat = [shared_response_manager, sat_solver, model] { 1008 sat_solver->NotifyThatModelIsUnsat(); 1009 shared_response_manager->NotifyThatImprovingProblemIsInfeasible( 1010 absl::StrCat(model->Name(), " [loading]")); 1011 }; 1012 1013 // We will add them all at once after model_proto is loaded. 1014 model->GetOrCreate<IntegerEncoder>()->DisableImplicationBetweenLiteral(); 1015 1016 auto* mapping = model->GetOrCreate<CpModelMapping>(); 1017 const SatParameters& parameters = *(model->GetOrCreate<SatParameters>()); 1018 const bool view_all_booleans_as_integers = 1019 (parameters.linearization_level() >= 2) || 1020 (parameters.search_branching() == SatParameters::FIXED_SEARCH && 1021 model_proto.search_strategy().empty()); 1022 LoadVariables(model_proto, view_all_booleans_as_integers, model); 1023 DetectOptionalVariables(model_proto, model); 1024 1025 // TODO(user): The core algo and symmetries seems to be problematic in some 1026 // cases. See for instance: neos-691058.mps.gz. This is probably because as 1027 // we modify the model, our symmetry might be wrong? investigate. 1028 if (!parameters.optimize_with_core() && parameters.symmetry_level() > 1 && 1029 !parameters.enumerate_all_solutions()) { 1030 LoadBooleanSymmetries(model_proto, model); 1031 } 1032 1033 ExtractEncoding(model_proto, model); 1034 ExtractElementEncoding(model_proto, model); 1035 PropagateEncodingFromEquivalenceRelations(model_proto, model); 1036 1037 // Check the model is still feasible before continuing. 1038 if (sat_solver->IsModelUnsat()) return unsat(); 1039 1040 // Fully encode variables as needed by the search strategy. 1041 AddFullEncodingFromSearchBranching(model_proto, model); 1042 1043 // Load the constraints. 1044 std::set<std::string> unsupported_types; 1045 int num_ignored_constraints = 0; 1046 for (const ConstraintProto& ct : model_proto.constraints()) { 1047 if (mapping->ConstraintIsAlreadyLoaded(&ct)) { 1048 ++num_ignored_constraints; 1049 continue; 1050 } 1051 1052 if (!LoadConstraint(ct, model)) { 1053 unsupported_types.insert(ConstraintCaseName(ct.constraint_case())); 1054 continue; 1055 } 1056 1057 // We propagate after each new Boolean constraint but not the integer 1058 // ones. So we call FinishPropagation() manually here. 1059 // 1060 // Note that we only do that in debug mode as this can be really slow on 1061 // certain types of problems with millions of constraints. 1062 if (DEBUG_MODE) { 1063 if (sat_solver->FinishPropagation()) { 1064 Trail* trail = model->GetOrCreate<Trail>(); 1065 const int old_num_fixed = trail->Index(); 1066 if (trail->Index() > old_num_fixed) { 1067 VLOG(3) << "Constraint fixed " << trail->Index() - old_num_fixed 1068 << " Boolean variable(s): " << ProtobufDebugString(ct); 1069 } 1070 } 1071 } 1072 if (sat_solver->IsModelUnsat()) { 1073 VLOG(2) << "UNSAT during extraction (after adding '" 1074 << ConstraintCaseName(ct.constraint_case()) << "'). " 1075 << ProtobufDebugString(ct); 1076 break; 1077 } 1078 } 1079 if (num_ignored_constraints > 0) { 1080 VLOG(3) << num_ignored_constraints << " constraints were skipped."; 1081 } 1082 if (!unsupported_types.empty()) { 1083 VLOG(1) << "There is unsupported constraints types in this model: "; 1084 for (const std::string& type : unsupported_types) { 1085 VLOG(1) << " - " << type; 1086 } 1087 return unsat(); 1088 } 1089 1090 model->GetOrCreate<IntegerEncoder>() 1091 ->AddAllImplicationsBetweenAssociatedLiterals(); 1092 if (!sat_solver->FinishPropagation()) return unsat(); 1093 } 1094 1095 void LoadFeasibilityPump(const CpModelProto& model_proto, Model* model) { 1096 LoadBaseModel(model_proto, model); 1097 1098 auto* mapping = model->GetOrCreate<CpModelMapping>(); 1099 const SatParameters& parameters = *(model->GetOrCreate<SatParameters>()); 1100 if (parameters.linearization_level() == 0) return; 1101 1102 // Add linear constraints to Feasibility Pump. 1103 const LinearRelaxation relaxation = 1104 ComputeLinearRelaxation(model_proto, model); 1105 const int num_lp_constraints = relaxation.linear_constraints.size(); 1106 if (num_lp_constraints == 0) return; 1107 auto* feasibility_pump = model->GetOrCreate<FeasibilityPump>(); 1108 for (int i = 0; i < num_lp_constraints; i++) { 1109 feasibility_pump->AddLinearConstraint(relaxation.linear_constraints[i]); 1110 } 1111 1112 if (model_proto.has_objective()) { 1113 for (int i = 0; i < model_proto.objective().coeffs_size(); ++i) { 1114 const IntegerVariable var = 1115 mapping->Integer(model_proto.objective().vars(i)); 1116 const int64_t coeff = model_proto.objective().coeffs(i); 1117 feasibility_pump->SetObjectiveCoefficient(var, IntegerValue(coeff)); 1118 } 1119 } 1120 } 1121 1122 // Loads a CpModelProto inside the given model. 1123 // This should only be called once on a given 'Model' class. 1124 // 1125 // TODO(user): move to cp_model_loader.h/.cc 1126 void LoadCpModel(const CpModelProto& model_proto, Model* model) { 1127 auto* shared_response_manager = model->GetOrCreate<SharedResponseManager>(); 1128 1129 LoadBaseModel(model_proto, model); 1130 1131 // Simple function for the few places where we do "return unsat()". 1132 auto* sat_solver = model->GetOrCreate<SatSolver>(); 1133 const auto unsat = [shared_response_manager, sat_solver, model] { 1134 sat_solver->NotifyThatModelIsUnsat(); 1135 shared_response_manager->NotifyThatImprovingProblemIsInfeasible( 1136 absl::StrCat(model->Name(), " [loading]")); 1137 }; 1138 1139 auto* mapping = model->GetOrCreate<CpModelMapping>(); 1140 const SatParameters& parameters = *(model->GetOrCreate<SatParameters>()); 1141 1142 // Auto detect "at least one of" constraints in the PrecedencesPropagator. 1143 // Note that we do that before we finish loading the problem (objective and 1144 // LP relaxation), because propagation will be faster at this point and it 1145 // should be enough for the purpose of this auto-detection. 1146 if (model->Mutable<PrecedencesPropagator>() != nullptr && 1147 parameters.auto_detect_greater_than_at_least_one_of()) { 1148 model->Mutable<PrecedencesPropagator>() 1149 ->AddGreaterThanAtLeastOneOfConstraints(model); 1150 if (!sat_solver->FinishPropagation()) return unsat(); 1151 } 1152 1153 // TODO(user): This should be done in the presolve instead. 1154 // TODO(user): We don't have a good deterministic time on all constraints, 1155 // so this might take more time than wanted. 1156 if (parameters.cp_model_probing_level() > 1) { 1157 Prober* prober = model->GetOrCreate<Prober>(); 1158 prober->ProbeBooleanVariables(/*deterministic_time_limit=*/1.0); 1159 if (model->GetOrCreate<SatSolver>()->IsModelUnsat()) { 1160 return unsat(); 1161 } 1162 if (!model->GetOrCreate<BinaryImplicationGraph>() 1163 ->ComputeTransitiveReduction()) { 1164 return unsat(); 1165 } 1166 } 1167 1168 // Create an objective variable and its associated linear constraint if 1169 // needed. 1170 IntegerVariable objective_var = kNoIntegerVariable; 1171 if (parameters.linearization_level() > 0) { 1172 // Linearize some part of the problem and register LP constraint(s). 1173 objective_var = AddLPConstraints(model_proto, model); 1174 } else if (model_proto.has_objective()) { 1175 const CpObjectiveProto& obj = model_proto.objective(); 1176 std::vector<std::pair<IntegerVariable, int64_t>> terms; 1177 terms.reserve(obj.vars_size()); 1178 for (int i = 0; i < obj.vars_size(); ++i) { 1179 terms.push_back( 1180 std::make_pair(mapping->Integer(obj.vars(i)), obj.coeffs(i))); 1181 } 1182 if (parameters.optimize_with_core()) { 1183 objective_var = GetOrCreateVariableWithTightBound(terms, model); 1184 } else { 1185 objective_var = GetOrCreateVariableGreaterOrEqualToSumOf(terms, model); 1186 } 1187 } 1188 1189 // Create the objective definition inside the Model so that it can be accessed 1190 // by the heuristics than needs it. 1191 if (objective_var != kNoIntegerVariable) { 1192 const CpObjectiveProto& objective_proto = model_proto.objective(); 1193 auto* objective_definition = model->GetOrCreate<ObjectiveDefinition>(); 1194 1195 objective_definition->scaling_factor = objective_proto.scaling_factor(); 1196 if (objective_definition->scaling_factor == 0.0) { 1197 objective_definition->scaling_factor = 1.0; 1198 } 1199 objective_definition->offset = objective_proto.offset(); 1200 objective_definition->objective_var = objective_var; 1201 1202 const int size = objective_proto.vars_size(); 1203 objective_definition->vars.resize(size); 1204 objective_definition->coeffs.resize(size); 1205 for (int i = 0; i < objective_proto.vars_size(); ++i) { 1206 // Note that if there is no mapping, then the variable will be 1207 // kNoIntegerVariable. 1208 objective_definition->vars[i] = mapping->Integer(objective_proto.vars(i)); 1209 objective_definition->coeffs[i] = IntegerValue(objective_proto.coeffs(i)); 1210 1211 // Fill the objective heuristics data. 1212 const int ref = objective_proto.vars(i); 1213 if (mapping->IsInteger(ref)) { 1214 const IntegerVariable var = mapping->Integer(objective_proto.vars(i)); 1215 objective_definition->objective_impacting_variables.insert( 1216 objective_proto.coeffs(i) > 0 ? var : NegationOf(var)); 1217 } 1218 } 1219 1220 // Register an objective special propagator. 1221 model->TakeOwnership( 1222 new LevelZeroEquality(objective_var, objective_definition->vars, 1223 objective_definition->coeffs, model)); 1224 } 1225 1226 // Intersect the objective domain with the given one if any. 1227 if (!model_proto.objective().domain().empty()) { 1228 const Domain user_domain = ReadDomainFromProto(model_proto.objective()); 1229 const Domain automatic_domain = 1230 model->GetOrCreate<IntegerTrail>()->InitialVariableDomain( 1231 objective_var); 1232 VLOG(3) << "Objective offset:" << model_proto.objective().offset() 1233 << " scaling_factor:" << model_proto.objective().scaling_factor(); 1234 VLOG(3) << "Automatic internal objective domain: " << automatic_domain; 1235 VLOG(3) << "User specified internal objective domain: " << user_domain; 1236 CHECK_NE(objective_var, kNoIntegerVariable); 1237 const bool ok = model->GetOrCreate<IntegerTrail>()->UpdateInitialDomain( 1238 objective_var, user_domain); 1239 if (!ok) { 1240 VLOG(2) << "UNSAT due to the objective domain."; 1241 return unsat(); 1242 } 1243 1244 // Make sure the sum take a value inside the objective domain by adding 1245 // the other side: objective <= sum terms. 1246 // 1247 // TODO(user): Use a better condition to detect when this is not useful. 1248 // Note also that for the core algorithm, we might need the other side too, 1249 // otherwise we could return feasible solution with an objective above the 1250 // user specified upper bound. 1251 if (!automatic_domain.IsIncludedIn(user_domain)) { 1252 std::vector<IntegerVariable> vars; 1253 std::vector<int64_t> coeffs; 1254 const CpObjectiveProto& obj = model_proto.objective(); 1255 for (int i = 0; i < obj.vars_size(); ++i) { 1256 vars.push_back(mapping->Integer(obj.vars(i))); 1257 coeffs.push_back(obj.coeffs(i)); 1258 } 1259 vars.push_back(objective_var); 1260 coeffs.push_back(-1); 1261 model->Add(WeightedSumGreaterOrEqual(vars, coeffs, 0)); 1262 } 1263 } 1264 1265 // Note that we do one last propagation at level zero once all the 1266 // constraints were added. 1267 SOLVER_LOG(model->GetOrCreate<SolverLogger>(), 1268 "Initial num_bool: ", sat_solver->NumVariables()); 1269 if (!sat_solver->FinishPropagation()) return unsat(); 1270 1271 if (model_proto.has_objective()) { 1272 // Report the initial objective variable bounds. 1273 auto* integer_trail = model->GetOrCreate<IntegerTrail>(); 1274 shared_response_manager->UpdateInnerObjectiveBounds( 1275 absl::StrCat(model->Name(), " initial_propagation"), 1276 integer_trail->LowerBound(objective_var), 1277 integer_trail->UpperBound(objective_var)); 1278 1279 // Watch improved objective best bounds. 1280 RegisterObjectiveBestBoundExport(objective_var, shared_response_manager, 1281 model); 1282 1283 // Import objective bounds. 1284 // TODO(user): Support objective bounds import in LNS and Core based 1285 // search. 1286 if (model->GetOrCreate<SatParameters>()->share_objective_bounds()) { 1287 RegisterObjectiveBoundsImport(shared_response_manager, model); 1288 } 1289 } 1290 1291 // Cache the links between model vars, IntegerVariables and lp constraints. 1292 // TODO(user): Cache this only if it is actually used. 1293 auto* integer_trail = model->GetOrCreate<IntegerTrail>(); 1294 auto* lp_dispatcher = model->GetOrCreate<LinearProgrammingDispatcher>(); 1295 auto* lp_vars = model->GetOrCreate<LPVariables>(); 1296 IntegerVariable size = integer_trail->NumIntegerVariables(); 1297 for (IntegerVariable positive_var(0); positive_var < size; 1298 positive_var += 2) { 1299 LPVariable lp_var; 1300 lp_var.positive_var = positive_var; 1301 lp_var.model_var = 1302 mapping->GetProtoVariableFromIntegerVariable(positive_var); 1303 lp_var.lp = gtl::FindWithDefault(*lp_dispatcher, positive_var, nullptr); 1304 1305 if (lp_var.model_var >= 0) { 1306 lp_vars->vars.push_back(lp_var); 1307 lp_vars->model_vars_size = 1308 std::max(lp_vars->model_vars_size, lp_var.model_var + 1); 1309 } 1310 } 1311 1312 // Initialize the fixed_search strategy. 1313 auto* search_heuristics = model->GetOrCreate<SearchHeuristics>(); 1314 search_heuristics->fixed_search = ConstructSearchStrategy( 1315 model_proto, mapping->GetVariableMapping(), objective_var, model); 1316 if (VLOG_IS_ON(3)) { 1317 search_heuristics->fixed_search = 1318 InstrumentSearchStrategy(model_proto, mapping->GetVariableMapping(), 1319 search_heuristics->fixed_search, model); 1320 } 1321 1322 // Initialize the "follow hint" strategy. 1323 std::vector<BooleanOrIntegerVariable> vars; 1324 std::vector<IntegerValue> values; 1325 for (int i = 0; i < model_proto.solution_hint().vars_size(); ++i) { 1326 const int ref = model_proto.solution_hint().vars(i); 1327 CHECK(RefIsPositive(ref)); 1328 BooleanOrIntegerVariable var; 1329 if (mapping->IsBoolean(ref)) { 1330 var.bool_var = mapping->Literal(ref).Variable(); 1331 } else { 1332 var.int_var = mapping->Integer(ref); 1333 } 1334 vars.push_back(var); 1335 values.push_back(IntegerValue(model_proto.solution_hint().values(i))); 1336 } 1337 search_heuristics->hint_search = FollowHint(vars, values, model); 1338 1339 // Create the CoreBasedOptimizer class if needed. 1340 if (parameters.optimize_with_core()) { 1341 // TODO(user): Remove code duplication with the solution_observer in 1342 // SolveLoadedCpModel(). 1343 const std::string solution_info = model->Name(); 1344 const auto solution_observer = [&model_proto, model, solution_info, 1345 shared_response_manager]() { 1346 CpSolverResponse response; 1347 FillSolutionInResponse(model_proto, *model, &response); 1348 response.set_solution_info(solution_info); 1349 shared_response_manager->NewSolution(response, model); 1350 }; 1351 1352 const auto& objective = *model->GetOrCreate<ObjectiveDefinition>(); 1353 CoreBasedOptimizer* core = 1354 new CoreBasedOptimizer(objective_var, objective.vars, objective.coeffs, 1355 solution_observer, model); 1356 model->Register<CoreBasedOptimizer>(core); 1357 model->TakeOwnership(core); 1358 } 1359 } 1360 1361 // Solves an already loaded cp_model_proto. 1362 // The final CpSolverResponse must be read from the shared_response_manager. 1363 // 1364 // TODO(user): This should be transformed so that it can be called many times 1365 // and resume from the last search state as if it wasn't interuped. That would 1366 // allow use to easily interleave different heuristics in the same thread. 1367 void SolveLoadedCpModel(const CpModelProto& model_proto, Model* model) { 1368 auto* shared_response_manager = model->GetOrCreate<SharedResponseManager>(); 1369 if (shared_response_manager->ProblemIsSolved()) return; 1370 1371 const std::string& solution_info = model->Name(); 1372 const auto solution_observer = [&model_proto, &model, &solution_info, 1373 &shared_response_manager]() { 1374 CpSolverResponse response; 1375 FillSolutionInResponse(model_proto, *model, &response); 1376 response.set_solution_info(solution_info); 1377 shared_response_manager->NewSolution(response, model); 1378 }; 1379 1380 // Reconfigure search heuristic if it was changed. 1381 ConfigureSearchHeuristics(model); 1382 1383 const auto& mapping = *model->GetOrCreate<CpModelMapping>(); 1384 SatSolver::Status status; 1385 const SatParameters& parameters = *model->GetOrCreate<SatParameters>(); 1386 if (parameters.use_probing_search()) { 1387 std::vector<BooleanVariable> bool_vars; 1388 std::vector<IntegerVariable> int_vars; 1389 IntegerTrail* integer_trail = model->GetOrCreate<IntegerTrail>(); 1390 absl::flat_hash_set<BooleanVariable> visited; 1391 for (int v = 0; v < model_proto.variables_size(); ++v) { 1392 if (mapping.IsBoolean(v)) { 1393 const BooleanVariable bool_var = mapping.Literal(v).Variable(); 1394 if (!visited.contains(bool_var)) { 1395 visited.insert(bool_var); 1396 bool_vars.push_back(bool_var); 1397 } 1398 } else { 1399 IntegerVariable var = mapping.Integer(v); 1400 if (integer_trail->IsFixed(var)) continue; 1401 int_vars.push_back(var); 1402 } 1403 } 1404 status = ContinuousProbing(bool_vars, int_vars, solution_observer, model); 1405 } else if (!model_proto.has_objective()) { 1406 while (true) { 1407 status = ResetAndSolveIntegerProblem( 1408 mapping.Literals(model_proto.assumptions()), model); 1409 if (status != SatSolver::Status::FEASIBLE) break; 1410 solution_observer(); 1411 if (!parameters.enumerate_all_solutions()) break; 1412 model->Add(ExcludeCurrentSolutionWithoutIgnoredVariableAndBacktrack()); 1413 } 1414 if (status == SatSolver::INFEASIBLE) { 1415 shared_response_manager->NotifyThatImprovingProblemIsInfeasible( 1416 solution_info); 1417 } 1418 if (status == SatSolver::ASSUMPTIONS_UNSAT) { 1419 shared_response_manager->NotifyThatImprovingProblemIsInfeasible( 1420 solution_info); 1421 1422 // Extract a good subset of assumptions and add it to the response. 1423 auto* time_limit = model->GetOrCreate<TimeLimit>(); 1424 auto* sat_solver = model->GetOrCreate<SatSolver>(); 1425 std::vector<Literal> core = sat_solver->GetLastIncompatibleDecisions(); 1426 MinimizeCoreWithPropagation(time_limit, sat_solver, &core); 1427 std::vector<int> core_in_proto_format; 1428 for (const Literal l : core) { 1429 core_in_proto_format.push_back( 1430 mapping.GetProtoVariableFromBooleanVariable(l.Variable())); 1431 if (!l.IsPositive()) { 1432 core_in_proto_format.back() = NegatedRef(core_in_proto_format.back()); 1433 } 1434 } 1435 shared_response_manager->AddUnsatCore(core_in_proto_format); 1436 } 1437 } else { 1438 // Optimization problem. 1439 const auto& objective = *model->GetOrCreate<ObjectiveDefinition>(); 1440 const IntegerVariable objective_var = objective.objective_var; 1441 CHECK_NE(objective_var, kNoIntegerVariable); 1442 1443 if (parameters.optimize_with_lb_tree_search()) { 1444 auto* search = model->GetOrCreate<LbTreeSearch>(); 1445 status = search->Search(solution_observer); 1446 } else if (parameters.optimize_with_core()) { 1447 // TODO(user): This doesn't work with splitting in chunk for now. It 1448 // shouldn't be too hard to fix. 1449 if (parameters.optimize_with_max_hs()) { 1450 status = MinimizeWithHittingSetAndLazyEncoding( 1451 objective, solution_observer, model); 1452 } else { 1453 status = model->Mutable<CoreBasedOptimizer>()->Optimize(); 1454 } 1455 } else { 1456 // TODO(user): This parameter break the splitting in chunk of a Solve(). 1457 // It should probably be moved into another SubSolver altogether. 1458 if (parameters.binary_search_num_conflicts() >= 0) { 1459 RestrictObjectiveDomainWithBinarySearch(objective_var, 1460 solution_observer, model); 1461 } 1462 status = MinimizeIntegerVariableWithLinearScanAndLazyEncoding( 1463 objective_var, solution_observer, model); 1464 } 1465 1466 // The search is done in both case. 1467 // 1468 // TODO(user): Remove the weird translation INFEASIBLE->FEASIBLE in the 1469 // function above? 1470 if (status == SatSolver::INFEASIBLE || status == SatSolver::FEASIBLE) { 1471 shared_response_manager->NotifyThatImprovingProblemIsInfeasible( 1472 solution_info); 1473 } 1474 } 1475 1476 // TODO(user): Remove from here when we split in chunk. We just want to 1477 // do that at the end. 1478 shared_response_manager->SetStatsFromModel(model); 1479 } 1480 1481 // Try to find a solution by following the hint and using a low conflict limit. 1482 // The CpModelProto must already be loaded in the Model. 1483 void QuickSolveWithHint(const CpModelProto& model_proto, Model* model) { 1484 if (!model_proto.has_solution_hint()) return; 1485 1486 auto* shared_response_manager = model->GetOrCreate<SharedResponseManager>(); 1487 if (shared_response_manager->ProblemIsSolved()) return; 1488 1489 // Temporarily change the parameters. 1490 auto* parameters = model->GetOrCreate<SatParameters>(); 1491 1492 // If the model was loaded with "optimize_with_core" then the objective 1493 // variable is not linked to its linear expression. Because of that, we can 1494 // return a solution that does not satisfy the objective domain. 1495 // 1496 // TODO(user): This is fixable, but then do we need the hint when optimizing 1497 // with core? 1498 if (parameters->optimize_with_core()) return; 1499 1500 const SatParameters saved_params = *parameters; 1501 parameters->set_max_number_of_conflicts(parameters->hint_conflict_limit()); 1502 parameters->set_search_branching(SatParameters::HINT_SEARCH); 1503 parameters->set_optimize_with_core(false); 1504 auto cleanup = ::absl::MakeCleanup( 1505 [parameters, saved_params]() { *parameters = saved_params; }); 1506 1507 // Solve decision problem. 1508 ConfigureSearchHeuristics(model); 1509 const auto& mapping = *model->GetOrCreate<CpModelMapping>(); 1510 const SatSolver::Status status = ResetAndSolveIntegerProblem( 1511 mapping.Literals(model_proto.assumptions()), model); 1512 1513 const std::string& solution_info = model->Name(); 1514 if (status == SatSolver::Status::FEASIBLE) { 1515 CpSolverResponse response; 1516 FillSolutionInResponse(model_proto, *model, &response); 1517 response.set_solution_info(absl::StrCat(solution_info, " [hint]")); 1518 shared_response_manager->NewSolution(response, model); 1519 1520 if (!model_proto.has_objective()) { 1521 if (parameters->enumerate_all_solutions()) { 1522 model->Add(ExcludeCurrentSolutionWithoutIgnoredVariableAndBacktrack()); 1523 } 1524 } else { 1525 // Restrict the objective. 1526 const IntegerVariable objective_var = 1527 model->GetOrCreate<ObjectiveDefinition>()->objective_var; 1528 model->GetOrCreate<SatSolver>()->Backtrack(0); 1529 IntegerTrail* integer_trail = model->GetOrCreate<IntegerTrail>(); 1530 if (!integer_trail->Enqueue( 1531 IntegerLiteral::LowerOrEqual( 1532 objective_var, 1533 shared_response_manager->GetInnerObjectiveUpperBound()), 1534 {}, {})) { 1535 shared_response_manager->NotifyThatImprovingProblemIsInfeasible( 1536 absl::StrCat(solution_info, " [hint]")); 1537 shared_response_manager->SetStatsFromModel(model); 1538 return; 1539 } 1540 } 1541 } 1542 1543 // This code is here to debug bad presolve during LNS that corrupt the hint. 1544 // Note that sometime the deterministic limit is hit before the hint can be 1545 // completed, so we don't report that has an error. 1546 if (parameters->debug_crash_on_bad_hint() && 1547 !model->GetOrCreate<TimeLimit>()->LimitReached() && 1548 status != SatSolver::Status::FEASIBLE) { 1549 LOG(FATAL) << "QuickSolveWithHint() didn't find a feasible solution. " 1550 << "The model name is '" << model_proto.name() << "'."; 1551 } 1552 } 1553 1554 // Solve a model with a different objective consisting of minimizing the L1 1555 // distance with the provided hint. Note that this method creates an in-memory 1556 // copy of the model and loads a local Model object from the copied model. 1557 void MinimizeL1DistanceWithHint(const CpModelProto& model_proto, Model* model) { 1558 Model local_model; 1559 1560 // Forward some shared class. 1561 local_model.Register<ModelSharedTimeLimit>( 1562 model->GetOrCreate<ModelSharedTimeLimit>()); 1563 local_model.Register<WallTimer>(model->GetOrCreate<WallTimer>()); 1564 1565 if (!model_proto.has_solution_hint()) return; 1566 1567 auto* shared_response_manager = model->GetOrCreate<SharedResponseManager>(); 1568 if (shared_response_manager->ProblemIsSolved()) return; 1569 1570 auto* parameters = local_model.GetOrCreate<SatParameters>(); 1571 // TODO(user): As of now the repair hint doesn't support when 1572 // enumerate_all_solutions is set since the solution is created on a different 1573 // model. 1574 if (parameters->enumerate_all_solutions()) return; 1575 1576 // Change the parameters. 1577 const SatParameters saved_params = *model->GetOrCreate<SatParameters>(); 1578 *parameters = saved_params; 1579 parameters->set_max_number_of_conflicts(parameters->hint_conflict_limit()); 1580 parameters->set_optimize_with_core(false); 1581 1582 // Update the model to introduce penalties to go away from hinted values. 1583 CpModelProto updated_model_proto = model_proto; 1584 updated_model_proto.clear_objective(); 1585 1586 // TODO(user): For boolean variables we can avoid creating new variables. 1587 for (int i = 0; i < model_proto.solution_hint().vars_size(); ++i) { 1588 const int var = model_proto.solution_hint().vars(i); 1589 const int64_t value = model_proto.solution_hint().values(i); 1590 1591 // Add a new var to represent the difference between var and value. 1592 const int new_var_index = updated_model_proto.variables_size(); 1593 IntegerVariableProto* var_proto = updated_model_proto.add_variables(); 1594 const int64_t min_domain = model_proto.variables(var).domain(0) - value; 1595 const int64_t max_domain = 1596 model_proto.variables(var).domain( 1597 model_proto.variables(var).domain_size() - 1) - 1598 value; 1599 var_proto->add_domain(min_domain); 1600 var_proto->add_domain(max_domain); 1601 1602 // new_var = var - value. 1603 ConstraintProto* const linear_constraint_proto = 1604 updated_model_proto.add_constraints(); 1605 LinearConstraintProto* linear = linear_constraint_proto->mutable_linear(); 1606 linear->add_vars(new_var_index); 1607 linear->add_coeffs(1); 1608 linear->add_vars(var); 1609 linear->add_coeffs(-1); 1610 linear->add_domain(-value); 1611 linear->add_domain(-value); 1612 1613 // abs_var = abs(new_var). 1614 const int abs_var_index = updated_model_proto.variables_size(); 1615 IntegerVariableProto* abs_var_proto = updated_model_proto.add_variables(); 1616 const int64_t abs_min_domain = 0; 1617 const int64_t abs_max_domain = 1618 std::max(std::abs(min_domain), std::abs(max_domain)); 1619 abs_var_proto->add_domain(abs_min_domain); 1620 abs_var_proto->add_domain(abs_max_domain); 1621 auto* abs_ct = updated_model_proto.add_constraints()->mutable_lin_max(); 1622 abs_ct->mutable_target()->add_vars(abs_var_index); 1623 abs_ct->mutable_target()->add_coeffs(1); 1624 LinearExpressionProto* left = abs_ct->add_exprs(); 1625 left->add_vars(new_var_index); 1626 left->add_coeffs(1); 1627 LinearExpressionProto* right = abs_ct->add_exprs(); 1628 right->add_vars(new_var_index); 1629 right->add_coeffs(-1); 1630 1631 updated_model_proto.mutable_objective()->add_vars(abs_var_index); 1632 updated_model_proto.mutable_objective()->add_coeffs(1); 1633 } 1634 1635 auto* local_response_manager = 1636 local_model.GetOrCreate<SharedResponseManager>(); 1637 local_response_manager->InitializeObjective(updated_model_proto); 1638 1639 // Solve optimization problem. 1640 LoadCpModel(updated_model_proto, &local_model); 1641 1642 ConfigureSearchHeuristics(&local_model); 1643 const auto& mapping = *local_model.GetOrCreate<CpModelMapping>(); 1644 const SatSolver::Status status = ResetAndSolveIntegerProblem( 1645 mapping.Literals(updated_model_proto.assumptions()), &local_model); 1646 1647 const std::string& solution_info = model->Name(); 1648 if (status == SatSolver::Status::FEASIBLE) { 1649 CpSolverResponse response; 1650 FillSolutionInResponse(model_proto, local_model, &response); 1651 if (DEBUG_MODE) { 1652 CpSolverResponse updated_response; 1653 FillSolutionInResponse(updated_model_proto, local_model, 1654 &updated_response); 1655 LOG(INFO) << "Found solution with repaired hint penalty = " 1656 << ComputeInnerObjective(updated_model_proto.objective(), 1657 updated_response); 1658 } 1659 response.set_solution_info(absl::StrCat(solution_info, " [repaired]")); 1660 shared_response_manager->NewSolution(response, &local_model); 1661 } 1662 } 1663 1664 // TODO(user): If this ever shows up in the profile, we could avoid copying 1665 // the mapping_proto if we are careful about how we modify the variable domain 1666 // before postsolving it. Note that 'num_variables_in_original_model' refers to 1667 // the model before presolve. 1668 void PostsolveResponseWithFullSolver(int num_variables_in_original_model, 1669 CpModelProto mapping_proto, 1670 const std::vector<int>& postsolve_mapping, 1671 std::vector<int64_t>* solution) { 1672 WallTimer wall_timer; 1673 wall_timer.Start(); 1674 1675 // Fix the correct variable in the mapping_proto. 1676 for (int i = 0; i < solution->size(); ++i) { 1677 auto* var_proto = mapping_proto.mutable_variables(postsolve_mapping[i]); 1678 var_proto->clear_domain(); 1679 var_proto->add_domain((*solution)[i]); 1680 var_proto->add_domain((*solution)[i]); 1681 } 1682 1683 // Postosolve parameters. 1684 // TODO(user): this problem is usually trivial, but we may still want to 1685 // impose a time limit or copy some of the parameters passed by the user. 1686 Model postsolve_model; 1687 postsolve_model.Register<WallTimer>(&wall_timer); 1688 { 1689 SatParameters& params = *postsolve_model.GetOrCreate<SatParameters>(); 1690 params.set_linearization_level(0); 1691 params.set_cp_model_probing_level(0); 1692 } 1693 1694 LoadCpModel(mapping_proto, &postsolve_model); 1695 SolveLoadedCpModel(mapping_proto, &postsolve_model); 1696 const CpSolverResponse postsolve_response = 1697 postsolve_model.GetOrCreate<SharedResponseManager>()->GetResponse(); 1698 CHECK(postsolve_response.status() == CpSolverStatus::FEASIBLE || 1699 postsolve_response.status() == CpSolverStatus::OPTIMAL) 1700 << CpSolverResponseStats(postsolve_response); 1701 1702 // We only copy the solution from the postsolve_response to the response. 1703 CHECK_LE(num_variables_in_original_model, 1704 postsolve_response.solution().size()); 1705 solution->assign( 1706 postsolve_response.solution().begin(), 1707 postsolve_response.solution().begin() + num_variables_in_original_model); 1708 } 1709 1710 void PostsolveResponseWrapper(const SatParameters& params, 1711 int num_variable_in_original_model, 1712 const CpModelProto& mapping_proto, 1713 const std::vector<int>& postsolve_mapping, 1714 std::vector<int64_t>* solution) { 1715 if (params.debug_postsolve_with_full_solver()) { 1716 PostsolveResponseWithFullSolver(num_variable_in_original_model, 1717 mapping_proto, postsolve_mapping, solution); 1718 } else { 1719 PostsolveResponse(num_variable_in_original_model, mapping_proto, 1720 postsolve_mapping, solution); 1721 } 1722 } 1723 1724 // TODO(user): Uniformize this function with the other one. 1725 CpSolverResponse SolvePureSatModel(const CpModelProto& model_proto, 1726 WallTimer* wall_timer, Model* model, 1727 SolverLogger* logger) { 1728 std::unique_ptr<SatSolver> solver(new SatSolver()); 1729 SatParameters parameters = *model->GetOrCreate<SatParameters>(); 1730 solver->SetParameters(parameters); 1731 model->GetOrCreate<TimeLimit>()->ResetLimitFromParameters(parameters); 1732 1733 // Create a DratProofHandler? 1734 std::unique_ptr<DratProofHandler> drat_proof_handler; 1735 #if !defined(__PORTABLE_PLATFORM__) 1736 if (!absl::GetFlag(FLAGS_drat_output).empty() || 1737 absl::GetFlag(FLAGS_drat_check)) { 1738 if (!absl::GetFlag(FLAGS_drat_output).empty()) { 1739 File* output; 1740 CHECK_OK(file::Open(absl::GetFlag(FLAGS_drat_output), "w", &output, 1741 file::Defaults())); 1742 drat_proof_handler = absl::make_unique<DratProofHandler>( 1743 /*in_binary_format=*/false, output, absl::GetFlag(FLAGS_drat_check)); 1744 } else { 1745 drat_proof_handler = absl::make_unique<DratProofHandler>(); 1746 } 1747 solver->SetDratProofHandler(drat_proof_handler.get()); 1748 } 1749 #endif // __PORTABLE_PLATFORM__ 1750 1751 auto get_literal = [](int ref) { 1752 if (ref >= 0) return Literal(BooleanVariable(ref), true); 1753 return Literal(BooleanVariable(NegatedRef(ref)), false); 1754 }; 1755 1756 std::vector<Literal> temp; 1757 const int num_variables = model_proto.variables_size(); 1758 solver->SetNumVariables(num_variables); 1759 if (drat_proof_handler != nullptr) { 1760 drat_proof_handler->SetNumVariables(num_variables); 1761 1762 // We load the model in the drat_proof_handler for the case where we want 1763 // to do in-memory checking. 1764 for (int ref = 0; ref < num_variables; ++ref) { 1765 const Domain domain = ReadDomainFromProto(model_proto.variables(ref)); 1766 if (domain.IsFixed()) { 1767 const Literal ref_literal = 1768 domain.Min() == 0 ? get_literal(ref).Negated() : get_literal(ref); 1769 drat_proof_handler->AddProblemClause({ref_literal}); 1770 } 1771 } 1772 for (const ConstraintProto& ct : model_proto.constraints()) { 1773 switch (ct.constraint_case()) { 1774 case ConstraintProto::ConstraintCase::kBoolAnd: { 1775 if (ct.enforcement_literal_size() == 0) { 1776 for (const int ref : ct.bool_and().literals()) { 1777 drat_proof_handler->AddProblemClause({get_literal(ref)}); 1778 } 1779 } else { 1780 // a => b 1781 const Literal not_a = 1782 get_literal(ct.enforcement_literal(0)).Negated(); 1783 for (const int ref : ct.bool_and().literals()) { 1784 drat_proof_handler->AddProblemClause({not_a, get_literal(ref)}); 1785 } 1786 } 1787 break; 1788 } 1789 case ConstraintProto::ConstraintCase::kBoolOr: 1790 temp.clear(); 1791 for (const int ref : ct.bool_or().literals()) { 1792 temp.push_back(get_literal(ref)); 1793 } 1794 for (const int ref : ct.enforcement_literal()) { 1795 temp.push_back(get_literal(ref).Negated()); 1796 } 1797 drat_proof_handler->AddProblemClause(temp); 1798 break; 1799 default: 1800 LOG(FATAL) << "Not supported"; 1801 } 1802 } 1803 } 1804 1805 for (const ConstraintProto& ct : model_proto.constraints()) { 1806 switch (ct.constraint_case()) { 1807 case ConstraintProto::ConstraintCase::kBoolAnd: { 1808 if (ct.enforcement_literal_size() == 0) { 1809 for (const int ref : ct.bool_and().literals()) { 1810 const Literal b = get_literal(ref); 1811 solver->AddUnitClause(b); 1812 } 1813 } else { 1814 // a => b 1815 const Literal not_a = 1816 get_literal(ct.enforcement_literal(0)).Negated(); 1817 for (const int ref : ct.bool_and().literals()) { 1818 const Literal b = get_literal(ref); 1819 solver->AddProblemClause({not_a, b}); 1820 } 1821 } 1822 break; 1823 } 1824 case ConstraintProto::ConstraintCase::kBoolOr: 1825 temp.clear(); 1826 for (const int ref : ct.bool_or().literals()) { 1827 temp.push_back(get_literal(ref)); 1828 } 1829 for (const int ref : ct.enforcement_literal()) { 1830 temp.push_back(get_literal(ref).Negated()); 1831 } 1832 solver->AddProblemClause(temp); 1833 break; 1834 default: 1835 LOG(FATAL) << "Not supported"; 1836 } 1837 } 1838 1839 // Deal with fixed variables. 1840 for (int ref = 0; ref < num_variables; ++ref) { 1841 const Domain domain = ReadDomainFromProto(model_proto.variables(ref)); 1842 if (domain.Min() == domain.Max()) { 1843 const Literal ref_literal = 1844 domain.Min() == 0 ? get_literal(ref).Negated() : get_literal(ref); 1845 solver->AddUnitClause(ref_literal); 1846 } 1847 } 1848 1849 SatSolver::Status status; 1850 CpSolverResponse response; 1851 if (parameters.cp_model_presolve()) { 1852 std::vector<bool> solution; 1853 status = SolveWithPresolve(&solver, model->GetOrCreate<TimeLimit>(), 1854 &solution, drat_proof_handler.get(), logger); 1855 if (status == SatSolver::FEASIBLE) { 1856 response.clear_solution(); 1857 for (int ref = 0; ref < num_variables; ++ref) { 1858 response.add_solution(solution[ref]); 1859 } 1860 } 1861 } else { 1862 status = solver->SolveWithTimeLimit(model->GetOrCreate<TimeLimit>()); 1863 if (status == SatSolver::FEASIBLE) { 1864 response.clear_solution(); 1865 for (int ref = 0; ref < num_variables; ++ref) { 1866 response.add_solution( 1867 solver->Assignment().LiteralIsTrue(get_literal(ref)) ? 1 : 0); 1868 } 1869 } 1870 } 1871 1872 // Tricky: the model local time limit is updated by the new functions, but 1873 // the old ones update time_limit directly. 1874 model->GetOrCreate<TimeLimit>()->AdvanceDeterministicTime( 1875 solver->model()->GetOrCreate<TimeLimit>()->GetElapsedDeterministicTime()); 1876 1877 switch (status) { 1878 case SatSolver::LIMIT_REACHED: { 1879 response.set_status(CpSolverStatus::UNKNOWN); 1880 break; 1881 } 1882 case SatSolver::FEASIBLE: { 1883 CHECK(SolutionIsFeasible( 1884 model_proto, std::vector<int64_t>(response.solution().begin(), 1885 response.solution().end()))); 1886 response.set_status(CpSolverStatus::OPTIMAL); 1887 break; 1888 } 1889 case SatSolver::INFEASIBLE: { 1890 response.set_status(CpSolverStatus::INFEASIBLE); 1891 break; 1892 } 1893 default: 1894 LOG(FATAL) << "Unexpected SatSolver::Status " << status; 1895 } 1896 response.set_num_booleans(solver->NumVariables()); 1897 response.set_num_branches(solver->num_branches()); 1898 response.set_num_conflicts(solver->num_failures()); 1899 response.set_num_binary_propagations(solver->num_propagations()); 1900 response.set_num_integer_propagations(0); 1901 response.set_wall_time(wall_timer->Get()); 1902 response.set_deterministic_time( 1903 model->Get<TimeLimit>()->GetElapsedDeterministicTime()); 1904 1905 if (status == SatSolver::INFEASIBLE && drat_proof_handler != nullptr) { 1906 WallTimer drat_timer; 1907 drat_timer.Start(); 1908 DratChecker::Status drat_status = drat_proof_handler->Check( 1909 absl::GetFlag(FLAGS_max_drat_time_in_seconds)); 1910 switch (drat_status) { 1911 case DratChecker::UNKNOWN: 1912 LOG(INFO) << "DRAT status: UNKNOWN"; 1913 break; 1914 case DratChecker::VALID: 1915 LOG(INFO) << "DRAT status: VALID"; 1916 break; 1917 case DratChecker::INVALID: 1918 LOG(ERROR) << "DRAT status: INVALID"; 1919 break; 1920 default: 1921 // Should not happen. 1922 break; 1923 } 1924 LOG(INFO) << "DRAT wall time: " << drat_timer.Get(); 1925 } else if (drat_proof_handler != nullptr) { 1926 // Always log a DRAT status to make it easier to extract it from a multirun 1927 // result with awk. 1928 LOG(INFO) << "DRAT status: NA"; 1929 LOG(INFO) << "DRAT wall time: NA"; 1930 LOG(INFO) << "DRAT user time: NA"; 1931 } 1932 return response; 1933 } 1934 1935 #if !defined(__PORTABLE_PLATFORM__) 1936 1937 // Small wrapper to simplify the constructions of the two SubSolver below. 1938 struct SharedClasses { 1939 CpModelProto const* model_proto; 1940 WallTimer* wall_timer; 1941 ModelSharedTimeLimit* time_limit; 1942 SharedBoundsManager* bounds; 1943 SharedResponseManager* response; 1944 SharedRelaxationSolutionRepository* relaxation_solutions; 1945 SharedLPSolutionRepository* lp_solutions; 1946 SharedIncompleteSolutionManager* incomplete_solutions; 1947 Model* global_model; 1948 1949 bool SearchIsDone() { 1950 if (response->ProblemIsSolved()) return true; 1951 if (time_limit->LimitReached()) return true; 1952 return false; 1953 } 1954 }; 1955 1956 // Encapsulate a full CP-SAT solve without presolve in the SubSolver API. 1957 class FullProblemSolver : public SubSolver { 1958 public: 1959 FullProblemSolver(const std::string& name, 1960 const SatParameters& local_parameters, bool split_in_chunks, 1961 SharedClasses* shared) 1962 : SubSolver(name), 1963 shared_(shared), 1964 split_in_chunks_(split_in_chunks), 1965 local_model_(absl::make_unique<Model>(name)) { 1966 // Setup the local model parameters and time limit. 1967 *(local_model_->GetOrCreate<SatParameters>()) = local_parameters; 1968 shared_->time_limit->UpdateLocalLimit( 1969 local_model_->GetOrCreate<TimeLimit>()); 1970 1971 if (shared->response != nullptr) { 1972 local_model_->Register<SharedResponseManager>(shared->response); 1973 } 1974 1975 if (shared->relaxation_solutions != nullptr) { 1976 local_model_->Register<SharedRelaxationSolutionRepository>( 1977 shared->relaxation_solutions); 1978 } 1979 1980 if (shared->lp_solutions != nullptr) { 1981 local_model_->Register<SharedLPSolutionRepository>(shared->lp_solutions); 1982 } 1983 1984 if (shared->incomplete_solutions != nullptr) { 1985 local_model_->Register<SharedIncompleteSolutionManager>( 1986 shared->incomplete_solutions); 1987 } 1988 } 1989 1990 bool TaskIsAvailable() override { 1991 if (shared_->SearchIsDone()) return false; 1992 1993 absl::MutexLock mutex_lock(&mutex_); 1994 return previous_task_is_completed_; 1995 } 1996 1997 std::function<void()> GenerateTask(int64_t task_id) override { 1998 { 1999 absl::MutexLock mutex_lock(&mutex_); 2000 previous_task_is_completed_ = false; 2001 } 2002 return [this]() { 2003 if (solving_first_chunk_) { 2004 LoadCpModel(*shared_->model_proto, local_model_.get()); 2005 2006 // Level zero variable bounds sharing. It is important to register 2007 // that after the probing that takes place in LoadCpModel() otherwise 2008 // we will have a mutex contention issue when all the thread probes 2009 // at the same time. 2010 if (shared_->bounds != nullptr) { 2011 RegisterVariableBoundsLevelZeroExport( 2012 *shared_->model_proto, shared_->bounds, local_model_.get()); 2013 RegisterVariableBoundsLevelZeroImport( 2014 *shared_->model_proto, shared_->bounds, local_model_.get()); 2015 } 2016 2017 if (local_model_->GetOrCreate<SatParameters>()->repair_hint()) { 2018 MinimizeL1DistanceWithHint(*shared_->model_proto, local_model_.get()); 2019 } else { 2020 QuickSolveWithHint(*shared_->model_proto, local_model_.get()); 2021 } 2022 2023 // No need for mutex since we only run one task at the time. 2024 solving_first_chunk_ = false; 2025 2026 if (split_in_chunks_) { 2027 // Abort first chunk and allow to schedule the next. 2028 absl::MutexLock mutex_lock(&mutex_); 2029 previous_task_is_completed_ = true; 2030 return; 2031 } 2032 } 2033 2034 auto* time_limit = local_model_->GetOrCreate<TimeLimit>(); 2035 if (split_in_chunks_) { 2036 // Configure time limit for chunk solving. Note that we do not want 2037 // to do that for the hint search for now. 2038 auto* params = local_model_->GetOrCreate<SatParameters>(); 2039 params->set_max_deterministic_time(1); 2040 time_limit->ResetLimitFromParameters(*params); 2041 shared_->time_limit->UpdateLocalLimit(time_limit); 2042 } 2043 2044 const double saved_dtime = time_limit->GetElapsedDeterministicTime(); 2045 SolveLoadedCpModel(*shared_->model_proto, local_model_.get()); 2046 { 2047 absl::MutexLock mutex_lock(&mutex_); 2048 deterministic_time_since_last_synchronize_ += 2049 time_limit->GetElapsedDeterministicTime() - saved_dtime; 2050 } 2051 2052 // Abort if the problem is solved. 2053 if (shared_->SearchIsDone()) { 2054 shared_->time_limit->Stop(); 2055 return; 2056 } 2057 2058 // In this mode, we allow to generate more task. 2059 if (split_in_chunks_) { 2060 absl::MutexLock mutex_lock(&mutex_); 2061 previous_task_is_completed_ = true; 2062 return; 2063 } 2064 2065 // Once a solver is done clear its memory and do not wait for the 2066 // destruction of the SubSolver. This is important because the full solve 2067 // might not be done at all, for instance this might have been configured 2068 // with stop_after_first_solution. 2069 local_model_.reset(); 2070 }; 2071 } 2072 2073 // TODO(user): A few of the information sharing we do between threads does not 2074 // happen here (bound sharing, RINS neighborhood, objective). Fix that so we 2075 // can have a deterministic parallel mode. 2076 void Synchronize() override { 2077 absl::MutexLock mutex_lock(&mutex_); 2078 deterministic_time_ += deterministic_time_since_last_synchronize_; 2079 shared_->time_limit->AdvanceDeterministicTime( 2080 deterministic_time_since_last_synchronize_); 2081 deterministic_time_since_last_synchronize_ = 0.0; 2082 } 2083 2084 std::string StatisticsString() const override { 2085 // The local model may have been deleted at the end of GenerateTask. 2086 // Do not crash in this case. 2087 // TODO(user): Revisit this case. 2088 if (local_model_ == nullptr) return std::string(); 2089 2090 const auto& lps = 2091 *local_model_->GetOrCreate<LinearProgrammingConstraintCollection>(); 2092 std::string lp_stats; 2093 if (!lps.empty() && 2094 local_model_->GetOrCreate<SatParameters>()->linearization_level() >= 2095 2) { 2096 for (const auto* lp : lps) { 2097 const std::string raw_statistics = lp->Statistics(); 2098 const std::vector<absl::string_view> lines = 2099 absl::StrSplit(raw_statistics, '\n', absl::SkipEmpty()); 2100 for (const absl::string_view& line : lines) { 2101 absl::StrAppend(&lp_stats, " ", line, "\n"); 2102 } 2103 } 2104 } 2105 return lp_stats; 2106 } 2107 2108 private: 2109 SharedClasses* shared_; 2110 const bool split_in_chunks_; 2111 std::unique_ptr<Model> local_model_; 2112 2113 // The first chunk is special. It is the one in which we load the model and 2114 // try to follow the hint. 2115 bool solving_first_chunk_ = true; 2116 2117 absl::Mutex mutex_; 2118 double deterministic_time_since_last_synchronize_ ABSL_GUARDED_BY(mutex_) = 2119 0.0; 2120 bool previous_task_is_completed_ ABSL_GUARDED_BY(mutex_) = true; 2121 }; 2122 2123 class FeasibilityPumpSolver : public SubSolver { 2124 public: 2125 FeasibilityPumpSolver(const SatParameters& local_parameters, 2126 SharedClasses* shared) 2127 : SubSolver("feasibility_pump"), 2128 shared_(shared), 2129 local_model_(absl::make_unique<Model>(name_)) { 2130 // Setup the local model parameters and time limit. 2131 *(local_model_->GetOrCreate<SatParameters>()) = local_parameters; 2132 shared_->time_limit->UpdateLocalLimit( 2133 local_model_->GetOrCreate<TimeLimit>()); 2134 2135 if (shared->response != nullptr) { 2136 local_model_->Register<SharedResponseManager>(shared->response); 2137 } 2138 2139 if (shared->relaxation_solutions != nullptr) { 2140 local_model_->Register<SharedRelaxationSolutionRepository>( 2141 shared->relaxation_solutions); 2142 } 2143 2144 if (shared->lp_solutions != nullptr) { 2145 local_model_->Register<SharedLPSolutionRepository>(shared->lp_solutions); 2146 } 2147 2148 if (shared->incomplete_solutions != nullptr) { 2149 local_model_->Register<SharedIncompleteSolutionManager>( 2150 shared->incomplete_solutions); 2151 } 2152 2153 // Level zero variable bounds sharing. 2154 if (shared_->bounds != nullptr) { 2155 RegisterVariableBoundsLevelZeroImport( 2156 *shared_->model_proto, shared_->bounds, local_model_.get()); 2157 } 2158 } 2159 2160 bool TaskIsAvailable() override { 2161 if (shared_->SearchIsDone()) return false; 2162 absl::MutexLock mutex_lock(&mutex_); 2163 return previous_task_is_completed_; 2164 } 2165 2166 std::function<void()> GenerateTask(int64_t task_id) override { 2167 return [this]() { 2168 { 2169 absl::MutexLock mutex_lock(&mutex_); 2170 if (!previous_task_is_completed_) return; 2171 previous_task_is_completed_ = false; 2172 } 2173 { 2174 absl::MutexLock mutex_lock(&mutex_); 2175 if (solving_first_chunk_) { 2176 LoadFeasibilityPump(*shared_->model_proto, local_model_.get()); 2177 // No new task will be scheduled for this worker if there is no 2178 // linear relaxation. 2179 if (local_model_->Get<FeasibilityPump>() == nullptr) return; 2180 solving_first_chunk_ = false; 2181 // Abort first chunk and allow to schedule the next. 2182 previous_task_is_completed_ = true; 2183 return; 2184 } 2185 } 2186 2187 auto* time_limit = local_model_->GetOrCreate<TimeLimit>(); 2188 const double saved_dtime = time_limit->GetElapsedDeterministicTime(); 2189 auto* feasibility_pump = local_model_->Mutable<FeasibilityPump>(); 2190 if (!feasibility_pump->Solve()) { 2191 shared_->response->NotifyThatImprovingProblemIsInfeasible(name_); 2192 } 2193 2194 { 2195 absl::MutexLock mutex_lock(&mutex_); 2196 deterministic_time_since_last_synchronize_ += 2197 time_limit->GetElapsedDeterministicTime() - saved_dtime; 2198 } 2199 2200 // Abort if the problem is solved. 2201 if (shared_->SearchIsDone()) { 2202 shared_->time_limit->Stop(); 2203 return; 2204 } 2205 2206 absl::MutexLock mutex_lock(&mutex_); 2207 previous_task_is_completed_ = true; 2208 }; 2209 } 2210 2211 void Synchronize() override { 2212 absl::MutexLock mutex_lock(&mutex_); 2213 deterministic_time_ += deterministic_time_since_last_synchronize_; 2214 shared_->time_limit->AdvanceDeterministicTime( 2215 deterministic_time_since_last_synchronize_); 2216 deterministic_time_since_last_synchronize_ = 0.0; 2217 } 2218 2219 // TODO(user, fdid): Display feasibility pump statistics. 2220 2221 private: 2222 SharedClasses* shared_; 2223 std::unique_ptr<Model> local_model_; 2224 2225 absl::Mutex mutex_; 2226 2227 // The first chunk is special. It is the one in which we load the linear 2228 // constraints. 2229 bool solving_first_chunk_ ABSL_GUARDED_BY(mutex_) = true; 2230 2231 double deterministic_time_since_last_synchronize_ ABSL_GUARDED_BY(mutex_) = 2232 0.0; 2233 bool previous_task_is_completed_ ABSL_GUARDED_BY(mutex_) = true; 2234 }; 2235 2236 // A Subsolver that generate LNS solve from a given neighborhood. 2237 class LnsSolver : public SubSolver { 2238 public: 2239 LnsSolver(std::unique_ptr<NeighborhoodGenerator> generator, 2240 const SatParameters& parameters, 2241 NeighborhoodGeneratorHelper* helper, SharedClasses* shared) 2242 : SubSolver(generator->name()), 2243 generator_(std::move(generator)), 2244 helper_(helper), 2245 parameters_(parameters), 2246 shared_(shared) {} 2247 2248 bool TaskIsAvailable() override { 2249 if (shared_->SearchIsDone()) return false; 2250 return generator_->ReadyToGenerate(); 2251 } 2252 2253 std::function<void()> GenerateTask(int64_t task_id) override { 2254 return [task_id, this]() { 2255 if (shared_->SearchIsDone()) return; 2256 2257 // Create a random number generator whose seed depends both on the task_id 2258 // and on the parameters_.random_seed() so that changing the later will 2259 // change the LNS behavior. 2260 const int32_t low = static_cast<int32_t>(task_id); 2261 const int32_t high = static_cast<int32_t>(task_id >> 32); 2262 std::seed_seq seed{low, high, parameters_.random_seed()}; 2263 random_engine_t random(seed); 2264 2265 NeighborhoodGenerator::SolveData data; 2266 data.difficulty = generator_->difficulty(); 2267 data.deterministic_limit = generator_->deterministic_limit(); 2268 2269 // Choose a base solution for this neighborhood. 2270 CpSolverResponse base_response; 2271 { 2272 const SharedSolutionRepository<int64_t>& repo = 2273 shared_->response->SolutionsRepository(); 2274 if (repo.NumSolutions() > 0) { 2275 base_response.set_status(CpSolverStatus::FEASIBLE); 2276 const SharedSolutionRepository<int64_t>::Solution solution = 2277 repo.GetRandomBiasedSolution(random); 2278 for (const int64_t value : solution.variable_values) { 2279 base_response.add_solution(value); 2280 } 2281 2282 // Note: We assume that the solution rank is the solution internal 2283 // objective. 2284 data.initial_best_objective = repo.GetSolution(0).rank; 2285 data.base_objective = solution.rank; 2286 } else { 2287 base_response.set_status(CpSolverStatus::UNKNOWN); 2288 2289 // If we do not have a solution, we use the current objective upper 2290 // bound so that our code that compute an "objective" improvement 2291 // works. 2292 // 2293 // TODO(user): this is non-deterministic. Fix. 2294 data.initial_best_objective = 2295 shared_->response->GetInnerObjectiveUpperBound(); 2296 data.base_objective = data.initial_best_objective; 2297 } 2298 } 2299 2300 Neighborhood neighborhood = 2301 generator_->Generate(base_response, data.difficulty, random); 2302 2303 if (!neighborhood.is_generated) return; 2304 2305 data.neighborhood_id = neighborhood.id; 2306 2307 const int64_t num_calls = std::max(int64_t{1}, generator_->num_calls()); 2308 const double fully_solved_proportion = 2309 static_cast<double>(generator_->num_fully_solved_calls()) / 2310 static_cast<double>(num_calls); 2311 std::string source_info = name(); 2312 if (!neighborhood.source_info.empty()) { 2313 absl::StrAppend(&source_info, "_", neighborhood.source_info); 2314 } 2315 const std::string solution_info = absl::StrFormat( 2316 "%s(d=%0.2f s=%i t=%0.2f p=%0.2f)", source_info, data.difficulty, 2317 task_id, data.deterministic_limit, fully_solved_proportion); 2318 2319 SatParameters local_params(parameters_); 2320 local_params.set_max_deterministic_time(data.deterministic_limit); 2321 local_params.set_stop_after_first_solution(false); 2322 local_params.set_log_search_progress(false); 2323 local_params.set_cp_model_probing_level(0); 2324 local_params.set_symmetry_level(0); 2325 local_params.set_solution_pool_size(0); 2326 2327 Model local_model(solution_info); 2328 *(local_model.GetOrCreate<SatParameters>()) = local_params; 2329 TimeLimit* local_time_limit = local_model.GetOrCreate<TimeLimit>(); 2330 local_time_limit->ResetLimitFromParameters(local_params); 2331 shared_->time_limit->UpdateLocalLimit(local_time_limit); 2332 2333 // Presolve and solve the LNS fragment. 2334 CpModelProto lns_fragment; 2335 CpModelProto mapping_proto; 2336 auto context = absl::make_unique<PresolveContext>( 2337 &local_model, &lns_fragment, &mapping_proto); 2338 2339 *lns_fragment.mutable_variables() = neighborhood.delta.variables(); 2340 { 2341 ModelCopy copier(context.get()); 2342 2343 // Copy and simplify the constraints from the initial model. 2344 if (!copier.ImportAndSimplifyConstraints( 2345 helper_->ModelProto(), neighborhood.constraints_to_ignore)) { 2346 return; 2347 } 2348 2349 // Copy and simplify the constraints from the delta model. 2350 if (!neighborhood.delta.constraints().empty() && 2351 !copier.ImportAndSimplifyConstraints(neighborhood.delta, {})) { 2352 return; 2353 } 2354 } 2355 2356 // Copy the rest of the model and overwrite the name. 2357 CopyEverythingExceptVariablesAndConstraintsFieldsIntoContext( 2358 helper_->ModelProto(), context.get()); 2359 lns_fragment.set_name(absl::StrCat("lns_", task_id)); 2360 2361 // Overwrite solution hinting. 2362 if (neighborhood.delta.has_solution_hint()) { 2363 *lns_fragment.mutable_solution_hint() = 2364 neighborhood.delta.solution_hint(); 2365 } 2366 2367 CpModelProto debug_copy; 2368 if (absl::GetFlag(FLAGS_cp_model_dump_problematic_lns)) { 2369 // We need to make a copy because the presolve is destructive. 2370 // It is why we do not do that by default. 2371 debug_copy = lns_fragment; 2372 } 2373 2374 #if !defined(__PORTABLE_PLATFORM__) 2375 #endif // __PORTABLE_PLATFORM__ 2376 2377 if (absl::GetFlag(FLAGS_cp_model_dump_lns)) { 2378 // TODO(user): export the delta too if needed. 2379 const std::string lns_name = 2380 absl::StrCat(absl::GetFlag(FLAGS_cp_model_dump_prefix), 2381 lns_fragment.name(), ".pb.txt"); 2382 LOG(INFO) << "Dumping LNS model to '" << lns_name << "'."; 2383 CHECK_OK(file::SetTextProto(lns_name, lns_fragment, file::Defaults())); 2384 } 2385 2386 std::vector<int> postsolve_mapping; 2387 const CpSolverStatus presolve_status = 2388 PresolveCpModel(context.get(), &postsolve_mapping); 2389 2390 // Release the context. 2391 context.reset(nullptr); 2392 neighborhood.delta.Clear(); 2393 2394 // TODO(user): Depending on the problem, we should probably use the 2395 // parameters that work bests (core, linearization_level, etc...) or 2396 // maybe we can just randomize them like for the base solution used. 2397 auto* local_response_manager = 2398 local_model.GetOrCreate<SharedResponseManager>(); 2399 local_response_manager->InitializeObjective(lns_fragment); 2400 2401 if (presolve_status == CpSolverStatus::UNKNOWN) { 2402 LoadCpModel(lns_fragment, &local_model); 2403 QuickSolveWithHint(lns_fragment, &local_model); 2404 SolveLoadedCpModel(lns_fragment, &local_model); 2405 } else { 2406 local_response_manager->MutableResponse()->set_status(presolve_status); 2407 } 2408 CpSolverResponse local_response = local_response_manager->GetResponse(); 2409 2410 // TODO(user): we actually do not need to postsolve if the solution is 2411 // not going to be used... 2412 if (local_params.cp_model_presolve() && 2413 (local_response.status() == CpSolverStatus::OPTIMAL || 2414 local_response.status() == CpSolverStatus::FEASIBLE)) { 2415 std::vector<int64_t> solution(local_response.solution().begin(), 2416 local_response.solution().end()); 2417 PostsolveResponseWrapper(local_params, 2418 helper_->ModelProto().variables_size(), 2419 mapping_proto, postsolve_mapping, &solution); 2420 local_response.mutable_solution()->Assign(solution.begin(), 2421 solution.end()); 2422 } 2423 2424 data.status = local_response.status(); 2425 data.deterministic_time = local_time_limit->GetElapsedDeterministicTime(); 2426 2427 bool new_solution = false; 2428 bool display_lns_info = false; 2429 if (generator_->IsRelaxationGenerator()) { 2430 bool has_feasible_solution = false; 2431 if (local_response.status() == CpSolverStatus::OPTIMAL || 2432 local_response.status() == CpSolverStatus::FEASIBLE) { 2433 has_feasible_solution = true; 2434 } 2435 2436 if (local_response.status() == CpSolverStatus::INFEASIBLE) { 2437 shared_->response->NotifyThatImprovingProblemIsInfeasible( 2438 local_response.solution_info()); 2439 } 2440 2441 if (shared_->model_proto->has_objective()) { 2442 // TODO(user): This is not deterministic since it is updated without 2443 // synchronization! So we shouldn't base the LNS score out of that. 2444 const IntegerValue current_obj_lb = 2445 shared_->response->GetInnerObjectiveLowerBound(); 2446 2447 const IntegerValue local_obj_lb = 2448 local_response_manager->GetInnerObjectiveLowerBound(); 2449 2450 const double scaled_local_obj_bound = ScaleObjectiveValue( 2451 lns_fragment.objective(), local_obj_lb.value()); 2452 2453 // Update the bound. 2454 const IntegerValue new_inner_obj_lb = IntegerValue( 2455 std::ceil(UnscaleObjectiveValue(shared_->model_proto->objective(), 2456 scaled_local_obj_bound) - 2457 1e-6)); 2458 data.new_objective_bound = new_inner_obj_lb; 2459 data.initial_best_objective_bound = current_obj_lb; 2460 if (new_inner_obj_lb > current_obj_lb) { 2461 shared_->response->UpdateInnerObjectiveBounds( 2462 solution_info, new_inner_obj_lb, kMaxIntegerValue); 2463 } 2464 } 2465 2466 // If we have a solution of the relaxed problem, we check if it is also 2467 // a valid solution of the non-relaxed one. 2468 if (has_feasible_solution) { 2469 if (SolutionIsFeasible( 2470 *shared_->model_proto, 2471 std::vector<int64_t>(local_response.solution().begin(), 2472 local_response.solution().end()))) { 2473 shared_->response->NewSolution(local_response, 2474 /*model=*/nullptr); 2475 2476 // Mark the solution optimal if the relaxation status is optimal. 2477 if (local_response.status() == CpSolverStatus::OPTIMAL) { 2478 shared_->response->NotifyThatImprovingProblemIsInfeasible( 2479 local_response.solution_info()); 2480 shared_->time_limit->Stop(); 2481 } 2482 } 2483 shared_->relaxation_solutions->NewRelaxationSolution(local_response); 2484 } 2485 } else { 2486 const std::vector<int64_t> solution(local_response.solution().begin(), 2487 local_response.solution().end()); 2488 2489 if (!solution.empty()) { 2490 // A solution that does not pass our validator indicates a bug. We 2491 // abort and dump the problematic model to facilitate debugging. 2492 // 2493 // TODO(user): In a production environment, we should probably just 2494 // ignore this fragment and continue. 2495 const bool feasible = 2496 SolutionIsFeasible(*shared_->model_proto, solution); 2497 if (!feasible) { 2498 if (absl::GetFlag(FLAGS_cp_model_dump_problematic_lns)) { 2499 const std::string name = 2500 absl::StrCat(absl::GetFlag(FLAGS_cp_model_dump_prefix), 2501 debug_copy.name(), ".pb.txt"); 2502 LOG(INFO) << "Dumping problematic LNS model to '" << name << "'."; 2503 CHECK_OK(file::SetTextProto(name, debug_copy, file::Defaults())); 2504 } 2505 LOG(FATAL) << "Infeasible LNS solution! " << solution_info 2506 << " solved with params " 2507 << local_params.ShortDebugString(); 2508 } 2509 } 2510 2511 // Special case if we solved a part of the full problem! 2512 // 2513 // TODO(user): This do not seem to work if they are symmetries loaded 2514 // into SAT. For now we just disable this if there is any symmetry. 2515 // See for instance spot5_1401.fzn. Be smarter about that: 2516 // - If there are connected compo in the inital model, we should only 2517 // compute generator that do not cross component? or are component 2518 // interchange useful? probably not. 2519 // - It should be fine if all our generator are fully or not at 2520 // all included in the variable we are fixing. So we can relax the 2521 // test here. Try on z26.mps or spot5_1401.fzn. 2522 if (local_response.status() == CpSolverStatus::OPTIMAL && 2523 !shared_->model_proto->has_symmetry() && !solution.empty() && 2524 neighborhood.is_simple && 2525 !neighborhood.variables_that_can_be_fixed_to_local_optimum 2526 .empty()) { 2527 display_lns_info = true; 2528 shared_->bounds->FixVariablesFromPartialSolution( 2529 solution, 2530 neighborhood.variables_that_can_be_fixed_to_local_optimum); 2531 } 2532 2533 // Finish to fill the SolveData now that the local solve is done. 2534 data.new_objective = data.base_objective; 2535 if (local_response.status() == CpSolverStatus::OPTIMAL || 2536 local_response.status() == CpSolverStatus::FEASIBLE) { 2537 data.new_objective = IntegerValue(ComputeInnerObjective( 2538 shared_->model_proto->objective(), local_response)); 2539 } 2540 2541 // Report any feasible solution we have. Optimization: We don't do that 2542 // if we just recovered the base solution. 2543 if (local_response.status() == CpSolverStatus::OPTIMAL || 2544 local_response.status() == CpSolverStatus::FEASIBLE) { 2545 const std::vector<int64_t> base_solution( 2546 base_response.solution().begin(), base_response.solution().end()); 2547 if (solution != base_solution) { 2548 new_solution = true; 2549 shared_->response->NewSolution(local_response, /*model=*/nullptr); 2550 } 2551 } 2552 if (!neighborhood.is_reduced && 2553 (local_response.status() == CpSolverStatus::OPTIMAL || 2554 local_response.status() == CpSolverStatus::INFEASIBLE)) { 2555 shared_->response->NotifyThatImprovingProblemIsInfeasible( 2556 local_response.solution_info()); 2557 shared_->time_limit->Stop(); 2558 } 2559 } 2560 2561 generator_->AddSolveData(data); 2562 2563 if (VLOG_IS_ON(2) || display_lns_info) { 2564 auto* logger = shared_->global_model->GetOrCreate<SolverLogger>(); 2565 std::string s = absl::StrCat(" LNS ", name(), ":"); 2566 if (new_solution) { 2567 const double base_obj = ScaleObjectiveValue( 2568 shared_->model_proto->objective(), 2569 ComputeInnerObjective(shared_->model_proto->objective(), 2570 base_response)); 2571 const double new_obj = ScaleObjectiveValue( 2572 shared_->model_proto->objective(), 2573 ComputeInnerObjective(shared_->model_proto->objective(), 2574 local_response)); 2575 absl::StrAppend(&s, " [new_sol:", base_obj, " -> ", new_obj, "]"); 2576 } 2577 if (neighborhood.is_simple) { 2578 absl::StrAppend( 2579 &s, " [", "relaxed:", neighborhood.num_relaxed_variables, 2580 " in_obj:", neighborhood.num_relaxed_variables_in_objective, 2581 " compo:", 2582 neighborhood.variables_that_can_be_fixed_to_local_optimum.size(), 2583 "]"); 2584 } 2585 SOLVER_LOG(logger, s, " [d:", data.difficulty, ", id:", task_id, 2586 ", dtime:", data.deterministic_time, "/", 2587 data.deterministic_limit, 2588 ", status:", ProtoEnumToString<CpSolverStatus>(data.status), 2589 ", #calls:", generator_->num_calls(), 2590 ", p:", fully_solved_proportion, "]"); 2591 } 2592 }; 2593 } 2594 2595 void Synchronize() override { 2596 generator_->Synchronize(); 2597 const double old = deterministic_time_; 2598 deterministic_time_ = generator_->deterministic_time(); 2599 shared_->time_limit->AdvanceDeterministicTime(deterministic_time_ - old); 2600 } 2601 2602 // TODO(user): Display LNS success rate. 2603 2604 private: 2605 std::unique_ptr<NeighborhoodGenerator> generator_; 2606 NeighborhoodGeneratorHelper* helper_; 2607 const SatParameters parameters_; 2608 SharedClasses* shared_; 2609 }; 2610 2611 void SolveCpModelParallel(const CpModelProto& model_proto, 2612 Model* global_model) { 2613 const SatParameters& parameters = *global_model->GetOrCreate<SatParameters>(); 2614 const int num_search_workers = parameters.num_search_workers(); 2615 CHECK(!parameters.enumerate_all_solutions()) 2616 << "Enumerating all solutions in parallel is not supported."; 2617 2618 std::unique_ptr<SharedBoundsManager> shared_bounds_manager; 2619 if (parameters.share_level_zero_bounds()) { 2620 shared_bounds_manager = absl::make_unique<SharedBoundsManager>(model_proto); 2621 } 2622 2623 std::unique_ptr<SharedRelaxationSolutionRepository> 2624 shared_relaxation_solutions; 2625 if (parameters.use_relaxation_lns()) { 2626 shared_relaxation_solutions = 2627 absl::make_unique<SharedRelaxationSolutionRepository>( 2628 /*num_solutions_to_keep=*/10); 2629 global_model->Register<SharedRelaxationSolutionRepository>( 2630 shared_relaxation_solutions.get()); 2631 } 2632 2633 auto shared_lp_solutions = absl::make_unique<SharedLPSolutionRepository>( 2634 /*num_solutions_to_keep=*/10); 2635 global_model->Register<SharedLPSolutionRepository>(shared_lp_solutions.get()); 2636 2637 // We currently only use the feasiblity pump if it is enabled and some other 2638 // parameters are not on. 2639 std::unique_ptr<SharedIncompleteSolutionManager> shared_incomplete_solutions; 2640 const bool use_feasibility_pump = parameters.use_feasibility_pump() && 2641 parameters.linearization_level() > 0 && 2642 !parameters.use_lns_only() && 2643 !parameters.interleave_search(); 2644 if (use_feasibility_pump) { 2645 shared_incomplete_solutions = 2646 absl::make_unique<SharedIncompleteSolutionManager>(); 2647 global_model->Register<SharedIncompleteSolutionManager>( 2648 shared_incomplete_solutions.get()); 2649 } 2650 2651 SharedClasses shared; 2652 shared.model_proto = &model_proto; 2653 shared.wall_timer = global_model->GetOrCreate<WallTimer>(); 2654 shared.time_limit = global_model->GetOrCreate<ModelSharedTimeLimit>(); 2655 shared.bounds = shared_bounds_manager.get(); 2656 shared.response = global_model->GetOrCreate<SharedResponseManager>(); 2657 shared.relaxation_solutions = shared_relaxation_solutions.get(); 2658 shared.lp_solutions = shared_lp_solutions.get(); 2659 shared.incomplete_solutions = shared_incomplete_solutions.get(); 2660 shared.global_model = global_model; 2661 2662 // The list of all the SubSolver that will be used in this parallel search. 2663 std::vector<std::unique_ptr<SubSolver>> subsolvers; 2664 2665 // Add a synchronization point for the shared classes. 2666 subsolvers.push_back(absl::make_unique<SynchronizationPoint>([&shared]() { 2667 shared.response->Synchronize(); 2668 shared.response->MutableSolutionsRepository()->Synchronize(); 2669 if (shared.bounds != nullptr) { 2670 shared.bounds->Synchronize(); 2671 } 2672 if (shared.relaxation_solutions != nullptr) { 2673 shared.relaxation_solutions->Synchronize(); 2674 } 2675 if (shared.lp_solutions != nullptr) { 2676 shared.lp_solutions->Synchronize(); 2677 } 2678 })); 2679 2680 if (parameters.use_lns_only()) { 2681 // Register something to find a first solution. Note that this is mainly 2682 // used for experimentation, and using no LP ususally result in a faster 2683 // first solution. 2684 SatParameters local_params = parameters; 2685 local_params.set_stop_after_first_solution(true); 2686 local_params.set_linearization_level(0); 2687 subsolvers.push_back(absl::make_unique<FullProblemSolver>( 2688 "first_solution", local_params, 2689 /*split_in_chunks=*/false, &shared)); 2690 } else { 2691 for (const SatParameters& local_params : GetDiverseSetOfParameters( 2692 parameters, model_proto, num_search_workers)) { 2693 // TODO(user): This is currently not supported here. 2694 if (parameters.optimize_with_max_hs()) continue; 2695 2696 subsolvers.push_back(absl::make_unique<FullProblemSolver>( 2697 local_params.name(), local_params, 2698 /*split_in_chunks=*/parameters.interleave_search(), &shared)); 2699 } 2700 } 2701 2702 // Add FeasibilityPumpSolver if enabled. 2703 if (use_feasibility_pump) { 2704 subsolvers.push_back( 2705 absl::make_unique<FeasibilityPumpSolver>(parameters, &shared)); 2706 } 2707 2708 // Add LNS SubSolver(s). 2709 2710 // Add the NeighborhoodGeneratorHelper as a special subsolver so that its 2711 // Synchronize() is called before any LNS neighborhood solvers. 2712 auto unique_helper = absl::make_unique<NeighborhoodGeneratorHelper>( 2713 &model_proto, ¶meters, shared.response, shared.time_limit, 2714 shared.bounds); 2715 NeighborhoodGeneratorHelper* helper = unique_helper.get(); 2716 subsolvers.push_back(std::move(unique_helper)); 2717 2718 // By default we use the user provided parameters. 2719 std::vector<SatParameters> lns_params = {parameters}; 2720 lns_params.back().set_name("default"); 2721 if (parameters.diversify_lns_params()) { 2722 std::vector<SatParameters> lns_params = 2723 GetDiverseSetOfParameters(parameters, model_proto, 6); 2724 } 2725 for (const SatParameters& local_params : lns_params) { 2726 // Only register following LNS SubSolver if there is an objective. 2727 if (model_proto.has_objective()) { 2728 // Enqueue all the possible LNS neighborhood subsolvers. 2729 // Each will have their own metrics. 2730 subsolvers.push_back(absl::make_unique<LnsSolver>( 2731 absl::make_unique<RelaxRandomVariablesGenerator>( 2732 helper, absl::StrCat("rnd_var_lns_", local_params.name())), 2733 local_params, helper, &shared)); 2734 subsolvers.push_back(absl::make_unique<LnsSolver>( 2735 absl::make_unique<RelaxRandomConstraintsGenerator>( 2736 helper, absl::StrCat("rnd_cst_lns_", local_params.name())), 2737 local_params, helper, &shared)); 2738 subsolvers.push_back(absl::make_unique<LnsSolver>( 2739 absl::make_unique<VariableGraphNeighborhoodGenerator>( 2740 helper, absl::StrCat("graph_var_lns_", local_params.name())), 2741 local_params, helper, &shared)); 2742 subsolvers.push_back(absl::make_unique<LnsSolver>( 2743 absl::make_unique<ConstraintGraphNeighborhoodGenerator>( 2744 helper, absl::StrCat("graph_cst_lns_", local_params.name())), 2745 local_params, helper, &shared)); 2746 2747 // TODO(user): If we have a model with scheduling + routing. We create 2748 // a lot of LNS generators. Investigate if we can reduce this number. 2749 if (!helper->TypeToConstraints(ConstraintProto::kNoOverlap).empty() || 2750 !helper->TypeToConstraints(ConstraintProto::kNoOverlap2D).empty() || 2751 !helper->TypeToConstraints(ConstraintProto::kCumulative).empty()) { 2752 subsolvers.push_back(absl::make_unique<LnsSolver>( 2753 absl::make_unique<SchedulingTimeWindowNeighborhoodGenerator>( 2754 helper, absl::StrCat("scheduling_time_window_lns_", 2755 local_params.name())), 2756 local_params, helper, &shared)); 2757 subsolvers.push_back(absl::make_unique<LnsSolver>( 2758 absl::make_unique<SchedulingNeighborhoodGenerator>( 2759 helper, 2760 absl::StrCat("scheduling_random_lns_", local_params.name())), 2761 local_params, helper, &shared)); 2762 } 2763 } 2764 2765 const int num_circuit = 2766 helper->TypeToConstraints(ConstraintProto::kCircuit).size(); 2767 const int num_routes = 2768 helper->TypeToConstraints(ConstraintProto::kRoutes).size(); 2769 if (num_circuit + num_routes > 0) { 2770 subsolvers.push_back(absl::make_unique<LnsSolver>( 2771 absl::make_unique<RoutingRandomNeighborhoodGenerator>( 2772 helper, absl::StrCat("routing_random_lns_", local_params.name())), 2773 local_params, helper, &shared)); 2774 2775 subsolvers.push_back(absl::make_unique<LnsSolver>( 2776 absl::make_unique<RoutingPathNeighborhoodGenerator>( 2777 helper, absl::StrCat("routing_path_lns_", local_params.name())), 2778 local_params, helper, &shared)); 2779 } 2780 if (num_routes > 0 || num_circuit > 1) { 2781 subsolvers.push_back(absl::make_unique<LnsSolver>( 2782 absl::make_unique<RoutingFullPathNeighborhoodGenerator>( 2783 helper, 2784 absl::StrCat("routing_full_path_lns_", local_params.name())), 2785 local_params, helper, &shared)); 2786 } 2787 2788 // TODO(user): for now this is not deterministic so we disable it on 2789 // interleave search. Fix. 2790 if (parameters.use_rins_lns() && !parameters.interleave_search()) { 2791 // Note that we always create the SharedLPSolutionRepository. This meets 2792 // the requirement of having at least one of 2793 // SharedRelaxationSolutionRepository or SharedLPSolutionRepository to 2794 // create RINS/RENS lns generators. 2795 2796 // RINS. 2797 subsolvers.push_back(absl::make_unique<LnsSolver>( 2798 absl::make_unique<RelaxationInducedNeighborhoodGenerator>( 2799 helper, shared.response, shared.relaxation_solutions, 2800 shared.lp_solutions, /*incomplete_solutions=*/nullptr, 2801 absl::StrCat("rins_lns_", local_params.name())), 2802 local_params, helper, &shared)); 2803 2804 // RENS. 2805 subsolvers.push_back(absl::make_unique<LnsSolver>( 2806 absl::make_unique<RelaxationInducedNeighborhoodGenerator>( 2807 helper, /*respons_manager=*/nullptr, shared.relaxation_solutions, 2808 shared.lp_solutions, shared.incomplete_solutions, 2809 absl::StrCat("rens_lns_", local_params.name())), 2810 local_params, helper, &shared)); 2811 } 2812 2813 if (parameters.use_relaxation_lns()) { 2814 subsolvers.push_back(absl::make_unique<LnsSolver>( 2815 absl::make_unique< 2816 ConsecutiveConstraintsRelaxationNeighborhoodGenerator>( 2817 helper, absl::StrCat("rnd_rel_lns_", local_params.name())), 2818 local_params, helper, &shared)); 2819 2820 subsolvers.push_back(absl::make_unique<LnsSolver>( 2821 absl::make_unique<WeightedRandomRelaxationNeighborhoodGenerator>( 2822 helper, absl::StrCat("wgt_rel_lns_", local_params.name())), 2823 local_params, helper, &shared)); 2824 } 2825 } 2826 2827 // Add a synchronization point for the gap integral that is executed last. 2828 // This way, after each batch, the proper deterministic time is updated and 2829 // then the function to integrate take the value of the new gap. 2830 subsolvers.push_back(absl::make_unique<SynchronizationPoint>( 2831 [&shared]() { shared.response->UpdateGapIntegral(); })); 2832 2833 // Log the name of all our SubSolvers. 2834 auto* logger = global_model->GetOrCreate<SolverLogger>(); 2835 if (logger->LoggingIsEnabled()) { 2836 std::vector<std::string> names; 2837 for (const auto& subsolver : subsolvers) { 2838 if (!subsolver->name().empty()) names.push_back(subsolver->name()); 2839 } 2840 SOLVER_LOG(logger, ""); 2841 SOLVER_LOG(logger, 2842 absl::StrFormat("Starting Search at %.2fs with %i " 2843 "workers and subsolvers: [ %s ]", 2844 shared.wall_timer->Get(), num_search_workers, 2845 absl::StrJoin(names, ", "))); 2846 } 2847 2848 // Launch the main search loop. 2849 if (parameters.interleave_search()) { 2850 DeterministicLoop(subsolvers, num_search_workers, 2851 parameters.interleave_batch_size()); 2852 } else { 2853 NonDeterministicLoop(subsolvers, num_search_workers); 2854 } 2855 2856 if (parameters.log_subsolver_statistics()) { 2857 SOLVER_LOG(logger, ""); 2858 SOLVER_LOG(logger, "Sub-solver search statistics:"); 2859 for (const auto& subsolver : subsolvers) { 2860 const std::string stats = subsolver->StatisticsString(); 2861 if (stats.empty()) continue; 2862 SOLVER_LOG(logger, absl::StrCat(" '", subsolver->name(), "':\n", stats)); 2863 } 2864 } 2865 } 2866 2867 #endif // __PORTABLE_PLATFORM__ 2868 2869 // If the option use_sat_inprocessing is true, then before postsolving a 2870 // solution, we need to make sure we add any new clause required for postsolving 2871 // to the mapping_model. 2872 void AddPostsolveClauses(const std::vector<int>& postsolve_mapping, 2873 Model* model, CpModelProto* mapping_proto) { 2874 auto* mapping = model->GetOrCreate<CpModelMapping>(); 2875 auto* postsolve = model->GetOrCreate<PostsolveClauses>(); 2876 for (const auto& clause : postsolve->clauses) { 2877 auto* ct = mapping_proto->add_constraints()->mutable_bool_or(); 2878 for (const Literal l : clause) { 2879 int var = mapping->GetProtoVariableFromBooleanVariable(l.Variable()); 2880 CHECK_NE(var, -1); 2881 var = postsolve_mapping[var]; 2882 ct->add_literals(l.IsPositive() ? var : NegatedRef(var)); 2883 } 2884 } 2885 postsolve->clauses.clear(); 2886 } 2887 2888 } // namespace 2889 2890 CpSolverResponse SolveCpModel(const CpModelProto& model_proto, Model* model) { 2891 auto* wall_timer = model->GetOrCreate<WallTimer>(); 2892 auto* user_timer = model->GetOrCreate<UserTimer>(); 2893 wall_timer->Start(); 2894 user_timer->Start(); 2895 2896 #if !defined(__PORTABLE_PLATFORM__) 2897 #endif // __PORTABLE_PLATFORM__ 2898 2899 #if !defined(__PORTABLE_PLATFORM__) 2900 // Dump initial model? 2901 if (absl::GetFlag(FLAGS_cp_model_dump_models)) { 2902 const std::string file = 2903 absl::StrCat(absl::GetFlag(FLAGS_cp_model_dump_prefix), "model.pb.txt"); 2904 LOG(INFO) << "Dumping cp model proto to '" << file << "'."; 2905 CHECK_OK(file::SetTextProto(file, model_proto, file::Defaults())); 2906 } 2907 #endif // __PORTABLE_PLATFORM__ 2908 2909 #if !defined(__PORTABLE_PLATFORM__) 2910 // Override parameters? 2911 if (!absl::GetFlag(FLAGS_cp_model_params).empty()) { 2912 SatParameters params = *model->GetOrCreate<SatParameters>(); 2913 SatParameters flag_params; 2914 CHECK(google::protobuf::TextFormat::ParseFromString( 2915 absl::GetFlag(FLAGS_cp_model_params), &flag_params)); 2916 params.MergeFrom(flag_params); 2917 *(model->GetOrCreate<SatParameters>()) = params; 2918 } 2919 #endif // __PORTABLE_PLATFORM__ 2920 2921 // Enable the logging component. 2922 const SatParameters& params = *model->GetOrCreate<SatParameters>(); 2923 SolverLogger* logger = model->GetOrCreate<SolverLogger>(); 2924 logger->EnableLogging(params.log_search_progress() || VLOG_IS_ON(1)); 2925 logger->SetLogToStdOut(params.log_to_stdout()); 2926 std::string log_string; 2927 if (params.log_to_response()) { 2928 logger->AddInfoLoggingCallback([&log_string](const std::string& message) { 2929 absl::StrAppend(&log_string, message, "\n"); 2930 }); 2931 } 2932 2933 auto* shared_response_manager = model->GetOrCreate<SharedResponseManager>(); 2934 shared_response_manager->set_dump_prefix( 2935 absl::GetFlag(FLAGS_cp_model_dump_prefix)); 2936 2937 #if !defined(__PORTABLE_PLATFORM__) 2938 // Note that the postprocessors are executed in reverse order, so this 2939 // will always dump the response just before it is returned since it is 2940 // the first one we register. 2941 if (absl::GetFlag(FLAGS_cp_model_dump_response)) { 2942 shared_response_manager->AddFinalResponsePostprocessor( 2943 [](CpSolverResponse* response) { 2944 const std::string file = absl::StrCat( 2945 absl::GetFlag(FLAGS_cp_model_dump_prefix), "response.pb.txt"); 2946 LOG(INFO) << "Dumping response proto to '" << file << "'."; 2947 CHECK_OK(file::SetTextProto(file, *response, file::Defaults())); 2948 }); 2949 } 2950 #endif // __PORTABLE_PLATFORM__ 2951 2952 // Always display the final response stats if requested. 2953 // This also copy the logs to the response if requested. 2954 shared_response_manager->AddFinalResponsePostprocessor( 2955 [logger, &model_proto, &log_string](CpSolverResponse* response) { 2956 SOLVER_LOG(logger, ""); 2957 SOLVER_LOG(logger, CpSolverResponseStats( 2958 *response, 2959 model_proto.has_objective() || 2960 model_proto.has_floating_point_objective())); 2961 if (!log_string.empty()) { 2962 response->set_solve_log(log_string); 2963 } 2964 }); 2965 2966 // Always add the timing information to a response. Note that it is important 2967 // to add this after the log/dump postprocessor since we execute them in 2968 // reverse order. 2969 auto* shared_time_limit = model->GetOrCreate<ModelSharedTimeLimit>(); 2970 shared_response_manager->AddResponsePostprocessor( 2971 [&wall_timer, &user_timer, 2972 &shared_time_limit](CpSolverResponse* response) { 2973 response->set_wall_time(wall_timer->Get()); 2974 response->set_user_time(user_timer->Get()); 2975 response->set_deterministic_time( 2976 shared_time_limit->GetElapsedDeterministicTime()); 2977 }); 2978 2979 // Validate parameters. 2980 // 2981 // Note that the few parameters we use before that are Booleans and thus 2982 // "safe". We need to delay the validation to return a proper response. 2983 { 2984 const std::string error = ValidateParameters(params); 2985 if (!error.empty()) { 2986 SOLVER_LOG(logger, "Invalid parameters: ", error); 2987 2988 // TODO(user): We currently reuse the MODEL_INVALID status even though it 2989 // is not the best name for this. Maybe we can add a PARAMETERS_INVALID 2990 // when it become needed. Or rename to INVALID_INPUT ? 2991 shared_response_manager->MutableResponse()->set_status( 2992 CpSolverStatus::MODEL_INVALID); 2993 shared_response_manager->MutableResponse()->set_solution_info(error); 2994 return shared_response_manager->GetResponse(); 2995 } 2996 } 2997 2998 // Initialize the time limit from the parameters. 2999 model->GetOrCreate<TimeLimit>()->ResetLimitFromParameters(params); 3000 3001 #if !defined(__PORTABLE_PLATFORM__) 3002 // Register SIGINT handler if requested by the parameters. 3003 if (params.catch_sigint_signal()) { 3004 model->GetOrCreate<SigintHandler>()->Register( 3005 [&shared_time_limit]() { shared_time_limit->Stop(); }); 3006 } 3007 #endif // __PORTABLE_PLATFORM__ 3008 3009 SOLVER_LOG(logger, ""); 3010 SOLVER_LOG(logger, "Starting CP-SAT solver v", OrToolsMajorVersion(), ".", 3011 OrToolsMinorVersion(), ".", OrToolsPatchVersion()); 3012 3013 // Initialize the number of workers if set to 0. 3014 if (params.num_search_workers() == 0) { 3015 #if !defined(__PORTABLE_PLATFORM__) 3016 // Sometimes, hardware_concurrency will return 0. So always default to 1. 3017 const int num_cores = 3018 params.enumerate_all_solutions() 3019 ? 1 3020 : std::max<int>(std::thread::hardware_concurrency(), 1); 3021 #else 3022 const int num_cores = 1; 3023 #endif 3024 SOLVER_LOG(logger, "Setting number of workers to ", num_cores); 3025 model->GetOrCreate<SatParameters>()->set_num_search_workers(num_cores); 3026 } 3027 3028 SOLVER_LOG(logger, "Parameters: ", params.ShortDebugString()); 3029 if (logger->LoggingIsEnabled() && params.use_absl_random()) { 3030 model->GetOrCreate<ModelRandomGenerator>()->LogSalt(); 3031 } 3032 3033 // Validate model_proto. 3034 // TODO(user): provide an option to skip this step for speed? 3035 { 3036 const std::string error = ValidateCpModel(model_proto); 3037 if (!error.empty()) { 3038 SOLVER_LOG(logger, "Invalid model: ", error); 3039 shared_response_manager->MutableResponse()->set_status( 3040 CpSolverStatus::MODEL_INVALID); 3041 shared_response_manager->MutableResponse()->set_solution_info(error); 3042 return shared_response_manager->GetResponse(); 3043 } 3044 } 3045 3046 SOLVER_LOG(logger, ""); 3047 SOLVER_LOG(logger, "Initial ", CpModelStats(model_proto)); 3048 3049 // Special case for pure-sat problem. 3050 // TODO(user): improve the normal presolver to do the same thing. 3051 // TODO(user): Support solution hint, but then the first TODO will make it 3052 // automatic. 3053 if (!params.use_sat_inprocessing() && !model_proto.has_objective() && 3054 !model_proto.has_floating_point_objective() && 3055 !model_proto.has_solution_hint() && !params.enumerate_all_solutions() && 3056 !params.use_lns_only() && params.num_search_workers() <= 1 && 3057 model_proto.assumptions().empty()) { 3058 bool is_pure_sat = true; 3059 for (const IntegerVariableProto& var : model_proto.variables()) { 3060 if (var.domain_size() != 2 || var.domain(0) < 0 || var.domain(1) > 1) { 3061 is_pure_sat = false; 3062 break; 3063 } 3064 } 3065 if (is_pure_sat) { 3066 for (const ConstraintProto& ct : model_proto.constraints()) { 3067 if (ct.constraint_case() != ConstraintProto::ConstraintCase::kBoolOr && 3068 ct.constraint_case() != ConstraintProto::ConstraintCase::kBoolAnd) { 3069 is_pure_sat = false; 3070 break; 3071 } 3072 } 3073 } 3074 if (is_pure_sat) { 3075 // TODO(user): All this duplication will go away when we are fast enough 3076 // on pure-sat model with the CpModel presolve... 3077 *shared_response_manager->MutableResponse() = 3078 SolvePureSatModel(model_proto, wall_timer, model, logger); 3079 if (params.fill_tightened_domains_in_response()) { 3080 *shared_response_manager->MutableResponse() 3081 ->mutable_tightened_variables() = model_proto.variables(); 3082 } 3083 return shared_response_manager->GetResponse(); 3084 } 3085 } 3086 3087 // Presolve and expansions. 3088 SOLVER_LOG(logger, ""); 3089 SOLVER_LOG(logger, 3090 absl::StrFormat("Starting presolve at %.2fs", wall_timer->Get())); 3091 CpModelProto new_cp_model_proto; 3092 CpModelProto mapping_proto; 3093 auto context = absl::make_unique<PresolveContext>(model, &new_cp_model_proto, 3094 &mapping_proto); 3095 3096 *context->working_model->mutable_variables() = model_proto.variables(); 3097 if (!ImportConstraintsWithBasicPresolveIntoContext(model_proto, 3098 context.get())) { 3099 VLOG(1) << "Model found infeasible during copy"; 3100 // TODO(user): At this point, the model is trivial, but we could exit 3101 // early. 3102 } 3103 CopyEverythingExceptVariablesAndConstraintsFieldsIntoContext(model_proto, 3104 context.get()); 3105 3106 // Checks for hints early in case they are forced to be hard constraints. 3107 if (params.fix_variables_to_their_hinted_value() && 3108 model_proto.has_solution_hint()) { 3109 SOLVER_LOG(context->logger(), "Fixing ", 3110 model_proto.solution_hint().vars().size(), 3111 " variables to their value in the solution hints."); 3112 for (int i = 0; i < model_proto.solution_hint().vars_size(); ++i) { 3113 const int var = model_proto.solution_hint().vars(i); 3114 const int64_t value = model_proto.solution_hint().values(i); 3115 if (!context->IntersectDomainWith(var, Domain(value))) { 3116 const IntegerVariableProto& var_proto = 3117 context->working_model->variables(var); 3118 const std::string var_name = var_proto.name().empty() 3119 ? absl::StrCat("var(", var, ")") 3120 : var_proto.name(); 3121 3122 const Domain var_domain = ReadDomainFromProto(var_proto); 3123 SOLVER_LOG(context->logger(), 3124 "Hint found infeasible when assigning variable '", var_name, 3125 "' with domain", var_domain.ToString(), " the value ", 3126 value); 3127 break; 3128 } 3129 } 3130 } 3131 3132 // If the hint is complete, we can use the solution checker to do more 3133 // validation. Note that after the model has been validated, we are sure there 3134 // are do duplicate variables in the solution hint, so we can just check the 3135 // size. 3136 // 3137 // TODO(user): Also check if the hint specifies all non-fixed variables. 3138 if (model_proto.has_solution_hint() && !context->ModelIsUnsat() && 3139 model_proto.solution_hint().vars().size() == 3140 model_proto.variables_size()) { 3141 std::vector<int64_t> solution(model_proto.variables_size(), 0); 3142 for (int i = 0; i < model_proto.solution_hint().vars_size(); ++i) { 3143 const int ref = model_proto.solution_hint().vars(i); 3144 const int64_t value = model_proto.solution_hint().values(i); 3145 solution[PositiveRef(ref)] = RefIsPositive(ref) ? value : -value; 3146 } 3147 if (SolutionIsFeasible(model_proto, solution)) { 3148 SOLVER_LOG(context->logger(), 3149 "The solution hint is complete and is feasible."); 3150 } else { 3151 // TODO(user): Change the code to make the solution checker more 3152 // informative by returning a message instead of just VLOGing it. 3153 SOLVER_LOG(context->logger(), 3154 "The solution hint is complete, but it is infeasible! we " 3155 "will try to repair it."); 3156 } 3157 } 3158 3159 // If the objective was a floating point one, do some postprocessing on the 3160 // final response. 3161 if (model_proto.has_floating_point_objective()) { 3162 shared_response_manager->AddFinalResponsePostprocessor( 3163 [¶ms, &model_proto, &mapping_proto, 3164 &logger](CpSolverResponse* response) { 3165 if (response->solution().empty()) return; 3166 3167 // Compute the true objective of the best returned solution. 3168 const auto& float_obj = model_proto.floating_point_objective(); 3169 double value = float_obj.offset(); 3170 const int num_terms = float_obj.vars().size(); 3171 for (int i = 0; i < num_terms; ++i) { 3172 value += float_obj.coeffs(i) * 3173 static_cast<double>(response->solution(float_obj.vars(i))); 3174 } 3175 response->set_objective_value(value); 3176 3177 // Also copy the scaled objective which must be in the mapping model. 3178 // This can be useful for some client, like if they want to do 3179 // multi-objective optimization in stages. 3180 if (!mapping_proto.has_objective()) return; 3181 const CpObjectiveProto& integer_obj = mapping_proto.objective(); 3182 *response->mutable_integer_objective() = integer_obj; 3183 3184 // If requested, compute a correct lb from the one on the integer 3185 // objective. We only do that if some error were introduced by the 3186 // scaling algorithm. 3187 if (params.mip_compute_true_objective_bound() && 3188 !integer_obj.scaling_was_exact()) { 3189 const int64_t integer_lb = response->inner_objective_lower_bound(); 3190 const double lb = ComputeTrueObjectiveLowerBound( 3191 model_proto, integer_obj, integer_lb); 3192 SOLVER_LOG(logger, "[Scaling] scaled_objective_bound: ", 3193 response->best_objective_bound(), 3194 " corrected_bound: ", lb, 3195 " delta: ", response->best_objective_bound() - lb); 3196 3197 // To avoid small errors that can be confusing, we take the 3198 // min/max with the objective value. 3199 if (float_obj.maximize()) { 3200 response->set_best_objective_bound( 3201 std::max(lb, response->objective_value())); 3202 } else { 3203 response->set_best_objective_bound( 3204 std::min(lb, response->objective_value())); 3205 } 3206 } 3207 3208 // Check the absolute gap, and display warning if needed. 3209 // TODO(user): Change status to IMPRECISE? 3210 if (response->status() == CpSolverStatus::OPTIMAL) { 3211 const double gap = std::abs(response->objective_value() - 3212 response->best_objective_bound()); 3213 if (gap > params.absolute_gap_limit()) { 3214 SOLVER_LOG(logger, 3215 "[Scaling] Warning: OPTIMAL was reported, yet the " 3216 "objective gap (", 3217 gap, ") is greater than requested absolute limit (", 3218 params.absolute_gap_limit(), ")."); 3219 } 3220 } 3221 }); 3222 } 3223 3224 if (params.num_search_workers() > 1 || model_proto.has_objective()) { 3225 // For the case where the assumptions are currently not supported, we just 3226 // assume they are fixed, and will always report all of them in the UNSAT 3227 // core if the problem turn out to be UNSAT. 3228 // 3229 // If the mode is not degraded, we will hopefully report a small subset 3230 // in case there is no feasible solution under these assumptions. 3231 shared_response_manager->AddFinalResponsePostprocessor( 3232 [&model_proto](CpSolverResponse* response) { 3233 if (response->status() != CpSolverStatus::INFEASIBLE) return; 3234 3235 // For now, just pass in all assumptions. 3236 *response->mutable_sufficient_assumptions_for_infeasibility() = 3237 model_proto.assumptions(); 3238 }); 3239 3240 context->InitializeNewDomains(); 3241 for (const int ref : model_proto.assumptions()) { 3242 if (!context->SetLiteralToTrue(ref)) { 3243 shared_response_manager->MutableResponse()->set_status( 3244 CpSolverStatus::INFEASIBLE); 3245 shared_response_manager->MutableResponse() 3246 ->add_sufficient_assumptions_for_infeasibility(ref); 3247 return shared_response_manager->GetResponse(); 3248 } 3249 } 3250 } 3251 3252 // Do the actual presolve. 3253 std::vector<int> postsolve_mapping; 3254 const CpSolverStatus presolve_status = 3255 PresolveCpModel(context.get(), &postsolve_mapping); 3256 if (presolve_status != CpSolverStatus::UNKNOWN) { 3257 SOLVER_LOG(logger, "Problem closed by presolve."); 3258 shared_response_manager->MutableResponse()->set_status(presolve_status); 3259 return shared_response_manager->GetResponse(); 3260 } 3261 3262 SOLVER_LOG(logger, ""); 3263 SOLVER_LOG(logger, "Presolved ", CpModelStats(new_cp_model_proto)); 3264 SOLVER_LOG(logger, ""); 3265 SOLVER_LOG(logger, "Preloading model."); 3266 3267 if (params.cp_model_presolve()) { 3268 shared_response_manager->AddSolutionPostprocessor( 3269 [&model_proto, ¶ms, &mapping_proto, &model, 3270 &postsolve_mapping](std::vector<int64_t>* solution) { 3271 AddPostsolveClauses(postsolve_mapping, model, &mapping_proto); 3272 PostsolveResponseWrapper(params, model_proto.variables_size(), 3273 mapping_proto, postsolve_mapping, solution); 3274 }); 3275 shared_response_manager->AddResponsePostprocessor( 3276 [&model_proto, ¶ms, &mapping_proto, &postsolve_mapping, wall_timer, 3277 model](CpSolverResponse* response) { 3278 // Map back the sufficient assumptions for infeasibility. 3279 for (int& ref : 3280 *(response 3281 ->mutable_sufficient_assumptions_for_infeasibility())) { 3282 ref = RefIsPositive(ref) 3283 ? postsolve_mapping[ref] 3284 : NegatedRef(postsolve_mapping[PositiveRef(ref)]); 3285 } 3286 if (!response->solution().empty()) { 3287 CHECK(SolutionIsFeasible( 3288 model_proto, 3289 std::vector<int64_t>(response->solution().begin(), 3290 response->solution().end()), 3291 &mapping_proto, &postsolve_mapping)) 3292 << "postsolved solution"; 3293 } 3294 if (params.fill_tightened_domains_in_response()) { 3295 // TODO(user): for now, we just use the domain infered during 3296 // presolve. 3297 if (mapping_proto.variables().size() >= 3298 model_proto.variables().size()) { 3299 for (int i = 0; i < model_proto.variables().size(); ++i) { 3300 *response->add_tightened_variables() = 3301 mapping_proto.variables(i); 3302 } 3303 } 3304 } 3305 }); 3306 } else { 3307 shared_response_manager->AddFinalResponsePostprocessor( 3308 [&model_proto](CpSolverResponse* response) { 3309 if (!response->solution().empty()) { 3310 CHECK(SolutionIsFeasible( 3311 model_proto, std::vector<int64_t>(response->solution().begin(), 3312 response->solution().end()))); 3313 } 3314 }); 3315 shared_response_manager->AddResponsePostprocessor( 3316 [&model_proto, ¶ms](CpSolverResponse* response) { 3317 // Truncate the solution in case model expansion added more variables. 3318 const int initial_size = model_proto.variables_size(); 3319 if (response->solution_size() > 0) { 3320 response->mutable_solution()->Truncate(initial_size); 3321 if (DEBUG_MODE || 3322 absl::GetFlag(FLAGS_cp_model_check_intermediate_solutions)) { 3323 CHECK(SolutionIsFeasible( 3324 model_proto, 3325 std::vector<int64_t>(response->solution().begin(), 3326 response->solution().end()))); 3327 } 3328 } 3329 if (params.fill_tightened_domains_in_response()) { 3330 *response->mutable_tightened_variables() = model_proto.variables(); 3331 } 3332 }); 3333 } 3334 3335 // Delete the context. 3336 context.reset(nullptr); 3337 3338 if (params.symmetry_level() > 1) { 3339 DetectAndAddSymmetryToProto(params, &new_cp_model_proto, logger); 3340 } 3341 3342 const auto& observers = model->GetOrCreate<SolutionObservers>()->observers; 3343 if (!observers.empty()) { 3344 shared_response_manager->AddSolutionCallback( 3345 [&observers](const CpSolverResponse& response) { 3346 for (const auto& observer : observers) { 3347 observer(response); 3348 } 3349 }); 3350 } 3351 3352 // If specified, we load the initial objective domain right away in the 3353 // response manager. Note that the presolve will always fill it with the 3354 // trivial min/max value if the user left it empty. This avoids to display 3355 // [-infinity, infinity] for the initial objective search space. 3356 if (new_cp_model_proto.has_objective()) { 3357 shared_response_manager->InitializeObjective(new_cp_model_proto); 3358 shared_response_manager->SetGapLimitsFromParameters(params); 3359 } 3360 3361 // Start counting the primal integral from the current determistic time and 3362 // initial objective domain gap that we just filled. 3363 shared_response_manager->UpdateGapIntegral(); 3364 3365 #if !defined(__PORTABLE_PLATFORM__) 3366 if (absl::GetFlag(FLAGS_cp_model_dump_models)) { 3367 const std::string presolved_file = absl::StrCat( 3368 absl::GetFlag(FLAGS_cp_model_dump_prefix), "presolved_model.pb.txt"); 3369 LOG(INFO) << "Dumping presolved cp model proto to '" << presolved_file 3370 << "'."; 3371 CHECK_OK(file::SetTextProto(presolved_file, new_cp_model_proto, 3372 file::Defaults())); 3373 3374 const std::string mapping_file = absl::StrCat( 3375 absl::GetFlag(FLAGS_cp_model_dump_prefix), "mapping_model.pb.txt"); 3376 LOG(INFO) << "Dumping mapping cp model proto to '" << mapping_file << "'."; 3377 CHECK_OK(file::SetTextProto(mapping_file, mapping_proto, file::Defaults())); 3378 } 3379 #endif // __PORTABLE_PLATFORM__ 3380 3381 if (params.stop_after_presolve() || shared_time_limit->LimitReached()) { 3382 int64_t num_terms = 0; 3383 for (const ConstraintProto& ct : new_cp_model_proto.constraints()) { 3384 num_terms += UsedVariables(ct).size(); 3385 } 3386 SOLVER_LOG( 3387 logger, "Stopped after presolve.", 3388 "\nPresolvedNumVariables: ", new_cp_model_proto.variables().size(), 3389 "\nPresolvedNumConstraints: ", new_cp_model_proto.constraints().size(), 3390 "\nPresolvedNumTerms: ", num_terms); 3391 3392 shared_response_manager->SetStatsFromModel(model); 3393 return shared_response_manager->GetResponse(); 3394 } 3395 3396 // Make sure everything stops when we have a first solution if requested. 3397 if (params.stop_after_first_solution()) { 3398 shared_response_manager->AddSolutionCallback( 3399 [shared_time_limit](const CpSolverResponse& response) { 3400 shared_time_limit->Stop(); 3401 }); 3402 } 3403 3404 #if defined(__PORTABLE_PLATFORM__) 3405 if (/* DISABLES CODE */ (false)) { 3406 // We ignore the multithreading parameter in this case. 3407 #else // __PORTABLE_PLATFORM__ 3408 if (params.num_search_workers() > 1 || params.interleave_search()) { 3409 SolveCpModelParallel(new_cp_model_proto, model); 3410 #endif // __PORTABLE_PLATFORM__ 3411 } else { 3412 SOLVER_LOG(logger, ""); 3413 SOLVER_LOG(logger, absl::StrFormat("Starting to load the model at %.2fs", 3414 wall_timer->Get())); 3415 shared_response_manager->SetUpdateGapIntegralOnEachChange(true); 3416 LoadCpModel(new_cp_model_proto, model); 3417 shared_response_manager->LoadDebugSolution(model); 3418 3419 SOLVER_LOG(logger, ""); 3420 SOLVER_LOG(logger, absl::StrFormat("Starting sequential search at %.2fs", 3421 wall_timer->Get())); 3422 if (params.repair_hint()) { 3423 MinimizeL1DistanceWithHint(new_cp_model_proto, model); 3424 } else { 3425 QuickSolveWithHint(new_cp_model_proto, model); 3426 } 3427 SolveLoadedCpModel(new_cp_model_proto, model); 3428 } 3429 3430 if (logger->LoggingIsEnabled()) { 3431 if (params.num_search_workers() <= 1) { 3432 const auto& lps = 3433 *model->GetOrCreate<LinearProgrammingConstraintCollection>(); 3434 if (!lps.empty()) { 3435 SOLVER_LOG(logger, ""); 3436 for (const auto* lp : lps) { 3437 SOLVER_LOG(logger, lp->Statistics()); 3438 } 3439 } 3440 } 3441 3442 if (params.num_search_workers() > 1) { 3443 SOLVER_LOG(logger, ""); 3444 shared_response_manager->DisplayImprovementStatistics(); 3445 } 3446 } 3447 3448 return shared_response_manager->GetResponse(); 3449 } 3450 3451 CpSolverResponse Solve(const CpModelProto& model_proto) { 3452 Model model; 3453 return SolveCpModel(model_proto, &model); 3454 } 3455 3456 CpSolverResponse SolveWithParameters(const CpModelProto& model_proto, 3457 const SatParameters& params) { 3458 Model model; 3459 model.Add(NewSatParameters(params)); 3460 return SolveCpModel(model_proto, &model); 3461 } 3462 3463 #if !defined(__PORTABLE_PLATFORM__) 3464 CpSolverResponse SolveWithParameters(const CpModelProto& model_proto, 3465 const std::string& params) { 3466 Model model; 3467 model.Add(NewSatParameters(params)); 3468 return SolveCpModel(model_proto, &model); 3469 } 3470 #endif // !__PORTABLE_PLATFORM__ 3471 3472 } // namespace sat 3473 } // namespace operations_research 3474