diff --git a/src/math/lp/CMakeLists.txt b/src/math/lp/CMakeLists.txt index 38aaa4706..3bd1104c1 100644 --- a/src/math/lp/CMakeLists.txt +++ b/src/math/lp/CMakeLists.txt @@ -14,6 +14,7 @@ z3_add_component(lp indexed_vector.cpp int_branch.cpp int_cube.cpp + int_gcd_test.cpp int_solver.cpp lar_solver.cpp lar_core_solver.cpp diff --git a/src/math/lp/gomory.cpp b/src/math/lp/gomory.cpp index d8d8b48c1..7cc8e8f32 100644 --- a/src/math/lp/gomory.cpp +++ b/src/math/lp/gomory.cpp @@ -378,9 +378,12 @@ bool gomory::is_gomory_cut_target(const row_strip& row) { int gomory::find_basic_var() { int result = -1; unsigned n = 0; - bool boxed = false; unsigned min_row_size = UINT_MAX; +#if 0 + bool boxed = false; mpq min_range; +#endif + // Prefer smaller row size // Prefer boxed to non-boxed diff --git a/src/math/lp/int_gcd_test.cpp b/src/math/lp/int_gcd_test.cpp new file mode 100644 index 000000000..b5a9e8c05 --- /dev/null +++ b/src/math/lp/int_gcd_test.cpp @@ -0,0 +1,192 @@ +/*++ +Copyright (c) 2020 Microsoft Corporation + +Module Name: + + int_gcd_test.cpp + +Abstract: + + Gcd_Test heuristic + +Author: + Nikolaj Bjorner (nbjorner) + Lev Nachmanson (levnach) + +Revision History: +--*/ + +#include "math/lp/int_solver.h" +#include "math/lp/lar_solver.h" +#include "math/lp/int_gcd_test.h" + +namespace lp { + + int_gcd_test::int_gcd_test(int_solver& lia): lia(lia), lra(lia.lra) {} + + lia_move int_gcd_test::operator()() { + lia.settings().stats().m_gcd_calls++; + TRACE("int_solver", tout << "gcd-test " << lia.settings().stats().m_gcd_calls << "\n";); + if (gcd_test()) { + return lia_move::undef; + } + else { + lia.settings().stats().m_gcd_conflicts++; + TRACE("gcd_test", tout << "gcd conflict\n";); + return lia_move::conflict; + } + } + + bool int_gcd_test::gcd_test() { + auto & A = lra.A_r(); // getting the matrix + for (unsigned i = 0; i < A.row_count(); i++) + if (!gcd_test_for_row(A, i)) + return false; + return true; + } + + static mpq get_denominators_lcm(const row_strip & row) { + mpq r(1); + for (auto & c : row) { + r = lcm(r, denominator(c.coeff())); + } + return r; + } + + bool int_gcd_test::gcd_test_for_row(static_matrix> & A, unsigned i) { + auto const& row = A.m_rows[i]; + auto & lcs = lra.m_mpq_lar_core_solver; + unsigned basic_var = lcs.m_r_basis[i]; + + if (!lia.column_is_int(basic_var) || lia.get_value(basic_var).is_int()) + return true; + mpq lcm_den = get_denominators_lcm(row); + mpq consts(0); + mpq gcds(0); + mpq least_coeff(0); + bool least_coeff_is_bounded = false; + unsigned j; + for (auto &c : A.m_rows[i]) { + j = c.var(); + const mpq& a = c.coeff(); + if (lra.column_is_fixed(j)) { + mpq aux = lcm_den * a; + consts += aux * lra.column_lower_bound(j).x; + } + else if (lra.column_is_real(j)) { + return true; + } + else if (gcds.is_zero()) { + gcds = abs(lcm_den * a); + least_coeff = gcds; + least_coeff_is_bounded = lra.column_is_bounded(j); + } + else { + mpq aux = abs(lcm_den * a); + gcds = gcd(gcds, aux); + if (aux < least_coeff) { + least_coeff = aux; + least_coeff_is_bounded = lra.column_is_bounded(j); + } + else if (least_coeff_is_bounded && aux == least_coeff) { + least_coeff_is_bounded = lra.column_is_bounded(j); + } + } + SASSERT(gcds.is_int()); + SASSERT(least_coeff.is_int()); + TRACE("gcd_test_bug", tout << "coeff: " << a << ", gcds: " << gcds + << " least_coeff: " << least_coeff << " consts: " << consts << "\n";); + + } + + if (gcds.is_zero()) { + // All variables are fixed. + // This theory guarantees that the assignment satisfies each row, and + // fixed integer variables are assigned to integer values. + return true; + } + + if (!(consts / gcds).is_int()) { + TRACE("gcd_test", tout << "row failed the GCD test:\n"; lia.display_row_info(tout, i);); + fill_explanation_from_fixed_columns(A.m_rows[i]); + return false; + } + + if (least_coeff.is_one() && !least_coeff_is_bounded) { + SASSERT(gcds.is_one()); + return true; + } + + if (least_coeff_is_bounded) { + return ext_gcd_test(A.m_rows[i], least_coeff, lcm_den, consts); + } + return true; + } + + bool int_gcd_test::ext_gcd_test(const row_strip & row, + mpq const & least_coeff, + mpq const & lcm_den, + mpq const & consts) { + TRACE("ext_gcd_test", tout << "row = "; lra.print_row(row, tout);); + mpq gcds(0); + mpq l(consts); + mpq u(consts); + + mpq a; + unsigned j; + for (const auto & c : row) { + j = c.var(); + TRACE("ext_gcd_test", tout << "col = "; lra.m_mpq_lar_core_solver.m_r_solver.print_column_bound_info(j, tout);); + const mpq & a = c.coeff(); + if (lra.column_is_fixed(j)) + continue; + SASSERT(!lra.column_is_real(j)); + mpq ncoeff = lcm_den * a; + SASSERT(ncoeff.is_int()); + mpq abs_ncoeff = abs(ncoeff); + if (abs_ncoeff == least_coeff) { + SASSERT(lra.column_is_bounded(j)); + if (ncoeff.is_pos()) { + // l += ncoeff * lra.column_lower_bound(j).x; + l.addmul(ncoeff, lra.column_lower_bound(j).x); + // u += ncoeff * lra.column_upper_bound(j).x; + u.addmul(ncoeff, lra.column_upper_bound(j).x); + } + else { + // l += ncoeff * upper_bound(j).get_rational(); + l.addmul(ncoeff, lra.column_upper_bound(j).x); + // u += ncoeff * lower_bound(j).get_rational(); + u.addmul(ncoeff, lra.column_lower_bound(j).x); + } + lia.add_to_explanation_from_fixed_or_boxed_column(j); + } + else if (gcds.is_zero()) { + gcds = abs_ncoeff; + } + else { + gcds = gcd(gcds, abs_ncoeff); + } + SASSERT(gcds.is_int()); + } + + if (gcds.is_zero()) { + return true; + } + + mpq l1 = ceil(l/gcds); + mpq u1 = floor(u/gcds); + + if (u1 < l1) { + fill_explanation_from_fixed_columns(row); + return false; + } + return true; + } + + void int_gcd_test::fill_explanation_from_fixed_columns(const row_strip & row) { + for (const auto & c : row) { + if (lra.column_is_fixed(c.var())) + lia.add_to_explanation_from_fixed_or_boxed_column(c.var()); + } + } +} diff --git a/src/math/lp/int_gcd_test.h b/src/math/lp/int_gcd_test.h new file mode 100644 index 000000000..00c97ad6f --- /dev/null +++ b/src/math/lp/int_gcd_test.h @@ -0,0 +1,42 @@ +/*++ +Copyright (c) 2020 Microsoft Corporation + +Module Name: + + int_gcd_test.h + +Abstract: + + Gcd_Test heuristic + +Author: + Nikolaj Bjorner (nbjorner) + Lev Nachmanson (levnach) + +Revision History: +--*/ +#pragma once + +#include "math/lp/lia_move.h" + +namespace lp { + class int_solver; + class lar_solver; + class int_gcd_test { + class int_solver& lia; + class lar_solver& lra; + + bool gcd_test(); + bool gcd_test_for_row(static_matrix> & A, unsigned i); + bool ext_gcd_test(const row_strip & row, + mpq const & least_coeff, + mpq const & lcm_den, + mpq const & consts); + void fill_explanation_from_fixed_columns(const row_strip & row); + + public: + int_gcd_test(int_solver& lia); + ~int_gcd_test() {} + lia_move operator()(); + }; +} diff --git a/src/math/lp/int_solver.cpp b/src/math/lp/int_solver.cpp index b61a4feb4..7aa6911ef 100644 --- a/src/math/lp/int_solver.cpp +++ b/src/math/lp/int_solver.cpp @@ -9,10 +9,60 @@ #include #include "math/lp/monic.h" #include "math/lp/gomory.h" -#include "math/lp/int_cube.h" #include "math/lp/int_branch.h" +#include "math/lp/int_gcd_test.h" +#include "math/lp/int_cube.h" namespace lp { +// this will allow to enable and disable tracking of the pivot rows +struct check_return_helper { + lar_solver& lra; + bool m_track_pivoted_rows; + check_return_helper(lar_solver& ls) : + lra(ls), + m_track_pivoted_rows(lra.get_track_pivoted_rows()) + { + TRACE("pivoted_rows", tout << "pivoted rows = " << lra.m_mpq_lar_core_solver.m_r_solver.m_pivoted_rows->size() << std::endl;); + lra.set_track_pivoted_rows(false); + } + ~check_return_helper() { + TRACE("pivoted_rows", tout << "pivoted rows = " << lra.m_mpq_lar_core_solver.m_r_solver.m_pivoted_rows->size() << std::endl;); + lra.set_track_pivoted_rows(m_track_pivoted_rows); + } +}; + +lia_move int_solver::check(lp::explanation * e) { + SASSERT(lra.ax_is_correct()); + if (!has_inf_int()) return lia_move::sat; + + m_t.clear(); + m_k.reset(); + m_ex = e; + m_ex->clear(); + m_upper = false; + lia_move r = lia_move::undef; + + gomory gc(*this); + int_cube cube(*this); + int_branch branch(*this); + int_gcd_test gcd(*this); + + if (should_run_gcd_test()) r = gcd(); + + check_return_helper pc(lra); + + if (settings().m_int_pivot_fixed_vars_from_basis) + lra.pivot_fixed_vars_from_basis(); + + + if (r == lia_move::undef) r = patch_nbasic_columns(); + ++m_number_of_calls; + if (r == lia_move::undef && should_find_cube()) r = cube(); + if (r == lia_move::undef && should_hnf_cut()) r = hnf_cut(); + if (r == lia_move::undef && should_gomory_cut()) r = gc(m_t, m_k, m_ex, m_upper); + if (r == lia_move::undef) r = branch(); + return r; +} std::ostream& int_solver::display_inf_rows(std::ostream& out) const { unsigned num = lra.A_r().column_count(); @@ -98,69 +148,19 @@ unsigned int_solver::column_count() const { return lra.column_count(); } -// this will allow to enable and disable tracking of the pivot rows -struct check_return_helper { - lar_solver& lra; - bool m_track_pivoted_rows; - check_return_helper(lar_solver& ls) : - lra(ls), - m_track_pivoted_rows(lra.get_track_pivoted_rows()) - { - TRACE("pivoted_rows", tout << "pivoted rows = " << lra.m_mpq_lar_core_solver.m_r_solver.m_pivoted_rows->size() << std::endl;); - lra.set_track_pivoted_rows(false); - } - ~check_return_helper() { - TRACE("pivoted_rows", tout << "pivoted rows = " << lra.m_mpq_lar_core_solver.m_r_solver.m_pivoted_rows->size() << std::endl;); - lra.set_track_pivoted_rows(m_track_pivoted_rows); - } -}; bool int_solver::should_find_cube() { return m_number_of_calls % settings().m_int_find_cube_period == 0; } -lia_move int_solver::find_cube() { - int_cube ic(*this); - if (should_find_cube()) { - return ic(); - } - else { - return lia_move::undef; - } -} - bool int_solver::should_run_gcd_test() { return settings().m_int_run_gcd_test; } -lia_move int_solver::run_gcd_test() { - if (should_run_gcd_test()) { - settings().stats().m_gcd_calls++; - TRACE("int_solver", tout << "gcd-test " << settings().stats().m_gcd_calls << "\n";); - if (!gcd_test()) { - settings().stats().m_gcd_conflicts++; - TRACE("gcd_test", tout << "gcd conflict\n";); - return lia_move::conflict; - } - } - return lia_move::undef; -} - bool int_solver::should_gomory_cut() { return m_number_of_calls % settings().m_int_gomory_cut_period == 0; } -lia_move int_solver::gomory_cut() { - TRACE("int_solver", tout << "gomory " << m_number_of_calls << "\n";); - gomory gc(*this); - if (should_gomory_cut()) { - return gc(m_t, m_k, m_ex, m_upper); - } - else { - return lia_move::undef; - } -} - void int_solver::try_add_term_to_A_for_hnf(unsigned i) { mpq rs; const lar_term* t = lra.terms()[i]; @@ -234,56 +234,16 @@ bool int_solver::should_hnf_cut() { } lia_move int_solver::hnf_cut() { - lia_move r = lia_move::undef; - if (should_hnf_cut()) { - r = make_hnf_cut(); - if (r == lia_move::undef) { - m_hnf_cut_period *= 2; - } - else { - m_hnf_cut_period = settings().hnf_cut_period(); - } + lia_move r = make_hnf_cut(); + if (r == lia_move::undef) { + m_hnf_cut_period *= 2; + } + else { + m_hnf_cut_period = settings().hnf_cut_period(); } return r; } -lia_move int_solver::check(lp::explanation * e) { - SASSERT(lra.ax_is_correct()); - if (!has_inf_int()) return lia_move::sat; - -#define CHECK_RET(fn) \ - r = fn; \ - if (r != lia_move::undef) { TRACE("int_solver", tout << #fn << "\n";); return r; } - - m_t.clear(); - m_k.reset(); - m_ex = e; - m_ex->clear(); - m_upper = false; - lia_move r; - - CHECK_RET(run_gcd_test()); - - check_return_helper pc(lra); - - if (settings().m_int_pivot_fixed_vars_from_basis) - lra.pivot_fixed_vars_from_basis(); - - CHECK_RET(patch_nbasic_columns()); - ++m_number_of_calls; - CHECK_RET(find_cube()); - CHECK_RET(hnf_cut()); - CHECK_RET(gomory_cut()); - CHECK_RET(branch_or_sat()); - return r; -} - -lia_move int_solver::branch_or_sat() { - int_branch b(*this); - return b(); -} - - void int_solver::set_value_for_nbasic_column_ignore_old_values(unsigned j, const impq & new_val) { lp_assert(!is_base(j)); auto & x = lra.m_mpq_lar_core_solver.m_r_x[j]; @@ -351,166 +311,14 @@ lia_move int_solver::patch_nbasic_columns() { return lia_move::undef; } -mpq get_denominators_lcm(const row_strip & row) { - mpq r(1); - for (auto & c : row) { - r = lcm(r, denominator(c.coeff())); - } - return r; -} - -bool int_solver::gcd_test_for_row(static_matrix> & A, unsigned i) { - auto const& row = A.m_rows[i]; - auto & lcs = lra.m_mpq_lar_core_solver; - unsigned basic_var = lcs.m_r_basis[i]; - - if (!column_is_int(basic_var) || get_value(basic_var).is_int()) - return true; - mpq lcm_den = get_denominators_lcm(row); - mpq consts(0); - mpq gcds(0); - mpq least_coeff(0); - bool least_coeff_is_bounded = false; - unsigned j; - for (auto &c : A.m_rows[i]) { - j = c.var(); - const mpq& a = c.coeff(); - if (lra.column_is_fixed(j)) { - mpq aux = lcm_den * a; - consts += aux * lra.column_lower_bound(j).x; - } - else if (lra.column_is_real(j)) { - return true; - } - else if (gcds.is_zero()) { - gcds = abs(lcm_den * a); - least_coeff = gcds; - least_coeff_is_bounded = lra.column_is_bounded(j); - } - else { - mpq aux = abs(lcm_den * a); - gcds = gcd(gcds, aux); - if (aux < least_coeff) { - least_coeff = aux; - least_coeff_is_bounded = lra.column_is_bounded(j); - } - else if (least_coeff_is_bounded && aux == least_coeff) { - least_coeff_is_bounded = lra.column_is_bounded(j); - } - } - SASSERT(gcds.is_int()); - SASSERT(least_coeff.is_int()); - TRACE("gcd_test_bug", tout << "coeff: " << a << ", gcds: " << gcds - << " least_coeff: " << least_coeff << " consts: " << consts << "\n";); - - } - - if (gcds.is_zero()) { - // All variables are fixed. - // This theory guarantees that the assignment satisfies each row, and - // fixed integer variables are assigned to integer values. - return true; - } - - if (!(consts / gcds).is_int()) { - TRACE("gcd_test", tout << "row failed the GCD test:\n"; display_row_info(tout, i);); - fill_explanation_from_fixed_columns(A.m_rows[i]); - return false; - } - - if (least_coeff.is_one() && !least_coeff_is_bounded) { - SASSERT(gcds.is_one()); - return true; - } - - if (least_coeff_is_bounded) { - return ext_gcd_test(A.m_rows[i], least_coeff, lcm_den, consts); - } - return true; -} - - +// TBD: move to gcd-test void int_solver::add_to_explanation_from_fixed_or_boxed_column(unsigned j) { constraint_index lc, uc; lra.get_bound_constraint_witnesses_for_column(j, lc, uc); m_ex->push_justification(lc); m_ex->push_justification(uc); } -void int_solver::fill_explanation_from_fixed_columns(const row_strip & row) { - for (const auto & c : row) { - if (!lra.column_is_fixed(c.var())) - continue; - add_to_explanation_from_fixed_or_boxed_column(c.var()); - } -} - -bool int_solver::gcd_test() { - auto & A = lra.A_r(); // getting the matrix - for (unsigned i = 0; i < A.row_count(); i++) - if (!gcd_test_for_row(A, i)) - return false; - return true; -} -bool int_solver::ext_gcd_test(const row_strip & row, - mpq const & least_coeff, - mpq const & lcm_den, - mpq const & consts) { - TRACE("ext_gcd_test", tout << "row = "; lra.print_row(row, tout);); - mpq gcds(0); - mpq l(consts); - mpq u(consts); - - mpq a; - unsigned j; - for (const auto & c : row) { - j = c.var(); - TRACE("ext_gcd_test", tout << "col = "; lra.m_mpq_lar_core_solver.m_r_solver.print_column_bound_info(j, tout);); - const mpq & a = c.coeff(); - if (lra.column_is_fixed(j)) - continue; - SASSERT(!lra.column_is_real(j)); - mpq ncoeff = lcm_den * a; - SASSERT(ncoeff.is_int()); - mpq abs_ncoeff = abs(ncoeff); - if (abs_ncoeff == least_coeff) { - SASSERT(lra.column_is_bounded(j)); - if (ncoeff.is_pos()) { - // l += ncoeff * lra.column_lower_bound(j).x; - l.addmul(ncoeff, lra.column_lower_bound(j).x); - // u += ncoeff * lra.column_upper_bound(j).x; - u.addmul(ncoeff, lra.column_upper_bound(j).x); - } - else { - // l += ncoeff * upper_bound(j).get_rational(); - l.addmul(ncoeff, lra.column_upper_bound(j).x); - // u += ncoeff * lower_bound(j).get_rational(); - u.addmul(ncoeff, lra.column_lower_bound(j).x); - } - add_to_explanation_from_fixed_or_boxed_column(j); - } - else if (gcds.is_zero()) { - gcds = abs_ncoeff; - } - else { - gcds = gcd(gcds, abs_ncoeff); - } - SASSERT(gcds.is_int()); - } - - if (gcds.is_zero()) { - return true; - } - - mpq l1 = ceil(l/gcds); - mpq u1 = floor(u/gcds); - - if (u1 < l1) { - fill_explanation_from_fixed_columns(row); - return false; - } - return true; -} int_solver::int_solver(lar_solver& lar_slv) : lra(lar_slv), diff --git a/src/math/lp/int_solver.h b/src/math/lp/int_solver.h index 42de3cc2e..1cfe4a74b 100644 --- a/src/math/lp/int_solver.h +++ b/src/math/lp/int_solver.h @@ -38,6 +38,7 @@ class int_solver { friend class gomory; friend class int_cube; friend class int_branch; + friend class int_gcd_test; public: // fields lar_solver& lra; @@ -107,7 +108,6 @@ private: bool should_gomory_cut(); bool should_hnf_cut(); - lia_move branch_or_sat(); lp_settings& settings(); const lp_settings& settings() const; bool at_bound(unsigned j) const; @@ -128,12 +128,10 @@ private: bool has_inf_int() const; public: bool is_term(unsigned j) const; - lia_move find_cube(); unsigned column_count() const; bool all_columns_are_bounded() const; void find_feasible_solution(); lia_move run_gcd_test(); - lia_move gomory_cut(); lia_move hnf_cut(); lia_move make_hnf_cut(); bool init_terms_for_hnf_cut();