diff --git a/cmake/flatzinc.cmake b/cmake/flatzinc.cmake index 6a74b799e8..55304f8935 100644 --- a/cmake/flatzinc.cmake +++ b/cmake/flatzinc.cmake @@ -70,6 +70,8 @@ add_library(flatzinc ortools/flatzinc/parser.yy.cc #ortools/flatzinc/parser_util.cc # Already #include in parser.tab.cc ortools/flatzinc/parser_util.h + ortools/flatzinc/presolve.cc + ortools/flatzinc/presolve.h ) ## Includes target_include_directories(flatzinc PUBLIC diff --git a/ortools/flatzinc/BUILD.bazel b/ortools/flatzinc/BUILD.bazel index 112bbd75bd..d3e8b22116 100644 --- a/ortools/flatzinc/BUILD.bazel +++ b/ortools/flatzinc/BUILD.bazel @@ -109,6 +109,21 @@ cc_library( ], ) +cc_library( + name = "presolve", + srcs = ["presolve.cc"], + hdrs = ["presolve.h"], + deps = [ + ":model", + "//ortools/base", + "//ortools/base:hash", + "//ortools/graph:cliques", + "//ortools/util:logging", + "//ortools/util:saturated_arithmetic", + "@com_google_absl//absl/strings", + ], +) + cc_library( name = "checker", srcs = ["checker.cc"], @@ -157,6 +172,7 @@ cc_binary( ":cp_model_fz_solver", ":model", ":parser_lib", + ":presolve", "//ortools/base", "//ortools/base:path", "//ortools/base:threadpool", diff --git a/ortools/flatzinc/checker.cc b/ortools/flatzinc/checker.cc index 6b179c311b..e21375ffe2 100644 --- a/ortools/flatzinc/checker.cc +++ b/ortools/flatzinc/checker.cc @@ -183,27 +183,6 @@ bool CheckArrayVarIntElement( return element == target; } -bool CheckOrtoolsArrayIntElement( - const Constraint& ct, const std::function& evaluator) { - const int64_t min_index = ct.arguments[1].values[0]; - const int64_t index = Eval(ct.arguments[0], evaluator) - min_index; - const int64_t element = EvalAt(ct.arguments[2], index, evaluator); - const int64_t target = Eval(ct.arguments[3], evaluator); - return element == target; -} - -bool CheckOrtoolsArrayIntElement2d( - const Constraint& ct, const std::function& evaluator) { - const int64_t min_index0 = ct.arguments[2].values[0]; - const int64_t min_index1 = ct.arguments[3].values[0]; - const int64_t span1 = ct.arguments[3].values[1] - min_index1 + 1; - const int64_t index0 = Eval(ct.arguments[0], evaluator) - min_index0; - const int64_t index1 = Eval(ct.arguments[1], evaluator) - min_index1; - const int64_t element = - EvalAt(ct.arguments[4], index0 * span1 + index1, evaluator); - const int64_t target = Eval(ct.arguments[5], evaluator); - return element == target; -} bool CheckAtMostInt(const Constraint& ct, const std::function& evaluator) { const int64_t expected = Eval(ct.arguments[0], evaluator); @@ -1205,6 +1184,7 @@ using CallMap = // They are created at compilation time when using the or-tools mzn library. CallMap CreateCallMap() { CallMap m; + m["fzn_all_different_int"] = CheckAllDifferentInt; m["alldifferent_except_0"] = CheckAlldifferentExcept0; m["among"] = CheckAmong; m["array_bool_and"] = CheckArrayBoolAnd; @@ -1213,137 +1193,130 @@ CallMap CreateCallMap() { m["array_bool_xor"] = CheckArrayBoolXor; m["array_int_element"] = CheckArrayIntElement; m["array_int_element_nonshifted"] = CheckArrayIntElementNonShifted; - m["array_int_maximum"] = CheckMaximumInt; - m["array_int_minimum"] = CheckMinimumInt; m["array_var_bool_element"] = CheckArrayVarIntElement; m["array_var_int_element"] = CheckArrayVarIntElement; m["at_most_int"] = CheckAtMostInt; m["bool_and"] = CheckBoolAnd; m["bool_clause"] = CheckBoolClause; + m["bool_eq"] = CheckIntEq; + m["bool2int"] = CheckIntEq; m["bool_eq_imp"] = CheckIntEqImp; m["bool_eq_reif"] = CheckIntEqReif; - m["bool_eq"] = CheckIntEq; + m["bool_ge"] = CheckIntGe; m["bool_ge_imp"] = CheckIntGeImp; m["bool_ge_reif"] = CheckIntGeReif; - m["bool_ge"] = CheckIntGe; + m["bool_gt"] = CheckIntGt; m["bool_gt_imp"] = CheckIntGtImp; m["bool_gt_reif"] = CheckIntGtReif; - m["bool_gt"] = CheckIntGt; + m["bool_le"] = CheckIntLe; m["bool_le_imp"] = CheckIntLeImp; m["bool_le_reif"] = CheckIntLeReif; - m["bool_le"] = CheckIntLe; m["bool_left_imp"] = CheckIntLe; m["bool_lin_eq"] = CheckIntLinEq; m["bool_lin_le"] = CheckIntLinLe; + m["bool_lt"] = CheckIntLt; m["bool_lt_imp"] = CheckIntLtImp; m["bool_lt_reif"] = CheckIntLtReif; - m["bool_lt"] = CheckIntLt; + m["bool_ne"] = CheckIntNe; m["bool_ne_imp"] = CheckIntNeImp; m["bool_ne_reif"] = CheckIntNeReif; - m["bool_ne"] = CheckIntNe; m["bool_not"] = CheckBoolNot; m["bool_or"] = CheckBoolOr; m["bool_right_imp"] = CheckIntGe; m["bool_xor"] = CheckBoolXor; - m["bool2int"] = CheckIntEq; + m["ortools_circuit"] = CheckCircuit; m["count_eq"] = CheckCountEq; + m["count"] = CheckCountEq; m["count_geq"] = CheckCountGeq; m["count_gt"] = CheckCountGt; m["count_leq"] = CheckCountLeq; m["count_lt"] = CheckCountLt; m["count_neq"] = CheckCountNeq; m["count_reif"] = CheckCountReif; - m["count"] = CheckCountEq; - m["diffn_k_with_sizes"] = CheckDiffnK; - m["diffn_nonstrict_k_with_sizes"] = CheckDiffnNonStrictK; - m["false_constraint"] = CheckFalseConstraint; - m["fixed_cumulative"] = CheckCumulative; - m["fzn_all_different_int"] = CheckAllDifferentInt; m["fzn_cumulative"] = CheckCumulative; - m["fzn_diffn_nonstrict"] = CheckDiffnNonStrict; + m["var_cumulative"] = CheckCumulative; + m["variable_cumulative"] = CheckCumulative; + m["fixed_cumulative"] = CheckCumulative; + m["ortools_cumulative_opt"] = CheckCumulativeOpt; m["fzn_diffn"] = CheckDiffn; - m["fzn_disjunctive_strict"] = CheckDisjunctiveStrict; + m["diffn_k_with_sizes"] = CheckDiffnK; + m["fzn_diffn_nonstrict"] = CheckDiffnNonStrict; + m["diffn_nonstrict_k_with_sizes"] = CheckDiffnNonStrictK; m["fzn_disjunctive"] = CheckDisjunctive; - m["global_cardinality_closed"] = CheckGlobalCardinalityClosed; - m["global_cardinality_low_up_closed"] = CheckGlobalCardinalityLowUpClosed; - m["global_cardinality_low_up"] = CheckGlobalCardinalityLowUp; - m["global_cardinality_old"] = CheckGlobalCardinalityOld; + m["fzn_disjunctive_strict"] = CheckDisjunctiveStrict; + m["ortools_disjunctive_strict_opt"] = CheckDisjunctiveStrictOpt; + m["false_constraint"] = CheckFalseConstraint; m["global_cardinality"] = CheckGlobalCardinality; + m["global_cardinality_closed"] = CheckGlobalCardinalityClosed; + m["global_cardinality_low_up"] = CheckGlobalCardinalityLowUp; + m["global_cardinality_low_up_closed"] = CheckGlobalCardinalityLowUpClosed; + m["global_cardinality_old"] = CheckGlobalCardinalityOld; m["int_abs"] = CheckIntAbs; m["int_div"] = CheckIntDiv; + m["int_eq"] = CheckIntEq; m["int_eq_imp"] = CheckIntEqImp; m["int_eq_reif"] = CheckIntEqReif; - m["int_eq"] = CheckIntEq; + m["int_ge"] = CheckIntGe; m["int_ge_imp"] = CheckIntGeImp; m["int_ge_reif"] = CheckIntGeReif; - m["int_ge"] = CheckIntGe; + m["int_gt"] = CheckIntGt; m["int_gt_imp"] = CheckIntGtImp; m["int_gt_reif"] = CheckIntGtReif; - m["int_gt"] = CheckIntGt; - m["int_in"] = CheckSetIn; + m["int_le"] = CheckIntLe; m["int_le_imp"] = CheckIntLeImp; m["int_le_reif"] = CheckIntLeReif; - m["int_le"] = CheckIntLe; + m["int_lin_eq"] = CheckIntLinEq; m["int_lin_eq_imp"] = CheckIntLinEqImp; m["int_lin_eq_reif"] = CheckIntLinEqReif; - m["int_lin_eq"] = CheckIntLinEq; + m["int_lin_ge"] = CheckIntLinGe; m["int_lin_ge_imp"] = CheckIntLinGeImp; m["int_lin_ge_reif"] = CheckIntLinGeReif; - m["int_lin_ge"] = CheckIntLinGe; + m["int_lin_le"] = CheckIntLinLe; m["int_lin_le_imp"] = CheckIntLinLeImp; m["int_lin_le_reif"] = CheckIntLinLeReif; - m["int_lin_le"] = CheckIntLinLe; + m["int_lin_ne"] = CheckIntLinNe; m["int_lin_ne_imp"] = CheckIntLinNeImp; m["int_lin_ne_reif"] = CheckIntLinNeReif; - m["int_lin_ne"] = CheckIntLinNe; + m["int_lt"] = CheckIntLt; m["int_lt_imp"] = CheckIntLtImp; m["int_lt_reif"] = CheckIntLtReif; - m["int_lt"] = CheckIntLt; m["int_max"] = CheckIntMax; m["int_min"] = CheckIntMin; m["int_minus"] = CheckIntMinus; m["int_mod"] = CheckIntMod; + m["int_ne"] = CheckIntNe; m["int_ne_imp"] = CheckIntNeImp; m["int_ne_reif"] = CheckIntNeReif; - m["int_ne"] = CheckIntNe; m["int_negate"] = CheckIntNegate; - m["int_not_in"] = CheckSetNotIn; m["int_plus"] = CheckIntPlus; m["int_times"] = CheckIntTimes; + m["ortools_inverse"] = CheckInverse; m["lex_less_bool"] = CheckLexLessInt; m["lex_less_int"] = CheckLexLessInt; m["lex_lesseq_bool"] = CheckLexLesseqInt; m["lex_lesseq_int"] = CheckLexLesseqInt; m["maximum_arg_int"] = CheckMaximumArgInt; m["maximum_int"] = CheckMaximumInt; + m["array_int_maximum"] = CheckMaximumInt; m["minimum_arg_int"] = CheckMinimumArgInt; m["minimum_int"] = CheckMinimumInt; - m["nvalue"] = CheckNvalue; - m["ortools_array_bool_element"] = CheckOrtoolsArrayIntElement; - m["ortools_array_int_element"] = CheckOrtoolsArrayIntElement; - m["ortools_array_var_bool_element"] = CheckOrtoolsArrayIntElement; - m["ortools_array_var_bool_element2d"] = CheckOrtoolsArrayIntElement2d; - m["ortools_array_var_int_element"] = CheckOrtoolsArrayIntElement; - m["ortools_array_var_int_element2d"] = CheckOrtoolsArrayIntElement2d; - m["ortools_circuit"] = CheckCircuit; - m["ortools_cumulative_opt"] = CheckCumulativeOpt; - m["ortools_disjunctive_strict_opt"] = CheckDisjunctiveStrictOpt; - m["ortools_inverse"] = CheckInverse; - m["ortools_network_flow_cost"] = CheckNetworkFlowCost; + m["array_int_minimum"] = CheckMinimumInt; m["ortools_network_flow"] = CheckNetworkFlow; + m["ortools_network_flow_cost"] = CheckNetworkFlowCost; + m["nvalue"] = CheckNvalue; m["ortools_regular"] = CheckRegular; - m["ortools_subcircuit"] = CheckSubCircuit; - m["ortools_table_bool"] = CheckTableInt; - m["ortools_table_int"] = CheckTableInt; m["regular_nfa"] = CheckRegularNfa; - m["set_in_reif"] = CheckSetInReif; m["set_in"] = CheckSetIn; + m["int_in"] = CheckSetIn; m["set_not_in"] = CheckSetNotIn; + m["int_not_in"] = CheckSetNotIn; + m["set_in_reif"] = CheckSetInReif; m["sliding_sum"] = CheckSlidingSum; m["sort"] = CheckSort; + m["ortools_subcircuit"] = CheckSubCircuit; m["symmetric_all_different"] = CheckSymmetricAllDifferent; - m["var_cumulative"] = CheckCumulative; - m["variable_cumulative"] = CheckCumulative; + m["ortools_table_bool"] = CheckTableInt; + m["ortools_table_int"] = CheckTableInt; return m; } diff --git a/ortools/flatzinc/cp_model_fz_solver.cc b/ortools/flatzinc/cp_model_fz_solver.cc index 25438654a6..d343d1b4b7 100644 --- a/ortools/flatzinc/cp_model_fz_solver.cc +++ b/ortools/flatzinc/cp_model_fz_solver.cc @@ -19,7 +19,6 @@ #include #include #include -#include #include #include "absl/container/flat_hash_map.h" @@ -77,9 +76,6 @@ struct CpModelProtoWithMapping { std::vector LookupVars(const fz::Argument& argument); std::vector LookupVarsOrValues(const fz::Argument& argument); - // Encoding literals. - int GetOrCreateVarEqValueLiteral(int var, int64_t value); - // Create and return the indices of the IntervalConstraint corresponding // to the flatzinc "interval" specified by a start var and a size var. // This method will cache intervals with the key . @@ -138,7 +134,6 @@ struct CpModelProtoWithMapping { absl::flat_hash_map, int> interval_key_to_index; absl::flat_hash_map var_to_lit_implies_greater_than_zero; - absl::flat_hash_map, int> value_encoding_literals; }; int CpModelProtoWithMapping::LookupConstant(int64_t value) { @@ -232,34 +227,6 @@ std::vector CpModelProtoWithMapping::LookupVarsOrValues( return result; } -int CpModelProtoWithMapping::GetOrCreateVarEqValueLiteral(int var, - int64_t value) { - const std::pair key = {var, value}; - const auto it = value_encoding_literals.find(key); - if (it != value_encoding_literals.end()) { - return it->second; - } - const int result = proto.variables_size(); - IntegerVariableProto* var_proto = proto.add_variables(); - var_proto->add_domain(0); - var_proto->add_domain(1); - value_encoding_literals[key] = result; - - ConstraintProto* pos_enforcement = AddEnforcedConstraint(result); - pos_enforcement->mutable_linear()->add_vars(var); - pos_enforcement->mutable_linear()->add_coeffs(1); - pos_enforcement->mutable_linear()->add_domain(value); - pos_enforcement->mutable_linear()->add_domain(value); - - ConstraintProto* neg_enforcement = AddEnforcedConstraint(NegatedRef(result)); - neg_enforcement->mutable_linear()->add_vars(var); - neg_enforcement->mutable_linear()->add_coeffs(1); - const Domain complement = Domain(value).Complement().IntersectionWith( - ReadDomainFromProto(proto.variables(var))); - FillDomainInProto(complement, neg_enforcement->mutable_linear()); - return result; -} - ConstraintProto* CpModelProtoWithMapping::AddEnforcedConstraint(int literal) { ConstraintProto* result = proto.add_constraints(); if (literal != kNoVar) { @@ -667,87 +634,46 @@ void CpModelProtoWithMapping::FillConstraint(const fz::Constraint& fz_ct, fz_ct.type == "array_var_int_element" || fz_ct.type == "array_var_bool_element" || fz_ct.type == "array_int_element_nonshifted") { - // Compatibility with the old format. - CHECK(fz_ct.arguments[0].type == fz::Argument::VAR_REF || - fz_ct.arguments[0].type == fz::Argument::INT_VALUE); - auto* arg = ct->mutable_element(); - arg->set_index(LookupVar(fz_ct.arguments[0])); - arg->set_target(LookupVar(fz_ct.arguments[2])); - - if (!absl::EndsWith(fz_ct.type, "_nonshifted")) { - // Add a dummy variable at position zero because flatzinc index start - // at 1. - // TODO(user): Make sure that zero is not in the index domain... - arg->add_vars(LookupConstant(0)); - } - for (const int var : LookupVars(fz_ct.arguments[1])) arg->add_vars(var); - } else if (fz_ct.type == "ortools_array_int_element" || - fz_ct.type == "ortools_array_bool_element" || - fz_ct.type == "ortools_array_var_int_element" || - fz_ct.type == "ortools_array_var_bool_element") { if (fz_ct.arguments[0].type == fz::Argument::VAR_REF || fz_ct.arguments[0].type == fz::Argument::INT_VALUE) { auto* arg = ct->mutable_element(); arg->set_index(LookupVar(fz_ct.arguments[0])); - arg->set_target(LookupVar(fz_ct.arguments[3])); - CHECK_EQ(fz_ct.arguments[1].type, fz::Argument::INT_INTERVAL); - const int64_t min_index = fz_ct.arguments[1].values.front(); - if (min_index > 0) { - const int zero_cst = LookupConstant(0); - for (int i = 0; i < min_index; ++i) { - arg->add_vars(zero_cst); - } + arg->set_target(LookupVar(fz_ct.arguments[2])); + + if (!absl::EndsWith(fz_ct.type, "_nonshifted")) { + // Add a dummy variable at position zero because flatzinc index start + // at 1. + // TODO(user): Make sure that zero is not in the index domain... + arg->add_vars(LookupConstant(0)); } - for (const int var : LookupVars(fz_ct.arguments[2])) arg->add_vars(var); - } - } else if (fz_ct.type == "ortools_array_var_int_element2d" || - fz_ct.type == "ortools_array_var_bool_element2d") { - const int index1 = LookupVar(fz_ct.arguments[0]); - const int index2 = LookupVar(fz_ct.arguments[1]); - const int target = LookupVar(fz_ct.arguments[5]); - - CHECK_EQ(fz_ct.arguments[2].type, fz::Argument::INT_INTERVAL); - CHECK_EQ(fz_ct.arguments[3].type, fz::Argument::INT_INTERVAL); - const int64_t min_1 = fz_ct.arguments[2].values[0]; - const int64_t max_1 = fz_ct.arguments[2].values[1]; - const int64_t min_2 = fz_ct.arguments[3].values[0]; - const int64_t max_2 = fz_ct.arguments[3].values[1]; - - if (fz_ct.arguments[4].type == fz::Argument::INT_LIST) { - // If the array is constant, we encode this as a table constraint. - auto* arg = ct->mutable_table(); - arg->add_vars(index1); - arg->add_vars(index2); - arg->add_vars(target); - - int i = 0; - for (int64_t val_1 = min_1; val_1 <= max_1; ++val_1) { - for (int64_t val_2 = min_2; val_2 <= max_2; ++val_2) { - arg->add_values(val_1); - arg->add_values(val_2); - arg->add_values(fz_ct.arguments[4].ValueAt(i++)); - } - } - CHECK_EQ(i, fz_ct.arguments[4].Size()); + for (const int var : LookupVars(fz_ct.arguments[1])) arg->add_vars(var); } else { - std::vector elems = LookupVars(fz_ct.arguments[4]); - int i = 0; - for (int64_t val_1 = min_1; val_1 <= max_1; ++val_1) { - const int lit1 = GetOrCreateVarEqValueLiteral(index1, val_1); - for (int64_t val_2 = min_2; val_2 <= max_2; ++val_2) { - const int lit2 = GetOrCreateVarEqValueLiteral(index2, val_2); - if (i != 0) ct = proto.add_constraints(); // new constraint. - ct->add_enforcement_literal(lit1); - ct->add_enforcement_literal(lit2); - ct->mutable_linear()->add_vars(target); - ct->mutable_linear()->add_coeffs(1); - ct->mutable_linear()->add_vars(elems[i++]); - ct->mutable_linear()->add_coeffs(-1); - ct->mutable_linear()->add_domain(0); - ct->mutable_linear()->add_domain(0); + // Special case added by the presolve or in flatzinc. We encode this + // as a table constraint. + CHECK(!absl::EndsWith(fz_ct.type, "_nonshifted")); + auto* arg = ct->mutable_table(); + + // the constraint is: + // values[coeff1 * vars[0] + coeff2 * vars[1] + offset] == target. + for (const int var : LookupVars(fz_ct.arguments[0])) arg->add_vars(var); + arg->add_vars(LookupVar(fz_ct.arguments[2])); // the target + + const std::vector& values = fz_ct.arguments[1].values; + const int64_t coeff1 = fz_ct.arguments[3].values[0]; + const int64_t coeff2 = fz_ct.arguments[3].values[1]; + const int64_t offset = fz_ct.arguments[4].values[0] - 1; + + for (const int64_t a : AllValuesInDomain(proto.variables(arg->vars(0)))) { + for (const int64_t b : + AllValuesInDomain(proto.variables(arg->vars(1)))) { + const int index = coeff1 * a + coeff2 * b + offset; + CHECK_GE(index, 0); + CHECK_LT(index, values.size()); + arg->add_values(a); + arg->add_values(b); + arg->add_values(values[index]); } } - CHECK_EQ(i, fz_ct.arguments[4].Size()); } } else if (fz_ct.type == "ortools_table_int") { auto* arg = ct->mutable_table(); diff --git a/ortools/flatzinc/fz.cc b/ortools/flatzinc/fz.cc index ee1ce9ef11..042102dae7 100644 --- a/ortools/flatzinc/fz.cc +++ b/ortools/flatzinc/fz.cc @@ -39,6 +39,7 @@ #include "ortools/flatzinc/cp_model_fz_solver.h" #include "ortools/flatzinc/model.h" #include "ortools/flatzinc/parser.h" +#include "ortools/flatzinc/presolve.h" #include "ortools/util/logging.h" ABSL_FLAG(double, time_limit, 0, "time limit in seconds."); @@ -49,6 +50,7 @@ ABSL_FLAG(bool, free_search, false, "If false, the solver must follow the defined search." "If true, other search are allowed."); ABSL_FLAG(int, threads, 0, "Number of threads the solver will use."); +ABSL_FLAG(bool, presolve, true, "Presolve the model to simplify it."); ABSL_FLAG(bool, statistics, false, "Print solver statistics after search."); ABSL_FLAG(bool, read_from_stdin, false, "Read the FlatZinc from stdin, not from a file."); @@ -151,6 +153,15 @@ Model ParseFlatzincModel(const std::string& input, bool input_is_filename, " parsed in ", timer.GetInMs(), " ms"); SOLVER_LOG(logger, ""); + // Presolve the model. + Presolver presolve(logger); + SOLVER_LOG(logger, "Presolve model"); + timer.Reset(); + timer.Start(); + presolve.Run(&model); + SOLVER_LOG(logger, " - done in ", timer.GetInMs(), " ms"); + SOLVER_LOG(logger); + // Print statistics. ModelStatistics stats(model, logger); stats.BuildStatistics(); diff --git a/ortools/flatzinc/mznlib/redefinitions-2.0.2.mzn b/ortools/flatzinc/mznlib/redefinitions-2.0.2.mzn deleted file mode 100644 index 74e28cd766..0000000000 --- a/ortools/flatzinc/mznlib/redefinitions-2.0.2.mzn +++ /dev/null @@ -1,26 +0,0 @@ -% Ignore. -predicate symmetry_breaking_constraint(var bool: b) = b; - -predicate redundant_constraint(var bool: b) = b; - -% array_var_bool_element_nonshifted. -predicate ortools_array_var_bool_element(var int: idx, - set of int: domain_of_x, - array [int] of var bool: x, - var bool: c); - -predicate array_var_bool_element_nonshifted(var int: idx, - array [int] of var bool: x, - var bool: c) = - ortools_array_var_bool_element(idx, index_set(x), x, c); - -% array_var_int_element_nonshifted. -predicate ortools_array_var_int_element(var int: idx, - set of int: domain_of_x, - array [int] of var int: x, - var int: c); - -predicate array_var_int_element_nonshifted(var int: idx, - array [int] of var int: x, - var int: c) = - ortools_array_var_int_element(idx, index_set(x), x, c); diff --git a/ortools/flatzinc/mznlib/redefinitions-2.5.2.mzn b/ortools/flatzinc/mznlib/redefinitions-2.5.2.mzn deleted file mode 100644 index e60bd8661f..0000000000 --- a/ortools/flatzinc/mznlib/redefinitions-2.5.2.mzn +++ /dev/null @@ -1,49 +0,0 @@ -% array_var_bool_element2d_nonshifted. -predicate ortools_array_var_bool_element2d(var int: idx1, - var int: idx2, - set of int: domain_of_x_1, - set of int: domain_of_x_2, - array[int] of var bool: x, - var bool: c); - -predicate array_var_bool_element2d_nonshifted(var int: idx1, - var int: idx2, - array[int,int] of var bool: x, - var bool: c) = - ortools_array_var_bool_element2d(idx1, - idx2, - index_set_1of2(x), - index_set_2of2(x), - array1d(x), - c); - -% array_var_int_element2d_nonshifted. -predicate ortools_array_var_int_element2d(var int: idx1, - var int: idx2, - set of int: domain_of_x_1, - set of int: domain_of_x_2, - array[int] of var int: x, - var int: c); - -predicate array_var_int_element2d_nonshifted(var int: idx1, - var int: idx2, - array[int,int] of var int: x, - var int: c) = - ortools_array_var_int_element2d(idx1, - idx2, - index_set_1of2(x), - index_set_2of2(x), - array1d(x), - c); - -predicate array_var_float_element2d_nonshifted(var int: idx1, var int: idx2, array[int,int] of var float: x, var float: c) = - let { - int: dim = card(index_set_2of2(x)); - int: min_flat = min(index_set_1of2(x))*dim+min(index_set_2of2(x))-1; - } in array_var_float_element_nonshifted((idx1*dim+idx2-min_flat)::domain, array1d(x), c); - -predicate array_var_set_element2d_nonshifted(var int: idx1, var int: idx2, array[int,int] of var set of int: x, var set of int: c) = - let { - int: dim = card(index_set_2of2(x)); - int: min_flat = min(index_set_1of2(x))*dim+min(index_set_2of2(x))-1; - } in array_var_set_element_nonshifted((idx1*dim+idx2-min_flat)::domain, array1d(x), c); diff --git a/ortools/flatzinc/parser_main.cc b/ortools/flatzinc/parser_main.cc index f571606c4f..eaebbedd00 100644 --- a/ortools/flatzinc/parser_main.cc +++ b/ortools/flatzinc/parser_main.cc @@ -25,6 +25,7 @@ #include "ortools/base/timer.h" #include "ortools/flatzinc/model.h" #include "ortools/flatzinc/parser.h" +#include "ortools/flatzinc/presolve.h" #include "ortools/util/logging.h" ABSL_FLAG(std::string, input, "", "Input file in the flatzinc format."); @@ -34,7 +35,7 @@ ABSL_FLAG(bool, statistics, false, "Print model statistics"); namespace operations_research { namespace fz { -void ParseFile(const std::string& filename) { +void ParseFile(const std::string& filename, bool presolve) { WallTimer timer; timer.Start(); @@ -57,6 +58,14 @@ void ParseFile(const std::string& filename) { Model model(problem_name); CHECK(ParseFlatzincFile(filename, &model)); + if (presolve) { + SOLVER_LOG(&logger, "Presolve model"); + timer.Reset(); + timer.Start(); + Presolver presolve(&logger); + presolve.Run(&model); + SOLVER_LOG(&logger, " - done in ", timer.GetInMs(), " ms"); + } if (absl::GetFlag(FLAGS_statistics)) { ModelStatistics stats(model, &logger); stats.BuildStatistics(); @@ -76,6 +85,7 @@ int main(int argc, char** argv) { absl::SetProgramUsageMessage(kUsage); absl::ParseCommandLine(argc, argv); google::InitGoogleLogging(argv[0]); - operations_research::fz::ParseFile(absl::GetFlag(FLAGS_input)); + operations_research::fz::ParseFile(absl::GetFlag(FLAGS_input), + absl::GetFlag(FLAGS_presolve)); return 0; } diff --git a/ortools/flatzinc/presolve.cc b/ortools/flatzinc/presolve.cc new file mode 100644 index 0000000000..2887f8db46 --- /dev/null +++ b/ortools/flatzinc/presolve.cc @@ -0,0 +1,495 @@ +// Copyright 2010-2024 Google LLC +// Licensed under the Apache License, Version 2.0 (the "License"); +// you may not use this file except in compliance with the License. +// You may obtain a copy of the License at +// +// http://www.apache.org/licenses/LICENSE-2.0 +// +// Unless required by applicable law or agreed to in writing, software +// distributed under the License is distributed on an "AS IS" BASIS, +// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +// See the License for the specific language governing permissions and +// limitations under the License. + +#include "ortools/flatzinc/presolve.h" + +#include +#include +#include +#include +#include + +#include "absl/container/flat_hash_set.h" +#include "absl/flags/flag.h" +#include "absl/log/check.h" +#include "absl/strings/string_view.h" +#include "absl/types/span.h" +#include "ortools/flatzinc/model.h" +#include "ortools/util/logging.h" + +ABSL_FLAG(bool, fz_floats_are_ints, false, + "Interpret floats as integers in all variables and constraints."); + +namespace operations_research { +namespace fz { +namespace { +enum PresolveState { ALWAYS_FALSE, ALWAYS_TRUE, UNDECIDED }; + +template +bool IsArrayBoolean(const std::vector& values) { + for (int i = 0; i < values.size(); ++i) { + if (values[i] != 0 && values[i] != 1) { + return false; + } + } + return true; +} + +template +bool AtMostOne0OrAtMostOne1(const std::vector& values) { + CHECK(IsArrayBoolean(values)); + int num_zero = 0; + int num_one = 0; + for (T val : values) { + if (val) { + num_one++; + } else { + num_zero++; + } + if (num_one > 1 && num_zero > 1) { + return false; + } + } + return true; +} + +template +void AppendIfNotInSet(T* value, absl::flat_hash_set* s, + std::vector* vec) { + if (s->insert(value).second) { + vec->push_back(value); + } + DCHECK_EQ(s->size(), vec->size()); +} + +} // namespace + +// Note on documentation +// +// In order to document presolve rules, we will use the following naming +// convention: +// - x, x1, xi, y, y1, yi denote integer variables +// - b, b1, bi denote boolean variables +// - c, c1, ci denote integer constants +// - t, t1, ti denote boolean constants +// - => x after a constraint denotes the target variable of this constraint. +// Arguments are listed in order. + +// Propagates cast constraint. +// Rule 1: +// Input: bool2int(b, c) or bool2int(t, x) +// Output: int_eq(...) +// +// Rule 2: +// Input: bool2int(b, x) +// Action: Replace all instances of x by b. +// Output: inactive constraint +void Presolver::PresolveBool2Int(Constraint* ct) { + DCHECK_EQ(ct->type, "bool2int"); + if (ct->arguments[0].HasOneValue() || ct->arguments[1].HasOneValue()) { + // Rule 1. + UpdateRuleStats("bool2int: rename to int_eq"); + ct->type = "int_eq"; + } else { + // Rule 2. + UpdateRuleStats("bool2int: merge boolean and integer variables."); + AddVariableSubstitution(ct->arguments[1].Var(), ct->arguments[0].Var()); + ct->MarkAsInactive(); + } +} + +// Propagates cast constraint. +// Rule 1: +// Input: int2float(x, y) +// Action: Replace all instances of y by x. +// Output: inactive constraint +void Presolver::PresolveInt2Float(Constraint* ct) { + DCHECK_EQ(ct->type, "int2float"); + // Rule 1. + UpdateRuleStats("int2float: merge integer and floating point variables."); + AddVariableSubstitution(ct->arguments[1].Var(), ct->arguments[0].Var()); + ct->MarkAsInactive(); +} + +void Presolver::PresolveStoreFlatteningMapping(Constraint* ct) { + CHECK_EQ(3, ct->arguments[1].variables.size()); + Variable* const var0 = ct->arguments[1].variables[0]; + Variable* const var1 = ct->arguments[1].variables[1]; + Variable* const var2 = ct->arguments[1].variables[2]; + const int64_t coeff0 = ct->arguments[0].values[0]; + const int64_t coeff1 = ct->arguments[0].values[1]; + const int64_t coeff2 = ct->arguments[0].values[2]; + const int64_t rhs = ct->arguments[2].Value(); + if (coeff0 == -1 && coeff2 == 1 && !array2d_index_map_.contains(var0)) { + array2d_index_map_[var0] = + Array2DIndexMapping(var1, coeff1, var2, -rhs, ct); + UpdateRuleStats("int_lin_eq: store 2d flattening mapping"); + } else if (coeff0 == -1 && coeff1 == 1 && + !array2d_index_map_.contains(var0)) { + array2d_index_map_[var0] = + Array2DIndexMapping(var2, coeff2, var1, -rhs, ct); + UpdateRuleStats("int_lin_eq: store 2d flattening mapping"); + } else if (coeff2 == -1 && coeff1 == 1 && + !array2d_index_map_.contains(var2)) { + array2d_index_map_[var2] = + Array2DIndexMapping(var0, coeff0, var1, -rhs, ct); + UpdateRuleStats("int_lin_eq: store 2d flattening mapping"); + } else if (coeff2 == -1 && coeff0 == 1 && + !array2d_index_map_.contains(var2)) { + array2d_index_map_[var2] = + Array2DIndexMapping(var1, coeff1, var0, -rhs, ct); + UpdateRuleStats("int_lin_eq: store 2d flattening mapping"); + } +} + +namespace { +bool IsIncreasingAndContiguous(absl::Span values) { + for (int i = 0; i < values.size() - 1; ++i) { + if (values[i + 1] != values[i] + 1) { + return false; + } + } + return true; +} + +bool AreOnesFollowedByMinusOne(absl::Span coeffs) { + CHECK(!coeffs.empty()); + for (int i = 0; i < coeffs.size() - 1; ++i) { + if (coeffs[i] != 1) { + return false; + } + } + return coeffs.back() == -1; +} + +template +bool IsStrictPrefix(const std::vector& v1, const std::vector& v2) { + if (v1.size() >= v2.size()) { + return false; + } + for (int i = 0; i < v1.size(); ++i) { + if (v1[i] != v2[i]) { + return false; + } + } + return true; +} +} // namespace + +// Rewrite array element: array_int_element: +// +// Rule 1: +// Input : array_int_element(x, [c1, .., cn], y) with x = a * x1 + x2 + b +// Output: array_int_element([x1, x2], [c_a1, .., c_am], b, [a, b]) +// to be interpreted by the extraction process. +// +// Rule 2: +// Input : array_int_element(x, [c1, .., cn], y) with x0 ci = c0 + i +// Output: int_lin_eq([-1, 1], [y, x], 1 - c) (e.g. y = x + c - 1) +void Presolver::PresolveSimplifyElement(Constraint* ct) { + if (ct->arguments[0].variables.size() != 1) return; + Variable* const index_var = ct->arguments[0].Var(); + + // Rule 1. + if (array2d_index_map_.contains(index_var)) { + UpdateRuleStats("array_int_element: rewrite as a 2d element"); + const Array2DIndexMapping& mapping = array2d_index_map_[index_var]; + // Rewrite constraint. + ct->arguments[0] = + Argument::VarRefArray({mapping.variable1, mapping.variable2}); + std::vector coefs; + coefs.push_back(mapping.coefficient); + coefs.push_back(1); + ct->arguments.push_back(Argument::IntegerList(coefs)); + ct->arguments.push_back(Argument::IntegerValue(mapping.offset)); + index_var->active = false; + mapping.constraint->MarkAsInactive(); + return; + } + + // Rule 2. + if (IsIncreasingAndContiguous(ct->arguments[1].values) && + ct->arguments[2].type == Argument::VAR_REF) { + const int64_t start = ct->arguments[1].values.front(); + Variable* const index = ct->arguments[0].Var(); + Variable* const target = ct->arguments[2].Var(); + UpdateRuleStats("array_int_element: rewrite as a linear constraint"); + + if (start == 1) { + ct->type = "int_eq"; + ct->RemoveArg(1); + } else { + // Rewrite constraint into a int_lin_eq + ct->type = "int_lin_eq"; + ct->arguments[0] = Argument::IntegerList({-1, 1}); + ct->arguments[1] = Argument::VarRefArray({target, index}); + ct->arguments[2] = Argument::IntegerValue(1 - start); + } + } +} + +void Presolver::Run(Model* model) { + // Should rewrite float constraints. + if (absl::GetFlag(FLAGS_fz_floats_are_ints)) { + // Treat float variables as int variables, convert constraints to int. + for (Constraint* const ct : model->constraints()) { + const std::string& id = ct->type; + if (id == "int2float") { + ct->type = "int_eq"; + } else if (id == "float_lin_le") { + ct->type = "int_lin_le"; + } else if (id == "float_lin_eq") { + ct->type = "int_lin_eq"; + } + } + } + + // Regroup increasing sequence of int_lin_eq([1,..,1,-1], [x1, ..., xn, yn]) + // into sequence of int_plus(x1, x2, y2), int_plus(y2, x3, y3)... + std::vector current_variables; + Variable* target_variable = nullptr; + Constraint* first_constraint = nullptr; + for (Constraint* const ct : model->constraints()) { + if (target_variable == nullptr) { + if (ct->type == "int_lin_eq" && ct->arguments[0].values.size() == 3 && + AreOnesFollowedByMinusOne(ct->arguments[0].values) && + ct->arguments[1].values.empty() && ct->arguments[2].Value() == 0) { + current_variables = ct->arguments[1].variables; + target_variable = current_variables.back(); + current_variables.pop_back(); + first_constraint = ct; + } + } else { + if (ct->type == "int_lin_eq" && + AreOnesFollowedByMinusOne(ct->arguments[0].values) && + ct->arguments[0].values.size() == current_variables.size() + 2 && + IsStrictPrefix(current_variables, ct->arguments[1].variables)) { + current_variables = ct->arguments[1].variables; + // Rewrite ct into int_plus. + ct->type = "int_plus"; + ct->arguments.clear(); + ct->arguments.push_back(Argument::VarRef(target_variable)); + ct->arguments.push_back( + Argument::VarRef(current_variables[current_variables.size() - 2])); + ct->arguments.push_back(Argument::VarRef(current_variables.back())); + target_variable = current_variables.back(); + current_variables.pop_back(); + + // We clean the first constraint too. + if (first_constraint != nullptr) { + first_constraint = nullptr; + } + } else { + current_variables.clear(); + target_variable = nullptr; + } + } + } + + // First pass. + for (Constraint* const ct : model->constraints()) { + if (ct->active && ct->type == "bool2int") { + PresolveBool2Int(ct); + } else if (ct->active && ct->type == "int2float") { + PresolveInt2Float(ct); + } else if (ct->active && ct->type == "int_lin_eq" && + ct->arguments[1].variables.size() == 3 && + ct->strong_propagation) { + PresolveStoreFlatteningMapping(ct); + } + } + if (!var_representative_map_.empty()) { + // Some new substitutions were introduced. Let's process them. + SubstituteEverywhere(model); + var_representative_map_.clear(); + var_representative_vector_.clear(); + } + + // Second pass. + for (Constraint* const ct : model->constraints()) { + if (ct->type == "array_int_element" || ct->type == "array_bool_element") { + PresolveSimplifyElement(ct); + } + } + + // Third pass: process objective with floating point coefficients. + Variable* float_objective_var = nullptr; + for (Variable* var : model->variables()) { + if (!var->active) continue; + if (var->domain.is_float) { + CHECK(float_objective_var == nullptr); + float_objective_var = var; + } + } + + Constraint* float_objective_ct = nullptr; + if (float_objective_var != nullptr) { + for (Constraint* ct : model->constraints()) { + if (!ct->active) continue; + if (ct->type == "float_lin_eq") { + CHECK(float_objective_ct == nullptr); + float_objective_ct = ct; + break; + } + } + } + + if (float_objective_ct != nullptr || float_objective_var != nullptr) { + CHECK(float_objective_ct != nullptr); + CHECK(float_objective_var != nullptr); + const int arity = float_objective_ct->arguments[0].Size(); + CHECK_EQ(float_objective_ct->arguments[1].variables[arity - 1], + float_objective_var); + CHECK_EQ(float_objective_ct->arguments[0].floats[arity - 1], -1.0); + for (int i = 0; i + 1 < arity; ++i) { + model->AddFloatingPointObjectiveTerm( + float_objective_ct->arguments[1].variables[i], + float_objective_ct->arguments[0].floats[i]); + } + model->SetFloatingPointObjectiveOffset( + -float_objective_ct->arguments[2].floats[0]); + model->ClearObjective(); + float_objective_var->active = false; + float_objective_ct->active = false; + } + + // Report presolve rules statistics. + if (!successful_rules_.empty()) { + for (const auto& rule : successful_rules_) { + if (rule.second == 1) { + SOLVER_LOG(logger_, " - rule '", rule.first, "' was applied 1 time"); + } else { + SOLVER_LOG(logger_, " - rule '", rule.first, "' was applied ", + rule.second, " times"); + } + } + } +} + +// ----- Substitution support ----- + +void Presolver::AddVariableSubstitution(Variable* from, Variable* to) { + CHECK(from != nullptr); + CHECK(to != nullptr); + // Apply the substitutions, if any. + from = FindRepresentativeOfVar(from); + to = FindRepresentativeOfVar(to); + if (to->temporary) { + // Let's switch to keep a non temporary as representative. + Variable* tmp = to; + to = from; + from = tmp; + } + if (from != to) { + CHECK(to->Merge(from->name, from->domain, from->temporary)); + from->active = false; + var_representative_map_[from] = to; + var_representative_vector_.push_back(from); + } +} + +Variable* Presolver::FindRepresentativeOfVar(Variable* var) { + if (var == nullptr) return nullptr; + Variable* start_var = var; + // First loop: find the top parent. + for (;;) { + const auto& it = var_representative_map_.find(var); + Variable* parent = it == var_representative_map_.end() ? var : it->second; + if (parent == var) break; + var = parent; + } + // Second loop: attach all the path to the top parent. + while (start_var != var) { + Variable* const parent = var_representative_map_[start_var]; + var_representative_map_[start_var] = var; + start_var = parent; + } + const auto& iter = var_representative_map_.find(var); + return iter == var_representative_map_.end() ? var : iter->second; +} + +void Presolver::SubstituteEverywhere(Model* model) { + // Rewrite the constraints. + for (Constraint* const ct : model->constraints()) { + if (ct != nullptr && ct->active) { + for (int i = 0; i < ct->arguments.size(); ++i) { + Argument& argument = ct->arguments[i]; + switch (argument.type) { + case Argument::VAR_REF: + case Argument::VAR_REF_ARRAY: { + for (int i = 0; i < argument.variables.size(); ++i) { + Variable* const old_var = argument.variables[i]; + Variable* const new_var = FindRepresentativeOfVar(old_var); + if (new_var != old_var) { + argument.variables[i] = new_var; + } + } + break; + } + default: { + } + } + } + } + } + // Rewrite the search. + for (Annotation* const ann : model->mutable_search_annotations()) { + SubstituteAnnotation(ann); + } + // Rewrite the output. + for (SolutionOutputSpecs* const output : model->mutable_output()) { + output->variable = FindRepresentativeOfVar(output->variable); + for (int i = 0; i < output->flat_variables.size(); ++i) { + output->flat_variables[i] = + FindRepresentativeOfVar(output->flat_variables[i]); + } + } + // Do not forget to merge domain that could have evolved asynchronously + // during presolve. + for (const auto& iter : var_representative_map_) { + iter.second->domain.IntersectWithDomain(iter.first->domain); + } + + // Change the objective variable. + Variable* const current_objective = model->objective(); + if (current_objective == nullptr) return; + Variable* const new_objective = FindRepresentativeOfVar(current_objective); + if (new_objective != current_objective) { + model->SetObjective(new_objective); + } +} + +void Presolver::SubstituteAnnotation(Annotation* ann) { + // TODO(user): Remove recursion. + switch (ann->type) { + case Annotation::ANNOTATION_LIST: + case Annotation::FUNCTION_CALL: { + for (int i = 0; i < ann->annotations.size(); ++i) { + SubstituteAnnotation(&ann->annotations[i]); + } + break; + } + case Annotation::VAR_REF: + case Annotation::VAR_REF_ARRAY: { + for (int i = 0; i < ann->variables.size(); ++i) { + ann->variables[i] = FindRepresentativeOfVar(ann->variables[i]); + } + break; + } + default: { + } + } +} + +} // namespace fz +} // namespace operations_research diff --git a/ortools/flatzinc/presolve.h b/ortools/flatzinc/presolve.h new file mode 100644 index 0000000000..38675968c7 --- /dev/null +++ b/ortools/flatzinc/presolve.h @@ -0,0 +1,115 @@ +// Copyright 2010-2024 Google LLC +// Licensed under the Apache License, Version 2.0 (the "License"); +// you may not use this file except in compliance with the License. +// You may obtain a copy of the License at +// +// http://www.apache.org/licenses/LICENSE-2.0 +// +// Unless required by applicable law or agreed to in writing, software +// distributed under the License is distributed on an "AS IS" BASIS, +// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +// See the License for the specific language governing permissions and +// limitations under the License. + +#ifndef OR_TOOLS_FLATZINC_PRESOLVE_H_ +#define OR_TOOLS_FLATZINC_PRESOLVE_H_ + +#include +#include +#include +#include +#include + +#include "absl/container/flat_hash_map.h" +#include "absl/container/flat_hash_set.h" +#include "absl/strings/match.h" +#include "ortools/base/hash.h" +#include "ortools/base/logging.h" +#include "ortools/base/types.h" +#include "ortools/flatzinc/model.h" +#include "ortools/util/logging.h" + +namespace operations_research { +namespace fz { +// The Presolver "pre-solves" a Model by applying some iterative +// transformations to it, which may simplify and/or shrink the model. +// +// TODO(user): Error reporting of unfeasible models. +class Presolver { + public: + explicit Presolver(SolverLogger* logger) : logger_(logger) {} + // Recursively apply all the pre-solve rules to the model, until exhaustion. + // The reduced model will: + // - Have some unused variables. + // - Have some unused constraints (marked as inactive). + // - Have some modified constraints (for example, they will no longer + // refer to unused variables). + void Run(Model* model); + + private: + // This struct stores the mapping of two index variables (of a 2D array; not + // included here) onto a single index variable (of the flattened 1D array). + // The original 2D array could be trimmed in the process; so we also need an + // offset. + // Eg. new_index_var = index_var1 * int_coeff + index_var2 + int_offset + struct Array2DIndexMapping { + Variable* variable1; + int64_t coefficient; + Variable* variable2; + int64_t offset; + Constraint* constraint; + + Array2DIndexMapping() + : variable1(nullptr), + coefficient(0), + variable2(nullptr), + offset(0), + constraint(nullptr) {} + Array2DIndexMapping(Variable* v1, int64_t c, Variable* v2, int64_t o, + Constraint* ct) + : variable1(v1), + coefficient(c), + variable2(v2), + offset(o), + constraint(ct) {} + }; + + // Substitution support. + void SubstituteEverywhere(Model* model); + void SubstituteAnnotation(Annotation* ann); + + // Presolve rules. + void PresolveBool2Int(Constraint* ct); + void PresolveInt2Float(Constraint* ct); + void PresolveStoreFlatteningMapping(Constraint* ct); + void PresolveSimplifyElement(Constraint* ct); + + // Helpers. + void UpdateRuleStats(const std::string& rule_name) { + successful_rules_[rule_name]++; + } + + // The presolver will discover some equivalence classes of variables [two + // variable are equivalent when replacing one by the other leads to the same + // logical model]. We will store them here, using a Union-find data structure. + // See http://en.wikipedia.org/wiki/Disjoint-set_data_structure. + // Note that the equivalence is directed. We prefer to replace all instances + // of 'from' with 'to', rather than the opposite. + void AddVariableSubstitution(Variable* from, Variable* to); + Variable* FindRepresentativeOfVar(Variable* var); + absl::flat_hash_map var_representative_map_; + std::vector var_representative_vector_; + + // Stores array2d_index_map_[z] = a * x + y + b. + absl::flat_hash_map array2d_index_map_; + + // Count applications of presolve rules. Use a sorted map for reporting + // purposes. + std::map successful_rules_; + + SolverLogger* logger_; +}; +} // namespace fz +} // namespace operations_research + +#endif // OR_TOOLS_FLATZINC_PRESOLVE_H_