diff --git a/src/math/lp/horner.cpp b/src/math/lp/horner.cpp index 2320f4091..43400c8eb 100644 --- a/src/math/lp/horner.cpp +++ b/src/math/lp/horner.cpp @@ -204,22 +204,92 @@ interv horner::interval_of_mul(const nex& e) { return a; } -bool horner::find_term_expr(rational& a, const lp::lar_term * &t, rational& b) const { +void horner::add_mul_to_vector(const nex& e, vector> &v) { + TRACE("nla_horner_details", tout << e << "\n";); + SASSERT(e.is_mul() && e.children().size() == 2); + rational r; + lpvar j = -1; + for (const nex & c : e.children()) { + switch (c.type()) { + case expr_type::SCALAR: + r = c.value(); + break; + case expr_type::VAR: + j = c.var(); + break; + default: + TRACE("nla_horner_details", tout << e.type() << "\n";); + SASSERT(false); + } + } + SASSERT((j + 1) && !(r.is_zero() || r.is_one())); + v.push_back(std::make_pair(r, j)); +} + +void horner::add_linear_to_vector(const nex& e, vector> &v) { + TRACE("nla_horner_details", tout << e << "\n";); + switch (e.type()) { + case expr_type::MUL: + add_mul_to_vector(e, v); + break; + case expr_type::VAR: + v.push_back(std::make_pair(rational(1), e.var())); + break; + default: + TRACE("nla_horner_details", tout << e.type() << "\n";); + UNREACHABLE(); + SASSERT(false); + } +} +// e = a * can_t + b, but we do not calculate the offset here +lp::lar_term horner::expression_to_canonical_form(nex& e, rational& a, rational& b) { + TRACE("nla_horner_details", tout << e << "\n";); + lpvar smallest_j; + vector> v; + b = rational(0); + unsigned a_index; + for (const nex& c : e.children()) { + if (c.is_scalar()) { + b += c.value(); + } else { + add_linear_to_vector(c, v); + if (v.size() == 1 || smallest_j < v.back().second) { + smallest_j = v.back().second; + a_index = v.size() - 1; + } + } + } + a = v[a_index].first; + lp::lar_term t; + + if (a.is_one()) { + for (unsigned k = 0; k < v.size(); k++) { + auto& p = v[k]; + t.add_coeff_var(p.first, p.second); + } + } else { + for (unsigned k = 0; k < v.size(); k++) { + auto& p = v[k]; + if (k != a_index) + t.add_coeff_var(p.first/a, p.second); + else + t.add_var(p.second); + } + } + TRACE("nla_horner_details", tout << a << "* ("; + lp::lar_solver::print_term_as_indices(t, tout) << ") + " << b << std::endl;); + return t; +} + +bool horner::find_term_expr(const nex& e, rational& a, const lp::lar_term* &t, rational& b) const { + nex n = e; + lp::lar_term can_t = expression_to_canonical_form(n, a, b); + TRACE("nla_horner_details", tout << "term = "; c().m_lar_solver.print_term(can_t, tout) << "\n";); + SASSERT(false); return false; } -interv horner::interval_of_sum(const nex& e) { - TRACE("nla_horner_details", tout << "e=" << e << "\n";); - SASSERT(e.is_sum()); - if (e.sum_is_linear()) { - const lp::lar_term * t; - rational a,b; - if (find_term_expr(a, t, b)) { - //todo create interval from a*t + b - SASSERT(false); - } - } - +interv horner::interval_of_sum_no_terms(const nex& e) { auto & es = e.children(); interv a = interval_of_expr(es[0]); if (m_intervals.is_inf(a)) { @@ -252,10 +322,42 @@ interv horner::interval_of_sum(const nex& e) { tout << " interv = "; m_intervals.display(tout, a);); return a; } + +bool horner::interval_from_term(const nex& e, interv & i) const { + rational a, b; + nex n = e; + lp::lar_term canonical_t = expression_to_canonical_form(n, a, b); + // TRACE("nla_horner_details", + SASSERT(false); + return false; +} + + +interv horner::interval_of_sum(const nex& e) { + TRACE("nla_horner_details", tout << "e=" << e << "\n";); + SASSERT(e.is_sum()); + interv i_e = interval_of_sum_no_terms(e); + if (e.sum_is_linear()) { + interv i_from_term ; + if (interval_from_term(e, i_from_term) + && + is_tighter(i_from_term, i_e)) + return i_from_term; + } + + return i_e; +} + // sets the dependencies also void horner::set_var_interval(lpvar v, interv& b) { m_intervals.set_var_interval_with_deps(v, b); TRACE("nla_horner_details_var", tout << "v = "; print_var(v, tout) << "\n"; m_intervals.display(tout, b);); } + +bool horner::is_tighter(const interv& a, const interv& b) const { + return m_intervals.is_tighter(a, b); +} + + } diff --git a/src/math/lp/horner.h b/src/math/lp/horner.h index d91ab8296..835cb2bcf 100644 --- a/src/math/lp/horner.h +++ b/src/math/lp/horner.h @@ -31,6 +31,7 @@ class horner : common { intervals m_intervals; public: typedef nla_expr nex; + typedef intervals::interval interv; horner(core *core); void horner_lemmas(); @@ -42,6 +43,7 @@ public: nex nexvar(lpvar j) const; intervals::interval interval_of_sum(const nex&); + intervals::interval interval_of_sum_no_terms(const nex&); intervals::interval interval_of_mul(const nex&); void set_interval_for_scalar(intervals::interval&, const rational&); void set_var_interval(lpvar j, intervals::interval&); @@ -49,6 +51,11 @@ public: template // T has an iterator of (coeff(), var()) bool row_has_monomial_to_refine(const T&) const; - bool find_term_expr(rational& a, const lp::lar_term * & t, rational& b) const; + bool find_term_expr(const nex& e, rational& a, const lp::lar_term * & t, rational& b) const; + static lp::lar_term expression_to_canonical_form(nex&, rational& a, rational & b); + static void add_linear_to_vector(const nex&, vector> &); + static void add_mul_to_vector(const nex&, vector> &); + bool is_tighter(const interv&, const interv&) const; + bool interval_from_term(const nex& e, interv&) const; }; // end of horner } diff --git a/src/math/lp/lar_solver.cpp b/src/math/lp/lar_solver.cpp index d45b4d4da..4c3ab8774 100644 --- a/src/math/lp/lar_solver.cpp +++ b/src/math/lp/lar_solver.cpp @@ -33,7 +33,8 @@ lar_solver::lar_solver() : m_status(lp_status::UNKNOWN), m_int_solver(nullptr), m_terms_start_index(1000000), m_var_register(0), - m_term_register(m_terms_start_index) + m_term_register(m_terms_start_index), + m_need_register_terms(true) // change to false {} void lar_solver::set_track_pivoted_rows(bool v) { @@ -1427,7 +1428,7 @@ std::ostream& lar_solver::print_term(lar_term const& term, std::ostream & out) c return out; } -std::ostream& lar_solver::print_term_as_indices(lar_term const& term, std::ostream & out) const { +std::ostream& lar_solver::print_term_as_indices(lar_term const& term, std::ostream & out) { print_linear_combination_of_column_indices_only(term.coeffs_as_vector(), out); return out; } @@ -1774,7 +1775,7 @@ void lar_solver::add_new_var_to_core_fields_for_mpq(bool register_in_basis) { var_index lar_solver::add_term_undecided(const vector> & coeffs) { - push_and_register_term(new lar_term(coeffs)); + push_term(new lar_term(coeffs)); return m_terms_start_index + m_terms.size() - 1; } @@ -1805,11 +1806,7 @@ bool lar_solver::term_coeffs_are_ok(const vector> & co return false; } #endif -void lar_solver::push_and_register_term(lar_term* t) { -#if Z3DEBUG_CHECK_UNIQUE_TERMS - lp_assert(m_set_of_terms.find(t) == m_set_of_terms.end()); - m_set_of_terms.insert(t); -#endif +void lar_solver::push_term(lar_term* t) { m_terms.push_back(t); } @@ -1829,8 +1826,8 @@ var_index lar_solver::add_term(const vector> & coeffs, lp_assert(all_vars_are_registered(coeffs)); if (strategy_is_undecided()) return add_term_undecided(coeffs); - - push_and_register_term(new lar_term(coeffs)); + lar_term * t = new lar_term(coeffs); + push_term(t); SASSERT(m_term_register.size() == m_terms.size()); unsigned adjusted_term_index = m_terms.size() - 1; var_index ret = m_terms_start_index + adjusted_term_index; @@ -1840,6 +1837,10 @@ var_index lar_solver::add_term(const vector> & coeffs, m_rows_with_changed_bounds.insert(A_r().row_count() - 1); } lp_assert(m_var_register.size() == A_r().column_count()); + if (m_need_register_terms) { + lar_term normalized_t = t->get_normalized_by_min_var(); + m_normalized_terms_to_columns[normalized_t] = A_r().column_count() - 1; + } return ret; } @@ -2255,10 +2256,6 @@ bool lar_solver::column_corresponds_to_term(unsigned j) const { return m_var_register.local_to_external(j) >= m_terms_start_index; } -var_index lar_solver::to_column(unsigned ext_j) const { - return m_var_register.external_to_local(ext_j); -} - bool lar_solver::tighten_term_bounds_by_delta(unsigned term_index, const impq& delta) { unsigned tj = term_index + m_terms_start_index; unsigned j; @@ -2374,6 +2371,14 @@ void lar_solver::set_cut_strategy(unsigned cut_frequency) { settings().set_hnf_cut_period(100000000); } } +void lar_solver::register_existing_terms() { + for (unsigned k = 0; k < m_terms.size(); k++) { + lar_term * t = m_terms[k]; + lar_term normalized_t = t->get_normalized_by_min_var(); + lpvar j = m_var_register.external_to_local(k + m_terms_start_index); + m_normalized_terms_to_columns[normalized_t] = j; + } +} } // namespace lp diff --git a/src/math/lp/lar_solver.h b/src/math/lp/lar_solver.h index e327a3553..196c87c86 100644 --- a/src/math/lp/lar_solver.h +++ b/src/math/lp/lar_solver.h @@ -47,7 +47,7 @@ namespace lp { class lar_solver : public column_namer { - + typedef unsigned lpvar; struct term_hasher { std::size_t operator()(const lar_term &t) const { @@ -73,7 +73,7 @@ class lar_solver : public column_namer { } }; - std::unordered_map m_normalized_terms_to_columns; + //////////////////// fields ////////////////////////// @@ -105,7 +105,8 @@ public: stacked_value m_term_count; vector m_terms; indexed_vector m_column_buffer; - // end of fields + bool m_need_register_terms; + std::unordered_map m_normalized_terms_to_columns; // end of fields unsigned terms_start_index() const { return m_terms_start_index; } const vector & terms() const { return m_terms; } @@ -181,8 +182,7 @@ public: var_index add_term_undecided(const vector> & coeffs); bool term_coeffs_are_ok(const vector> & coeffs); - void push_and_register_term(lar_term* t); - + void push_term(lar_term* t); void add_row_for_term(const lar_term * term, unsigned term_ext_index); void add_row_from_term_no_constraint(const lar_term * term, unsigned term_ext_index); @@ -505,7 +505,7 @@ public: std::ostream& print_term(lar_term const& term, std::ostream & out) const; - std::ostream& print_term_as_indices(lar_term const& term, std::ostream & out) const; + static std::ostream& print_term_as_indices(lar_term const& term, std::ostream & out); std::ostream& print_constraint(const lar_base_constraint * c, std::ostream & out) const; std::ostream& print_constraint_indices_only(const lar_base_constraint * c, std::ostream & out) const; @@ -628,7 +628,6 @@ public: lar_core_solver & get_core_solver() { return m_mpq_lar_core_solver; } bool column_corresponds_to_term(unsigned) const; void catch_up_in_updating_int_solver(); - var_index to_column(unsigned ext_j) const; bool tighten_term_bounds_by_delta(unsigned, const impq&); void round_to_integer_solution(); void update_delta_for_terms(const impq & delta, unsigned j, const vector&); @@ -644,5 +643,6 @@ public: void fix_Ax_b_on_rounded_rows(); void fix_Ax_b_on_rounded_row(unsigned); void collect_rounded_rows_to_fix(); + void register_existing_terms(); }; } diff --git a/src/math/lp/lar_term.h b/src/math/lp/lar_term.h index cefbf0ce6..469e34f22 100644 --- a/src/math/lp/lar_term.h +++ b/src/math/lp/lar_term.h @@ -39,7 +39,7 @@ public: } } - void add_var(unsigned j) { + void add_var(lpvar j) { rational c(1); add_monomial(c, j); } @@ -65,8 +65,8 @@ public: // some terms get used in add constraint // it is the same as the offset in the m_constraints - vector> coeffs_as_vector() const { - vector> ret; + vector> coeffs_as_vector() const { + vector> ret; for (const auto & p : m_coeffs) { ret.push_back(std::make_pair(p.m_value, p.m_key)); } @@ -93,7 +93,7 @@ public: m_coeffs.insert(k, b); } - bool contains(unsigned j) const { + bool contains(lpvar j) const { return m_coeffs.contains(j); } @@ -116,11 +116,11 @@ public: } struct ival { - unsigned m_var; + lpvar m_var; const mpq & m_coeff; - ival(unsigned var, const mpq & val) : m_var(var), m_coeff(val) { + ival(lpvar var, const mpq & val) : m_var(var), m_coeff(val) { } - unsigned var() const { return m_var;} + lpvar var() const { return m_var;} const mpq & coeff() const { return m_coeff; } }; diff --git a/src/math/lp/nla_intervals.h b/src/math/lp/nla_intervals.h index 4761b92d5..3b3a52d27 100644 --- a/src/math/lp/nla_intervals.h +++ b/src/math/lp/nla_intervals.h @@ -185,7 +185,46 @@ public: m_config.add_deps(a, b, deps, i); } + bool is_tighter_on_lower(const interval& a, const interval& b) const { + if (lower_is_inf(a)) + return false; + if (lower_is_inf(b)) + return true; + if (rational(lower(a)) < rational(lower(b))) + return true; + if (lower(a) > lower(b)) + return false; + + if (!a.m_lower_open) + return false; + if (b.m_lower_open) + return false; + + return true; + } + + bool is_tighter_on_upper(const interval& a, const interval& b) const { + if (upper_is_inf(a)) + return false; + if (upper_is_inf(b)) + return true; + if (rational(upper(a)) > rational(upper(b))) + return true; + if (rational(upper(a)) < rational(upper(b))) + return false; + + if (!a.m_upper_open) + return false; + if (b.m_upper_open) + return false; + + return true; + } + bool is_tighter(const interval& a, const interval& b) const { + return (is_tighter_on_lower(a, b) && !is_tighter_on_upper(b, a)) || + (is_tighter_on_upper(a, b) && is_tighter_on_lower(b, a)); + } bool upper_is_inf(const interval& a) const { return m_config.upper_is_inf(a); } bool lower_is_inf(const interval& a) const { return m_config.lower_is_inf(a); } diff --git a/src/smt/theory_lra.cpp b/src/smt/theory_lra.cpp index b9291f857..da82c350f 100644 --- a/src/smt/theory_lra.cpp +++ b/src/smt/theory_lra.cpp @@ -649,6 +649,10 @@ class theory_lra::imp { } TRACE("arith", tout << "v" << v << " := " << mk_pp(t, m) << " " << _has_var << "\n";); if (!_has_var) { + if (m_solver->m_need_register_terms == false) { + m_solver->m_need_register_terms = true; + m_solver->register_existing_terms(); + } m_switcher.add_monomial(register_theory_var_in_lar_solver(v), vars.size(), vars.c_ptr()); } }