mirrors/z3
3
0
Fork 0
mirror of https://github.com/Z3Prover/z3 synced 2025-04-15 13:28:47 +00:00
Code Activity

no calling cut_solver when there are bounded columns

Browse source 
Signed-off-by: Lev Nachmanson <levnach@hotmail.com>

use special bounds inf find_cube for x+y, x-y

Signed-off-by: Lev Nachmanson <levnach@hotmail.com>

bug fixes in column patching, add stats to patching, restructure int_solver

Signed-off-by: Lev Nachmanson <levnach@hotmail.com>

comment out m_old_values from int_solver

Signed-off-by: Lev Nachmanson <levnach@hotmail.com>

avoid calling pivot_fixed_vars_from_basis() in int_solver

Signed-off-by: Lev Nachmanson <levnach@hotmail.com>

fix the return value from path_nbasic_columns

Signed-off-by: Lev Nachmanson <levnach@hotmail.com>

fix the return value from path_nbasic_columns

Signed-off-by: Lev Nachmanson <levnach@hotmail.com>

work in patch_columns

Signed-off-by: Lev Nachmanson <levnach@hotmail.com>

work on int_solver check()

Signed-off-by: Lev Nachmanson <levnach@hotmail.com>

exit from find_free_interval() when l >= u

Signed-off-by: Lev Nachmanson <levnach@hotmail.com>

experiment with branching on nbasic columns

Signed-off-by: Lev Nachmanson <levnach@hotmail.com>

remove m_old_values

Signed-off-by: Lev Nachmanson <levnach@hotmail.com>

add rounding to patch_columns

Signed-off-by: Lev Nachmanson <levnach@hotmail.com>

qflia

Signed-off-by: Lev Nachmanson <levnach@hotmail.com>

patch all columns, round non-patched, branch or basic columns

Signed-off-by: Lev Nachmanson <levnach@hotmail.com>

refactor int_solver::check()

Signed-off-by: Lev Nachmanson <levnach@hotmail.com>

restore  move_non_basic_columns_to_bounds() after a failure in find_cube()

Signed-off-by: Lev Nachmanson <levnach@hotmail.com>

optimize gomory cuts search

Signed-off-by: Lev Nachmanson <levnach@hotmail.com>

produce gomory cuts without moving columns to bounds

Signed-off-by: Lev Nachmanson <levnach@hotmail.com>

call find_feasible_solution() after moving columns

Signed-off-by: Lev Nachmanson <levnach@hotmail.com>

alway move colums to bounds before gomory cut

Signed-off-by: Lev Nachmanson <levnach@hotmail.com>

merge from best branch

Signed-off-by: Lev Nachmanson <levnach@hotmail.com>
This commit is contained in:
Lev Nachmanson 2018-04-03 17:29:58 -07:00
parent 7e82ab595e
commit c04bcb411d
13 changed files with 341 additions and 276 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

2
src/smt/theory_arith.h
View file

@ -42,7 +42,7 @@ Revision History:
namespace smt { namespace smt {
struct theory_arith_stats { struct theory_arith_stats {
unsigned m_conflicts, m_add_rows, m_pivots, m_diseq_cs, m_gomory_cuts, m_branches, m_gcd_tests; unsigned m_conflicts, m_add_rows, m_pivots, m_diseq_cs, m_gomory_cuts, m_branches, m_gcd_tests, m_patches, m_patches_succ;
unsigned m_assert_lower, m_assert_upper, m_assert_diseq, m_core2th_eqs, m_core2th_diseqs; unsigned m_assert_lower, m_assert_upper, m_assert_diseq, m_core2th_eqs, m_core2th_diseqs;
unsigned m_th2core_eqs, m_th2core_diseqs, m_bound_props, m_offset_eqs, m_fixed_eqs, m_offline_eqs; unsigned m_th2core_eqs, m_th2core_diseqs, m_bound_props, m_offset_eqs, m_fixed_eqs, m_offline_eqs;
unsigned m_max_min; unsigned m_max_min;

3
src/smt/theory_arith_int.h
View file

@ -1335,7 +1335,7 @@ namespace smt {
} }
} }
}); });
m_stats.m_patches++;
patch_int_infeasible_vars(); patch_int_infeasible_vars();
fix_non_base_vars(); fix_non_base_vars();
@ -1368,6 +1368,7 @@ namespace smt {
theory_var int_var = find_infeasible_int_base_var(); theory_var int_var = find_infeasible_int_base_var();
if (int_var == null_theory_var) { if (int_var == null_theory_var) {
m_stats.m_patches_succ++;
TRACE("arith_int_incomp", tout << "FC_DONE 2...\n"; display(tout);); TRACE("arith_int_incomp", tout << "FC_DONE 2...\n"; display(tout););
return m_liberal_final_check || !m_changed_assignment ? FC_DONE : FC_CONTINUE; return m_liberal_final_check || !m_changed_assignment ? FC_DONE : FC_CONTINUE;
} }

2
src/smt/theory_arith_pp.h
View file

@ -38,6 +38,8 @@ namespace smt {
st.update("arith gcd tests", m_stats.m_gcd_tests); st.update("arith gcd tests", m_stats.m_gcd_tests);
st.update("arith ineq splits", m_stats.m_branches); st.update("arith ineq splits", m_stats.m_branches);
st.update("arith gomory cuts", m_stats.m_gomory_cuts); st.update("arith gomory cuts", m_stats.m_gomory_cuts);
st.update("arith patches", m_stats.m_patches);
st.update("arith patches_succ", m_stats.m_patches_succ);
st.update("arith max-min", m_stats.m_max_min); st.update("arith max-min", m_stats.m_max_min);
st.update("arith grobner", m_stats.m_gb_compute_basis); st.update("arith grobner", m_stats.m_gb_compute_basis);
st.update("arith pseudo nonlinear", m_stats.m_nl_linear); st.update("arith pseudo nonlinear", m_stats.m_nl_linear);

16
src/smt/theory_lra.cpp
View file

@ -306,9 +306,11 @@ class theory_lra::imp {
reset_variable_values(); reset_variable_values();
m_solver->settings().bound_propagation() = BP_NONE != propagation_mode(); m_solver->settings().bound_propagation() = BP_NONE != propagation_mode();
m_solver->set_track_pivoted_rows(lp.bprop_on_pivoted_rows()); m_solver->set_track_pivoted_rows(lp.bprop_on_pivoted_rows());
m_solver->settings().m_int_branch_cut_gomory_threshold = ctx().get_fparams().m_arith_branch_cut_ratio; m_solver->settings().m_int_gomory_cut_period = ctx().get_fparams().m_arith_branch_cut_ratio;
m_solver->settings().m_int_branch_cut_solver = std::max(8u, ctx().get_fparams().m_arith_branch_cut_ratio); m_solver->settings().m_int_cuts_etc_period = ctx().get_fparams().m_arith_branch_cut_ratio;
m_solver->settings().m_run_gcd_test = ctx().get_fparams().m_arith_gcd_test; m_solver->settings().m_int_cut_solver_period = std::max(8u, ctx().get_fparams().m_arith_branch_cut_ratio);
m_solver->settings().m_int_run_gcd_test = ctx().get_fparams().m_arith_gcd_test;
m_solver->settings().set_random_seed(ctx().get_fparams().m_random_seed); m_solver->settings().set_random_seed(ctx().get_fparams().m_random_seed);
//m_solver->settings().set_ostream(0); //m_solver->settings().set_ostream(0);
m_lia = alloc(lp::int_solver, m_solver.get()); m_lia = alloc(lp::int_solver, m_solver.get());
@ -1286,7 +1288,7 @@ public:
lp::explanation ex; // TBD, this should be streamlined accross different explanations lp::explanation ex; // TBD, this should be streamlined accross different explanations
bool upper; bool upper;
switch(m_lia->check(term, k, ex, upper)) { switch(m_lia->check(term, k, ex, upper)) {
case lp::lia_move::ok: case lp::lia_move::sat:
return l_true; return l_true;
case lp::lia_move::branch: { case lp::lia_move::branch: {
app_ref b = mk_bound(term, k, !upper); app_ref b = mk_bound(term, k, !upper);
@ -1321,8 +1323,8 @@ public:
m_explanation = ex.m_explanation; m_explanation = ex.m_explanation;
set_conflict1(); set_conflict1();
return l_false; return l_false;
case lp::lia_move::give_up: case lp::lia_move::undef:
TRACE("arith", tout << "lia giveup\n";); TRACE("arith", tout << "lia undef\n";);
return l_undef; return l_undef;
default: default:
UNREACHABLE(); UNREACHABLE();
@ -2818,6 +2820,8 @@ public:
st.update("gcd-conflict", m_solver->settings().st().m_gcd_conflicts); st.update("gcd-conflict", m_solver->settings().st().m_gcd_conflicts);
st.update("cube-calls", m_solver->settings().st().m_cube_calls); st.update("cube-calls", m_solver->settings().st().m_cube_calls);
st.update("cube-success", m_solver->settings().st().m_cube_success); st.update("cube-success", m_solver->settings().st().m_cube_success);
st.update("arith-patches", m_solver->settings().st().m_patches);
st.update("arith-patches-success", m_solver->settings().st().m_patches_success);
} }
}; };

14
src/util/lp/bound_analyzer_on_row.h
View file

@ -116,7 +116,7 @@ public :
unsigned j; unsigned j;
void analyze() { void analyze() {
for (auto c : m_row) { for (const auto & c : m_row) {
if ((m_column_of_l == -2) && (m_column_of_u == -2)) if ((m_column_of_l == -2) && (m_column_of_u == -2))
break; break;
analyze_bound_on_var_on_coeff(c.var(), c.coeff()); analyze_bound_on_var_on_coeff(c.var(), c.coeff());
@ -225,7 +225,7 @@ public :
int strict = 0; int strict = 0;
mpq total; mpq total;
lp_assert(is_zero(total)); lp_assert(is_zero(total));
for (auto p : m_row) { for (const auto& p : m_row) {
bool str; bool str;
total -= monoid_min(p.coeff(), p.var(), str); total -= monoid_min(p.coeff(), p.var(), str);
if (str) if (str)
@ -233,7 +233,7 @@ public :
} }
for (auto p : m_row) { for (const auto &p : m_row) {
bool str; bool str;
bool a_is_pos = is_pos(p.coeff()); bool a_is_pos = is_pos(p.coeff());
mpq bound = total / p.coeff() + monoid_min_no_mult(a_is_pos, p.var(), str); mpq bound = total / p.coeff() + monoid_min_no_mult(a_is_pos, p.var(), str);
@ -250,14 +250,14 @@ public :
int strict = 0; int strict = 0;
mpq total; mpq total;
lp_assert(is_zero(total)); lp_assert(is_zero(total));
for (auto p : m_row) { for (const auto &p : m_row) {
bool str; bool str;
total -= monoid_max(p.coeff(), p.var(), str); total -= monoid_max(p.coeff(), p.var(), str);
if (str) if (str)
strict++; strict++;
} }
for (auto p : m_row) { for (const auto& p : m_row) {
bool str; bool str;
bool a_is_pos = is_pos(p.coeff()); bool a_is_pos = is_pos(p.coeff());
mpq bound = total / p.coeff() + monoid_max_no_mult(a_is_pos, p.var(), str); mpq bound = total / p.coeff() + monoid_max_no_mult(a_is_pos, p.var(), str);
@ -279,7 +279,7 @@ public :
unsigned j; unsigned j;
mpq bound = -m_rs.x; mpq bound = -m_rs.x;
bool strict = false; bool strict = false;
for (auto p : m_row) { for (const auto& p : m_row) {
j = p.var(); j = p.var();
if (j == static_cast<unsigned>(m_column_of_u)) { if (j == static_cast<unsigned>(m_column_of_u)) {
u_coeff = p.coeff(); u_coeff = p.coeff();
@ -308,7 +308,7 @@ public :
unsigned j; unsigned j;
mpq bound = -m_rs.x; mpq bound = -m_rs.x;
bool strict = false; bool strict = false;
for (auto p : m_row) { for (const auto &p : m_row) {
j = p.var(); j = p.var();
if (j == static_cast<unsigned>(m_column_of_l)) { if (j == static_cast<unsigned>(m_column_of_l)) {
l_coeff = p.coeff(); l_coeff = p.coeff();

2
src/util/lp/column_namer.h
View file

@ -27,7 +27,7 @@ public:
template <typename T> template <typename T>
void print_row(const row_strip<T> & row, std::ostream & out) const { void print_row(const row_strip<T> & row, std::ostream & out) const {
vector<std::pair<T, unsigned>> coeff; vector<std::pair<T, unsigned>> coeff;
for (auto p : row) { for (auto & p : row) {
coeff.push_back(std::make_pair(p.coeff(), p.var())); coeff.push_back(std::make_pair(p.coeff(), p.var()));
} }
print_linear_combination_of_column_indices(coeff, out); print_linear_combination_of_column_indices(coeff, out);

2
src/util/lp/cut_solver.h
View file

@ -1715,7 +1715,7 @@ public:
m_cancelled = true; m_cancelled = true;
return true; return true;
} }
unsigned bound = m_asserts.size() * 200 / (1 + m_settings.m_int_branch_cut_solver); unsigned bound = m_asserts.size() * 200 / (1 + m_settings.m_int_cut_solver_period);
if (m_trail.size() > bound || m_number_of_conflicts > bound) { if (m_trail.size() > bound || m_number_of_conflicts > bound) {
m_cancelled = true; m_cancelled = true;
return true; return true;

457
src/util/lp/int_solver.cpp
View file

@ -10,18 +10,6 @@
#include <utility> #include <utility>
namespace lp { namespace lp {
void int_solver::failed() {
auto & lcs = m_lar_solver->m_mpq_lar_core_solver;
for (unsigned j : m_old_values_set.m_index) {
lcs.m_r_x[j] = m_old_values_data[j];
lp_assert(lcs.m_r_solver.column_is_feasible(j));
lcs.m_r_solver.remove_column_from_inf_set(j);
}
lp_assert(lcs.m_r_solver.calc_current_x_is_feasible_include_non_basis());
lp_assert(lcs.m_r_solver.current_x_is_feasible());
m_old_values_set.clear();
}
void int_solver::trace_inf_rows() const { void int_solver::trace_inf_rows() const {
unsigned num = m_lar_solver->A_r().column_count(); unsigned num = m_lar_solver->A_r().column_count();
@ -44,6 +32,13 @@ void int_solver::trace_inf_rows() const {
); );
} }
bool int_solver::all_columns_are_bounded() const {
for (unsigned j = 0; j < m_lar_solver->column_count(); j++)
if (m_lar_solver->column_is_bounded(j) == false)
return false;
return true;
}
bool int_solver::has_inf_int() const { bool int_solver::has_inf_int() const {
return m_lar_solver->has_inf_int(); return m_lar_solver->has_inf_int();
} }
@ -71,6 +66,14 @@ int int_solver::get_kth_inf_int(unsigned k) const {
return -1; return -1;
} }
int int_solver::find_inf_int_nbasis_column() const {
for (unsigned j : m_lar_solver->r_nbasis())
if (! column_is_int_inf(j) )
return j;
return -1;
}
int int_solver::find_inf_int_boxed_base_column_with_smallest_range(unsigned & inf_int_count) { int int_solver::find_inf_int_boxed_base_column_with_smallest_range(unsigned & inf_int_count) {
inf_int_count = 0; inf_int_count = 0;
int result = -1; int result = -1;
@ -117,9 +120,11 @@ int int_solver::find_inf_int_boxed_base_column_with_smallest_range(unsigned & in
bool int_solver::is_gomory_cut_target(const row_strip<mpq>& row) { bool int_solver::is_gomory_cut_target(const row_strip<mpq>& row) {
// All non base variables must be at their bounds and assigned to rationals (that is, infinitesimals are not allowed). // All non base variables must be at their bounds and assigned to rationals (that is, infinitesimals are not allowed).
unsigned j; unsigned j;
for (auto p : row) { for (const auto & p : row) {
j = p.var(); j = p.var();
if (is_base(j)) continue; if (is_base(j)) continue;
if (!at_bound(j))
return false;
if (!is_zero(get_value(j).y)) { if (!is_zero(get_value(j).y)) {
TRACE("gomory_cut", tout << "row is not gomory cut target:\n"; TRACE("gomory_cut", tout << "row is not gomory cut target:\n";
display_column(tout, j); display_column(tout, j);
@ -131,7 +136,7 @@ bool int_solver::is_gomory_cut_target(const row_strip<mpq>& row) {
} }
void int_solver::real_case_in_gomory_cut(const mpq & a, unsigned x_j, mpq & k, lar_term& pol, explanation & expl, const mpq& f_0, const mpq& one_minus_f_0) { void int_solver::real_case_in_gomory_cut(const mpq & a, unsigned x_j, const mpq& f_0, const mpq& one_minus_f_0) {
TRACE("gomory_cut_detail_real", tout << "real\n";); TRACE("gomory_cut_detail_real", tout << "real\n";);
mpq new_a; mpq new_a;
if (at_low(x_j)) { if (at_low(x_j)) {
@ -142,9 +147,9 @@ void int_solver::real_case_in_gomory_cut(const mpq & a, unsigned x_j, mpq & k, l
new_a = a / f_0; new_a = a / f_0;
new_a.neg(); new_a.neg();
} }
k.addmul(new_a, lower_bound(x_j).x); // is it a faster operation than m_k->addmul(new_a, lower_bound(x_j).x); // is it a faster operation than
// k += lower_bound(x_j).x * new_a; // k += lower_bound(x_j).x * new_a;
expl.push_justification(column_lower_bound_constraint(x_j), new_a); m_ex->push_justification(column_lower_bound_constraint(x_j), new_a);
} }
else { else {
lp_assert(at_upper(x_j)); lp_assert(at_upper(x_j));
@ -155,11 +160,11 @@ void int_solver::real_case_in_gomory_cut(const mpq & a, unsigned x_j, mpq & k, l
else { else {
new_a = a / one_minus_f_0; new_a = a / one_minus_f_0;
} }
k.addmul(new_a, upper_bound(x_j).x); // k += upper_bound(x_j).x * new_a; m_k->addmul(new_a, upper_bound(x_j).x); // k += upper_bound(x_j).x * new_a;
expl.push_justification(column_upper_bound_constraint(x_j), new_a); m_ex->push_justification(column_upper_bound_constraint(x_j), new_a);
} }
TRACE("gomory_cut_detail_real", tout << a << "*v" << x_j << " k: " << k << "\n";); TRACE("gomory_cut_detail_real", tout << a << "*v" << x_j << " k: " << *m_k << "\n";);
pol.add_monomial(new_a, x_j); m_t->add_monomial(new_a, x_j);
} }
constraint_index int_solver::column_upper_bound_constraint(unsigned j) const { constraint_index int_solver::column_upper_bound_constraint(unsigned j) const {
@ -171,12 +176,13 @@ constraint_index int_solver::column_lower_bound_constraint(unsigned j) const {
} }
void int_solver::int_case_in_gomory_cut(const mpq & a, unsigned x_j, mpq & k, lar_term & t, explanation& expl, mpq & lcm_den, const mpq& f_0, const mpq& one_minus_f_0) { void int_solver::int_case_in_gomory_cut(const mpq & a, unsigned x_j,
mpq & lcm_den, const mpq& f_0, const mpq& one_minus_f_0) {
lp_assert(is_int(x_j)); lp_assert(is_int(x_j));
lp_assert(!a.is_int()); lp_assert(!a.is_int());
mpq f_j = fractional_part(a); mpq f_j = fractional_part(a);
TRACE("gomory_cut_detail", TRACE("gomory_cut_detail",
tout << a << " x_j" << x_j << " k = " << k << "\n"; tout << a << " x_j" << x_j << " k = " << *m_k << "\n";
tout << "f_j: " << f_j << "\n"; tout << "f_j: " << f_j << "\n";
tout << "f_0: " << f_0 << "\n"; tout << "f_0: " << f_0 << "\n";
tout << "1 - f_0: " << 1 - f_0 << "\n"; tout << "1 - f_0: " << 1 - f_0 << "\n";
@ -191,8 +197,8 @@ void int_solver::int_case_in_gomory_cut(const mpq & a, unsigned x_j, mpq & k, la
else { else {
new_a = (1 - f_j) / f_0; new_a = (1 - f_j) / f_0;
} }
k.addmul(new_a, lower_bound(x_j).x); m_k->addmul(new_a, lower_bound(x_j).x);
expl.push_justification(column_lower_bound_constraint(x_j), new_a); m_ex->push_justification(column_lower_bound_constraint(x_j), new_a);
} }
else { else {
lp_assert(at_upper(x_j)); lp_assert(at_upper(x_j));
@ -203,19 +209,19 @@ void int_solver::int_case_in_gomory_cut(const mpq & a, unsigned x_j, mpq & k, la
new_a = (mpq(1) - f_j) / one_minus_f_0; new_a = (mpq(1) - f_j) / one_minus_f_0;
} }
new_a.neg(); // the upper terms are inverted new_a.neg(); // the upper terms are inverted
k.addmul(new_a, upper_bound(x_j).x); m_k->addmul(new_a, upper_bound(x_j).x);
expl.push_justification(column_upper_bound_constraint(x_j), new_a); m_ex->push_justification(column_upper_bound_constraint(x_j), new_a);
} }
TRACE("gomory_cut_detail", tout << "new_a: " << new_a << " k: " << k << "\n";); TRACE("gomory_cut_detail", tout << "new_a: " << new_a << " k: " << *m_k << "\n";);
t.add_monomial(new_a, x_j); m_t->add_monomial(new_a, x_j);
lcm_den = lcm(lcm_den, denominator(new_a)); lcm_den = lcm(lcm_den, denominator(new_a));
} }
lia_move int_solver::report_conflict_from_gomory_cut(mpq & k) { lia_move int_solver::report_conflict_from_gomory_cut() {
TRACE("empty_pol",); TRACE("empty_pol",);
lp_assert(k.is_pos()); lp_assert(m_k->is_pos());
// conflict 0 >= k where k is positive // conflict 0 >= k where k is positive
k.neg(); // returning 0 <= -k m_k->neg(); // returning 0 <= -k
return lia_move::conflict; return lia_move::conflict;
} }
@ -263,43 +269,43 @@ void int_solver::gomory_cut_adjust_t_and_k(vector<std::pair<mpq, unsigned>> & po
} }
} }
bool int_solver::current_solution_is_inf_on_cut(const lar_term& t, const mpq& k) const { bool int_solver::current_solution_is_inf_on_cut() const {
const auto & x = m_lar_solver->m_mpq_lar_core_solver.m_r_x; const auto & x = m_lar_solver->m_mpq_lar_core_solver.m_r_x;
impq v = t.apply(x); impq v = m_t->apply(x);
TRACE( TRACE(
"current_solution_is_inf_on_cut", tout << "v = " << v << " k = " << k << std::endl; "current_solution_is_inf_on_cut", tout << "v = " << v << " k = " << (*m_k) << std::endl;
if (v <=k) { if (v <=(*m_k)) {
tout << "v <= k - it should not happen!\n"; tout << "v <= k - it should not happen!\n";
} }
); );
return v > k; return v > (*m_k);
} }
void int_solver::adjust_term_and_k_for_some_ints_case_gomory(lar_term& t, mpq& k, mpq &lcm_den) { void int_solver::adjust_term_and_k_for_some_ints_case_gomory(mpq &lcm_den) {
lp_assert(!t.is_empty()); lp_assert(!m_t->is_empty());
auto pol = t.coeffs_as_vector(); auto pol = m_t->coeffs_as_vector();
t.clear(); m_t->clear();
if (pol.size() == 1) { if (pol.size() == 1) {
TRACE("gomory_cut_detail", tout << "pol.size() is 1" << std::endl;); TRACE("gomory_cut_detail", tout << "pol.size() is 1" << std::endl;);
unsigned v = pol[0].second; unsigned v = pol[0].second;
lp_assert(is_int(v)); lp_assert(is_int(v));
const mpq& a = pol[0].first; const mpq& a = pol[0].first;
k /= a; (*m_k) /= a;
if (a.is_pos()) { // we have av >= k if (a.is_pos()) { // we have av >= k
if (!k.is_int()) if (!(*m_k).is_int())
k = ceil(k); (*m_k) = ceil((*m_k));
// switch size // switch size
t.add_monomial(- mpq(1), v); m_t->add_monomial(- mpq(1), v);
k.neg(); (*m_k).neg();
} else { } else {
if (!k.is_int()) if (!(*m_k).is_int())
k = floor(k); (*m_k) = floor((*m_k));
t.add_monomial(mpq(1), v); m_t->add_monomial(mpq(1), v);
} }
} else { } else {
TRACE("gomory_cut_detail", tout << "pol.size() > 1" << std::endl;); TRACE("gomory_cut_detail", tout << "pol.size() > 1" << std::endl;);
lcm_den = lcm(lcm_den, denominator(k)); lcm_den = lcm(lcm_den, denominator((*m_k)));
lp_assert(lcm_den.is_pos()); lp_assert(lcm_den.is_pos());
if (!lcm_den.is_one()) { if (!lcm_den.is_one()) {
// normalize coefficients of integer parameters to be integers. // normalize coefficients of integer parameters to be integers.
@ -307,41 +313,41 @@ void int_solver::adjust_term_and_k_for_some_ints_case_gomory(lar_term& t, mpq& k
pi.first *= lcm_den; pi.first *= lcm_den;
SASSERT(!is_int(pi.second) || pi.first.is_int()); SASSERT(!is_int(pi.second) || pi.first.is_int());
} }
k *= lcm_den; (*m_k) *= lcm_den;
} }
// negate everything to return -pol <= -k // negate everything to return -pol <= -(*m_k)
for (const auto & pi: pol) for (const auto & pi: pol)
t.add_monomial(-pi.first, pi.second); m_t->add_monomial(-pi.first, pi.second);
k.neg(); (*m_k).neg();
} }
TRACE("gomory_cut_detail", tout << "k = " << k << std::endl;); TRACE("gomory_cut_detail", tout << "k = " << (*m_k) << std::endl;);
lp_assert(k.is_int()); lp_assert((*m_k).is_int());
} }
lia_move int_solver::mk_gomory_cut(lar_term& t, mpq& k, explanation & expl, unsigned inf_col, const row_strip<mpq> & row) { lia_move int_solver::mk_gomory_cut( unsigned inf_col, const row_strip<mpq> & row) {
lp_assert(column_is_int_inf(inf_col)); lp_assert(column_is_int_inf(inf_col));
TRACE("gomory_cut", TRACE("gomory_cut",
tout << "applying cut at:\n"; m_lar_solver->print_row(row, tout); tout << std::endl; tout << "applying cut at:\n"; m_lar_solver->print_row(row, tout); tout << std::endl;
for (auto p : row) { for (auto & p : row) {
m_lar_solver->m_mpq_lar_core_solver.m_r_solver.print_column_info(p.var(), tout); m_lar_solver->m_mpq_lar_core_solver.m_r_solver.print_column_info(p.var(), tout);
} }
tout << "inf_col = " << inf_col << std::endl; tout << "inf_col = " << inf_col << std::endl;
); );
// gomory will be t <= k and the current solution has a property t > k // gomory will be t <= k and the current solution has a property t > k
k = 1; *m_k = 1;
mpq lcm_den(1); mpq lcm_den(1);
unsigned x_j; unsigned x_j;
mpq a; mpq a;
bool some_int_columns = false; bool some_int_columns = false;
mpq f_0 = int_solver::fractional_part(get_value(inf_col)); mpq f_0 = int_solver::fractional_part(get_value(inf_col));
mpq one_min_f_0 = 1 - f_0; mpq one_min_f_0 = 1 - f_0;
for (auto p : row) { for (auto & p : row) {
x_j = p.var(); x_j = p.var();
if (x_j == inf_col) if (x_j == inf_col)
continue; continue;
@ -349,34 +355,28 @@ lia_move int_solver::mk_gomory_cut(lar_term& t, mpq& k, explanation & expl, unsi
a = p.coeff(); a = p.coeff();
a.neg(); a.neg();
if (is_real(x_j)) if (is_real(x_j))
real_case_in_gomory_cut(a, x_j, k, t, expl, f_0, one_min_f_0); real_case_in_gomory_cut(a, x_j, f_0, one_min_f_0);
else { else {
if (a.is_int()) continue; // f_j will be zero and no monomial will be added if (a.is_int()) continue; // f_j will be zero and no monomial will be added
some_int_columns = true; some_int_columns = true;
int_case_in_gomory_cut(a, x_j, k, t, expl, lcm_den, f_0, one_min_f_0); int_case_in_gomory_cut(a, x_j, lcm_den, f_0, one_min_f_0);
} }
} }
if (t.is_empty()) if (m_t->is_empty())
return report_conflict_from_gomory_cut(k); return report_conflict_from_gomory_cut();
if (some_int_columns) if (some_int_columns)
adjust_term_and_k_for_some_ints_case_gomory(t, k, lcm_den); adjust_term_and_k_for_some_ints_case_gomory(lcm_den);
lp_assert(current_solution_is_inf_on_cut(t, k)); lp_assert(current_solution_is_inf_on_cut());
m_lar_solver->subs_term_columns(t); m_lar_solver->subs_term_columns(*m_t);
TRACE("gomory_cut", tout<<"precut:"; m_lar_solver->print_term(t, tout); tout << " <= " << k << std::endl;); TRACE("gomory_cut", tout<<"precut:"; m_lar_solver->print_term(*m_t, tout); tout << " <= " << *m_k << std::endl;);
return lia_move::cut; return lia_move::cut;
} }
void int_solver::init_check_data() {
unsigned n = m_lar_solver->A_r().column_count();
m_old_values_set.resize(n);
m_old_values_data.resize(n);
}
int int_solver::find_free_var_in_gomory_row(const row_strip<mpq>& row) { int int_solver::find_free_var_in_gomory_row(const row_strip<mpq>& row) {
unsigned j; unsigned j;
for (auto p : row) { for (const auto & p : row) {
j = p.var(); j = p.var();
if (!is_base(j) && is_free(j)) if (!is_base(j) && is_free(j))
return static_cast<int>(j); return static_cast<int>(j);
@ -384,21 +384,17 @@ int int_solver::find_free_var_in_gomory_row(const row_strip<mpq>& row) {
return -1; return -1;
} }
lia_move int_solver::proceed_with_gomory_cut(lar_term& t, mpq& k, explanation& ex, unsigned j, bool & upper) { lia_move int_solver::proceed_with_gomory_cut(unsigned j) {
lia_move ret;
const row_strip<mpq>& row = m_lar_solver->get_row(row_of_basic_column(j)); const row_strip<mpq>& row = m_lar_solver->get_row(row_of_basic_column(j));
int free_j = find_free_var_in_gomory_row(row); int free_j = find_free_var_in_gomory_row(row);
if (free_j != -1) { if (free_j != -1)
ret = create_branch_on_column(j, t, k, true, upper); return lia_move::undef;
} else if (!is_gomory_cut_target(row)) { if (!is_gomory_cut_target(row)) {
bool upper; return lia_move::undef;
ret = create_branch_on_column(j, t, k, false, upper);
} else {
upper = false;
ret = mk_gomory_cut(t, k, ex, j, row);
} }
return ret;
*m_upper = false;
return mk_gomory_cut(j, row);
} }
@ -460,14 +456,14 @@ struct pivoted_rows_tracking_control {
} }
}; };
void int_solver::copy_explanations_from_cut_solver(explanation &ex) { void int_solver::copy_explanations_from_cut_solver() {
TRACE("propagate_and_backjump_step_int", TRACE("propagate_and_backjump_step_int",
for (unsigned j: m_cut_solver.m_explanation) for (unsigned j: m_cut_solver.m_explanation)
m_lar_solver->print_constraint(m_lar_solver->constraints()[j], tout);); m_lar_solver->print_constraint(m_lar_solver->constraints()[j], tout););
for (unsigned j : m_cut_solver.m_explanation) { for (unsigned j : m_cut_solver.m_explanation) {
ex.push_justification(j); m_ex->push_justification(j);
} }
m_cut_solver.m_explanation.clear(); m_cut_solver.m_explanation.clear();
} }
@ -490,18 +486,43 @@ void int_solver::catch_up_in_adding_constraints_to_cut_solver() {
} }
} }
impq int_solver::get_cube_delta_for_term(const lar_term& t) const {
if (t.size() == 2) {
bool seen_minus = false;
bool seen_plus = false;
for(const auto & p : t) {
const mpq & c = p.coeff();
if (c == one_of_type<mpq>()) {
seen_plus = true;
} else if (c == -one_of_type<mpq>()) {
seen_minus = true;
} else {
goto usual_delta;
}
}
if (seen_minus && seen_plus)
return zero_of_type<impq>();
return impq(0, 1);
}
usual_delta:
mpq delta = zero_of_type<mpq>();
for (const auto & p : t) {
delta += abs(p.coeff());
}
delta *= mpq(1, 2);
return impq(delta);
}
bool int_solver::tighten_term_for_cube(unsigned i) { bool int_solver::tighten_term_for_cube(unsigned i) {
unsigned ti = i + m_lar_solver->terms_start_index(); unsigned ti = i + m_lar_solver->terms_start_index();
if (!m_lar_solver->term_is_used_as_row(ti)) if (!m_lar_solver->term_is_used_as_row(ti))
return true; return true;
const lar_term* t = m_lar_solver->terms()[i]; const lar_term* t = m_lar_solver->terms()[i];
mpq delta = zero_of_type<mpq>();
for (const auto & p : *t) { impq delta = get_cube_delta_for_term(*t);
delta += abs(p.coeff());
}
delta *= mpq(1, 2);
TRACE("cube", m_lar_solver->print_term_as_indices(*t, tout); tout << ", delta = " << delta;); TRACE("cube", m_lar_solver->print_term_as_indices(*t, tout); tout << ", delta = " << delta;);
if (is_zero(delta))
return true;
return m_lar_solver->tighten_term_bounds_by_delta(i, delta); return m_lar_solver->tighten_term_bounds_by_delta(i, delta);
} }
@ -515,7 +536,7 @@ bool int_solver::tighten_terms_for_cube() {
} }
bool int_solver::find_cube() { bool int_solver::find_cube() {
if (m_branch_cut_counter % settings().m_int_branch_find_cube != 0) if (m_branch_cut_counter % settings().m_int_find_cube_period != 0)
return false; return false;
settings().st().m_cube_calls++; settings().st().m_cube_calls++;
@ -535,95 +556,120 @@ bool int_solver::find_cube() {
TRACE("cube", tout << "cannot find a feasiblie solution";); TRACE("cube", tout << "cannot find a feasiblie solution";);
m_lar_solver->pop(); m_lar_solver->pop();
move_non_basic_columns_to_bounds(); move_non_basic_columns_to_bounds();
m_lar_solver->find_feasible_solution(); find_feasible_solution();
lp_assert(m_lar_solver->get_status() == lp_status::OPTIMAL); lp_assert(m_cut_solver.cancel() || is_feasible());
// it can happen that we found an integer solution here // it can happen that we found an integer solution here
return !m_lar_solver->r_basis_has_inf_int(); return !m_lar_solver->r_basis_has_inf_int();
} }
m_lar_solver->round_to_integer_solution();
m_lar_solver->pop(); m_lar_solver->pop();
m_lar_solver->round_to_integer_solution();
lp_assert(m_cut_solver.cancel() || is_feasible());
return true; return true;
} }
lia_move int_solver::check(lar_term& t, mpq& k, explanation& ex, bool & upper) { void int_solver::find_feasible_solution() {
init_check_data(); m_lar_solver->find_feasible_solution();
// it is a reimplementation of lp_assert(lp_status::OPTIMAL == m_lar_solver->get_status() || lp_status::FEASIBLE == m_lar_solver->get_status());
// final_check_status theory_arith<Ext>::check_int_feasibility() }
// from theory_arith_int.h with the addition of cut_solver
if (!has_inf_int()) lia_move int_solver::run_gcd_test() {
return lia_move::ok; if (settings().m_int_run_gcd_test) {
if (settings().m_run_gcd_test) {
settings().st().m_gcd_calls++; settings().st().m_gcd_calls++;
if (!gcd_test(ex)) { if (!gcd_test()) {
TRACE("gcd_test", tout << "conflict";);
settings().st().m_gcd_conflicts++; settings().st().m_gcd_conflicts++;
return lia_move::conflict; return lia_move::conflict;
} }
} else {
TRACE("gcd_test", tout << "no test";);
} }
return lia_move::undef;
}
lia_move int_solver::call_cut_solver() {
if ((m_branch_cut_counter) % settings().m_int_cut_solver_period != 0 || !all_columns_are_bounded())
return lia_move::undef;
TRACE("check_main_int", tout<<"cut_solver";);
catch_up_in_adding_constraints_to_cut_solver();
auto check_res = m_cut_solver.check();
settings().st().m_cut_solver_calls++;
switch (check_res) {
case cut_solver::lbool::l_false:
copy_explanations_from_cut_solver();
settings().st().m_cut_solver_false++;
return lia_move::conflict;
case cut_solver::lbool::l_true:
settings().st().m_cut_solver_true++;
copy_values_from_cut_solver();
lp_assert(m_lar_solver->all_constraints_hold());
return lia_move::sat;
case cut_solver::lbool::l_undef:
settings().st().m_cut_solver_undef++;
if (m_cut_solver.try_getting_cut(*m_t, *m_k, m_lar_solver->m_mpq_lar_core_solver.m_r_x)) {
m_lar_solver->subs_term_columns(*m_t);
TRACE("cut_solver_cuts",
tout<<"precut from cut_solver:"; m_lar_solver->print_term(*m_t, tout); tout << " <= " << *m_k << std::endl;);
return lia_move::cut;
}
default:
return lia_move::undef;
}
}
lia_move int_solver::check(lar_term& t, mpq& k, explanation& ex, bool & upper) {
if (!has_inf_int())
return lia_move::sat;
m_t = &t; m_k = &k; m_ex = &ex; m_upper = &upper;
if (run_gcd_test() == lia_move::conflict)
return lia_move::conflict;
pivoted_rows_tracking_control pc(m_lar_solver); pivoted_rows_tracking_control pc(m_lar_solver);
/* if (m_params.m_arith_euclidean_solver) apply_euclidean_solver(); */ if(settings().m_int_pivot_fixed_vars_from_basis)
//m_lar_solver->pivot_fixed_vars_from_basis(); m_lar_solver->pivot_fixed_vars_from_basis();
patch_int_infeasible_nbasic_columns();
if (!has_inf_int()) if (patch_nbasic_columns() == lia_move::sat)
return lia_move::ok; return lia_move::sat;
++m_branch_cut_counter; ++m_branch_cut_counter;
if (find_cube()){ if (find_cube()){
settings().st().m_cube_success++; settings().st().m_cube_success++;
return lia_move::ok; return lia_move::sat;
} }
TRACE("cube", tout << "cube did not succeed";);
lia_move r = call_cut_solver();
if (r != lia_move::undef)
return r;
if ((m_branch_cut_counter) % settings().m_int_branch_cut_solver == 0) { if ((m_branch_cut_counter) % settings().m_int_gomory_cut_period == 0) {
TRACE("check_main_int", tout<<"cut_solver";);
catch_up_in_adding_constraints_to_cut_solver();
auto check_res = m_cut_solver.check();
settings().st().m_cut_solver_calls++;
switch (check_res) {
case cut_solver::lbool::l_false:
copy_explanations_from_cut_solver(ex);
settings().st().m_cut_solver_false++;
return lia_move::conflict;
case cut_solver::lbool::l_true:
settings().st().m_cut_solver_true++;
copy_values_from_cut_solver();
lp_assert(m_lar_solver->all_constraints_hold());
return lia_move::ok;
case cut_solver::lbool::l_undef:
settings().st().m_cut_solver_undef++;
if (m_cut_solver.try_getting_cut(t, k, m_lar_solver->m_mpq_lar_core_solver.m_r_x)) {
m_lar_solver->subs_term_columns(t);
TRACE("cut_solver_cuts",
tout<<"precut from cut_solver:"; m_lar_solver->print_term(t, tout); tout << " <= " << k << std::endl;);
return lia_move::cut;
}
break;
default:
return lia_move::give_up;
}
}
if ((m_branch_cut_counter) % settings().m_int_branch_cut_gomory_threshold == 0) {
TRACE("check_main_int", tout << "gomory";); TRACE("check_main_int", tout << "gomory";);
if (move_non_basic_columns_to_bounds()) { if (move_non_basic_columns_to_bounds()) {
lp_status st = m_lar_solver->find_feasible_solution(); lp_status st = m_lar_solver->find_feasible_solution();
lp_assert(non_basic_columns_are_at_bounds()); lp_assert(non_basic_columns_are_at_bounds());
if (st != lp_status::FEASIBLE && st != lp_status::OPTIMAL) { if (st != lp_status::FEASIBLE && st != lp_status::OPTIMAL) {
TRACE("arith_int", tout << "give_up\n";); TRACE("arith_int", tout << "give_up\n";);
return lia_move::give_up; return lia_move::undef;
} }
} }
int j = find_inf_int_base_column(); int j = find_inf_int_base_column();
if (j == -1) return lia_move::ok; if (j == -1) {
j = find_inf_int_nbasis_column();
return j == -1? lia_move::sat : create_branch_on_column(j);
}
TRACE("arith_int", tout << "j = " << j << " does not have an integer assignment: " << get_value(j) << "\n";); TRACE("arith_int", tout << "j = " << j << " does not have an integer assignment: " << get_value(j) << "\n";);
return proceed_with_gomory_cut(t, k, ex, j, upper); r = proceed_with_gomory_cut(j);
if (r != lia_move::undef)
return r;
return create_branch_on_column(j);
} }
TRACE("check_main_int", tout << "branch"; ); TRACE("check_main_int", tout << "branch"; );
return create_branch_on_column(find_inf_int_base_column(), t, k, false, upper); int j = find_inf_int_base_column();
if (j == -1) {
j = find_inf_int_nbasis_column();
if (j == -1)
return lia_move::sat;
}
return create_branch_on_column(j);
} }
bool int_solver::move_non_basic_column_to_bounds(unsigned j) { bool int_solver::move_non_basic_column_to_bounds(unsigned j) {
@ -668,6 +714,9 @@ bool int_solver::move_non_basic_columns_to_bounds() {
if (move_non_basic_column_to_bounds(j)) if (move_non_basic_column_to_bounds(j))
change = true; change = true;
} }
if (settings().simplex_strategy() == simplex_strategy_enum::tableau_costs)
m_lar_solver->update_x_and_inf_costs_for_columns_with_changed_bounds_tableau();
return change; return change;
} }
@ -683,27 +732,26 @@ void int_solver::set_value_for_nbasic_column_ignore_old_values(unsigned j, const
void int_solver::set_value_for_nbasic_column(unsigned j, const impq & new_val) { void int_solver::set_value_for_nbasic_column(unsigned j, const impq & new_val) {
lp_assert(!is_base(j)); lp_assert(!is_base(j));
auto & x = m_lar_solver->m_mpq_lar_core_solver.m_r_x[j]; auto & x = m_lar_solver->m_mpq_lar_core_solver.m_r_x[j];
if (m_lar_solver->has_int_var() && !m_old_values_set.contains(j)) {
m_old_values_set.insert(j);
m_old_values_data[j] = x;
}
auto delta = new_val - x; auto delta = new_val - x;
x = new_val; x = new_val;
m_lar_solver->change_basic_columns_dependend_on_a_given_nb_column(j, delta); m_lar_solver->change_basic_columns_dependend_on_a_given_nb_column(j, delta);
} }
void int_solver::patch_int_infeasible_non_basic_column(unsigned j) { void int_solver::patch_nbasic_column(unsigned j) {
if (!is_int(j)) return;
bool inf_l, inf_u;
impq l, u;
mpq m;
if (!get_value(j).is_int() || !get_freedom_interval_for_column(j, inf_l, l, inf_u, u, m)) {
move_non_basic_column_to_bounds(j);
return;
}
auto & lcs = m_lar_solver->m_mpq_lar_core_solver; auto & lcs = m_lar_solver->m_mpq_lar_core_solver;
impq & val = lcs.m_r_x[j]; impq & val = lcs.m_r_x[j];
bool val_is_int = val.is_int(); bool val_is_int = val.is_int();
if (settings().m_int_patch_only_integer_values) {
if (!val_is_int)
return;
}
bool inf_l, inf_u;
impq l, u;
mpq m;
if (!get_freedom_interval_for_column(j, inf_l, l, inf_u, u, m)) {
return;
}
bool m_is_one = m.is_one(); bool m_is_one = m.is_one();
if (m.is_one() && val_is_int) if (m.is_one() && val_is_int)
return; return;
@ -722,7 +770,6 @@ void int_solver::patch_int_infeasible_non_basic_column(unsigned j) {
if (inf_u || l <= u) { if (inf_u || l <= u) {
TRACE("patch_int", TRACE("patch_int",
tout << "patching with l: " << l << '\n';); tout << "patching with l: " << l << '\n';);
set_value_for_nbasic_column(j, l); set_value_for_nbasic_column(j, l);
} }
else { else {
@ -742,28 +789,29 @@ void int_solver::patch_int_infeasible_non_basic_column(unsigned j) {
tout << "patching with 0\n";); tout << "patching with 0\n";);
} }
} }
void int_solver::patch_int_infeasible_nbasic_columns() { lia_move int_solver::patch_nbasic_columns() {
settings().st().m_patches++;
lp_assert(is_feasible()); lp_assert(is_feasible());
for (unsigned j : m_lar_solver->m_mpq_lar_core_solver.m_r_nbasis) { for (unsigned j : m_lar_solver->m_mpq_lar_core_solver.m_r_nbasis) {
patch_int_infeasible_non_basic_column(j); patch_nbasic_column(j);
if (!is_feasible())
break;
} }
if (!is_feasible()) { lp_assert(is_feasible());
move_non_basic_columns_to_bounds(); if (!has_inf_int()) {
m_lar_solver->find_feasible_solution(); settings().st().m_patches_success++;
return lia_move::sat;
} }
return lia_move::undef;
} }
mpq get_denominators_lcm(const row_strip<mpq> & row) { mpq get_denominators_lcm(const row_strip<mpq> & row) {
mpq r(1); mpq r(1);
for (auto c : row) { for (auto & c : row) {
r = lcm(r, denominator(c.coeff())); r = lcm(r, denominator(c.coeff()));
} }
return r; return r;
} }
bool int_solver::gcd_test_for_row(static_matrix<mpq, numeric_pair<mpq>> & A, unsigned i, explanation & ex) { bool int_solver::gcd_test_for_row(static_matrix<mpq, numeric_pair<mpq>> & A, unsigned i) {
mpq lcm_den = get_denominators_lcm(A.m_rows[i]); mpq lcm_den = get_denominators_lcm(A.m_rows[i]);
mpq consts(0); mpq consts(0);
mpq gcds(0); mpq gcds(0);
@ -812,7 +860,7 @@ bool int_solver::gcd_test_for_row(static_matrix<mpq, numeric_pair<mpq>> & A, uns
if (!(consts / gcds).is_int()) { if (!(consts / gcds).is_int()) {
TRACE("gcd_test", tout << "row failed the GCD test:\n"; display_row_info(tout, i);); TRACE("gcd_test", tout << "row failed the GCD test:\n"; display_row_info(tout, i););
fill_explanation_from_fixed_columns(A.m_rows[i], ex); fill_explanation_from_fixed_columns(A.m_rows[i]);
return false; return false;
} }
@ -822,29 +870,29 @@ bool int_solver::gcd_test_for_row(static_matrix<mpq, numeric_pair<mpq>> & A, uns
} }
if (least_coeff_is_bounded) { if (least_coeff_is_bounded) {
return ext_gcd_test(A.m_rows[i], least_coeff, lcm_den, consts, ex); return ext_gcd_test(A.m_rows[i], least_coeff, lcm_den, consts);
} }
return true; return true;
} }
void int_solver::add_to_explanation_from_fixed_or_boxed_column(unsigned j, explanation & ex) { void int_solver::add_to_explanation_from_fixed_or_boxed_column(unsigned j) {
constraint_index lc, uc; constraint_index lc, uc;
m_lar_solver->get_bound_constraint_witnesses_for_column(j, lc, uc); m_lar_solver->get_bound_constraint_witnesses_for_column(j, lc, uc);
ex.m_explanation.push_back(std::make_pair(mpq(1), lc)); m_ex->m_explanation.push_back(std::make_pair(mpq(1), lc));
ex.m_explanation.push_back(std::make_pair(mpq(1), uc)); m_ex->m_explanation.push_back(std::make_pair(mpq(1), uc));
} }
void int_solver::fill_explanation_from_fixed_columns(const row_strip<mpq> & row, explanation & ex) { void int_solver::fill_explanation_from_fixed_columns(const row_strip<mpq> & row) {
for (const auto & c : row) { for (const auto & c : row) {
if (!m_lar_solver->column_is_fixed(c.var())) if (!m_lar_solver->column_is_fixed(c.var()))
continue; continue;
add_to_explanation_from_fixed_or_boxed_column(c.var(), ex); add_to_explanation_from_fixed_or_boxed_column(c.var());
} }
} }
bool int_solver::gcd_test(explanation & ex) { bool int_solver::gcd_test() {
auto & A = m_lar_solver->A_r(); // getting the matrix auto & A = m_lar_solver->A_r(); // getting the matrix
for (unsigned i = 0; i < A.row_count(); i++) for (unsigned i = 0; i < A.row_count(); i++)
if (!gcd_test_for_row(A, i, ex)) { if (!gcd_test_for_row(A, i)) {
return false; return false;
} }
@ -854,7 +902,7 @@ bool int_solver::gcd_test(explanation & ex) {
bool int_solver::ext_gcd_test(const row_strip<mpq> & row, bool int_solver::ext_gcd_test(const row_strip<mpq> & row,
mpq const & least_coeff, mpq const & least_coeff,
mpq const & lcm_den, mpq const & lcm_den,
mpq const & consts, explanation& ex) { mpq const & consts) {
mpq gcds(0); mpq gcds(0);
mpq l(consts); mpq l(consts);
mpq u(consts); mpq u(consts);
@ -884,7 +932,7 @@ bool int_solver::ext_gcd_test(const row_strip<mpq> & row,
// u += ncoeff * lower_bound(j).get_rational(); // u += ncoeff * lower_bound(j).get_rational();
u.addmul(ncoeff, m_lar_solver->column_lower_bound(j).x); u.addmul(ncoeff, m_lar_solver->column_lower_bound(j).x);
} }
add_to_explanation_from_fixed_or_boxed_column(j, ex); add_to_explanation_from_fixed_or_boxed_column(j);
} }
else if (gcds.is_zero()) { else if (gcds.is_zero()) {
gcds = abs_ncoeff; gcds = abs_ncoeff;
@ -903,7 +951,7 @@ bool int_solver::ext_gcd_test(const row_strip<mpq> & row,
mpq u1 = floor(u/gcds); mpq u1 = floor(u/gcds);
if (u1 < l1) { if (u1 < l1) {
fill_explanation_from_fixed_columns(row, ex); fill_explanation_from_fixed_columns(row);
return false; return false;
} }
@ -926,9 +974,6 @@ int_solver::int_solver(lar_solver* lar_slv) :
[this]() {return m_lar_solver->A_r().column_count();}, [this]() {return m_lar_solver->A_r().column_count();},
[this](unsigned j) {return get_value(j);}, [this](unsigned j) {return get_value(j);},
settings()) { settings()) {
lp_assert(m_old_values_set.size() == 0);
m_old_values_set.resize(lar_slv->A_r().column_count());
m_old_values_data.resize(lar_slv->A_r().column_count(), zero_of_type<impq>());
m_lar_solver->set_int_solver(this); m_lar_solver->set_int_solver(this);
} }
@ -996,7 +1041,7 @@ bool int_solver::get_freedom_interval_for_column(unsigned j, bool & inf_l, impq
unsigned row_index; unsigned row_index;
lp_assert(settings().use_tableau()); lp_assert(settings().use_tableau());
const auto & A = m_lar_solver->A_r(); const auto & A = m_lar_solver->A_r();
for (auto c : A.column(j)) { for (const auto &c : A.column(j)) {
row_index = c.var(); row_index = c.var();
const mpq & a = c.coeff(); const mpq & a = c.coeff();
@ -1017,7 +1062,7 @@ bool int_solver::get_freedom_interval_for_column(unsigned j, bool & inf_l, impq
if (has_low(i)) if (has_low(i))
set_upper(u, inf_u, xj + (xi - lcs.m_r_lower_bounds()[i]) / a); set_upper(u, inf_u, xj + (xi - lcs.m_r_lower_bounds()[i]) / a);
} }
if (!inf_l && !inf_u && l == u) break;; if (!inf_l && !inf_u && l >= u) break;
} }
TRACE("freedom_interval", TRACE("freedom_interval",
@ -1029,11 +1074,9 @@ bool int_solver::get_freedom_interval_for_column(unsigned j, bool & inf_l, impq
if (inf_u) tout << "oo"; else tout << u; if (inf_u) tout << "oo"; else tout << u;
tout << "]\n"; tout << "]\n";
tout << "val = " << get_value(j) << "\n"; tout << "val = " << get_value(j) << "\n";
tout << "return " << (inf_l || inf_u || l <= u);
); );
lp_assert(inf_l || l <= get_value(j)); return (inf_l || inf_u || l <= u);
lp_assert(inf_u || u >= get_value(j));
return true;
} }
bool int_solver::is_int(unsigned j) const { bool int_solver::is_int(unsigned j) const {
@ -1238,22 +1281,22 @@ const impq& int_solver::lower_bound(unsigned j) const {
return m_lar_solver->column_lower_bound(j); return m_lar_solver->column_lower_bound(j);
} }
lia_move int_solver::create_branch_on_column(int j, lar_term& t, mpq& k, bool free_column, bool & upper) { lia_move int_solver::create_branch_on_column(int j) {
TRACE("check_main_int", tout << "branching" << std::endl;); TRACE("check_main_int", tout << "branching" << std::endl;);
lp_assert(t.is_empty()); lp_assert(m_t->is_empty());
lp_assert(j != -1); lp_assert(j != -1);
t.add_monomial(mpq(1), m_lar_solver->adjust_column_index_to_term_index(j)); m_t->add_monomial(mpq(1), m_lar_solver->adjust_column_index_to_term_index(j));
if (free_column) { if (is_free(j)) {
upper = true; *m_upper = true;
k = mpq(0); *m_k = mpq(0);
} else { } else {
upper = left_branch_is_more_narrow_than_right(j); *m_upper = left_branch_is_more_narrow_than_right(j);
k = upper? floor(get_value(j)) : ceil(get_value(j)); *m_k = *m_upper? floor(get_value(j)) : ceil(get_value(j));
} }
TRACE("arith_int", tout << "branching v" << j << " = " << get_value(j) << "\n"; TRACE("arith_int", tout << "branching v" << j << " = " << get_value(j) << "\n";
display_column(tout, j); display_column(tout, j);
tout << "k = " << k << std::endl; tout << "k = " << *m_k << std::endl;
); );
return lia_move::branch; return lia_move::branch;

76
src/util/lp/int_solver.h
View file

@ -30,14 +30,14 @@ class lar_solver;
template <typename T, typename X> template <typename T, typename X>
struct lp_constraint; struct lp_constraint;
enum class lia_move { enum class lia_move {
ok, sat,
branch, branch,
cut, cut,
conflict, conflict,
continue_with_check, continue_with_check,
give_up, undef,
unsat unsat
}; };
struct explanation { struct explanation {
vector<std::pair<mpq, constraint_index>> m_explanation; vector<std::pair<mpq, constraint_index>> m_explanation;
@ -52,11 +52,13 @@ struct explanation {
class int_solver { class int_solver {
public: public:
// fields // fields
lar_solver *m_lar_solver; lar_solver * m_lar_solver;
int_set m_old_values_set; unsigned m_branch_cut_counter;
vector<impq> m_old_values_data; cut_solver m_cut_solver;
unsigned m_branch_cut_counter; lar_term* m_t; // the term to return in the cut
cut_solver m_cut_solver; mpq *m_k; // the right side of the cut
explanation *m_ex; // the conflict explanation
bool *m_upper; // we have a cut m_t*x <= k if m_upper is true nad m_t*x >= k otherwise
// methods // methods
int_solver(lar_solver* lp); int_solver(lar_solver* lp);
@ -69,7 +71,7 @@ private:
// how to tighten bounds for integer variables. // how to tighten bounds for integer variables.
bool gcd_test_for_row(static_matrix<mpq, numeric_pair<mpq>> & A, unsigned i, explanation &); bool gcd_test_for_row(static_matrix<mpq, numeric_pair<mpq>> & A, unsigned i);
// gcd test // gcd test
// 5*x + 3*y + 6*z = 5 // 5*x + 3*y + 6*z = 5
@ -79,25 +81,17 @@ private:
// this is unsolvable because 5/3 is not an integer. // this is unsolvable because 5/3 is not an integer.
// so we create a lemma that rules out this condition. // so we create a lemma that rules out this condition.
// //
bool gcd_test(explanation & ); // returns false in case of failure. Creates a theory lemma in case of failure. bool gcd_test(); // returns false in case of failure. Creates a theory lemma in case of failure.
// create goromy cuts
// either creates a conflict or a bound.
// branch and bound:
// decide what to branch and bound on
// creates a fresh inequality.
bool branch(const lp_constraint<mpq, mpq> & new_inequality); bool branch(const lp_constraint<mpq, mpq> & new_inequality);
bool ext_gcd_test(const row_strip<mpq>& row, bool ext_gcd_test(const row_strip<mpq>& row,
mpq const & least_coeff, mpq const & least_coeff,
mpq const & lcm_den, mpq const & lcm_den,
mpq const & consts, mpq const & consts);
explanation & ex); void fill_explanation_from_fixed_columns(const row_strip<mpq> & row);
void fill_explanation_from_fixed_columns(const row_strip<mpq> & row, explanation &); void add_to_explanation_from_fixed_or_boxed_column(unsigned j);
void add_to_explanation_from_fixed_or_boxed_column(unsigned j, explanation &); void patch_nbasic_column(unsigned j);
void patch_int_infeasible_non_basic_column(unsigned j); lia_move patch_nbasic_columns();
void patch_int_infeasible_nbasic_columns();
bool get_freedom_interval_for_column(unsigned j, bool & inf_l, impq & l, bool & inf_u, impq & u, mpq & m); bool get_freedom_interval_for_column(unsigned j, bool & inf_l, impq & l, bool & inf_u, impq & u, mpq & m);
const impq & lower_bound(unsigned j) const; const impq & lower_bound(unsigned j) const;
const impq & upper_bound(unsigned j) const; const impq & upper_bound(unsigned j) const;
@ -111,7 +105,6 @@ private:
void set_value_for_nbasic_column(unsigned j, const impq & new_val); void set_value_for_nbasic_column(unsigned j, const impq & new_val);
void set_value_for_nbasic_column_ignore_old_values(unsigned j, const impq & new_val); void set_value_for_nbasic_column_ignore_old_values(unsigned j, const impq & new_val);
bool non_basic_columns_are_at_bounds() const; bool non_basic_columns_are_at_bounds() const;
void failed();
bool is_feasible() const; bool is_feasible() const;
const impq & get_value(unsigned j) const; const impq & get_value(unsigned j) const;
bool column_is_int_inf(unsigned j) const; bool column_is_int_inf(unsigned j) const;
@ -122,11 +115,10 @@ private:
lp_settings& settings(); lp_settings& settings();
bool move_non_basic_columns_to_bounds(); bool move_non_basic_columns_to_bounds();
void branch_infeasible_int_var(unsigned); void branch_infeasible_int_var(unsigned);
lia_move mk_gomory_cut(lar_term& t, mpq& k,explanation & ex, unsigned inf_col, const row_strip<mpq>& row); lia_move mk_gomory_cut(unsigned inf_col, const row_strip<mpq>& row);
lia_move report_conflict_from_gomory_cut(mpq & k); lia_move report_conflict_from_gomory_cut();
void adjust_term_and_k_for_some_ints_case_gomory(lar_term& t, mpq& k, mpq& lcm_den); void adjust_term_and_k_for_some_ints_case_gomory(mpq& lcm_den);
void init_check_data(); lia_move proceed_with_gomory_cut(unsigned j);
lia_move proceed_with_gomory_cut(lar_term& t, mpq& k, explanation& ex, unsigned j, bool & upper);
int find_free_var_in_gomory_row(const row_strip<mpq>& ); int find_free_var_in_gomory_row(const row_strip<mpq>& );
bool is_gomory_cut_target(const row_strip<mpq>&); bool is_gomory_cut_target(const row_strip<mpq>&);
bool at_bound(unsigned j) const; bool at_bound(unsigned j) const;
@ -147,19 +139,19 @@ public:
return n.x - floor(n.x); return n.x - floor(n.x);
} }
private: private:
void real_case_in_gomory_cut(const mpq & a, unsigned x_j, mpq & k, lar_term& t, explanation & ex, const mpq& f_0, const mpq& one_minus_f_0); void real_case_in_gomory_cut(const mpq & a, unsigned x_j, const mpq& f_0, const mpq& one_minus_f_0);
void int_case_in_gomory_cut(const mpq & a, unsigned x_j, mpq & k, lar_term& t, explanation& ex, mpq & lcm_den, const mpq& f_0, const mpq& one_minus_f_0); void int_case_in_gomory_cut(const mpq & a, unsigned x_j, mpq & lcm_den, const mpq& f_0, const mpq& one_minus_f_0);
constraint_index column_upper_bound_constraint(unsigned j) const; constraint_index column_upper_bound_constraint(unsigned j) const;
constraint_index column_lower_bound_constraint(unsigned j) const; constraint_index column_lower_bound_constraint(unsigned j) const;
void display_row_info(std::ostream & out, unsigned row_index) const; void display_row_info(std::ostream & out, unsigned row_index) const;
void gomory_cut_adjust_t_and_k(vector<std::pair<mpq, unsigned>> & pol, lar_term & t, mpq &k, bool num_ints, mpq &lcm_den); void gomory_cut_adjust_t_and_k(vector<std::pair<mpq, unsigned>> & pol, lar_term & t, mpq &k, bool num_ints, mpq &lcm_den);
bool current_solution_is_inf_on_cut(const lar_term& t, const mpq& k) const; bool current_solution_is_inf_on_cut() const;
public: public:
bool shift_var(unsigned j, unsigned range); bool shift_var(unsigned j, unsigned range);
private: private:
unsigned random(); unsigned random();
bool has_inf_int() const; bool has_inf_int() const;
lia_move create_branch_on_column(int j, lar_term& t, mpq& k, bool free_column, bool & upper); lia_move create_branch_on_column(int j);
void catch_up_in_adding_constraints_to_cut_solver(); void catch_up_in_adding_constraints_to_cut_solver();
public: public:
template <typename T> template <typename T>
@ -168,7 +160,7 @@ public:
void get_int_coeffs_from_constraint(const lar_base_constraint* c, vector<cut_solver::monomial>& coeff, T & rs); void get_int_coeffs_from_constraint(const lar_base_constraint* c, vector<cut_solver::monomial>& coeff, T & rs);
bool is_term(unsigned j) const; bool is_term(unsigned j) const;
void add_constraint_to_cut_solver(unsigned,const lar_base_constraint*); void add_constraint_to_cut_solver(unsigned,const lar_base_constraint*);
void copy_explanations_from_cut_solver(explanation &); void copy_explanations_from_cut_solver();
void pop(unsigned); void pop(unsigned);
void push(); void push();
void copy_values_from_cut_solver(); void copy_values_from_cut_solver();
@ -177,5 +169,11 @@ public:
bool tighten_terms_for_cube(); bool tighten_terms_for_cube();
bool tighten_term_for_cube(unsigned); bool tighten_term_for_cube(unsigned);
unsigned column_count() const; unsigned column_count() const;
bool all_columns_are_bounded() const;
impq get_cube_delta_for_term(const lar_term&) const;
void find_feasible_solution();
int find_inf_int_nbasis_column() const;
lia_move run_gcd_test();
lia_move call_cut_solver();
}; };
} }

1
src/util/lp/lar_core_solver.h
View file

@ -831,5 +831,6 @@ public:
} }
const vector<unsigned>& r_basis() const { return m_r_basis; } const vector<unsigned>& r_basis() const { return m_r_basis; }
const vector<unsigned>& r_nbasis() const { return m_r_nbasis; }
}; };
} }

14
src/util/lp/lar_solver.cpp
View file

@ -2131,7 +2131,7 @@ var_index lar_solver:: to_var_index(unsigned ext_j) const {
return it->second.internal_j(); return it->second.internal_j();
} }
bool lar_solver::tighten_term_bounds_by_delta(unsigned term_index, const mpq& delta) { bool lar_solver::tighten_term_bounds_by_delta(unsigned term_index, const impq& delta) {
unsigned tj = term_index + m_terms_start_index; unsigned tj = term_index + m_terms_start_index;
auto it = m_ext_vars_to_columns.find(tj); auto it = m_ext_vars_to_columns.find(tj);
if (it == m_ext_vars_to_columns.end()) if (it == m_ext_vars_to_columns.end())
@ -2141,17 +2141,23 @@ bool lar_solver::tighten_term_bounds_by_delta(unsigned term_index, const mpq& de
TRACE("cube", tout << "delta = " << delta << std::endl; TRACE("cube", tout << "delta = " << delta << std::endl;
m_int_solver->display_column(tout, j); ); m_int_solver->display_column(tout, j); );
if (slv.column_has_upper_bound(j) && slv.column_has_lower_bound(j)) { if (slv.column_has_upper_bound(j) && slv.column_has_lower_bound(j)) {
if (slv.m_upper_bounds[j].x - delta < slv.m_lower_bounds[j].x + delta) { if (slv.m_upper_bounds[j] - delta < slv.m_lower_bounds[j] + delta) {
TRACE("cube", tout << "cannot tighten, delta = " << delta;); TRACE("cube", tout << "cannot tighten, delta = " << delta;);
return false; return false;
} }
} }
TRACE("cube", tout << "can tighten";); TRACE("cube", tout << "can tighten";);
if (slv.column_has_upper_bound(j)) { if (slv.column_has_upper_bound(j)) {
add_var_bound(tj, lconstraint_kind::LE, slv.m_upper_bounds[j].x - delta); if (!is_zero(delta.y))
add_var_bound(tj, lconstraint_kind::LT, slv.m_upper_bounds[j].x - delta.x);
else
add_var_bound(tj, lconstraint_kind::LE, slv.m_upper_bounds[j].x - delta.x);
} }
if (slv.column_has_lower_bound(j)) { if (slv.column_has_lower_bound(j)) {
add_var_bound(tj, lconstraint_kind::GE, slv.m_lower_bounds[j].x + delta); if (!is_zero(delta.y))
add_var_bound(tj, lconstraint_kind::GT, slv.m_lower_bounds[j].x + delta.x);
else
add_var_bound(tj, lconstraint_kind::GE, slv.m_lower_bounds[j].x + delta.x);
} }
return true; return true;
} }

3
src/util/lp/lar_solver.h
View file

@ -572,11 +572,12 @@ public:
bool column_corresponds_to_term(unsigned) const; bool column_corresponds_to_term(unsigned) const;
void catch_up_in_updating_int_solver(); void catch_up_in_updating_int_solver();
var_index to_var_index(unsigned ext_j) const; var_index to_var_index(unsigned ext_j) const;
bool tighten_term_bounds_by_delta(unsigned, const mpq&); bool tighten_term_bounds_by_delta(unsigned, const impq&);
void round_to_integer_solution(); void round_to_integer_solution();
void update_delta_for_terms(const impq & delta, unsigned j, const vector<unsigned>&); void update_delta_for_terms(const impq & delta, unsigned j, const vector<unsigned>&);
void fill_vars_to_terms(vector<vector<unsigned>> & vars_to_terms); void fill_vars_to_terms(vector<vector<unsigned>> & vars_to_terms);
unsigned column_count() const { return A_r().column_count(); } unsigned column_count() const { return A_r().column_count(); }
const vector<unsigned> & r_basis() const { return m_mpq_lar_core_solver.r_basis(); } const vector<unsigned> & r_basis() const { return m_mpq_lar_core_solver.r_basis(); }
const vector<unsigned> & r_nbasis() const { return m_mpq_lar_core_solver.r_nbasis(); }
}; };
} }

25
src/util/lp/lp_settings.h
View file

@ -110,6 +110,8 @@ struct stats {
unsigned m_gcd_conflicts; unsigned m_gcd_conflicts;
unsigned m_cube_calls; unsigned m_cube_calls;
unsigned m_cube_success; unsigned m_cube_success;
unsigned m_patches;
unsigned m_patches_success;
stats() { reset(); } stats() { reset(); }
void reset() { memset(this, 0, sizeof(*this)); } void reset() { memset(this, 0, sizeof(*this)); }
}; };
@ -229,10 +231,14 @@ public:
max_row_length_for_bound_propagation(300), max_row_length_for_bound_propagation(300),
backup_costs(true), backup_costs(true),
column_number_threshold_for_using_lu_in_lar_solver(4000), column_number_threshold_for_using_lu_in_lar_solver(4000),
m_int_branch_cut_gomory_threshold(4), m_int_gomory_cut_period(4),
m_int_branch_cut_solver(8), m_int_cut_solver_period(8),
m_run_gcd_test(true), m_int_find_cube_period(4),
m_cut_solver_cycle_on_var(10) m_int_cuts_etc_period(4),
m_int_run_gcd_test(true),
m_cut_solver_cycle_on_var(10),
m_int_pivot_fixed_vars_from_basis(false),
m_int_patch_only_integer_values(true)
{} {}
void set_resource_limit(lp_resource_limit& lim) { m_resource_limit = &lim; } void set_resource_limit(lp_resource_limit& lim) { m_resource_limit = &lim; }
@ -339,11 +345,14 @@ public:
unsigned max_row_length_for_bound_propagation; unsigned max_row_length_for_bound_propagation;
bool backup_costs; bool backup_costs;
unsigned column_number_threshold_for_using_lu_in_lar_solver; unsigned column_number_threshold_for_using_lu_in_lar_solver;
unsigned m_int_branch_cut_gomory_threshold; unsigned m_int_gomory_cut_period;
unsigned m_int_branch_cut_solver; unsigned m_int_cut_solver_period;
unsigned m_int_branch_find_cube; unsigned m_int_find_cube_period;
bool m_run_gcd_test; unsigned m_int_cuts_etc_period;
bool m_int_run_gcd_test;
unsigned m_cut_solver_cycle_on_var; unsigned m_cut_solver_cycle_on_var;
bool m_int_pivot_fixed_vars_from_basis;
bool m_int_patch_only_integer_values;
}; // end of lp_settings class }; // end of lp_settings class