From 739515263250b5cbba94b02f2379250c90d6368e Mon Sep 17 00:00:00 2001 From: Nikolaj Bjorner Date: Sun, 23 Nov 2025 08:59:55 -0800 Subject: [PATCH] factor out coi, use polynomial elaboration for nlsat solver (#8039) * factor out coi, use polynomial elaboration for nlsat solver Signed-off-by: Nikolaj Bjorner * remove unused functionality Signed-off-by: Nikolaj Bjorner --------- Signed-off-by: Nikolaj Bjorner --- src/ast/sls/sls_arith_base.cpp | 4 +- src/math/lp/CMakeLists.txt | 1 + src/math/lp/nla_coi.cpp | 88 +++++++++ src/math/lp/nla_coi.h | 43 +++++ src/math/lp/nla_core.h | 6 + src/math/lp/nra_solver.cpp | 341 ++++++++++++++++++++------------- src/math/lp/nra_solver.h | 2 + 7 files changed, 349 insertions(+), 136 deletions(-) create mode 100644 src/math/lp/nla_coi.cpp create mode 100644 src/math/lp/nla_coi.h diff --git a/src/ast/sls/sls_arith_base.cpp b/src/ast/sls/sls_arith_base.cpp index 5dc768206..df40b0251 100644 --- a/src/ast/sls/sls_arith_base.cpp +++ b/src/ast/sls/sls_arith_base.cpp @@ -2620,8 +2620,10 @@ namespace sls { display(out, ad) << "\n"; } }; - for (var_t v = 0; v < m_vars.size(); ++v) { + for (var_t v = 0; v < m_vars.size(); ++v) { if (!eval_is_correct(v)) { +// if (m.rlimit().is_canceled()) +// return; report_error(verbose_stream(), v); TRACE(arith, report_error(tout, v)); UNREACHABLE(); diff --git a/src/math/lp/CMakeLists.txt b/src/math/lp/CMakeLists.txt index 5c156d38a..729b47782 100644 --- a/src/math/lp/CMakeLists.txt +++ b/src/math/lp/CMakeLists.txt @@ -24,6 +24,7 @@ z3_add_component(lp monomial_bounds.cpp nex_creator.cpp nla_basics_lemmas.cpp + nla_coi.cpp nla_common.cpp nla_core.cpp nla_divisions.cpp diff --git a/src/math/lp/nla_coi.cpp b/src/math/lp/nla_coi.cpp new file mode 100644 index 000000000..fcab22021 --- /dev/null +++ b/src/math/lp/nla_coi.cpp @@ -0,0 +1,88 @@ + +/*++ + Copyright (c) 2025 Microsoft Corporation + + Author: + Lev Nachmanson (levnach) + Nikolaj Bjorner (nbjorner) + --*/ + +#include "math/lp/nla_core.h" +#include "math/lp/nla_coi.h" + +namespace nla { + + void coi::init() { + indexed_uint_set visited; + unsigned_vector todo; + vector var2occurs; + m_term_set.reset(); + m_mon_set.reset(); + m_constraint_set.reset(); + m_var_set.reset(); + auto& lra = c.lra_solver(); + + for (auto ci : lra.constraints().indices()) { + auto const& c = lra.constraints()[ci]; + if (c.is_auxiliary()) + continue; + for (auto const& [coeff, v] : c.coeffs()) { + var2occurs.reserve(v + 1); + var2occurs[v].constraints.push_back(ci); + } + } + + for (auto const& m : c.emons()) { + for (auto v : m.vars()) { + var2occurs.reserve(v + 1); + var2occurs[v].monics.push_back(m.var()); + } + } + + for (const auto *t : lra.terms() ) { + for (auto const iv : *t) { + auto v = iv.j(); + var2occurs.reserve(v + 1); + var2occurs[v].terms.push_back(t->j()); + } + } + + for (auto const& m : c.to_refine()) + todo.push_back(m); + + for (unsigned i = 0; i < todo.size(); ++i) { + auto v = todo[i]; + if (visited.contains(v)) + continue; + visited.insert(v); + m_var_set.insert(v); + var2occurs.reserve(v + 1); + for (auto ci : var2occurs[v].constraints) { + m_constraint_set.insert(ci); + auto const& c = lra.constraints()[ci]; + for (auto const& [coeff, w] : c.coeffs()) + todo.push_back(w); + } + for (auto w : var2occurs[v].monics) + todo.push_back(w); + + for (auto ti : var2occurs[v].terms) { + for (auto iv : lra.get_term(ti)) + todo.push_back(iv.j()); + todo.push_back(ti); + } + + if (lra.column_has_term(v)) { + m_term_set.insert(v); + for (auto kv : lra.get_term(v)) + todo.push_back(kv.j()); + } + + if (c.is_monic_var(v)) { + m_mon_set.insert(v); + for (auto w : c.emons()[v]) + todo.push_back(w); + } + } + } +} \ No newline at end of file diff --git a/src/math/lp/nla_coi.h b/src/math/lp/nla_coi.h new file mode 100644 index 000000000..d05f08fbd --- /dev/null +++ b/src/math/lp/nla_coi.h @@ -0,0 +1,43 @@ + +/*++ + Copyright (c) 2025 Microsoft Corporation + + Abstract: + Class for computing the cone of influence for NL constraints. + It includes variables that come from monomials that have incorrect evaluation and + transitively all constraints and variables that are connected. + + Author: + Lev Nachmanson (levnach) + Nikolaj Bjorner (nbjorner) + --*/ + +#pragma once + +namespace nla { + + class core; + + class coi { + core& c; + indexed_uint_set m_mon_set, m_constraint_set; + indexed_uint_set m_term_set, m_var_set; + + struct occurs { + unsigned_vector constraints; + unsigned_vector monics; + unsigned_vector terms; + }; + + public: + coi(core& c) : c(c) {} + + void init(); + + indexed_uint_set const& mons() const { return m_mon_set; } + indexed_uint_set const& constraints() const { return m_constraint_set; } + indexed_uint_set& terms() { return m_term_set; } + indexed_uint_set const &vars() { return m_var_set; } + + }; +} \ No newline at end of file diff --git a/src/math/lp/nla_core.h b/src/math/lp/nla_core.h index 6121a79a7..baacbc8e8 100644 --- a/src/math/lp/nla_core.h +++ b/src/math/lp/nla_core.h @@ -450,6 +450,12 @@ public: nla_throttle& throttle() { return m_throttle; } const nla_throttle& throttle() const { return m_throttle; } + lp::lar_solver& lra_solver() { return lra; } + + indexed_uint_set const& to_refine() const { + return m_to_refine; + } + }; // end of core struct pp_mon { diff --git a/src/math/lp/nra_solver.cpp b/src/math/lp/nra_solver.cpp index fec10fe0e..bcb33c5b7 100644 --- a/src/math/lp/nra_solver.cpp +++ b/src/math/lp/nra_solver.cpp @@ -9,6 +9,7 @@ #include #include "math/lp/lar_solver.h" #include "math/lp/nra_solver.h" +#include "math/lp/nla_coi.h" #include "nlsat/nlsat_solver.h" #include "math/polynomial/polynomial.h" #include "math/polynomial/algebraic_numbers.h" @@ -25,114 +26,143 @@ typedef nla::mon_eq mon_eq; typedef nla::variable_map_type variable_map_type; struct solver::imp { + lp::lar_solver& lra; reslimit& m_limit; params_ref m_params; u_map m_lp2nl; // map from lar_solver variables to nlsat::solver variables - indexed_uint_set m_term_set; scoped_ptr m_nlsat; scoped_ptr m_values; // values provided by LRA solver scoped_ptr m_tmp1, m_tmp2; + nla::coi m_coi; nla::core& m_nla_core; imp(lp::lar_solver& s, reslimit& lim, params_ref const& p, nla::core& nla_core): lra(s), m_limit(lim), m_params(p), + m_coi(nla_core), m_nla_core(nla_core) {} bool need_check() { return m_nla_core.m_to_refine.size() != 0; } - indexed_uint_set m_mon_set, m_constraint_set; - - struct occurs { - unsigned_vector constraints; - unsigned_vector monics; - unsigned_vector terms; - }; - - void init_cone_of_influence() { - indexed_uint_set visited; - unsigned_vector todo; - vector var2occurs; - m_term_set.reset(); - m_mon_set.reset(); - m_constraint_set.reset(); - - for (auto ci : lra.constraints().indices()) { - auto const& c = lra.constraints()[ci]; - if (c.is_auxiliary()) - continue; - for (auto const& [coeff, v] : c.coeffs()) { - var2occurs.reserve(v + 1); - var2occurs[v].constraints.push_back(ci); - } - } - - for (auto const& m : m_nla_core.emons()) { - for (auto v : m.vars()) { - var2occurs.reserve(v + 1); - var2occurs[v].monics.push_back(m.var()); - } - } - - for (const auto *t : lra.terms() ) { - for (auto const iv : *t) { - auto v = iv.j(); - var2occurs.reserve(v + 1); - var2occurs[v].terms.push_back(t->j()); - } - } - - for (auto const& m : m_nla_core.m_to_refine) - todo.push_back(m); - - for (unsigned i = 0; i < todo.size(); ++i) { - auto v = todo[i]; - if (visited.contains(v)) - continue; - visited.insert(v); - var2occurs.reserve(v + 1); - for (auto ci : var2occurs[v].constraints) { - m_constraint_set.insert(ci); - auto const& c = lra.constraints()[ci]; - for (auto const& [coeff, w] : c.coeffs()) - todo.push_back(w); - } - for (auto w : var2occurs[v].monics) - todo.push_back(w); - - for (auto ti : var2occurs[v].terms) { - for (auto iv : lra.get_term(ti)) - todo.push_back(iv.j()); - todo.push_back(ti); - } - - if (lra.column_has_term(v)) { - m_term_set.insert(v); - for (auto kv : lra.get_term(v)) - todo.push_back(kv.j()); - } - - if (m_nla_core.is_monic_var(v)) { - m_mon_set.insert(v); - for (auto w : m_nla_core.emons()[v]) - todo.push_back(w); - } - } - } - void reset() { m_values = nullptr; m_tmp1 = nullptr; m_tmp2 = nullptr; m_nlsat = alloc(nlsat::solver, m_limit, m_params, false); m_values = alloc(scoped_anum_vector, am()); - m_term_set.reset(); m_lp2nl.reset(); } + // Create polynomial definition for variable v used in setup_assignment_solver. + // Side-effects: updates m_vars2mon when v is a monic variable. + void mk_definition(unsigned v, polynomial_ref_vector &definitions, vector& denominators) { + auto &pm = m_nlsat->pm(); + polynomial::polynomial_ref p(pm); + rational den(1); + if (m_nla_core.emons().is_monic_var(v)) { + auto const &m = m_nla_core.emons()[v]; + for (auto v2 : m.vars()) { + polynomial_ref pw(definitions.get(v2), m_nlsat->pm()); + if (!p) + p = pw; + else + p = p * pw; + } + } + else if (lra.column_has_term(v)) { + for (auto const &[w, coeff] : lra.get_term(v)) { + den = lcm(denominator(coeff), den); + } + for (auto const &[w, coeff] : lra.get_term(v)) { + auto coeff1 = den * coeff; + polynomial_ref pw(definitions.get(w), m_nlsat->pm()); + if (!p) + p = constant(coeff1) * pw; + else + p = p + (constant(coeff1) * pw); + } + } + else { + p = pm.mk_polynomial(lp2nl(v)); // nlsat var index equals v (verified above when created) + } + definitions.push_back(p); + denominators.push_back(den); + } + + void setup_solver_poly() { + m_coi.init(); + auto &pm = m_nlsat->pm(); + polynomial_ref_vector definitions(pm); + vector denominators; + for (unsigned v = 0; v < lra.number_of_vars(); ++v) { + if (m_coi.vars().contains(v)) { + auto j = m_nlsat->mk_var(lra.var_is_int(v)); + m_lp2nl.insert(v, j); // we don't really need this. It is going to be the identify map. + mk_definition(v, definitions, denominators); + } + else { + definitions.push_back(nullptr); + denominators.push_back(rational(0)); + } + } + + // we rely on that all information encoded into the tableau is present as a constraint. + for (auto ci : m_coi.constraints()) { + auto &c = lra.constraints()[ci]; + auto &pm = m_nlsat->pm(); + auto k = c.kind(); + auto rhs = c.rhs(); + auto lhs = c.coeffs(); + rational den = denominator(rhs); + for (auto [coeff, v] : lhs) + den = lcm(lcm(den, denominator(coeff)), denominators[v]); + polynomial::polynomial_ref p(pm); + p = pm.mk_const(-den * rhs); + + for (auto [coeff, v] : lhs) { + polynomial_ref poly(pm); + poly = definitions.get(v); + poly = poly * constant(den * coeff / denominators[v]); + p = p + poly; + } + auto lit = add_constraint(p, ci, k); + } + definitions.reset(); + } + + void setup_solver_terms() { + m_coi.init(); + // add linear inequalities from lra_solver + for (auto ci : m_coi.constraints()) + add_constraint(ci); + + // add polynomial definitions. + for (auto const &m : m_coi.mons()) + add_monic_eq(m_nla_core.emons()[m]); + + // add term definitions. + for (unsigned i : m_coi.terms()) + add_term(i); + } + + + + polynomial::polynomial_ref sub(polynomial::polynomial *a, polynomial::polynomial *b) { + return polynomial_ref(m_nlsat->pm().sub(a, b), m_nlsat->pm()); + } + polynomial::polynomial_ref mul(polynomial::polynomial *a, polynomial::polynomial *b) { + return polynomial_ref(m_nlsat->pm().mul(a, b), m_nlsat->pm()); + } + polynomial::polynomial_ref var(lp::lpvar v) { + return polynomial_ref(m_nlsat->pm().mk_polynomial(lp2nl(v)), m_nlsat->pm()); + } + polynomial::polynomial_ref constant(rational const& r) { + return polynomial_ref(m_nlsat->pm().mk_const(r), m_nlsat->pm()); + } + /** \brief one-shot nlsat check. A one shot checker is the least functionality that can @@ -147,24 +177,14 @@ struct solver::imp { lbool check() { SASSERT(need_check()); reset(); - vector core; - - init_cone_of_influence(); - // add linear inequalities from lra_solver - for (auto ci : m_constraint_set) - add_constraint(ci); + vector core; - // add polynomial definitions. - for (auto const& m : m_mon_set) - add_monic_eq(m_nla_core.emons()[m]); + smt_params_helper p(m_params); - // add term definitions. - for (unsigned i : m_term_set) - add_term(i); + setup_solver_poly(); TRACE(nra, m_nlsat->display(tout)); - smt_params_helper p(m_params); if (p.arith_nl_log()) { static unsigned id = 0; std::stringstream strm; @@ -196,7 +216,8 @@ struct solver::imp { } } m_nlsat->collect_statistics(st); - TRACE(nra, + TRACE(nra, tout << "nra result " << r << "\n"); + CTRACE(nra, false, m_nlsat->display(tout << r << "\n"); display(tout); for (auto [j, x] : m_lp2nl) tout << "j" << j << " := x" << x << "\n";); @@ -223,14 +244,15 @@ struct solver::imp { case l_false: { lp::explanation ex; m_nlsat->get_core(core); - for (auto c : core) { - unsigned idx = static_cast(static_cast(c) - this); - ex.push_back(idx); - TRACE(nra, lra.display_constraint(tout << "ex: " << idx << ": ", idx) << "\n";); - } nla::lemma_builder lemma(m_nla_core, __FUNCTION__); + for (auto c : core) { + unsigned idx = static_cast(static_cast(c) - this); + ex.push_back(idx); + } + lemma &= ex; m_nla_core.set_use_nra_model(true); + TRACE(nra, tout << lemma << "\n"); break; } case l_undef: @@ -272,12 +294,25 @@ struct solver::imp { coeffs.push_back(mpz(1)); coeffs.push_back(mpz(-1)); polynomial::polynomial_ref p(pm.mk_polynomial(2, coeffs.data(), mls), pm); - polynomial::polynomial* ps[1] = { p }; - bool even[1] = { false }; - nlsat::literal lit = m_nlsat->mk_ineq_literal(nlsat::atom::kind::EQ, 1, ps, even); + auto lit = mk_literal(p.get(), lp::lconstraint_kind::EQ); + nlsat::assumption a = nullptr; m_nlsat->mk_clause(1, &lit, nullptr); } + nlsat::literal mk_literal(polynomial::polynomial *p, lp::lconstraint_kind k) { + polynomial::polynomial *ps[1] = { p }; + bool is_even[1] = { false }; + switch (k) { + case lp::lconstraint_kind::LE: return ~m_nlsat->mk_ineq_literal(nlsat::atom::kind::GT, 1, ps, is_even); + case lp::lconstraint_kind::GE: return ~m_nlsat->mk_ineq_literal(nlsat::atom::kind::LT, 1, ps, is_even); + case lp::lconstraint_kind::LT: return m_nlsat->mk_ineq_literal(nlsat::atom::kind::LT, 1, ps, is_even); + case lp::lconstraint_kind::GT: return m_nlsat->mk_ineq_literal(nlsat::atom::kind::GT, 1, ps, is_even); + case lp::lconstraint_kind::EQ: return m_nlsat->mk_ineq_literal(nlsat::atom::kind::EQ, 1, ps, is_even); + default: UNREACHABLE(); // unreachable + } + throw default_exception("uexpected operator"); + } + void add_constraint(unsigned idx) { auto& c = lra.constraints()[idx]; auto& pm = m_nlsat->pm(); @@ -297,30 +332,26 @@ struct solver::imp { } rhs *= den; polynomial::polynomial_ref p(pm.mk_linear(sz, coeffs.data(), vars.data(), -rhs), pm); - polynomial::polynomial* ps[1] = { p }; - bool is_even[1] = { false }; + nlsat::literal lit = mk_literal(p.get(), k); + nlsat::assumption a = this + idx; + m_nlsat->mk_clause(1, &lit, a); + } + + nlsat::literal add_constraint(polynomial::polynomial *p, unsigned idx, lp::lconstraint_kind k) { + polynomial::polynomial *ps[1] = {p}; + bool is_even[1] = {false}; nlsat::literal lit; nlsat::assumption a = this + idx; switch (k) { - case lp::lconstraint_kind::LE: - lit = ~m_nlsat->mk_ineq_literal(nlsat::atom::kind::GT, 1, ps, is_even); - break; - case lp::lconstraint_kind::GE: - lit = ~m_nlsat->mk_ineq_literal(nlsat::atom::kind::LT, 1, ps, is_even); - break; - case lp::lconstraint_kind::LT: - lit = m_nlsat->mk_ineq_literal(nlsat::atom::kind::LT, 1, ps, is_even); - break; - case lp::lconstraint_kind::GT: - lit = m_nlsat->mk_ineq_literal(nlsat::atom::kind::GT, 1, ps, is_even); - break; - case lp::lconstraint_kind::EQ: - lit = m_nlsat->mk_ineq_literal(nlsat::atom::kind::EQ, 1, ps, is_even); - break; - default: - UNREACHABLE(); // unreachable + case lp::lconstraint_kind::LE: lit = ~m_nlsat->mk_ineq_literal(nlsat::atom::kind::GT, 1, ps, is_even); break; + case lp::lconstraint_kind::GE: lit = ~m_nlsat->mk_ineq_literal(nlsat::atom::kind::LT, 1, ps, is_even); break; + case lp::lconstraint_kind::LT: lit = m_nlsat->mk_ineq_literal(nlsat::atom::kind::LT, 1, ps, is_even); break; + case lp::lconstraint_kind::GT: lit = m_nlsat->mk_ineq_literal(nlsat::atom::kind::GT, 1, ps, is_even); break; + case lp::lconstraint_kind::EQ: lit = m_nlsat->mk_ineq_literal(nlsat::atom::kind::EQ, 1, ps, is_even); break; + default: UNREACHABLE(); // unreachable } m_nlsat->mk_clause(1, &lit, a); + return lit; } bool check_monic(mon_eq const& m) { @@ -370,7 +401,7 @@ struct solver::imp { for (auto const& m : m_nla_core.emons()) if (any_of(m.vars(), [&](lp::lpvar v) { return m_lp2nl.contains(v); })) add_monic_eq_bound(m); - for (unsigned i : m_term_set) + for (unsigned i : m_coi.terms()) add_term(i); for (auto const& [v, w] : m_lp2nl) { if (lra.column_has_lower_bound(v)) @@ -418,6 +449,7 @@ struct solver::imp { ex.push_back(ci); nla::lemma_builder lemma(m_nla_core, __FUNCTION__); lemma &= ex; + TRACE(nra, tout << lemma << "\n"); break; } case l_undef: @@ -554,8 +586,8 @@ struct solver::imp { if (!m_lp2nl.find(v, r)) { r = m_nlsat->mk_var(is_int(v)); m_lp2nl.insert(v, r); - if (!m_term_set.contains(v) && lra.column_has_term(v)) { - m_term_set.insert(v); + if (!m_coi.terms().contains(v) && lra.column_has_term(v)) { + m_coi.terms().insert(v); } } return r; @@ -586,20 +618,55 @@ struct solver::imp { m_nlsat->mk_clause(1, &lit, nullptr); } - nlsat::anum const& value(lp::lpvar v) { - polynomial::var pv; - if (m_lp2nl.find(v, pv)) - return m_nlsat->value(pv); - else { - for (unsigned w = m_values->size(); w <= v; ++w) { - scoped_anum a(am()); - am().set(a, m_nla_core.val(w).to_mpq()); + nlsat::anum const &value(lp::lpvar v) { + init_values(v + 1); + return (*m_values)[v]; + } + + void init_values(unsigned sz) { + if (m_values->size() >= sz) + return; + unsigned w; + scoped_anum a(am()); + for (unsigned v = m_values->size(); v < sz; ++v) { + if (m_nla_core.emons().is_monic_var(v)) { + am().set(a, 1); + auto &m = m_nla_core.emon(v); + for (auto x : m.vars()) + am().mul(a, (*m_values)[x], a); m_values->push_back(a); - } - return (*m_values)[v]; + } + else if (lra.column_has_term(v)) { + scoped_anum b(am()); + am().set(a, 0); + for (auto const &[w, coeff] : lra.get_term(v)) { + am().set(b, coeff.to_mpq()); + am().mul(b, (*m_values)[w], b); + am().add(a, b, a); + } + m_values->push_back(a); + } + else if (m_lp2nl.find(v, w)) { + m_values->push_back(m_nlsat->value(w)); + } + else { + am().set(a, m_nla_core.val(v).to_mpq()); + m_values->push_back(a); + } } } + + void set_value(lp::lpvar v, rational const& value) { + if (!m_values) + m_values = alloc(scoped_anum_vector, am()); + scoped_anum a(am()); + am().set(a, value.to_mpq()); + while (m_values->size() <= v) + m_values->push_back(a); + am().set((*m_values)[v], a); + } + nlsat::anum_manager& am() { return m_nlsat->am(); } @@ -680,4 +747,8 @@ void solver::updt_params(params_ref& p) { m_imp->updt_params(p); } +void solver::set_value(lp::lpvar v, rational const& value) { + m_imp->set_value(v, value); +} + } diff --git a/src/math/lp/nra_solver.h b/src/math/lp/nra_solver.h index 90f022ba6..b009b3c12 100644 --- a/src/math/lp/nra_solver.h +++ b/src/math/lp/nra_solver.h @@ -59,6 +59,8 @@ namespace nra { nlsat::anum_manager& am(); + void set_value(lp::lpvar v, rational const &value); + scoped_anum& tmp1(); scoped_anum& tmp2();