diff --git a/src/math/dd/dd_pdd.cpp b/src/math/dd/dd_pdd.cpp index 857665583..455fc4039 100644 --- a/src/math/dd/dd_pdd.cpp +++ b/src/math/dd/dd_pdd.cpp @@ -1308,7 +1308,7 @@ namespace dd { else if (l == lo() && h == hi()) return *this; else - return m.mk_var(v)*h + l; + return m.mk_var(var())*h + l; } diff --git a/src/math/grobner/pdd_solver.cpp b/src/math/grobner/pdd_solver.cpp index dd79b506b..51aed555c 100644 --- a/src/math/grobner/pdd_solver.cpp +++ b/src/math/grobner/pdd_solver.cpp @@ -342,6 +342,7 @@ namespace dd { for (equation* e : m_solved) dealloc(e); for (equation* e : m_to_simplify) dealloc(e); for (equation* e : m_processed) dealloc(e); + m_subst.reset(); m_solved.reset(); m_processed.reset(); m_to_simplify.reset(); @@ -354,16 +355,39 @@ namespace dd { void solver::add(pdd const& p, u_dependency * dep) { if (p.is_zero()) return; equation * eq = alloc(equation, p, dep); - if (check_conflict(*eq)) { + if (check_conflict(*eq)) return; - } push_equation(to_simplify, eq); - if (!m_var2level.empty()) { + if (!m_var2level.empty()) m_levelp1 = std::max(m_var2level[p.var()]+1, m_levelp1); - } update_stats_max_degree_and_size(*eq); - } + } + + void solver::add_subst(unsigned v, pdd const& p, u_dependency* dep) { + SASSERT(m_processed.empty()); + SASSERT(m_solved.empty()); + + m_subst.push_back({v, p, dep}); + + for (auto* e : m_to_simplify) { + auto r = e->poly().subst_pdd(v, p); + if (r == e->poly()) + continue; + *e = m_dep_manager.mk_join(dep, e->dep()); + *e = r; + } + } + + void solver::simplify(pdd& p, u_dependency*& d) { + for (auto const& [v, q, d2] : m_subst) { + pdd r = p.subst_pdd(v, q); + if (r != p) { + p = r; + d = m_dep_manager.mk_join(d, d2); + } + } + } bool solver::canceled() { return m_limit.is_canceled(); @@ -446,9 +470,24 @@ namespace dd { } std::ostream& solver::display(std::ostream& out) const { - out << "solved\n"; for (auto e : m_solved) display(out, *e); - out << "processed\n"; for (auto e : m_processed) display(out, *e); - out << "to_simplify\n"; for (auto e : m_to_simplify) display(out, *e); + if (!m_solved.empty()) { + out << "solved\n"; for (auto e : m_solved) display(out, *e); + } + if (!m_processed.empty()) { + out << "processed\n"; for (auto e : m_processed) display(out, *e); + } + if (!m_to_simplify.empty()) { + out << "to_simplify\n"; for (auto e : m_to_simplify) display(out, *e); + } + if (!m_subst.empty()) { + out << "subst\n"; + for (auto const& [v, p, d] : m_subst) { + out << "v" << v << " := " << p; + if (m_print_dep) + m_print_dep(d, out); + out << "\n"; + } + } return display_statistics(out); } diff --git a/src/math/grobner/pdd_solver.h b/src/math/grobner/pdd_solver.h index 150a70df8..3f3c4e73e 100644 --- a/src/math/grobner/pdd_solver.h +++ b/src/math/grobner/pdd_solver.h @@ -118,6 +118,7 @@ private: equation_vector m_solved; // equations with solved variables, triangular equation_vector m_processed; equation_vector m_to_simplify; + vector> m_subst; mutable u_dependency_manager m_dep_manager; equation_vector m_all_eqs; equation* m_conflict; @@ -136,6 +137,9 @@ public: void add(pdd const& p) { add(p, nullptr); } void add(pdd const& p, u_dependency * dep); + void simplify(pdd& p, u_dependency*& dep); + void add_subst(unsigned v, pdd const& p, u_dependency* dep); + void simplify(); void saturate(); diff --git a/src/math/lp/nla_core.cpp b/src/math/lp/nla_core.cpp index e6eb8ad26..9e3e168e9 100644 --- a/src/math/lp/nla_core.cpp +++ b/src/math/lp/nla_core.cpp @@ -1716,7 +1716,7 @@ void core::configure_grobner() { for (lpvar k : emons()[j].vars()) r *= pdd_expr(rational::one(), k, dep); r -= val_of_fixed_var_with_deps(j, dep); - m_pdd_grobner.add(r, dep); + add_eq_to_grobner(r, dep); } } } @@ -1880,6 +1880,76 @@ dd::pdd core::pdd_expr(const rational& c, lpvar j, u_dependency*& dep) { return r; } +/** + \brief convert p == 0 into a solved form v == r, such that + v has bounds [lo, oo) iff r has bounds [lo', oo) + v has bounds (oo,hi] iff r has bounds (oo,hi'] + + The solved form allows the Grobner solver identify more bounds conflicts. + A bad leading term can miss bounds conflicts. + For example for x + y + z == 0 where x, y : [0, oo) and z : (oo,0] + we prefer to solve z == -x - y instead of x == -z - y + because the solution -z - y has neither an upper, nor a lower bound. + */ +bool core::is_solved(dd::pdd const& p, unsigned& v, dd::pdd& r) { + if (!p.is_linear()) + return false; + r = p; + unsigned num_lo = 0, num_hi = 0; + unsigned lo = 0, hi = 0; + rational lc, hc, c; + while (!r.is_val()) { + SASSERT(r.hi().is_val()); + v = r.var(); + rational val = r.hi().val(); + switch (m_lar_solver.get_column_type(v)) { + case lp::column_type::lower_bound: + if (val > 0) num_lo++, lo = v, lc = val; else num_hi++, hi = v, hc = val; + break; + case lp::column_type::upper_bound: + if (val < 0) num_lo++, lo = v, lc = val; else num_hi++, hi = v, hc = val; + break; + case lp::column_type::fixed: + case lp::column_type::boxed: + break; + default: + return false; + } + if (num_lo > 1 && num_hi > 1) + return false; + r = r.lo(); + } + if (num_lo == 1 && num_hi > 1) { + v = lo; + c = lc; + } + else if (num_hi == 1 && num_lo > 1) { + v = hi; + c = hc; + } + else + return false; + + r = c*m_pdd_manager.mk_var(v) - p; + if (c != 1) + r = r * (1/c); + return true; +} + +/** + \brief add an equality to grobner solver, convert it to solved form if available. +*/ +void core::add_eq_to_grobner(dd::pdd& p, u_dependency* dep) { + unsigned v; + dd::pdd q(m_pdd_manager); + m_pdd_grobner.simplify(p, dep); + if (is_solved(p, v, q)) + m_pdd_grobner.add_subst(v, q, dep); + else + m_pdd_grobner.add(p, dep); +} + + void core::add_row_to_grobner(const vector> & row) { u_dependency *dep = nullptr; rational val; @@ -1887,7 +1957,7 @@ void core::add_row_to_grobner(const vector> & row) { for (const auto &p : row) sum += pdd_expr(p.coeff(), p.var(), dep); TRACE("grobner", print_row(row, tout) << " " << sum << "\n"); - m_pdd_grobner.add(sum, dep); + add_eq_to_grobner(sum, dep); } @@ -1907,7 +1977,6 @@ void core::find_nl_cluster() { TRACE("grobner", tout << "vars in cluster: "; for (lpvar j : active_var_set()) tout << "j" << j << " "; tout << "\n"; display_matrix_of_m_rows(tout); - /*emons().display(tout << "emons\n");*/ ); } diff --git a/src/math/lp/nla_core.h b/src/math/lp/nla_core.h index e50075564..ee07c7db1 100644 --- a/src/math/lp/nla_core.h +++ b/src/math/lp/nla_core.h @@ -466,7 +466,9 @@ public: void display_matrix_of_m_rows(std::ostream & out) const; void set_active_vars_weights(nex_creator&); unsigned get_var_weight(lpvar) const; - void add_row_to_grobner(const vector> & row); + void add_row_to_grobner(const vector> & row); + bool is_solved(dd::pdd const& p, unsigned& v, dd::pdd& r); + void add_eq_to_grobner(dd::pdd& p, u_dependency* dep); bool check_pdd_eq(const dd::solver::equation*); const rational& val_of_fixed_var_with_deps(lpvar j, u_dependency*& dep); dd::pdd pdd_expr(const rational& c, lpvar j, u_dependency*&);