From 453a5d2de1242b2ad7620c6538c540807e0e7081 Mon Sep 17 00:00:00 2001 From: Lev Nachmanson Date: Wed, 5 Jun 2019 15:57:33 -0700 Subject: [PATCH] hook up nla_solver it lp bound propagation Signed-off-by: Lev Nachmanson --- src/math/lp/nla_intervals.cpp | 107 +++++++++++++++++++++++ src/math/lp/nla_intervals.h | 158 ++++++++++++++++++++++++++++++++++ 2 files changed, 265 insertions(+) create mode 100644 src/math/lp/nla_intervals.cpp create mode 100644 src/math/lp/nla_intervals.h diff --git a/src/math/lp/nla_intervals.cpp b/src/math/lp/nla_intervals.cpp new file mode 100644 index 000000000..14d683fb5 --- /dev/null +++ b/src/math/lp/nla_intervals.cpp @@ -0,0 +1,107 @@ + + +#include "math/lp/nla_core.h" +#include "math/interval/interval_def.h" +#include "math/lp/nla_intervals.h" +namespace nla { + + bool intervals::check() { + m_region.reset(); + for (auto const& m : c().emons()) { + if (!check(m)) { + return false; + } + } + for (auto const& t : m_solver.terms()) { + if (!check(*t)) { + return false; + } + } + return true; + } + + bool intervals::check(monomial const& m) { + interval a, b, c, d; + m_imanager.set(a, rational(1).to_mpq()); + set_interval(m.var(), d); + if (m_imanager.lower_is_inf(d) && m_imanager.upper_is_inf(d)) { + return true; + } + for (lpvar v : m.vars()) { + // TBD allow for division to get range of a + // m = a*b*c, where m and b*c are bounded, then interval for a is m/b*c + if (m_imanager.lower_is_inf(a) && m_imanager.upper_is_inf(a)) { + return true; + } + // TBD: deal with powers v^n interval instead of multiplying v*v .. * v + set_interval(v, b); + interval_deps deps; + m_imanager.mul(a, b, c, deps); + m_imanager.set(a, c); + m_config.set_deps(a, b, deps, a); + } + if (m_imanager.before(a, d)) { + svector cs; + m_dep_manager.linearize(a.m_upper_dep, cs); + m_dep_manager.linearize(d.m_lower_dep, cs); + for (auto ci : cs) { + //expl.push_justification(ci); + } + // TBD conflict + return false; + } + if (m_imanager.before(d, a)) { + svector cs; + m_dep_manager.linearize(a.m_lower_dep, cs); + m_dep_manager.linearize(d.m_upper_dep, cs); + for (auto ci : cs) { + //expl.push_justification(ci); + } + // TBD conflict + return false; + } + // could also perform bounds propagation: + // a has tighter lower/upper bound than m.var(), + // -> transfer bound to m.var() + // all but one variable has bound + // -> transfer bound to that variable using division + return true; + } + + void intervals::set_interval(lpvar v, interval& b) { + lp::constraint_index ci; + rational val; + bool is_strict; + if (m_solver.has_lower_bound(v, ci, val, is_strict)) { + m_config.set_lower(b, val); + m_config.set_lower_is_open(b, is_strict); + m_config.set_lower_is_inf(b, false); + b.m_lower_dep = mk_dep(ci); + } + else { + m_config.set_lower_is_open(b, true); + m_config.set_lower_is_inf(b, true); + b.m_lower_dep = nullptr; + } + if (m_solver.has_upper_bound(v, ci, val, is_strict)) { + m_config.set_upper(b, val); + m_config.set_upper_is_open(b, is_strict); + m_config.set_upper_is_inf(b, false); + b.m_upper_dep = mk_dep(ci); + } + else { + m_config.set_upper_is_open(b, true); + m_config.set_upper_is_inf(b, true); + b.m_upper_dep = nullptr; + } + } + + intervals::ci_dependency *intervals::mk_dep(lp::constraint_index ci) { + return m_dep_manager.mk_leaf(ci); + } + + bool intervals::check(lp::lar_term const& t) { + // convert term into factors for improved precision + return true; + } +} diff --git a/src/math/lp/nla_intervals.h b/src/math/lp/nla_intervals.h new file mode 100644 index 000000000..7213adba8 --- /dev/null +++ b/src/math/lp/nla_intervals.h @@ -0,0 +1,158 @@ +/*++ + Copyright (c) 2017 Microsoft Corporation + + Module Name: + + + + Abstract: + + + + Author: + Nikolaj Bjorner (nbjorner) + Lev Nachmanson (levnach) + + Revision History: + + + --*/ +#pragma once +#include "util/dependency.h" +#include "math/lp/nla_common.h" +#include "math/lp/lar_solver.h" +#include "math/interval/interval.h" + + +namespace nla { + class core; + + class intervals : common { + + class ci_value_manager { + public: + void inc_ref(lp::constraint_index const & v) { + } + + void dec_ref(lp::constraint_index const & v) { + } + }; + + struct ci_dependency_config { + typedef ci_value_manager value_manager; + typedef small_object_allocator allocator; + static const bool ref_count = false; + typedef lp::constraint_index value; + }; + + typedef dependency_manager ci_dependency_manager; + + typedef ci_dependency_manager::dependency ci_dependency; + + class im_config { + unsynch_mpq_manager& m_manager; + ci_dependency_manager& m_dep_manager; + + public: + typedef unsynch_mpq_manager numeral_manager; + typedef mpq numeral; + + // Every configuration object must provide an interval type. + // The actual fields are irrelevant, the interval manager + // accesses interval data using the following API. + + struct interval { + interval(): + m_lower(), m_upper(), + m_lower_open(1), m_upper_open(1), + m_lower_inf(1), m_upper_inf(1), + m_lower_dep(nullptr), m_upper_dep(nullptr) {} + numeral m_lower; + numeral m_upper; + unsigned m_lower_open:1; + unsigned m_upper_open:1; + unsigned m_lower_inf:1; + unsigned m_upper_inf:1; + ci_dependency * m_lower_dep; // justification for the lower bound + ci_dependency * m_upper_dep; // justification for the upper bound + }; + + + void set_deps(interval const& a, interval const& b, interval_deps const& deps, interval& i) { + ci_dependency* lo = mk_dependency(a, b, deps.m_lower_deps); + ci_dependency* hi = mk_dependency(a, b, deps.m_upper_deps); + i.m_lower_dep = lo; + i.m_upper_dep = hi; + } + + // Should be NOOPs for precise numeral types. + // For imprecise types (e.g., floats) it should set the rounding mode. + void round_to_minus_inf() {} + void round_to_plus_inf() {} + void set_rounding(bool to_plus_inf) {} + + // Getters + numeral const & lower(interval const & a) const { return a.m_lower; } + numeral const & upper(interval const & a) const { return a.m_upper; } + numeral & lower(interval & a) { return a.m_lower; } + numeral & upper(interval & a) { return a.m_upper; } + bool lower_is_open(interval const & a) const { return a.m_lower_open; } + bool upper_is_open(interval const & a) const { return a.m_upper_open; } + bool lower_is_inf(interval const & a) const { return a.m_lower_inf; } + bool upper_is_inf(interval const & a) const { return a.m_upper_inf; } + + // Setters + void set_lower(interval & a, numeral const & n) { m_manager.set(a.m_lower, n); } + void set_upper(interval & a, numeral const & n) { m_manager.set(a.m_upper, n); } + void set_lower(interval & a, rational const & n) { set_lower(a, n.to_mpq()); } + void set_upper(interval & a, rational const & n) { set_upper(a, n.to_mpq()); } + void set_lower_is_open(interval & a, bool v) { a.m_lower_open = v; } + void set_upper_is_open(interval & a, bool v) { a.m_upper_open = v; } + void set_lower_is_inf(interval & a, bool v) { a.m_lower_inf = v; } + void set_upper_is_inf(interval & a, bool v) { a.m_upper_inf = v; } + + // Reference to numeral manager + numeral_manager & m() const { return m_manager; } + + im_config(numeral_manager & m, ci_dependency_manager& d):m_manager(m), m_dep_manager(d) {} + private: + ci_dependency* mk_dependency(interval const& a, interval const& b, bound_deps bd) { + ci_dependency* dep = nullptr; + if (dep_in_lower1(bd)) { + dep = m_dep_manager.mk_join(dep, a.m_lower_dep); + } + if (dep_in_lower2(bd)) { + dep = m_dep_manager.mk_join(dep, b.m_lower_dep); + } + if (dep_in_upper1(bd)) { + dep = m_dep_manager.mk_join(dep, a.m_upper_dep); + } + if (dep_in_upper2(bd)) { + dep = m_dep_manager.mk_join(dep, b.m_upper_dep); + } + return dep; + } + }; + + ci_dependency_manager m_dep_manager; + im_config m_config; + interval_manager m_imanager; + region m_region; + lp::lar_solver& m_solver; + + typedef interval_manager::interval interval; + + bool check(monomial const& m); + + void set_interval(lpvar v, interval & b); + + ci_dependency* mk_dep(lp::constraint_index ci); + + bool check(lp::lar_term const& t); + + public: + intervals(core* c); + bool check(); + }; + +} // end of namespace nla