diff --git a/src/util/lp/general_matrix.h b/src/util/lp/general_matrix.h index ce510eb6e..1c643161a 100644 --- a/src/util/lp/general_matrix.h +++ b/src/util/lp/general_matrix.h @@ -105,12 +105,12 @@ public: } template - void init_row_from_container(int i, const T & c, std::function column_fix) { + void init_row_from_container(int i, const T & c, std::function column_fix, const mpq& sign) { auto & row = m_data[adjust_row(i)]; lp_assert(row_is_initialized_correctly(row)); for (const auto & p : c) { unsigned j = adjust_column(column_fix(p.var())); - row[j] = p.coeff(); + row[j] = sign * p.coeff(); } } diff --git a/src/util/lp/hnf_cutter.h b/src/util/lp/hnf_cutter.h index dff0b8ef8..061b3a130 100644 --- a/src/util/lp/hnf_cutter.h +++ b/src/util/lp/hnf_cutter.h @@ -29,6 +29,7 @@ class hnf_cutter { var_register m_var_register; general_matrix m_A; vector m_terms; + vector m_terms_upper; svector m_constraints_for_explanation; vector m_right_sides; lp_settings & m_settings; @@ -54,15 +55,22 @@ public: // m_A will be filled from scratch in init_matrix_A m_var_register.clear(); m_terms.clear(); + m_terms_upper.clear(); m_constraints_for_explanation.clear(); m_right_sides.clear(); m_abs_max = zero_of_type(); m_overflow = false; } - void add_term(const lar_term* t, const mpq &rs, constraint_index ci) { + void add_term(const lar_term* t, const mpq &rs, constraint_index ci, bool upper_bound) { m_terms.push_back(t); - m_right_sides.push_back(rs); + m_terms_upper.push_back(upper_bound); + if (upper_bound) + m_right_sides.push_back(rs); + else + m_right_sides.push_back(-rs); + m_constraints_for_explanation.push_back(ci); + for (const auto &p : *t) { m_var_register.add_var(p.var()); mpq t = abs(ceil(p.coeff())); @@ -76,7 +84,8 @@ public: } void initialize_row(unsigned i) { - m_A.init_row_from_container(i, * m_terms[i], [this](unsigned j) { return m_var_register.add_var(j);}); + mpq sign = m_terms_upper[i]? one_of_type(): - one_of_type(); + m_A.init_row_from_container(i, * m_terms[i], [this](unsigned j) { return m_var_register.add_var(j);}, sign); } void init_matrix_A() { diff --git a/src/util/lp/int_solver.cpp b/src/util/lp/int_solver.cpp index f49ecabfd..e58abb4fd 100644 --- a/src/util/lp/int_solver.cpp +++ b/src/util/lp/int_solver.cpp @@ -537,8 +537,9 @@ void int_solver::try_add_term_to_A_for_hnf(unsigned i) { mpq rs; const lar_term* t = m_lar_solver->terms()[i]; constraint_index ci; - if (!hnf_cutter_is_full() && m_lar_solver->get_equality_and_right_side_for_term_on_current_x(i, rs, ci)) { - m_hnf_cutter.add_term(t, rs, ci); + bool upper_bound; + if (!hnf_cutter_is_full() && m_lar_solver->get_equality_and_right_side_for_term_on_current_x(i, rs, ci, upper_bound)) { + m_hnf_cutter.add_term(t, rs, ci, upper_bound); } } diff --git a/src/util/lp/lar_solver.cpp b/src/util/lp/lar_solver.cpp index 99a0c5883..f8d16ab8c 100644 --- a/src/util/lp/lar_solver.cpp +++ b/src/util/lp/lar_solver.cpp @@ -2227,8 +2227,19 @@ void lar_solver::round_to_integer_solution() { update_delta_for_terms(del, j, vars_to_terms[j]); } } +// return true if all y coords are zeroes +bool lar_solver::sum_first_coords(const lar_term& t, mpq & val) const { + val = zero_of_type(); + for (const auto & c : t) { + const auto & x = m_mpq_lar_core_solver.m_r_x[c.var()]; + if (!is_zero(x.y)) + return false; + val += x.x * c.coeff(); + } + return true; +} -bool lar_solver::get_equality_and_right_side_for_term_on_current_x(unsigned term_index, mpq & rs, constraint_index& ci) const { +bool lar_solver::get_equality_and_right_side_for_term_on_current_x(unsigned term_index, mpq & rs, constraint_index& ci, bool &upper_bound) const { unsigned tj = term_index + m_terms_start_index; unsigned j; bool is_int; @@ -2236,26 +2247,29 @@ bool lar_solver::get_equality_and_right_side_for_term_on_current_x(unsigned term return false; // the term does not have a bound because it does not correspond to a column if (!is_int) // todo - allow for the next version of hnf return false; - impq term_val; - bool term_val_ready = false; + bool rs_is_calculated = false; mpq b; bool is_strict; + const lar_term& t = *terms()[term_index]; if (has_upper_bound(j, ci, b, is_strict) && !is_strict) { lp_assert(b.is_int()); - term_val = terms()[term_index]->apply(m_mpq_lar_core_solver.m_r_x); - term_val_ready = true; - if (term_val.x == b) { - rs = b; + if (!sum_first_coords(t, rs)) + return false; + rs_is_calculated = true; + if (rs == b) { + upper_bound = true; return true; } } if (has_lower_bound(j, ci, b, is_strict) && !is_strict) { - if (!term_val_ready) - term_val = terms()[term_index]->apply(m_mpq_lar_core_solver.m_r_x); + if (!rs_is_calculated){ + if (!sum_first_coords(t, rs)) + return false; + } lp_assert(b.is_int()); - if (term_val.x == b) { - rs = b; + if (rs == b) { + upper_bound = false; return true; } } diff --git a/src/util/lp/lar_solver.h b/src/util/lp/lar_solver.h index f17aa4a0d..72705dfb3 100644 --- a/src/util/lp/lar_solver.h +++ b/src/util/lp/lar_solver.h @@ -573,9 +573,10 @@ public: unsigned column_count() const { return A_r().column_count(); } const vector & r_basis() const { return m_mpq_lar_core_solver.r_basis(); } const vector & r_nbasis() const { return m_mpq_lar_core_solver.r_nbasis(); } - bool get_equality_and_right_side_for_term_on_current_x(unsigned i, mpq &rs, constraint_index& ci) const; + bool get_equality_and_right_side_for_term_on_current_x(unsigned i, mpq &rs, constraint_index& ci, bool &upper_bound) const; bool remove_from_basis(unsigned); lar_term get_term_to_maximize(unsigned ext_j) const; void set_cut_strategy(unsigned cut_frequency); + bool sum_first_coords(const lar_term& t, mpq & val) const; }; }