From e31cecf5dbe7b6d0542434451e8ed64497bd3075 Mon Sep 17 00:00:00 2001 From: Lev Nachmanson Date: Thu, 21 Sep 2023 11:27:53 -0700 Subject: [PATCH] transfer propagate monomial bounds to nla_solver --- src/math/lp/lp_bound_propagator.h | 165 ----------------------- src/math/lp/nla_core.cpp | 213 +++++++++++++++++++++++++++++- src/math/lp/nla_core.h | 22 ++- src/math/lp/nla_solver.cpp | 11 +- src/math/lp/nla_solver.h | 4 +- src/sat/smt/arith_internalize.cpp | 2 +- src/smt/theory_lra.cpp | 9 +- 7 files changed, 243 insertions(+), 183 deletions(-) diff --git a/src/math/lp/lp_bound_propagator.h b/src/math/lp/lp_bound_propagator.h index 99755f606..7ed872aa2 100644 --- a/src/math/lp/lp_bound_propagator.h +++ b/src/math/lp/lp_bound_propagator.h @@ -120,171 +120,6 @@ public: m_column_types = &lp().get_column_types(); } - bool is_linear(const svector& m, lpvar& zero_var, lpvar& non_fixed) { - zero_var = non_fixed = null_lpvar; - unsigned n_of_non_fixed = 0; - for (lpvar v : m) { - if (!this->column_is_fixed(v)) { - n_of_non_fixed++; - non_fixed = v; - continue; - } - const auto & b = get_lower_bound(v).x; - if (b.is_zero()) { - zero_var = v; - return true; - } - - } - return n_of_non_fixed <= 1; - } - - void add_bounds_for_zero_var(lpvar monic_var, lpvar zero_var) { - auto& lps = lp(); - auto lambda = [zero_var,&lps]() { - return lps.get_bound_constraint_witnesses_for_column(zero_var); - }; - TRACE("add_bound", lp().print_column_info(zero_var, tout) << std::endl;); - add_lower_bound_monic(monic_var, mpq(0), false, lambda); - add_upper_bound_monic(monic_var, mpq(0), false, lambda); - } - - void add_lower_bound_monic(lpvar j, const mpq& v, bool is_strict, std::function explain_dep) { - TRACE("add_bound", lp().print_column_info(j, tout) << std::endl;); - j = lp().column_to_reported_index(j); - unsigned k; - if (!m_improved_lower_bounds.find(j, k)) { - m_improved_lower_bounds.insert(j,static_cast(m_ibounds.size())); - m_ibounds.push_back(implied_bound(v, j, true, is_strict, explain_dep)); - } - else { - auto& found_bound = m_ibounds[k]; - if (v > found_bound.m_bound || (v == found_bound.m_bound && !found_bound.m_strict && is_strict)) { - found_bound = implied_bound(v, j, true, is_strict, explain_dep); - TRACE("add_bound", lp().print_implied_bound(found_bound, tout);); - } - } - } - - void add_upper_bound_monic(lpvar j, const mpq& bound_val, bool is_strict, std::function explain_bound) { - j = lp().column_to_reported_index(j); - unsigned k; - if (!m_improved_upper_bounds.find(j, k)) { - m_improved_upper_bounds.insert(j, static_cast(m_ibounds.size())); - m_ibounds.push_back(implied_bound(bound_val, j, false, is_strict, explain_bound)); - } - else { - auto& found_bound = m_ibounds[k]; - if (bound_val > found_bound.m_bound || (bound_val == found_bound.m_bound && !found_bound.m_strict && is_strict)) { - found_bound = implied_bound(bound_val, j, false, is_strict, explain_bound); - TRACE("add_bound", lp().print_implied_bound(found_bound, tout);); - } - } - } - - void propagate_monic(lpvar monic_var, const svector& vars) { - lpvar non_fixed, zero_var; - if (!is_linear(vars, zero_var, non_fixed)) - return; - - if (zero_var != null_lpvar) - add_bounds_for_zero_var(monic_var, zero_var); - else { - rational k = rational(1); - for (auto v : vars) - if (v != non_fixed) { - k *= lp().get_column_value(v).x; - if (k.is_big()) return; - } - - if (non_fixed != null_lpvar) - propagate_monic_with_non_fixed(monic_var, vars, non_fixed, k); - else // all variables are fixed - propagate_monic_with_all_fixed(monic_var, vars, k); - } - } - - void propagate_monic_with_non_fixed(lpvar monic_var, const svector& vars, lpvar non_fixed, const rational& k) { - lp::impq bound_value; - bool is_strict; - auto& lps = lp(); - - if (lower_bound_is_available(non_fixed)) { - bound_value = lp().column_lower_bound(non_fixed); - is_strict = !bound_value.y.is_zero(); - auto lambda = [vars, non_fixed,&lps]() { - u_dependency* dep = lps.get_column_lower_bound_witness(non_fixed); - for (auto v : vars) - if (v != non_fixed) - dep = lps.join_deps(dep, lps.get_bound_constraint_witnesses_for_column(v)); - return dep; - }; - if (k.is_pos()) - add_lower_bound_monic(monic_var, k * bound_value.x, is_strict, lambda); - else - add_upper_bound_monic(monic_var, k * bound_value.x, is_strict, lambda); - } - - if (upper_bound_is_available(non_fixed)) { - bound_value = lp().column_upper_bound(non_fixed); - is_strict = !bound_value.y.is_zero(); - auto lambda = [vars, non_fixed,&lps]() { - u_dependency* dep = lps.get_column_upper_bound_witness(non_fixed); - for (auto v : vars) - if (v != non_fixed) - dep = lps.join_deps(dep, lps.get_bound_constraint_witnesses_for_column(v)); - return dep; - }; - if (k.is_neg()) - add_lower_bound_monic(monic_var, k * bound_value.x, is_strict, lambda); - else - add_upper_bound_monic(monic_var, k * bound_value.x, is_strict, lambda); - } - - if (lower_bound_is_available(monic_var)) { - auto lambda = [vars, monic_var, non_fixed,&lps]() { - u_dependency* dep = lps.get_column_lower_bound_witness(monic_var); - for (auto v : vars) { - if (v != non_fixed) { - dep = lps.join_deps(dep, lps.get_bound_constraint_witnesses_for_column(v)); - } - } - return dep; - }; - bound_value = lp().column_lower_bound(monic_var); - is_strict = !bound_value.y.is_zero(); - if (k.is_pos()) - add_lower_bound_monic(non_fixed, bound_value.x / k, is_strict, lambda); - else - add_upper_bound_monic(non_fixed, bound_value.x / k, is_strict, lambda); - } - - if (upper_bound_is_available(monic_var)) { - bound_value = lp().column_upper_bound(monic_var); - is_strict = !bound_value.y.is_zero(); - auto lambda = [vars, monic_var, non_fixed,&lps]() { - u_dependency* dep = lps.get_column_upper_bound_witness(monic_var); - for (auto v : vars) { - if (v != non_fixed) { - dep = lps.join_deps(dep, lps.get_bound_constraint_witnesses_for_column(v)); - } - } - return dep; - }; - if (k.is_neg()) - add_lower_bound_monic(non_fixed, bound_value.x / k, is_strict, lambda); - else - add_upper_bound_monic(non_fixed, bound_value.x / k, is_strict, lambda); - } - } - - void propagate_monic_with_all_fixed(lpvar monic_var, const svector& vars, const rational& k) { - auto& lps = lp(); - auto lambda = [vars,&lps]() { return lps.get_bound_constraint_witnesses_for_columns(vars); }; - add_lower_bound_monic(monic_var, k, false, lambda); - add_upper_bound_monic(monic_var, k, false, lambda); - } - column_type get_column_type(unsigned j) const { return (*m_column_types)[j]; } diff --git a/src/math/lp/nla_core.cpp b/src/math/lp/nla_core.cpp index 63dd29d9b..4e16793ec 100644 --- a/src/math/lp/nla_core.cpp +++ b/src/math/lp/nla_core.cpp @@ -17,11 +17,12 @@ Author: #include "math/grobner/pdd_solver.h" #include "math/dd/pdd_interval.h" #include "math/dd/pdd_eval.h" +#include "nla_core.h" namespace nla { typedef lp::lar_term term; -core::core(lp::lar_solver& s, params_ref const& p, reslimit& lim) : m_evars(), +core::core(lp::lar_solver& s, params_ref const& p, reslimit& lim, std_vector& implied_bounds) : m_evars(), lra(s), m_reslim(lim), m_params(p), @@ -37,7 +38,8 @@ core::core(lp::lar_solver& s, params_ref const& p, reslimit& lim) : m_evars(), m_grobner(this), m_emons(m_evars), m_use_nra_model(false), - m_nra(s, m_nra_lim, *this) { + m_nra(s, m_nra_lim, *this), + m_implied_bounds(implied_bounds) { m_nlsat_delay = lp_settings().nlsat_delay(); lra.m_find_monics_with_changed_bounds_func = [&](const indexed_uint_set& columns_with_changed_bounds) { for (const auto& m : m_emons) { @@ -1837,5 +1839,210 @@ bool core::improve_bounds() { return bounds_improved; } -} // end of nla +bool core::is_linear(const svector& m, lpvar& zero_var, lpvar& non_fixed) +{ + zero_var = non_fixed = null_lpvar; + unsigned n_of_non_fixed = 0; + for (lpvar v : m) { + if (!this->var_is_fixed(v)) { + n_of_non_fixed++; + non_fixed = v; + continue; + } + const auto& b = get_lower_bound(v); + if (b.is_zero()) { + zero_var = v; + return true; + } + } + return n_of_non_fixed <= 1; +} + +void core::add_lower_bound_monic(lpvar j, const lp::mpq& v, bool is_strict, std::function explain_dep) +{ + TRACE("add_bound", lra.print_column_info(j, tout) << std::endl;); + j = lra.column_to_reported_index(j); + unsigned k; + if (!m_improved_lower_bounds.find(j, k)) { + m_improved_lower_bounds.insert(j, static_cast(m_implied_bounds.size())); + m_implied_bounds.push_back(lp::implied_bound(v, j, true, is_strict, explain_dep)); + } + else { + auto& found_bound = m_implied_bounds[k]; + if (v > found_bound.m_bound || (v == found_bound.m_bound && !found_bound.m_strict && is_strict)) { + found_bound = lp::implied_bound(v, j, true, is_strict, explain_dep); + TRACE("add_bound", lra.print_implied_bound(found_bound, tout);); + } + } +} + +void core::add_upper_bound_monic(lpvar j, const lp::mpq& bound_val, bool is_strict, std::function explain_dep) +{ + j = lra.column_to_reported_index(j); + unsigned k; + if (!m_improved_upper_bounds.find(j, k)) { + m_improved_upper_bounds.insert(j, static_cast(m_implied_bounds.size())); + m_implied_bounds.push_back(lp::implied_bound(bound_val, j, false, is_strict, explain_dep)); + } + else { + auto& found_bound = m_implied_bounds[k]; + if (bound_val > found_bound.m_bound || (bound_val == found_bound.m_bound && !found_bound.m_strict && is_strict)) { + found_bound = lp::implied_bound(bound_val, j, false, is_strict, explain_dep); + TRACE("add_bound", lra.print_implied_bound(found_bound, tout);); + } + } +} + +bool core::upper_bound_is_available(unsigned j) const +{ + switch (get_column_type(j)) { + case lp::column_type::fixed: + case lp::column_type::boxed: + case lp::column_type::upper_bound: + return true; + default: + return false; + } +} + +bool core::lower_bound_is_available(unsigned j) const +{ + switch (get_column_type(j)) { + case lp::column_type::fixed: + case lp::column_type::boxed: + case lp::column_type::lower_bound: + return true; + default: + return false; + } +} + +void core::propagate_monic_with_non_fixed(lpvar monic_var, const svector& vars, lpvar non_fixed, const rational& k) +{ + lp::impq bound_value; + bool is_strict; + auto& lps = lra; + + if (lower_bound_is_available(non_fixed)) { + bound_value = lra.column_lower_bound(non_fixed); + is_strict = !bound_value.y.is_zero(); + auto lambda = [vars, non_fixed, &lps]() { + u_dependency* dep = lps.get_column_lower_bound_witness(non_fixed); + for (auto v : vars) + if (v != non_fixed) + dep = lps.join_deps(dep, lps.get_bound_constraint_witnesses_for_column(v)); + return dep; + }; + if (k.is_pos()) + add_lower_bound_monic(monic_var, k * bound_value.x, is_strict, lambda); + else + add_upper_bound_monic(monic_var, k * bound_value.x, is_strict, lambda); + } + + if (upper_bound_is_available(non_fixed)) { + bound_value = lra.column_upper_bound(non_fixed); + is_strict = !bound_value.y.is_zero(); + auto lambda = [vars, non_fixed, &lps]() { + u_dependency* dep = lps.get_column_upper_bound_witness(non_fixed); + for (auto v : vars) + if (v != non_fixed) + dep = lps.join_deps(dep, lps.get_bound_constraint_witnesses_for_column(v)); + return dep; + }; + if (k.is_neg()) + add_lower_bound_monic(monic_var, k * bound_value.x, is_strict, lambda); + else + add_upper_bound_monic(monic_var, k * bound_value.x, is_strict, lambda); + } + + if (lower_bound_is_available(monic_var)) { + auto lambda = [vars, monic_var, non_fixed, &lps]() { + u_dependency* dep = lps.get_column_lower_bound_witness(monic_var); + for (auto v : vars) { + if (v != non_fixed) { + dep = lps.join_deps(dep, lps.get_bound_constraint_witnesses_for_column(v)); + } + } + return dep; + }; + bound_value = lra.column_lower_bound(monic_var); + is_strict = !bound_value.y.is_zero(); + if (k.is_pos()) + add_lower_bound_monic(non_fixed, bound_value.x / k, is_strict, lambda); + else + add_upper_bound_monic(non_fixed, bound_value.x / k, is_strict, lambda); + } + + if (upper_bound_is_available(monic_var)) { + bound_value = lra.column_upper_bound(monic_var); + is_strict = !bound_value.y.is_zero(); + auto lambda = [vars, monic_var, non_fixed, &lps]() { + u_dependency* dep = lps.get_column_upper_bound_witness(monic_var); + for (auto v : vars) { + if (v != non_fixed) { + dep = lps.join_deps(dep, lps.get_bound_constraint_witnesses_for_column(v)); + } + } + return dep; + }; + if (k.is_neg()) + add_lower_bound_monic(non_fixed, bound_value.x / k, is_strict, lambda); + else + add_upper_bound_monic(non_fixed, bound_value.x / k, is_strict, lambda); + } +} + +void core::propagate_monic_with_all_fixed(lpvar monic_var, const svector& vars, const rational& k) +{ + auto* lps = &lra; + auto lambda = [vars, lps]() { return lps->get_bound_constraint_witnesses_for_columns(vars); }; + add_lower_bound_monic(monic_var, k, false, lambda); + add_upper_bound_monic(monic_var, k, false, lambda); +} + +void core::add_bounds_for_zero_var(lpvar monic_var, lpvar zero_var) +{ + auto* lps = &lra; + auto lambda = [zero_var, lps]() { + return lps->get_bound_constraint_witnesses_for_column(zero_var); + }; + TRACE("add_bound", lra.print_column_info(zero_var, tout) << std::endl;); + add_lower_bound_monic(monic_var, lp::mpq(0), false, lambda); + add_upper_bound_monic(monic_var, lp::mpq(0), false, lambda); +} + +void core::calculate_implied_bounds_for_monic(lp::lpvar monic_var) +{ + lpvar non_fixed, zero_var; + const auto& vars = m_emons[monic_var].vars(); + if (!is_linear(vars, zero_var, non_fixed)) + return; + + if (zero_var != null_lpvar) + add_bounds_for_zero_var(monic_var, zero_var); + else { + rational k = rational(1); + for (auto v : vars) + if (v != non_fixed) { + k *= val(v); + if (k.is_big()) return; + } + + if (non_fixed != null_lpvar) + propagate_monic_with_non_fixed(monic_var, vars, non_fixed, k); + else // all variables are fixed + propagate_monic_with_all_fixed(monic_var, vars, k); + } +} + +void core::init_bound_propagation() +{ + this->m_implied_bounds.clear(); + this->m_improved_lower_bounds.reset(); + this->m_improved_upper_bounds.reset(); + this->m_column_types = &lra.get_column_types(); +} + +} // namespace nla + diff --git a/src/math/lp/nla_core.h b/src/math/lp/nla_core.h index 3b888f8ef..54c7e1df1 100644 --- a/src/math/lp/nla_core.h +++ b/src/math/lp/nla_core.h @@ -103,7 +103,10 @@ class core { emonics m_emons; svector m_add_buffer; mutable indexed_uint_set m_active_var_set; - + // these maps map a column index to the corresponding index in ibounds + u_map m_improved_lower_bounds; + u_map m_improved_upper_bounds; + const vector* m_column_types; reslimit m_nra_lim; bool m_use_nra_model = false; @@ -114,12 +117,13 @@ class core { void check_weighted(unsigned sz, std::pair>* checks); void add_bounds(); + std_vector & m_implied_bounds; // try to improve bounds for variables in monomials. bool improve_bounds(); public: // constructor - core(lp::lar_solver& s, params_ref const& p, reslimit&); + core(lp::lar_solver& s, params_ref const& p, reslimit&, std_vector & implied_bounds); const auto& monics_with_changed_bounds() const { return m_monics_with_changed_bounds; } void reset_monics_with_changed_bounds() { m_monics_with_changed_bounds.reset(); } void insert_to_refine(lpvar j); @@ -431,15 +435,23 @@ public: void set_use_nra_model(bool m); bool use_nra_model() const { return m_use_nra_model; } void collect_statistics(::statistics&); + bool is_linear(const svector& m, lpvar& zero_var, lpvar& non_fixed); + void add_bounds_for_zero_var(lpvar monic_var, lpvar zero_var); + void propagate_monic_with_non_fixed(lpvar monic_var, const svector& vars, lpvar non_fixed, const rational& k); + void propagate_monic_with_all_fixed(lpvar monic_var, const svector& vars, const rational& k); + void add_lower_bound_monic(lpvar j, const lp::mpq& v, bool is_strict, std::function explain_dep); + void add_upper_bound_monic(lpvar j, const lp::mpq& v, bool is_strict, std::function explain_dep); + bool upper_bound_is_available(unsigned j) const; + bool lower_bound_is_available(unsigned j) const; private: - void restore_patched_values(); + lp::column_type get_column_type(unsigned j) const { return (*m_column_types)[j]; } void constrain_nl_in_tableau(); bool solve_tableau(); void restore_tableau(); void save_tableau(); bool integrality_holds(); - - + void calculate_implied_bounds_for_monic(lp::lpvar v); + void init_bound_propagation(); }; // end of core struct pp_mon { diff --git a/src/math/lp/nla_solver.cpp b/src/math/lp/nla_solver.cpp index 8be918f01..3475a3509 100644 --- a/src/math/lp/nla_solver.cpp +++ b/src/math/lp/nla_solver.cpp @@ -54,8 +54,8 @@ namespace nla { m_core->pop(n); } - solver::solver(lp::lar_solver& s, params_ref const& p, reslimit& limit): - m_core(alloc(core, s, p, limit)) { + solver::solver(lp::lar_solver& s, params_ref const& p, reslimit& limit, std_vector & implied_bounds): + m_core(alloc(core, s, p, limit, implied_bounds)) { } bool solver::influences_nl_var(lpvar j) const { @@ -88,6 +88,9 @@ namespace nla { m_core->collect_statistics(st); } + void solver::calculate_implied_bounds_for_monic(lp::lpvar v) { + m_core->calculate_implied_bounds_for_monic(v); + } // ensure r = x^y, add abstraction/refinement lemmas lbool solver::check_power(lpvar r, lpvar x, lpvar y, vector& lemmas) { return m_core->check_power(r, x, y, lemmas); @@ -97,4 +100,8 @@ namespace nla { m_core->check_bounded_divisions(lemmas); } + void solver::init_bound_propagation() { + m_core->init_bound_propagation(); + } + } diff --git a/src/math/lp/nla_solver.h b/src/math/lp/nla_solver.h index 07bf095a6..a4de90320 100644 --- a/src/math/lp/nla_solver.h +++ b/src/math/lp/nla_solver.h @@ -24,7 +24,7 @@ namespace nla { core* m_core; public: - solver(lp::lar_solver& s, params_ref const& p, reslimit& limit); + solver(lp::lar_solver& s, params_ref const& p, reslimit& limit, std_vector & implied_bounds); ~solver(); const auto& monics_with_changed_bounds() const { return m_core->monics_with_changed_bounds(); } void reset_monics_with_changed_bounds() { m_core->reset_monics_with_changed_bounds(); } @@ -48,5 +48,7 @@ namespace nla { nlsat::anum_manager& am(); nlsat::anum const& am_value(lp::var_index v) const; void collect_statistics(::statistics & st); + void calculate_implied_bounds_for_monic(lp::lpvar v); + void init_bound_propagation(); }; } diff --git a/src/sat/smt/arith_internalize.cpp b/src/sat/smt/arith_internalize.cpp index 3174ad775..5893c8520 100644 --- a/src/sat/smt/arith_internalize.cpp +++ b/src/sat/smt/arith_internalize.cpp @@ -61,7 +61,7 @@ namespace arith { void solver::ensure_nla() { if (!m_nla) { - m_nla = alloc(nla::solver, *m_solver.get(), s().params(), m.limit()); + m_nla = alloc(nla::solver, *m_solver.get(), s().params(), m.limit(), m_implied_bounds); for (auto const& _s : m_scopes) { (void)_s; m_nla->push(); diff --git a/src/smt/theory_lra.cpp b/src/smt/theory_lra.cpp index aeb138d3c..ef7ab7091 100644 --- a/src/smt/theory_lra.cpp +++ b/src/smt/theory_lra.cpp @@ -264,7 +264,7 @@ class theory_lra::imp { void ensure_nla() { if (!m_nla) { - m_nla = alloc(nla::solver, *m_solver.get(), ctx().get_params(), m.limit()); + m_nla = alloc(nla::solver, *m_solver.get(), ctx().get_params(), m.limit(), m_implied_bounds); for (auto const& _s : m_scopes) { (void)_s; m_nla->push(); @@ -2198,13 +2198,10 @@ public: finish_bound_propagation(); } - void calculate_implied_bounds_for_monic(lpvar monic_var, const svector& vars) { - m_bp.propagate_monic(monic_var, vars); - } - void propagate_bounds_for_touched_monomials() { + m_nla->init_bound_propagation(); for (unsigned v : m_nla->monics_with_changed_bounds()) { - calculate_implied_bounds_for_monic(v, m_nla->get_core().emons()[v].vars()); + m_nla->calculate_implied_bounds_for_monic(v); } m_nla->reset_monics_with_changed_bounds(); }