diff --git a/ortools/sat/BUILD.bazel b/ortools/sat/BUILD.bazel index 3f40cda1c5..a9c7234d86 100644 --- a/ortools/sat/BUILD.bazel +++ b/ortools/sat/BUILD.bazel @@ -196,6 +196,7 @@ cc_library( ":lb_tree_search", ":linear_programming_constraint", ":linear_relaxation", + ":max_hs", ":model", ":optimization", ":parameters_validation", @@ -1169,6 +1170,38 @@ cc_library( ], ) +cc_library( + name = "max_hs", + srcs = ["max_hs.cc"], + hdrs = ["max_hs.h"], + deps = [ + ":boolean_problem", + ":boolean_problem_cc_proto", + ":cp_model_mapping", + ":cp_model_utils", + ":encoding", + ":integer", + ":integer_expr", + ":integer_search", + ":linear_relaxation", + ":model", + ":optimization", + ":pb_constraint", + ":sat_base", + ":sat_parameters_cc_proto", + ":sat_solver", + ":util", + "//ortools/base", + "//ortools/base:map_util", + "//ortools/linear_solver:linear_solver_cc_proto", + "//ortools/port:proto_utils", + "//ortools/util:strong_integers", + "//ortools/util:time_limit", + "@com_google_absl//absl/strings", + "@com_google_absl//absl/strings:str_format", + ], +) + cc_library( name = "util", srcs = ["util.cc"], diff --git a/ortools/sat/cp_model_search.cc b/ortools/sat/cp_model_search.cc index c6c7dd3220..645a778ddf 100644 --- a/ortools/sat/cp_model_search.cc +++ b/ortools/sat/cp_model_search.cc @@ -23,6 +23,10 @@ #include "ortools/sat/integer.h" #include "ortools/sat/util.h" +// TODO(user): remove this when the code is stable and does not use SCIP +// anymore. +ABSL_FLAG(bool, cp_model_use_max_hs, false, "Use max_hs in search portfolio."); + namespace operations_research { namespace sat { @@ -450,7 +454,6 @@ std::vector GetDiverseSetOfParameters( new_params.set_search_branching(SatParameters::AUTOMATIC_SEARCH); new_params.set_optimize_with_core(true); new_params.set_optimize_with_max_hs(true); - new_params.set_find_multiple_cores(false); strategies["max_hs"] = new_params; } @@ -562,6 +565,9 @@ std::vector GetDiverseSetOfParameters( names.push_back("quick_restart_no_lp"); names.push_back("lb_tree_search"); names.push_back("probing"); +#if !defined(__PORTABLE_PLATFORM__) && defined(USE_SCIP) + if (absl::GetFlag(FLAGS_cp_model_use_max_hs)) names.push_back("max_hs"); +#endif // !defined(__PORTABLE_PLATFORM__) && defined(USE_SCIP) } else { for (const std::string& name : base_params.subsolvers()) { names.push_back(name); diff --git a/ortools/sat/cp_model_solver.cc b/ortools/sat/cp_model_solver.cc index 316afeaede..4cd07a3f1e 100644 --- a/ortools/sat/cp_model_solver.cc +++ b/ortools/sat/cp_model_solver.cc @@ -78,6 +78,7 @@ #include "ortools/sat/linear_programming_constraint.h" #include "ortools/sat/linear_relaxation.h" #include "ortools/sat/lp_utils.h" +#include "ortools/sat/max_hs.h" #include "ortools/sat/optimization.h" #include "ortools/sat/parameters_validation.h" #include "ortools/sat/precedences.h" @@ -1080,7 +1081,8 @@ void LoadBaseModel(const CpModelProto& model_proto, Model* model) { const bool view_all_booleans_as_integers = (parameters.linearization_level() >= 2) || (parameters.search_branching() == SatParameters::FIXED_SEARCH && - model_proto.search_strategy().empty()); + model_proto.search_strategy().empty()) || + parameters.optimize_with_max_hs(); LoadVariables(model_proto, view_all_booleans_as_integers, model); DetectOptionalVariables(model_proto, model); @@ -1412,11 +1414,18 @@ void LoadCpModel(const CpModelProto& model_proto, Model* model) { }; const auto& objective = *model->GetOrCreate(); - CoreBasedOptimizer* core = - new CoreBasedOptimizer(objective_var, objective.vars, objective.coeffs, - solution_observer, model); - model->Register(core); - model->TakeOwnership(core); + if (parameters.optimize_with_max_hs()) { + HittingSetOptimizer* max_hs = new HittingSetOptimizer( + model_proto, objective, solution_observer, model); + model->Register(max_hs); + model->TakeOwnership(max_hs); + } else { + CoreBasedOptimizer* core = + new CoreBasedOptimizer(objective_var, objective.vars, + objective.coeffs, solution_observer, model); + model->Register(core); + model->TakeOwnership(core); + } } } @@ -1506,8 +1515,7 @@ void SolveLoadedCpModel(const CpModelProto& model_proto, Model* model) { // TODO(user): This doesn't work with splitting in chunk for now. It // shouldn't be too hard to fix. if (parameters.optimize_with_max_hs()) { - status = MinimizeWithHittingSetAndLazyEncoding( - objective, solution_observer, model); + status = model->Mutable()->Optimize(); } else { status = model->Mutable()->Optimize(); } diff --git a/ortools/sat/encoding.cc b/ortools/sat/encoding.cc index bf29ba1ba2..fe8e0a058f 100644 --- a/ortools/sat/encoding.cc +++ b/ortools/sat/encoding.cc @@ -377,9 +377,10 @@ std::vector ReduceNodesAndExtractAssumptions( } // Fix the nodes right-most variables that are above the gap. + // If we closed the problem, we abort and return and empty vector. if (upper_bound != kCoefficientMax) { const Coefficient gap = upper_bound - *lower_bound; - if (gap <= 0) return {}; + if (gap < 0) return {}; for (EncodingNode* n : *nodes) { n->ApplyUpperBound((gap / n->weight()).value(), solver); } diff --git a/ortools/sat/integer_expr.cc b/ortools/sat/integer_expr.cc index 735f1ae8d0..6478ec0212 100644 --- a/ortools/sat/integer_expr.cc +++ b/ortools/sat/integer_expr.cc @@ -1404,7 +1404,7 @@ bool FixedModuloPropagator::PropagateOuterBounds() { bool FixedModuloPropagator::PropagateBoundsWhenExprIsPositive( AffineExpression expr, AffineExpression target) { const IntegerValue min_target = integer_trail_->LowerBound(target); - DCHECK_GE(min_target, 0); + DCHECK_GE(min_target, 0) << target.DebugString(); const IntegerValue max_target = integer_trail_->UpperBound(target); // The propagation rules below will not be triggered if the domain of target diff --git a/ortools/sat/optimization.cc b/ortools/sat/optimization.cc index ee410e24b5..ce592f1108 100644 --- a/ortools/sat/optimization.cc +++ b/ortools/sat/optimization.cc @@ -35,11 +35,6 @@ #include "ortools/base/map_util.h" #include "ortools/base/stl_util.h" #include "ortools/base/timer.h" -#include "ortools/util/strong_integers.h" -#if !defined(__PORTABLE_PLATFORM__) && defined(USE_SCIP) -#include "ortools/linear_solver/linear_solver.h" -#include "ortools/linear_solver/linear_solver.pb.h" -#endif // __PORTABLE_PLATFORM__ #include "ortools/port/proto_utils.h" #include "ortools/sat/boolean_problem.h" #include "ortools/sat/encoding.h" @@ -47,6 +42,7 @@ #include "ortools/sat/pb_constraint.h" #include "ortools/sat/sat_parameters.pb.h" #include "ortools/sat/util.h" +#include "ortools/util/strong_integers.h" #include "ortools/util/time_limit.h" namespace operations_research { @@ -258,7 +254,17 @@ void MinimizeCoreWithPropagation(TimeLimit* limit, SatSolver* solver, solver->SetAssumptionLevel(0); if (candidate.size() < core->size()) { VLOG(1) << "minimization " << core->size() << " -> " << candidate.size(); - core->assign(candidate.begin(), candidate.end()); + + // We want to preserve the order of literal in the response. + absl::flat_hash_set set; + for (const Literal l : candidate) set.insert(l.Index()); + int new_size = 0; + for (const Literal l : *core) { + if (set.contains(l.Index())) { + (*core)[new_size++] = l; + } + } + core->resize(new_size); } } @@ -993,6 +999,9 @@ SatSolver::Status SolveWithCardinalityEncodingAndCore( int max_depth = 0; std::string previous_core_info = ""; for (int iter = 0;; ++iter) { + // TODO(user): We are suboptimal here because we use for upper bound the + // current best objective, not best_obj - 1. This code is not really used + // but we should still fix it. const std::vector assumptions = ReduceNodesAndExtractAssumptions( upper_bound, stratified_lower_bound, &lower_bound, &nodes, solver); if (assumptions.empty()) return SatSolver::FEASIBLE; @@ -1266,63 +1275,6 @@ SatSolver::Status FindCores(std::vector assumptions, return SatSolver::ASSUMPTIONS_UNSAT; } -// Slightly different algo than FindCores() which aim to extract more cores, but -// not necessarily non-overlaping ones. -SatSolver::Status FindMultipleCoresForMaxHs( - std::vector assumptions, Model* model, - std::vector>* cores) { - cores->clear(); - SatSolver* sat_solver = model->GetOrCreate(); - TimeLimit* limit = model->GetOrCreate(); - const double saved_dlimit = limit->GetDeterministicLimit(); - auto cleanup = ::absl::MakeCleanup([limit, saved_dlimit]() { - limit->ChangeDeterministicLimit(saved_dlimit); - }); - - bool first_loop = true; - do { - if (limit->LimitReached()) return SatSolver::LIMIT_REACHED; - - // The order of assumptions do not matter. - // Randomizing it should improve diversity. - auto* random = model->GetOrCreate(); - std::shuffle(assumptions.begin(), assumptions.end(), *random); - - const SatSolver::Status result = - ResetAndSolveIntegerProblem(assumptions, model); - if (result != SatSolver::ASSUMPTIONS_UNSAT) return result; - std::vector core = sat_solver->GetLastIncompatibleDecisions(); - if (sat_solver->parameters().minimize_core()) { - MinimizeCoreWithPropagation(limit, sat_solver, &core); - } - CHECK(!core.empty()); - cores->push_back(core); - if (!sat_solver->parameters().find_multiple_cores()) break; - - // Pick a random literal from the core and remove it from the set of - // assumptions. - CHECK(!core.empty()); - const Literal random_literal = - core[absl::Uniform(*random, 0, core.size())]; - for (int i = 0; i < assumptions.size(); ++i) { - if (assumptions[i] == random_literal) { - std::swap(assumptions[i], assumptions.back()); - assumptions.pop_back(); - break; - } - } - - // Once we found at least one core, we impose a time limit to not spend - // too much time finding more. - if (first_loop) { - limit->ChangeDeterministicLimit( - std::min(saved_dlimit, limit->GetElapsedDeterministicTime() + 1.0)); - first_loop = false; - } - } while (!assumptions.empty()); - return SatSolver::ASSUMPTIONS_UNSAT; -} - } // namespace CoreBasedOptimizer::CoreBasedOptimizer( @@ -1555,7 +1507,139 @@ bool CoreBasedOptimizer::CoverOptimization() { return PropagateObjectiveBounds(); } +SatSolver::Status CoreBasedOptimizer::OptimizeWithSatEncoding( + const std::vector& literals, + const std::vector& coefficients, IntegerValue offset) { + // Create one initial nodes per variables with cost. + // TODO(user): We could create EncodingNode out of IntegerVariable. + std::deque repository; + Coefficient unused = 0; + std::vector nodes = + CreateInitialEncodingNodes(literals, coefficients, &unused, &repository); + CHECK_EQ(unused, 0); + + // Initialize the bounds. + // This is in term of number of variables not at their minimal value. + Coefficient lower_bound(0); + + // This is used by the "stratified" approach. + // TODO(user): Take into account parameters. + Coefficient stratified_lower_bound(0); + for (EncodingNode* n : nodes) { + stratified_lower_bound = std::max(stratified_lower_bound, n->weight()); + } + + // Start the algorithm. + int max_depth = 0; + std::string previous_core_info = ""; + for (int iter = 0;;) { + if (time_limit_->LimitReached()) return SatSolver::LIMIT_REACHED; + if (!sat_solver_->ResetToLevelZero()) return SatSolver::INFEASIBLE; + + const Coefficient upper_bound( + (integer_trail_->UpperBound(objective_var_) - offset).value()); + const std::vector assumptions = ReduceNodesAndExtractAssumptions( + upper_bound, stratified_lower_bound, &lower_bound, &nodes, sat_solver_); + if (assumptions.empty()) { + stratified_lower_bound = + MaxNodeWeightSmallerThan(nodes, stratified_lower_bound); + if (stratified_lower_bound > 0) continue; + + // We do not have any assumptions anymore, but we still need to see + // if the problem is feasible or not! + } + const IntegerValue new_obj_lb(lower_bound.value() + offset.value()); + if (new_obj_lb > integer_trail_->LowerBound(objective_var_)) { + if (!integer_trail_->Enqueue( + IntegerLiteral::GreaterOrEqual(objective_var_, new_obj_lb), {}, + {})) { + return SatSolver::INFEASIBLE; + } + + // Report the improvement. + // Note that we have a callback that will do the same, but doing it + // earlier allow us to add more information. + const int num_bools = sat_solver_->NumVariables(); + const int num_fixed = sat_solver_->NumFixedVariables(); + model_->GetOrCreate()->UpdateInnerObjectiveBounds( + absl::StrFormat("BoolCore num_cores:%d [%s] assumptions:%u " + "depth:%d fixed_bools:%d/%d", + iter, previous_core_info, nodes.size(), max_depth, + num_fixed, num_bools), + new_obj_lb, integer_trail_->LevelZeroUpperBound(objective_var_)); + } + + // Solve under the assumptions. + // + // TODO(user): Find multiple core like in the "main" algorithm. + // this is just trying to solve with assumptions not involving the newly + // found core. + const SatSolver::Status result = + ResetAndSolveIntegerProblem(assumptions, model_); + if (result == SatSolver::FEASIBLE) { + if (!ProcessSolution()) return SatSolver::INFEASIBLE; + if (stop_) return SatSolver::LIMIT_REACHED; + + // If not all assumptions were taken, continue with a lower stratified + // bound. Otherwise we have an optimal solution. + stratified_lower_bound = + MaxNodeWeightSmallerThan(nodes, stratified_lower_bound); + if (stratified_lower_bound > 0) continue; + return SatSolver::INFEASIBLE; + } + if (result != SatSolver::ASSUMPTIONS_UNSAT) return result; + + // We have a new core. + std::vector core = sat_solver_->GetLastIncompatibleDecisions(); + if (parameters_->minimize_core()) { + MinimizeCoreWithPropagation(time_limit_, sat_solver_, &core); + } + + // Compute the min weight of all the nodes in the core. + // The lower bound will be increased by that much. + const Coefficient min_weight = ComputeCoreMinWeight(nodes, core); + previous_core_info = + absl::StrFormat("core:%u mw:%d", core.size(), min_weight.value()); + + // We only count an iter when we found a core. + ++iter; + ProcessCore(core, min_weight, &repository, &nodes, sat_solver_); + max_depth = std::max(max_depth, nodes.back()->depth()); + } + + return SatSolver::FEASIBLE; // shouldn't reach here. +} + SatSolver::Status CoreBasedOptimizer::Optimize() { + // Hack: If the objective is fully Boolean, we use the + // OptimizeWithSatEncoding() version as it seems to be better. + // + // TODO(user): Try to understand exactly why and merge both code path. + if (!parameters_->interleave_search()) { + IntegerValue offset(0); + std::vector literals; + std::vector coefficients; + bool all_booleans = true; + for (const ObjectiveTerm& term : terms_) { + const IntegerVariable var = term.var; + const IntegerValue coeff = term.weight; + const IntegerValue lb = integer_trail_->LowerBound(var); + const IntegerValue ub = integer_trail_->UpperBound(var); + if (ub - lb == 1) { + offset += lb * coeff; + literals.push_back(integer_encoder_->GetOrCreateAssociatedLiteral( + IntegerLiteral::GreaterOrEqual(var, ub))); + coefficients.push_back(Coefficient(coeff.value())); + } else { + all_booleans = false; + break; + } + } + if (all_booleans) { + return OptimizeWithSatEncoding(literals, coefficients, offset); + } + } + // TODO(user): The core is returned in the same order as the assumptions, // so we don't really need this map, we could just do a linear scan to // recover which node are part of the core. This however needs to be properly @@ -1765,247 +1849,10 @@ SatSolver::Status CoreBasedOptimizer::Optimize() { } // Abort if we reached the time limit. Note that we still add any cores we - // found in case the solve is splitted in "chunk". + // found in case the solve is split in "chunk". if (result == SatSolver::LIMIT_REACHED) return result; } } -// TODO(user): take the MPModelRequest or MPModelProto directly, so that we can -// have initial constraints! -// -// TODO(user): remove code duplication with MinimizeWithCoreAndLazyEncoding(); -SatSolver::Status MinimizeWithHittingSetAndLazyEncoding( - const ObjectiveDefinition& objective_definition, - const std::function& feasible_solution_observer, Model* model) { -#if !defined(__PORTABLE_PLATFORM__) && defined(USE_SCIP) - - IntegerVariable objective_var = objective_definition.objective_var; - std::vector variables = objective_definition.vars; - std::vector coefficients = objective_definition.coeffs; - - auto* sat_solver = model->GetOrCreate(); - auto* integer_trail = model->GetOrCreate(); - auto* integer_encoder = model->GetOrCreate(); - auto* time_limit = model->GetOrCreate(); - - // This will be called each time a feasible solution is found. - const auto process_solution = [&]() { - // We don't assume that objective_var is linked with its linear term, so - // we recompute the objective here. - IntegerValue objective(0); - for (int i = 0; i < variables.size(); ++i) { - objective += - coefficients[i] * IntegerValue(model->Get(Value(variables[i]))); - } - if (objective > integer_trail->UpperBound(objective_var)) return true; - - if (feasible_solution_observer != nullptr) { - feasible_solution_observer(); - } - - // Constrain objective_var. This has a better result when objective_var is - // used in an LP relaxation for instance. - sat_solver->Backtrack(0); - sat_solver->SetAssumptionLevel(0); - if (!integer_trail->Enqueue( - IntegerLiteral::LowerOrEqual(objective_var, objective - 1), {}, - {})) { - return false; - } - return true; - }; - - // This is the "generalized" hitting set problem we will solve. Each time - // we find a core, a new constraint will be added to this problem. - MPModelRequest request; - request.set_solver_specific_parameters("limits/gap = 0"); - request.set_solver_type(MPModelRequest::SCIP_MIXED_INTEGER_PROGRAMMING); - - MPModelProto& hs_model = *request.mutable_model(); - const int num_variables_in_objective = variables.size(); - for (int i = 0; i < num_variables_in_objective; ++i) { - if (coefficients[i] < 0) { - variables[i] = NegationOf(variables[i]); - coefficients[i] = -coefficients[i]; - } - const IntegerVariable var = variables[i]; - MPVariableProto* var_proto = hs_model.add_variable(); - var_proto->set_lower_bound(integer_trail->LowerBound(var).value()); - var_proto->set_upper_bound(integer_trail->UpperBound(var).value()); - var_proto->set_objective_coefficient(coefficients[i].value()); - var_proto->set_is_integer(true); - } - - MPSolutionResponse response; - - // This is used by the "stratified" approach. We will only consider terms with - // a weight not lower than this threshold. The threshold will decrease as the - // algorithm progress. - IntegerValue stratified_threshold = kMaxIntegerValue; - - // TODO(user): The core is returned in the same order as the assumptions, - // so we don't really need this map, we could just do a linear scan to - // recover which node are part of the core. - std::map> assumption_to_indices; - - // New Booleans variable in the MIP model to represent X >= cte. - std::map, int> created_var; - - const SatParameters& parameters = *(model->GetOrCreate()); - - // Start the algorithm. - SatSolver::Status result; - for (int iter = 0;; ++iter) { - // TODO(user): Even though we keep the same solver, currently the solve is - // not really done incrementally. It might be hard to improve though. - // - // TODO(user): C^c is broken when using SCIP. - MPSolver::SolveWithProto(request, &response); - if (response.status() != MPSolverResponseStatus::MPSOLVER_OPTIMAL) { - // We currently abort if we have a non-optimal result. - // This is correct if we had a limit reached, but not in the other cases. - // - // TODO(user): It is actually easy to use a FEASIBLE result. If when - // passing it to SAT it is no feasbile, we can still create cores. If it - // is feasible, we have a solution, but we cannot increase the lower - // bound. - return SatSolver::LIMIT_REACHED; - } - CHECK_EQ(response.status(), MPSolverResponseStatus::MPSOLVER_OPTIMAL); - - const IntegerValue mip_objective( - static_cast(std::round(response.objective_value()))); - VLOG(1) << "constraints: " << hs_model.constraint_size() - << " variables: " << hs_model.variable_size() << " hs_lower_bound: " - << objective_definition.ScaleIntegerObjective(mip_objective) - << " strat: " << stratified_threshold; - - // Update the objective lower bound with our current bound. - // - // Note(user): This is not needed for correctness, but it might cause - // more propagation and is nice to have for reporting/logging purpose. - if (!integer_trail->Enqueue( - IntegerLiteral::GreaterOrEqual(objective_var, mip_objective), {}, - {})) { - result = SatSolver::INFEASIBLE; - break; - } - - sat_solver->Backtrack(0); - sat_solver->SetAssumptionLevel(0); - std::vector assumptions; - assumption_to_indices.clear(); - IntegerValue next_stratified_threshold(0); - for (int i = 0; i < num_variables_in_objective; ++i) { - const IntegerValue hs_value( - static_cast(response.variable_value(i))); - if (hs_value == integer_trail->UpperBound(variables[i])) continue; - - // Only consider the terms above the threshold. - if (coefficients[i] < stratified_threshold) { - next_stratified_threshold = - std::max(next_stratified_threshold, coefficients[i]); - } else { - // It is possible that different variables have the same associated - // literal. So we do need to consider this case. - assumptions.push_back(integer_encoder->GetOrCreateAssociatedLiteral( - IntegerLiteral::LowerOrEqual(variables[i], hs_value))); - assumption_to_indices[assumptions.back().Index()].push_back(i); - } - } - - // No assumptions with the current stratified_threshold? use the new one. - if (assumptions.empty() && next_stratified_threshold > 0) { - CHECK_LT(next_stratified_threshold, stratified_threshold); - stratified_threshold = next_stratified_threshold; - --iter; // "false" iteration, the lower bound does not increase. - continue; - } - - // TODO(user): we could also randomly shuffle the assumptions to find more - // cores for only one MIP solve. - // - // TODO(user): Use the real weights and exploit the extra cores. - std::vector> cores; - result = FindMultipleCoresForMaxHs(assumptions, model, &cores); - if (result == SatSolver::FEASIBLE) { - if (!process_solution()) return SatSolver::INFEASIBLE; - if (parameters.stop_after_first_solution()) { - return SatSolver::LIMIT_REACHED; - } - if (cores.empty()) { - // If not all assumptions were taken, continue with a lower stratified - // bound. Otherwise we have an optimal solution. - stratified_threshold = next_stratified_threshold; - if (stratified_threshold == 0) break; - --iter; // "false" iteration, the lower bound does not increase. - continue; - } - } else if (result == SatSolver::LIMIT_REACHED) { - // Hack: we use a local limit internally that we restore at the end. - // However we still return LIMIT_REACHED in this case... - if (time_limit->LimitReached()) break; - } else if (result != SatSolver::ASSUMPTIONS_UNSAT) { - break; - } - - sat_solver->Backtrack(0); - sat_solver->SetAssumptionLevel(0); - for (const std::vector& core : cores) { - if (core.size() == 1) { - for (const int index : - gtl::FindOrDie(assumption_to_indices, core.front().Index())) { - hs_model.mutable_variable(index)->set_lower_bound( - integer_trail->LowerBound(variables[index]).value()); - } - continue; - } - - // Add the corresponding constraint to hs_model. - MPConstraintProto* ct = hs_model.add_constraint(); - ct->set_lower_bound(1.0); - for (const Literal lit : core) { - for (const int index : - gtl::FindOrDie(assumption_to_indices, lit.Index())) { - const double lb = integer_trail->LowerBound(variables[index]).value(); - const double hs_value = response.variable_value(index); - if (hs_value == lb) { - ct->add_var_index(index); - ct->add_coefficient(1.0); - ct->set_lower_bound(ct->lower_bound() + lb); - } else { - const std::pair key = {index, hs_value}; - if (!gtl::ContainsKey(created_var, key)) { - const int new_var_index = hs_model.variable_size(); - created_var[key] = new_var_index; - - MPVariableProto* new_var = hs_model.add_variable(); - new_var->set_lower_bound(0); - new_var->set_upper_bound(1); - new_var->set_is_integer(true); - - // (new_var == 1) => x > hs_value. - // (x - lb) - (hs_value - lb + 1) * new_var >= 0. - MPConstraintProto* implication = hs_model.add_constraint(); - implication->set_lower_bound(lb); - implication->add_var_index(index); - implication->add_coefficient(1.0); - implication->add_var_index(new_var_index); - implication->add_coefficient(lb - hs_value - 1); - } - ct->add_var_index(gtl::FindOrDieNoPrint(created_var, key)); - ct->add_coefficient(1.0); - } - } - } - } - } - - return result; -#else // !__PORTABLE_PLATFORM__ && USE_SCIP - LOG(FATAL) << "Not supported."; -#endif // !__PORTABLE_PLATFORM__ && USE_SCIP -} - } // namespace sat } // namespace operations_research diff --git a/ortools/sat/optimization.h b/ortools/sat/optimization.h index 99699dba57..58d160744e 100644 --- a/ortools/sat/optimization.h +++ b/ortools/sat/optimization.h @@ -34,7 +34,7 @@ namespace sat { // be inferred by propagation by any subset of the other literal, it will be // removed. // -// Note that this function doest NOT preserve the order of Literal in the core. +// Note that the literal of the minimized core will stay in the same order. // // TODO(user): Avoid spending too much time trying to minimize a core. void MinimizeCoreWithPropagation(TimeLimit* limit, SatSolver* solver, @@ -154,6 +154,19 @@ class CoreBasedOptimizer { // some of the work already done, so it might just never find anything. SatSolver::Status Optimize(); + // A different way to encode the objective as core are found. This one do + // not introduce IntegerVariable and encode everything in Boolean. + // + // It seems to be more powerful, but it isn't completely implemented yet. + // TODO(user): + // - Make it work for integer variable in the objective. + // - Only call it if the objective domain is not too large? + // - Support resuming for interleaved search. + // - Implement all core heurisitics. + SatSolver::Status OptimizeWithSatEncoding( + const std::vector& literals, + const std::vector& coefficients, IntegerValue offset); + private: CoreBasedOptimizer(const CoreBasedOptimizer&) = delete; CoreBasedOptimizer& operator=(const CoreBasedOptimizer&) = delete; @@ -210,24 +223,6 @@ class CoreBasedOptimizer { bool stop_ = false; }; -// Generalization of the max-HS algorithm (HS stands for Hitting Set). This is -// similar to MinimizeWithCoreAndLazyEncoding() but it uses a hybrid approach -// with a MIP solver to handle the discovered infeasibility cores. -// -// See, Jessica Davies and Fahiem Bacchus, "Solving MAXSAT by Solving a -// Sequence of Simpler SAT Instances", -// http://www.cs.toronto.edu/~jdavies/daviesCP11.pdf" -// -// Note that an implementation of this approach won the 2016 max-SAT competition -// on the industrial category, see -// http://maxsat.ia.udl.cat/results/#wpms-industrial -// -// TODO(user): This function requires linking with the SCIP MIP solver which is -// quite big, maybe we should find a way not to do that. -SatSolver::Status MinimizeWithHittingSetAndLazyEncoding( - const ObjectiveDefinition& objective_definition, - const std::function& feasible_solution_observer, Model* model); - } // namespace sat } // namespace operations_research