From 9c18ede68724132b148f6c75ea92242f6647aecf Mon Sep 17 00:00:00 2001 From: Lev Nachmanson Date: Wed, 5 Jun 2019 15:26:11 -0700 Subject: [PATCH] hook up nla_solver it lp bound propagation Signed-off-by: Lev Nachmanson --- src/math/lp/bound_analyzer_on_row.h | 89 ++++++------------- src/math/lp/emonomials.cpp | 1 + src/math/lp/lar_solver.cpp | 14 +-- src/math/lp/lar_solver.h | 16 ++-- src/math/lp/lp_bound_propagator.cpp | 76 ++++++++++++++-- ...und_propagator.h => lp_bound_propagator.h} | 18 +++- src/math/lp/nla_core.cpp | 15 ++++ src/math/lp/nla_core.h | 6 +- src/math/lp/nla_solver.cpp | 20 +++++ src/math/lp/nla_solver.h | 7 +- src/smt/theory_lra.cpp | 8 +- 11 files changed, 177 insertions(+), 93 deletions(-) rename src/math/lp/{bound_propagator.h => lp_bound_propagator.h} (66%) diff --git a/src/math/lp/bound_analyzer_on_row.h b/src/math/lp/bound_analyzer_on_row.h index 22c40a753..528980f74 100644 --- a/src/math/lp/bound_analyzer_on_row.h +++ b/src/math/lp/bound_analyzer_on_row.h @@ -21,7 +21,7 @@ Revision History: #include "util/vector.h" #include "implied_bound.h" #include "test_bound_analyzer.h" -#include "math/lp/bound_propagator.h" +#include "math/lp/lp_bound_propagator.h" // We have an equality : sum by j of row[j]*x[j] = rs // We try to pin a var by pushing the total by using the variable bounds // In a loop we drive the partial sum down, denoting the variables of this process by _u. @@ -30,7 +30,7 @@ namespace lp { template // C plays a role of a container class bound_analyzer_on_row { const C& m_row; - bound_propagator & m_bp; + lp_bound_propagator & m_bp; unsigned m_row_or_term_index; int m_column_of_u; // index of an unlimited from above monoid // -1 means that such a value is not found, -2 means that at least two of such monoids were found @@ -44,7 +44,7 @@ public : unsigned bj, // basis column for the row const numeric_pair& rs, unsigned row_or_term_index, - bound_propagator & bp + lp_bound_propagator & bp ) : m_row(it), @@ -55,8 +55,6 @@ public : m_rs(rs) {} - - unsigned j; void analyze() { for (const auto & c : m_row) { if ((m_column_of_l == -2) && (m_column_of_u == -2)) @@ -80,86 +78,53 @@ public : } bool upper_bound_is_available(unsigned j) const { - switch (m_bp.get_column_type(j)) - { - case column_type::fixed: - case column_type::boxed: - case column_type::upper_bound: - return true; - default: - return false; - } + return m_bp.upper_bound_is_available(j); } bool lower_bound_is_available(unsigned j) const { - switch (m_bp.get_column_type(j)) - { - case column_type::fixed: - case column_type::boxed: - case column_type::lower_bound: - return true; - default: - return false; - } + return m_bp.lower_bound_is_available(j); } - const impq & ub(unsigned j) const { + impq ub(unsigned j) const { lp_assert(upper_bound_is_available(j)); return m_bp.get_upper_bound(j); } - const impq & lb(unsigned j) const { + impq lb(unsigned j) const { lp_assert(lower_bound_is_available(j)); return m_bp.get_lower_bound(j); } + mpq u_strict(unsigned j, bool & strict) const { + const impq u = ub(j); + strict = !is_zero(u.y); + return u.x; + } - const mpq & monoid_max_no_mult(bool a_is_pos, unsigned j, bool & strict) const { - if (a_is_pos) { - strict = !is_zero(ub(j).y); - return ub(j).x; - } - strict = !is_zero(lb(j).y); - return lb(j).x; + mpq l_strict(unsigned j, bool & strict) const { + const impq l = lb(j); + strict = !is_zero(l.y); + return l.x; + } + + mpq monoid_max_no_mult(bool a_is_pos, unsigned j, bool & strict) const { + return a_is_pos? u_strict(j, strict) : l_strict(j, strict); } mpq monoid_max(const mpq & a, unsigned j) const { - if (is_pos(a)) { - return a * ub(j).x; - } - return a * lb(j).x; + return a * (is_pos(a) ? ub(j).x : lb(j).x); } mpq monoid_max(const mpq & a, unsigned j, bool & strict) const { - if (is_pos(a)) { - strict = !is_zero(ub(j).y); - return a * ub(j).x; - } - strict = !is_zero(lb(j).y); - return a * lb(j).x; + return a * (is_pos(a)? u_strict(j, strict) : l_strict(j, strict)); } - const mpq & monoid_min_no_mult(bool a_is_pos, unsigned j, bool & strict) const { - if (!a_is_pos) { - strict = !is_zero(ub(j).y); - return ub(j).x; - } - strict = !is_zero(lb(j).y); - return lb(j).x; + mpq monoid_min_no_mult(bool a_is_pos, unsigned j, bool & strict) const { + return a_is_pos? l_strict(j, strict) : u_strict(j, strict); } mpq monoid_min(const mpq & a, unsigned j, bool& strict) const { - if (is_neg(a)) { - strict = !is_zero(ub(j).y); - return a * ub(j).x; - } - - strict = !is_zero(lb(j).y); - return a * lb(j).x; + return a * (is_neg(a)? u_strict(j, strict) : l_strict(j, strict)); } mpq monoid_min(const mpq & a, unsigned j) const { - if (is_neg(a)) { - return a * ub(j).x; - } - - return a * lb(j).x; + return a * (is_neg(a)? ub(j).x : lb(j).x ); } @@ -338,7 +303,7 @@ public : unsigned bj, // basis column for the row const numeric_pair& rs, unsigned row_or_term_index, - bound_propagator & bp + lp_bound_propagator & bp ) { bound_analyzer_on_row a(row, bj, rs, row_or_term_index, bp); a.analyze(); diff --git a/src/math/lp/emonomials.cpp b/src/math/lp/emonomials.cpp index dd73cc171..3f93a7309 100644 --- a/src/math/lp/emonomials.cpp +++ b/src/math/lp/emonomials.cpp @@ -282,6 +282,7 @@ bool emonomials::is_visited(monomial const& m) const { class of equal up-to var_eqs monomials. */ void emonomials::add(lpvar v, unsigned sz, lpvar const* vs) { + TRACE("nla_solver", tout << "v = " << v << "\n";); unsigned idx = m_monomials.size(); m_monomials.push_back(monomial(v, sz, vs, idx)); lpvar last_var = UINT_MAX; diff --git a/src/math/lp/lar_solver.cpp b/src/math/lp/lar_solver.cpp index 53ef1f84f..92613579a 100644 --- a/src/math/lp/lar_solver.cpp +++ b/src/math/lp/lar_solver.cpp @@ -158,7 +158,7 @@ bool lar_solver::implied_bound_is_correctly_explained(implied_bound const & be, void lar_solver::analyze_new_bounds_on_row( unsigned row_index, - bound_propagator & bp) { + lp_bound_propagator & bp) { lp_assert(!use_tableau()); unsigned j = m_mpq_lar_core_solver.m_r_basis[row_index]; // basis column for the row bound_analyzer_on_row> @@ -173,7 +173,7 @@ void lar_solver::analyze_new_bounds_on_row( void lar_solver::analyze_new_bounds_on_row_tableau( unsigned row_index, - bound_propagator & bp ) { + lp_bound_propagator & bp ) { if (A_r().m_rows[row_index].size() > settings().max_row_length_for_bound_propagation) return; @@ -199,7 +199,7 @@ void lar_solver::substitute_basis_var_in_terms_for_row(unsigned i) { } } -void lar_solver::calculate_implied_bounds_for_row(unsigned i, bound_propagator & bp) { +void lar_solver::calculate_implied_bounds_for_row(unsigned i, lp_bound_propagator & bp) { if(use_tableau()) { analyze_new_bounds_on_row_tableau(i, bp); } else { @@ -219,12 +219,12 @@ unsigned lar_solver::map_term_index_to_column_index(unsigned j) const { return m_var_register.external_to_local(j); } -void lar_solver::propagate_bounds_on_a_term(const lar_term& t, bound_propagator & bp, unsigned term_offset) { +void lar_solver::propagate_bounds_on_a_term(const lar_term& t, lp_bound_propagator & bp, unsigned term_offset) { lp_assert(false); // not implemented } -void lar_solver::explain_implied_bound(implied_bound & ib, bound_propagator & bp) { +void lar_solver::explain_implied_bound(implied_bound & ib, lp_bound_propagator & bp) { unsigned i = ib.m_row_or_term_index; int bound_sign = ib.m_is_lower_bound? 1: -1; int j_sign = (ib.m_coeff_before_j_is_pos ? 1 :-1) * bound_sign; @@ -252,7 +252,7 @@ bool lar_solver::term_is_used_as_row(unsigned term) const { return m_var_register.external_is_used(term); } -void lar_solver::propagate_bounds_on_terms(bound_propagator & bp) { +void lar_solver::propagate_bounds_on_terms(lp_bound_propagator & bp) { for (unsigned i = 0; i < m_terms.size(); i++) { if (term_is_used_as_row(i + m_terms_start_index)) continue; // this term is used a left side of a constraint, @@ -263,7 +263,7 @@ void lar_solver::propagate_bounds_on_terms(bound_propagator & bp) { // goes over touched rows and tries to induce bounds -void lar_solver::propagate_bounds_for_touched_rows(bound_propagator & bp) { +void lar_solver::propagate_bounds_for_touched_rows(lp_bound_propagator & bp) { if (!use_tableau()) return; // todo: consider to remove the restriction diff --git a/src/math/lp/lar_solver.h b/src/math/lp/lar_solver.h index b3aa2ab7b..74516db64 100644 --- a/src/math/lp/lar_solver.h +++ b/src/math/lp/lar_solver.h @@ -41,7 +41,7 @@ Revision History: #include "math/lp/conversion_helper.h" #include "math/lp/int_solver.h" #include "math/lp/nra_solver.h" -#include "math/lp/bound_propagator.h" +#include "math/lp/lp_bound_propagator.h" namespace lp { @@ -264,16 +264,16 @@ public: void analyze_new_bounds_on_row( unsigned row_index, - bound_propagator & bp); + lp_bound_propagator & bp); void analyze_new_bounds_on_row_tableau( unsigned row_index, - bound_propagator & bp); + lp_bound_propagator & bp); void substitute_basis_var_in_terms_for_row(unsigned i); - void calculate_implied_bounds_for_row(unsigned i, bound_propagator & bp); + void calculate_implied_bounds_for_row(unsigned i, lp_bound_propagator & bp); unsigned adjust_column_index_to_term_index(unsigned j) const; @@ -313,19 +313,19 @@ public: } - void propagate_bounds_on_a_term(const lar_term& t, bound_propagator & bp, unsigned term_offset); + void propagate_bounds_on_a_term(const lar_term& t, lp_bound_propagator & bp, unsigned term_offset); - void explain_implied_bound(implied_bound & ib, bound_propagator & bp); + void explain_implied_bound(implied_bound & ib, lp_bound_propagator & bp); bool term_is_used_as_row(unsigned term) const; - void propagate_bounds_on_terms(bound_propagator & bp); + void propagate_bounds_on_terms(lp_bound_propagator & bp); // goes over touched rows and tries to induce bounds - void propagate_bounds_for_touched_rows(bound_propagator & bp); + void propagate_bounds_for_touched_rows(lp_bound_propagator & bp); lp_status get_status() const; diff --git a/src/math/lp/lp_bound_propagator.cpp b/src/math/lp/lp_bound_propagator.cpp index ee25d8dd0..afe9d4371 100644 --- a/src/math/lp/lp_bound_propagator.cpp +++ b/src/math/lp/lp_bound_propagator.cpp @@ -3,19 +3,79 @@ Author: Lev Nachmanson */ #include "math/lp/lar_solver.h" +#include "math/lp/nla_solver.h" namespace lp { -bound_propagator::bound_propagator(lar_solver & ls): - m_lar_solver(ls) {} -column_type bound_propagator::get_column_type(unsigned j) const { + +lp_bound_propagator::lp_bound_propagator(lar_solver & ls, nla::solver* nla): + m_lar_solver(ls), m_nla_solver(nla) {} + +column_type lp_bound_propagator::get_column_type(unsigned j) const { return m_lar_solver.m_mpq_lar_core_solver.m_column_types()[j]; } -const impq & bound_propagator::get_lower_bound(unsigned j) const { - return m_lar_solver.m_mpq_lar_core_solver.m_r_lower_bounds()[j]; + +bool lp_bound_propagator::lower_bound_is_available(unsigned j) const { + if (lower_bound_is_available_for_column(j)) + return true; + if (!m_nla_solver->is_monomial_var(j)) + return false; + return m_nla_solver->monomial_has_lower_bound(j); } -const impq & bound_propagator::get_upper_bound(unsigned j) const { - return m_lar_solver.m_mpq_lar_core_solver.m_r_upper_bounds()[j]; + +bool lp_bound_propagator::upper_bound_is_available_for_column(unsigned j) const { + switch (get_column_type(j)) { + case column_type::fixed: + case column_type::boxed: + case column_type::upper_bound: + return true; + default: + return false; + } } -void bound_propagator::try_add_bound(mpq v, unsigned j, bool is_low, bool coeff_before_j_is_pos, unsigned row_or_term_index, bool strict) { + +bool lp_bound_propagator::lower_bound_is_available_for_column(unsigned j) const { + switch (get_column_type(j)) { + case column_type::fixed: + case column_type::boxed: + case column_type::lower_bound: + return true; + default: + return false; + } +} +bool lp_bound_propagator::upper_bound_is_available(unsigned j) const { + if (upper_bound_is_available_for_column(j)) + return true; + if (!m_nla_solver->is_monomial_var(j)) + return false; + return m_nla_solver->monomial_has_upper_bound(j); +} + +impq lp_bound_propagator::get_lower_bound(unsigned j) const { + lp_assert(lower_bound_is_available(j)); + if (lower_bound_is_available_for_column(j)) { + const impq& l = m_lar_solver.m_mpq_lar_core_solver.m_r_lower_bounds()[j]; + if (!(m_nla_solver && m_nla_solver->is_monomial_var(j))) + return l; + if (! m_nla_solver->monomial_has_lower_bound(j)) + return std::max(l, m_nla_solver->get_lower_bound(j)); + } + return m_nla_solver->get_lower_bound(j); +} + +impq lp_bound_propagator::get_upper_bound(unsigned j) const { + lp_assert(upper_bound_is_available(j)); + if (upper_bound_is_available_for_column(j)) { + const impq& l = m_lar_solver.m_mpq_lar_core_solver.m_r_upper_bounds()[j]; + if (!(m_nla_solver && m_nla_solver->is_monomial_var(j))) + return l; + if (! m_nla_solver->monomial_has_upper_bound(j)) + return std::min(l, m_nla_solver->get_upper_bound(j)); + } + return m_nla_solver->get_upper_bound(j); +} + +void lp_bound_propagator::try_add_bound(mpq v, unsigned j, bool is_low, bool coeff_before_j_is_pos, unsigned row_or_term_index, bool strict) { + CTRACE("nla_solver", m_nla_solver && m_nla_solver->is_monomial_var(j), tout << "mon_var = " << j << "\n"; ); j = m_lar_solver.adjust_column_index_to_term_index(j); lconstraint_kind kind = is_low? GE : LE; diff --git a/src/math/lp/bound_propagator.h b/src/math/lp/lp_bound_propagator.h similarity index 66% rename from src/math/lp/bound_propagator.h rename to src/math/lp/lp_bound_propagator.h index a75cabf6d..a3f4c9802 100644 --- a/src/math/lp/bound_propagator.h +++ b/src/math/lp/lp_bound_propagator.h @@ -4,19 +4,29 @@ */ #pragma once #include "math/lp/lp_settings.h" +namespace nla { +// forward definition +class solver; +} namespace lp { class lar_solver; -class bound_propagator { +class lp_bound_propagator { std::unordered_map m_improved_lower_bounds; // these maps map a column index to the corresponding index in ibounds std::unordered_map m_improved_upper_bounds; lar_solver & m_lar_solver; + nla::solver * m_nla_solver; public: vector m_ibounds; public: - bound_propagator(lar_solver & ls); + lp_bound_propagator(lar_solver & ls, nla::solver * nla); column_type get_column_type(unsigned) const; - const impq & get_lower_bound(unsigned) const; - const impq & get_upper_bound(unsigned) const; + + bool lower_bound_is_available_for_column(unsigned) const; + bool upper_bound_is_available_for_column(unsigned) const; + bool lower_bound_is_available(unsigned) const; + bool upper_bound_is_available(unsigned) const; + impq get_lower_bound(unsigned) const; + impq get_upper_bound(unsigned) const; void try_add_bound(mpq v, unsigned j, bool is_low, bool coeff_before_j_is_pos, unsigned row_or_term_index, bool strict); virtual bool bound_is_interesting(unsigned vi, lp::lconstraint_kind kind, diff --git a/src/math/lp/nla_core.cpp b/src/math/lp/nla_core.cpp index 16afe48ad..3a3fb8b6e 100644 --- a/src/math/lp/nla_core.cpp +++ b/src/math/lp/nla_core.cpp @@ -1331,6 +1331,21 @@ lbool core::test_check( return check(l); } +lp::impq core::get_upper_bound_of_monomial(lpvar j) const { + SASSERT(false); +} + +lp::impq core::get_lower_bound_of_monomial(lpvar j) const { + SASSERT(false); +} + +bool core::monomial_has_lower_bound(lpvar j) const { + SASSERT(false); +} + +bool core::monomial_has_upper_bound(lpvar j) const { + SASSERT(false); +} } // end of nla diff --git a/src/math/lp/nla_core.h b/src/math/lp/nla_core.h index 466df1213..06f4871b7 100644 --- a/src/math/lp/nla_core.h +++ b/src/math/lp/nla_core.h @@ -96,7 +96,7 @@ public: lp::lar_term subs_terms_to_columns(const lp::lar_term& t) const; bool ineq_holds(const ineq& n) const; bool lemma_holds(const lemma& l) const; - + bool is_monomial_var(lpvar j) const { return m_emons.is_monomial_var(j); } rational val(lpvar j) const { return m_lar_solver.get_column_value_rational(j); } rational val(const monomial& m) const { return m_lar_solver.get_column_value_rational(m.var()); } @@ -338,6 +338,10 @@ public: lbool test_check(vector& l); lpvar map_to_root(lpvar) const; + lp::impq get_upper_bound_of_monomial(lpvar j) const; + lp::impq get_lower_bound_of_monomial(lpvar j) const; + bool monomial_has_lower_bound(lpvar j) const; + bool monomial_has_upper_bound(lpvar j) const; }; // end of core struct pp_mon { diff --git a/src/math/lp/nla_solver.cpp b/src/math/lp/nla_solver.cpp index 7df1fee72..8665d1612 100644 --- a/src/math/lp/nla_solver.cpp +++ b/src/math/lp/nla_solver.cpp @@ -31,6 +31,10 @@ void solver::add_monomial(lpvar v, unsigned sz, lpvar const* vs) { m_core->add(v, sz, vs); } +bool solver::is_monomial_var(lpvar v) const { + return m_core->is_monomial_var(v); +} + bool solver::need_check() { return true; } lbool solver::check(vector& l) { @@ -57,5 +61,21 @@ solver::solver(lp::lar_solver& s) { solver::~solver() { dealloc(m_core); } +lp::impq solver::get_upper_bound(lpvar j) const { + SASSERT(is_monomial_var(j) && m_core->monomial_has_upper_bound(j)); + return m_core->get_upper_bound_of_monomial(j); +} + +lp::impq solver::get_lower_bound(lpvar j) const { + SASSERT(is_monomial_var(j) && m_core->monomial_has_lower_bound(j)); + return m_core->get_lower_bound_of_monomial(j); +} + +bool solver::monomial_has_lower_bound(lpvar j) const { + return m_core->monomial_has_lower_bound(j); +} +bool solver::monomial_has_upper_bound(lpvar j) const { + return m_core->monomial_has_upper_bound(j); +} } diff --git a/src/math/lp/nla_solver.h b/src/math/lp/nla_solver.h index 44f57908c..292298742 100644 --- a/src/math/lp/nla_solver.h +++ b/src/math/lp/nla_solver.h @@ -32,7 +32,7 @@ namespace nla { class solver { core* m_core; public: - void add_monomial(lp::var_index v, unsigned sz, lp::var_index const* vs); + void add_monomial(lpvar v, unsigned sz, lpvar const* vs); solver(lp::lar_solver& s); ~solver(); @@ -42,5 +42,10 @@ public: bool need_check(); lbool check(vector&); std::ostream& display(std::ostream& out); + bool is_monomial_var(lpvar) const; + lp::impq get_lower_bound(lpvar j) const; + lp::impq get_upper_bound(lpvar j) const; + bool monomial_has_lower_bound(lpvar j) const; + bool monomial_has_upper_bound(lpvar j) const; }; } diff --git a/src/smt/theory_lra.cpp b/src/smt/theory_lra.cpp index e31f6d02e..f00d53085 100644 --- a/src/smt/theory_lra.cpp +++ b/src/smt/theory_lra.cpp @@ -2321,9 +2321,13 @@ public: return false; } - struct local_bound_propagator: public lp::bound_propagator { + nla::solver* get_nla_solver() { + return m_switcher.m_nla ? m_switcher.m_nla->get() : nullptr; + } + + struct local_bound_propagator: public lp::lp_bound_propagator { imp & m_imp; - local_bound_propagator(imp& i) : bound_propagator(*i.m_solver), m_imp(i) {} + local_bound_propagator(imp& i) : lp_bound_propagator(*i.m_solver, i.get_nla_solver()), m_imp(i) {} bool bound_is_interesting(unsigned j, lp::lconstraint_kind kind, const rational & v) override { return m_imp.bound_is_interesting(j, kind, v);