diff --git a/src/math/dd/pdd_interval.h b/src/math/dd/pdd_interval.h index 70e6c9a9b..55d2c9b21 100644 --- a/src/math/dd/pdd_interval.h +++ b/src/math/dd/pdd_interval.h @@ -49,12 +49,8 @@ public: get_interval(p.hi(), hi); get_interval(p.lo(), lo); if (deps) { - interval_deps_combine_rule combine_rule; - m_dep_intervals.mul(hi, a, t, combine_rule); - m_dep_intervals.combine_deps(hi, a, combine_rule, t); - combine_rule.reset(); - m_dep_intervals.add(t, lo, ret, combine_rule); - m_dep_intervals.combine_deps(t, lo, combine_rule, ret); + m_dep_intervals.mul(hi, a, t); + m_dep_intervals.add(t, lo, ret); } else { m_dep_intervals.mul(hi, a, t); m_dep_intervals.add(t, lo, ret); diff --git a/src/math/interval/dep_intervals.h b/src/math/interval/dep_intervals.h index f2fca6131..b0b752e36 100644 --- a/src/math/interval/dep_intervals.h +++ b/src/math/interval/dep_intervals.h @@ -134,13 +134,15 @@ private: } }; + public: + typedef interval_manager::interval interval; + mutable unsynch_mpq_manager m_num_manager; mutable u_dependency_manager m_dep_manager; im_config m_config; mutable interval_manager m_imanager; - typedef interval_manager::interval interval; unsynch_mpq_manager& num_manager() { return m_num_manager; } const unsynch_mpq_manager& num_manager() const { return m_num_manager; } @@ -148,34 +150,6 @@ public: u_dependency* mk_leaf(unsigned d) { return m_dep_manager.mk_leaf(d); } u_dependency* mk_join(u_dependency* a, u_dependency* b) { return m_dep_manager.mk_join(a, b); } - template - void mul(const rational& r, const interval& a, interval& b) const { - if (r.is_zero()) return; - m_imanager.mul(r.to_mpq(), a, b); - if (wd == with_deps) { - if (r.is_pos()) { - b.m_lower_dep = a.m_lower_dep; - b.m_upper_dep = a.m_upper_dep; - } - else { - b.m_upper_dep = a.m_lower_dep; - b.m_lower_dep = a.m_upper_dep; - } - } - } - - void mul(const interval& a, const interval& b, interval& c, interval_deps_combine_rule& deps) { m_imanager.mul(a, b, c, deps); } - void add(const interval& a, const interval& b, interval& c, interval_deps_combine_rule& deps) { m_imanager.add(a, b, c, deps); } - - void combine_deps(interval const& a, interval const& b, interval_deps_combine_rule const& deps, interval& i) const { - SASSERT(&a != &i && &b != &i); - m_config.add_deps(a, b, deps, i); - } - - void combine_deps(interval const& a, interval_deps_combine_rule const& deps, interval& i) const { - SASSERT(&a != &i); - m_config.add_deps(a, deps, i); - } public: u_dependency_manager& dep_manager() { return m_dep_manager; } @@ -188,12 +162,20 @@ public: std::ostream& display(std::ostream& out, const interval& i) const; void set_lower(interval& a, rational const& n) const { m_config.set_lower(a, n.to_mpq()); } void set_upper(interval& a, rational const& n) const { m_config.set_upper(a, n.to_mpq()); } - void set_lower_is_open(interval& a, bool strict) { m_config.set_lower_is_open(a, strict); } - void set_lower_is_inf(interval& a, bool inf) { m_config.set_lower_is_inf(a, inf); } - void set_upper_is_open(interval& a, bool strict) { m_config.set_upper_is_open(a, strict); } - void set_upper_is_inf(interval& a, bool inf) { m_config.set_upper_is_inf(a, inf); } + void set_lower_is_open(interval& a, bool strict) const { m_config.set_lower_is_open(a, strict); } + void set_lower_is_inf(interval& a, bool inf) const { m_config.set_lower_is_inf(a, inf); } + void set_upper_is_open(interval& a, bool strict) const { m_config.set_upper_is_open(a, strict); } + void set_upper_is_inf(interval& a, bool inf) const { m_config.set_upper_is_inf(a, inf); } void set_lower_dep(interval& a, u_dependency* d) const { m_config.set_lower_dep(a, d); } void set_upper_dep(interval& a, u_dependency* d) const { m_config.set_upper_dep(a, d); } + void set_value(interval& a, rational const& n) const { + set_lower(a, n); + set_upper(a, n); + set_lower_is_open(a, false); + set_upper_is_open(a, false); + set_lower_is_inf(a, false); + set_upper_is_inf(a, false); + } bool is_zero(const interval& a) const { return m_config.is_zero(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); } @@ -220,6 +202,8 @@ public: template void power(const interval& a, unsigned n, interval& b) { + if (n == 1 && &a == &b) + return; if (with_deps == wd) { interval_deps_combine_rule combine_rule; m_imanager.power(a, n, b, combine_rule); @@ -232,6 +216,58 @@ public: display(tout, b) << "\n"; ); } + template + void mul(const rational& r, const interval& a, interval& b) const { + if (r.is_zero()) return; + m_imanager.mul(r.to_mpq(), a, b); + if (wd == with_deps) { + if (r.is_pos()) { + b.m_lower_dep = a.m_lower_dep; + b.m_upper_dep = a.m_upper_dep; + } + else { + b.m_upper_dep = a.m_lower_dep; + b.m_lower_dep = a.m_upper_dep; + } + } + } + + template + void mul(const interval& a, const interval& b, interval& c) { + if (wd == with_deps) { + interval_deps_combine_rule comb_rule; + mul(a, b, c, comb_rule); + combine_deps(a, b, comb_rule, c); + } + else { + mul(a, b, c); + } + } + + template + void div(const interval& a, const interval& b, interval& c) { + if (wd == with_deps) { + interval_deps_combine_rule comb_rule; + m_imanager.div(a, b, c, comb_rule); + combine_deps(a, b, comb_rule, c); + } + else { + m_imanager.div(a, b, c); + } + } + + template + void add(const interval& a, const interval& b, interval& c) { + if (wd == with_deps) { + interval_deps_combine_rule comb_rule; + add(a, b, c, comb_rule); + combine_deps(a, b, comb_rule, c); + } + else { + add(a, b, c); + } + } + template void copy_upper_bound(const interval& a, interval& i) const { SASSERT(a.m_upper_inf == false); @@ -370,6 +406,21 @@ public: copy_upper_bound(b, i); } +private: + void mul(const interval& a, const interval& b, interval& c, interval_deps_combine_rule& deps) { m_imanager.mul(a, b, c, deps); } + void add(const interval& a, const interval& b, interval& c, interval_deps_combine_rule& deps) { m_imanager.add(a, b, c, deps); } + + void combine_deps(interval const& a, interval const& b, interval_deps_combine_rule const& deps, interval& i) const { + SASSERT(&a != &i && &b != &i); + m_config.add_deps(a, b, deps, i); + } + + void combine_deps(interval const& a, interval_deps_combine_rule const& deps, interval& i) const { + SASSERT(&a != &i); + m_config.add_deps(a, deps, i); + } + }; typedef _scoped_interval scoped_dep_interval; +typedef dep_intervals::interval dep_interval; diff --git a/src/math/lp/monomial_bounds.cpp b/src/math/lp/monomial_bounds.cpp index 71c0595a0..d5eb03c0f 100644 --- a/src/math/lp/monomial_bounds.cpp +++ b/src/math/lp/monomial_bounds.cpp @@ -20,40 +20,50 @@ namespace nla { bool propagated = false; for (lpvar v : c().m_to_refine) { monic const& m = c().emons()[v]; - if (propagate_up(m)) + if (propagate(m)) propagated = true; } return propagated; } - bool monomial_bounds::propagate_up(monic const& m) { + void monomial_bounds::compute_product(unsigned start, monic const& m, scoped_dep_interval& product) { auto & intervals = c().m_intervals; auto & dep_intervals = intervals.get_dep_intervals(); - scoped_dep_interval product(dep_intervals), di(dep_intervals); - lpvar v = m.vars()[0]; - var2interval(v, product); - for (unsigned i = 1; i < m.size(); ++i) { - var2interval(m.vars()[i], di); - dep_intervals.mul(product, di, product); + scoped_dep_interval vi(dep_intervals); + for (unsigned i = start; i < m.size(); ) { + lpvar v = m.vars()[i]; + unsigned power = 1; + var2interval(v, vi); + ++i; + for (; i < m.size() && m.vars()[i] == v; ++i) { + ++power; + } + dep_intervals.power(vi, power, vi); + dep_intervals.mul(product, vi, product); } - if (dep_intervals.is_below(product, c().val(m.var()))) { + } + + bool monomial_bounds::propagate_value(dep_interval& range, lpvar v) { + auto & intervals = c().m_intervals; + auto & dep_intervals = intervals.get_dep_intervals(); + if (dep_intervals.is_below(range, c().val(v))) { lp::explanation ex; - dep_intervals.get_upper_dep(product, ex); - new_lemma lemma(c(), "propagate up - upper bound of product is below value"); + dep_intervals.get_upper_dep(range, ex); + new_lemma lemma(c(), "propagate up - upper bound of range is below value"); lemma &= ex; - auto const& upper = dep_intervals.upper(product); - auto cmp = dep_intervals.upper_is_open(product) ? llc::LT : llc::LE; - lemma |= ineq(m.var(), cmp, upper); + auto const& upper = dep_intervals.upper(range); + auto cmp = dep_intervals.upper_is_open(range) ? llc::LT : llc::LE; + lemma |= ineq(v, cmp, upper); return true; } - else if (dep_intervals.is_above(product, c().val(m.var()))) { + else if (dep_intervals.is_above(range, c().val(v))) { lp::explanation ex; - dep_intervals.get_lower_dep(product, ex); - new_lemma lemma(c(), "propagate up - lower bound of product is above value"); + dep_intervals.get_lower_dep(range, ex); + new_lemma lemma(c(), "propagate up - lower bound of range is above value"); lemma &= ex; - auto const& lower = dep_intervals.upper(product); - auto cmp = dep_intervals.lower_is_open(product) ? llc::GT : llc::GE; - lemma |= ineq(m.var(), cmp, lower); + auto const& lower = dep_intervals.upper(range); + auto cmp = dep_intervals.lower_is_open(range) ? llc::GT : llc::GE; + lemma |= ineq(v, cmp, lower); return true; } else { @@ -72,11 +82,58 @@ namespace nla { dep_intervals.set_lower(i, bound); dep_intervals.set_lower_dep(i, dep_intervals.mk_leaf(ci)); } + else { + dep_intervals.set_lower_is_inf(i, true); + } if (c().has_upper_bound(v, ci, bound, is_strict)) { dep_intervals.set_upper_is_open(i, is_strict); dep_intervals.set_upper(i, bound); dep_intervals.set_upper_dep(i, dep_intervals.mk_leaf(ci)); } + else { + dep_intervals.set_upper_is_inf(i, true); + } } + + bool monomial_bounds::propagate(monic const& m) { + auto & intervals = c().m_intervals; + auto & dep_intervals = intervals.get_dep_intervals(); + scoped_dep_interval product(dep_intervals); + scoped_dep_interval vi(dep_intervals), mi(dep_intervals); + scoped_dep_interval other_product(dep_intervals); + var2interval(m.var(), mi); + if (dep_intervals.lower_is_inf(mi) && dep_intervals.upper_is_inf(mi)) + return false; + dep_intervals.set_value(product, rational::one()); + for (unsigned i = 0; i < m.size(); ) { + lpvar v = m.vars()[i]; + ++i; + unsigned power = 1; + for (; i < m.size() && v == m.vars()[i]; ++i) + ++power; + var2interval(v, vi); + dep_intervals.power(vi, power, vi); + if (power == 1) { + dep_intervals.set(other_product, product); + compute_product(i, m, other_product); + if (propagate_down(m, mi, v, other_product)) + return true; + } + dep_intervals.mul(product, vi, product); + } + return propagate_value(product, m.var()); + } + + bool monomial_bounds::propagate_down(monic const& m, dep_interval& mi, lpvar v, dep_interval& product) { + auto & intervals = c().m_intervals; + auto & dep_intervals = intervals.get_dep_intervals(); + if (!dep_intervals.separated_from_zero(product)) + return false; + scoped_dep_interval range(dep_intervals); + dep_intervals.set(range, mi); + dep_intervals.div(range, product, range); + return propagate_value(range, v); + } + } diff --git a/src/math/lp/monomial_bounds.h b/src/math/lp/monomial_bounds.h index 146514f92..714bcbd51 100644 --- a/src/math/lp/monomial_bounds.h +++ b/src/math/lp/monomial_bounds.h @@ -19,7 +19,11 @@ namespace nla { class monomial_bounds : common { void var2interval(lpvar v, scoped_dep_interval& i); - bool propagate_up(monic const& m); + bool propagate_down(monic const& m, lpvar u); + bool propagate_value(dep_interval& range, lpvar v); + void compute_product(unsigned start, monic const& m, scoped_dep_interval& i); + bool propagate(monic const& m); + bool propagate_down(monic const& m, dep_interval& mi, lpvar v, dep_interval& product); public: monomial_bounds(core* core); bool operator()(); diff --git a/src/math/lp/nla_intervals.cpp b/src/math/lp/nla_intervals.cpp index a25c8a30b..ad8c62be0 100644 --- a/src/math/lp/nla_intervals.cpp +++ b/src/math/lp/nla_intervals.cpp @@ -327,14 +327,7 @@ bool intervals::interval_of_sum_no_term(const nex_sum& e, scoped_dep_interval & scoped_dep_interval c(get_dep_intervals()); TRACE("nla_intervals_details", tout << "sdi = "; display(tout, sdi) << "\nb = "; display(tout, b) << "\n";); - if (wd == e_with_deps::with_deps) { - interval_deps_combine_rule combine_rule; - m_dep_intervals.add(sdi, b, c, combine_rule); - m_dep_intervals.combine_deps(sdi, b, combine_rule, c); - } - else { - m_dep_intervals.add(sdi, b, c); - } + m_dep_intervals.add(sdi, b, c); m_dep_intervals.set(sdi, c); TRACE("nla_intervals_details", tout << *e[k] << ", "; display(tout, sdi); tout << "\n";); @@ -426,14 +419,7 @@ bool intervals::interval_of_mul(const nex_mul& e, scoped_dep_interval& a, const return false; TRACE("nla_intervals_details", tout << "ep = " << ep << ", "; display(tout, b); ); scoped_dep_interval c(get_dep_intervals()); - if (wd == e_with_deps::with_deps) { - interval_deps_combine_rule comb_rule; - m_dep_intervals.mul(a, b, c, comb_rule); - TRACE("nla_intervals_details", tout << "c before combine_deps() "; display(tout, c);); - m_dep_intervals.combine_deps(a, b, comb_rule, c); - } else { - m_dep_intervals.mul(a, b, c); - } + m_dep_intervals.mul(a, b, c); TRACE("nla_intervals_details", tout << "a "; display(tout, a);); TRACE("nla_intervals_details", tout << "c "; display(tout, c);); m_dep_intervals.set(a, c);