diff --git a/src/util/lp/CMakeLists.txt b/src/util/lp/CMakeLists.txt index c3dd8df1f..61617089e 100644 --- a/src/util/lp/CMakeLists.txt +++ b/src/util/lp/CMakeLists.txt @@ -25,7 +25,9 @@ z3_add_component(lp matrix.cpp mon_eq.cpp nla_basics_lemmas.cpp + nla_common.cpp nla_core.cpp + nla_order_lemmas.cpp nla_solver.cpp nra_solver.cpp permutation_matrix.cpp diff --git a/src/util/lp/nla_basics_lemmas.cpp b/src/util/lp/nla_basics_lemmas.cpp index 1e18b71f9..5bf0c7c2e 100644 --- a/src/util/lp/nla_basics_lemmas.cpp +++ b/src/util/lp/nla_basics_lemmas.cpp @@ -21,11 +21,8 @@ #include "util/lp/nla_core.h" #include "util/lp/factorization_factory_imp.h" namespace nla { -template rational basics::vvr(T const& t) const { return m_core->vvr(t); } -rational basics::vvr(lpvar t) const { return m_core->vvr(t); } -template lpvar basics::var(T const& t) const { return m_core->var(t); } -basics::basics(core * c) : m_core(c) {} +basics::basics(core * c) : common(c) {} // Monomials m and n vars have the same values, up to "sign" // Generate a lemma if values of m.var() and n.var() are not the same up to sign bool basics::basic_sign_lemma_on_two_monomials(const monomial& m, const monomial& n, const rational& sign) { @@ -225,7 +222,6 @@ void basics::add_fixed_zero_lemma(const monomial& m, lpvar j) { c().mk_ineq(m.var(), llc::EQ); TRACE("nla_solver", c().print_lemma(tout);); } -void basics::add_empty_lemma() { c().add_empty_lemma(); } void basics::negate_strict_sign(lpvar j) { TRACE("nla_solver_details", c().print_var(j, tout);); if (!vvr(j).is_zero()) { @@ -845,11 +841,5 @@ void basics::basic_lemma_for_mon_non_zero_model_based(const rooted_mon& rm, cons basic_lemma_for_mon_non_zero_model_based_mf(f); } -bool basics::done() const { return c().done(); } - -template void basics::explain(const T& t) { - c().explain(t, c().current_expl()); -} -template void basics::explain(const monomial& t); } diff --git a/src/util/lp/nla_basics_lemmas.h b/src/util/lp/nla_basics_lemmas.h index 658beb753..d370411d4 100644 --- a/src/util/lp/nla_basics_lemmas.h +++ b/src/util/lp/nla_basics_lemmas.h @@ -21,11 +21,12 @@ #include "util/lp/monomial.h" #include "util/lp/rooted_mons.h" #include "util/lp/factorization.h" +#include "util/lp/nla_common.h" namespace nla { struct core; -struct basics { +struct basics: common { core* m_core; core& c() { return *m_core; } const core& c() const { return *m_core; } @@ -85,9 +86,6 @@ struct basics { void basic_lemma_for_mon(const rooted_mon& rm, bool derived); // use basic multiplication properties to create a lemma bool basic_lemma(bool derived); - template rational vvr(T const& t) const; - rational vvr(lpvar) const; - template lpvar var(T const& t) const; void generate_sign_lemma(const monomial& m, const monomial& n, const rational& sign); void generate_zero_lemmas(const monomial& m); lpvar find_best_zero(const monomial& m, unsigned_vector & fixed_zeros) const; @@ -97,14 +95,11 @@ struct basics { void generate_strict_case_zero_lemma(const monomial& m, unsigned zero_j, int sign_of_zj); void add_fixed_zero_lemma(const monomial& m, lpvar j); - void add_empty_lemma(); void negate_strict_sign(lpvar j); - bool done() const; // x != 0 or y = 0 => |xy| >= |y| void proportion_lemma_model_based(const rooted_mon& rm, const factorization& factorization); // x != 0 or y = 0 => |xy| >= |y| bool proportion_lemma_derived(const rooted_mon& rm, const factorization& factorization); - template void explain(const T&); // if there are no zero factors then |m| >= |m[factor_index]| void generate_pl_on_mon(const monomial& m, unsigned factor_index); diff --git a/src/util/lp/nla_common.cpp b/src/util/lp/nla_common.cpp new file mode 100644 index 000000000..798500e58 --- /dev/null +++ b/src/util/lp/nla_common.cpp @@ -0,0 +1,127 @@ +/*++ + Copyright (c) 2017 Microsoft Corporation + + Module Name: + + + + Abstract: + + + + Author: + Nikolaj Bjorner (nbjorner) + Lev Nachmanson (levnach) + + Revision History: + + + --*/ +#include "util/lp/nla_common.h" +#include "util/lp/nla_core.h" + +namespace nla { +bool common::done() const { return m_core->done(); } + +template void common::explain(const T& t) { + c().explain(t, c().current_expl()); +} +template void common::explain(const monomial& t); +template void common::explain(const factor& t); +template void common::explain(const rooted_mon& t); +template void common::explain(const factorization& t); + +void common::explain(lpvar j) { c().explain(j, c().current_expl()); } + +template rational common::vvr(T const& t) const { return m_core->vvr(t); } +template rational common::vvr(monomial const& t) const; +template rational common::vvr(rooted_mon const& t) const; +template rational common::vvr(factor const& t) const; +rational common::vvr(lpvar t) const { return m_core->vvr(t); } +template lpvar common::var(T const& t) const { return m_core->var(t); } +template lpvar common::var(factor const& t) const; +template lpvar common::var(rooted_mon const& t) const; +void common::add_empty_lemma() { c().add_empty_lemma(); } +template rational common::canonize_sign(const T& t) const { + return c().canonize_sign(t); +} +template rational common::canonize_sign(const rooted_mon& t) const; +template rational common::canonize_sign(const factor& t) const; +rational common::canonize_sign(lpvar j) const { + return c().canonize_sign_of_var(j); +} +void common::mk_ineq(lp::lar_term& t, llc cmp, const rational& rs){ + c().mk_ineq(t, cmp, rs); +} +void common::mk_ineq(const rational& a, lpvar j, const rational& b, lpvar k, llc cmp, const rational& rs){ + c().mk_ineq(a, j, b, j, cmp, rs); +} + +void common::mk_ineq(lpvar j, const rational& b, lpvar k, llc cmp, const rational& rs){ + c().mk_ineq(j, b, k, cmp, rs); +} + +void common::mk_ineq(lpvar j, const rational& b, lpvar k, llc cmp){ + c().mk_ineq(j, b, k, cmp); +} + +void common::mk_ineq(const rational& a, lpvar j, const rational& b, lpvar k, llc cmp) { + c().mk_ineq(a, j, b, k, cmp); +} + +void common::mk_ineq(const rational& a ,lpvar j, lpvar k, llc cmp, const rational& rs) { + c().mk_ineq(a, j, k, cmp, rs); +} + +void common::mk_ineq(lpvar j, lpvar k, llc cmp, const rational& rs) { + c().mk_ineq(j, k, cmp, rs);} + +void common::mk_ineq(lpvar j, llc cmp, const rational& rs){ + c().mk_ineq(j, cmp, rs);} + +void common::mk_ineq(const rational& a, lpvar j, llc cmp, const rational& rs) { + c().mk_ineq(a, j, cmp, rs); +} +void common::mk_ineq(const rational& a, lpvar j, llc cmp){ + c().mk_ineq(a, j, cmp); +} + +void common::mk_ineq(lpvar j, llc cmp){ + c().mk_ineq(j, cmp); +} +std::ostream& common::print_lemma(std::ostream& out) const { + return c().print_lemma(out); +} +template +std::ostream& common::print_product(const T & m, std::ostream& out) const { + return c().print_product(m, out); +} +template +std::ostream& common::print_product(const monomial & m, std::ostream& out) const; + +template std::ostream& common::print_product(const rooted_mon & m, std::ostream& out) const; + + +std::ostream& common::print_monomial(const monomial & m, std::ostream& out) const { + return c().print_monomial(m, out); +} +std::ostream& common::print_rooted_monomial(const rooted_mon& rm, std::ostream& out) const { + return c().print_rooted_monomial(rm, out); +} +std::ostream& common::print_rooted_monomial_with_vars(const rooted_mon& rm, std::ostream& out) const { + return c().print_rooted_monomial_with_vars(rm, out); +} +std::ostream& common::print_factor(const factor & f, std::ostream& out) const { + return c().print_factor(f, out); +} + +std::ostream& common::print_var(lpvar j, std::ostream& out) const { + return c().print_var(j, out); +} + +bool common::check_monomial(const monomial& m) const { + return c().check_monomial(m); +} + +} + diff --git a/src/util/lp/nla_common.h b/src/util/lp/nla_common.h new file mode 100644 index 000000000..57facac55 --- /dev/null +++ b/src/util/lp/nla_common.h @@ -0,0 +1,90 @@ +/*++ + Copyright (c) 2017 Microsoft Corporation + + Module Name: + + + + Abstract: + + + + Author: + Nikolaj Bjorner (nbjorner) + Lev Nachmanson (levnach) + + Revision History: + + + --*/ +#pragma once +#include "util/rational.h" +#include "util/lp/nla_defs.h" +#include "util/lp/lar_term.h" +#include "util/lp/monomial.h" +#include "util/lp/factorization.h" +#include "util/lp/rooted_mons.h" +namespace nla { +inline llc negate(llc cmp) { + switch(cmp) { + case llc::LE: return llc::GT; + case llc::LT: return llc::GE; + case llc::GE: return llc::LT; + case llc::GT: return llc::LE; + case llc::EQ: return llc::NE; + case llc::NE: return llc::EQ; + default: SASSERT(false); + }; + return cmp; // not reachable +} + +struct core; +struct common { + core* m_core; + common(core* c): m_core(c) {} + core& c() { return *m_core; } + const core& c() const { return *m_core; } + template rational vvr(T const& t) const; + rational vvr(lpvar) const; + template lpvar var(T const& t) const; + bool done() const; + template void explain(const T&); + void explain(lpvar); + void add_empty_lemma(); + template rational canonize_sign(const T&) const; + rational canonize_sign(lpvar) const; + void mk_ineq(lp::lar_term& t, llc cmp, const rational& rs); + void mk_ineq(const rational& a, lpvar j, const rational& b, lpvar k, llc cmp, const rational& rs); + + void mk_ineq(lpvar j, const rational& b, lpvar k, llc cmp, const rational& rs); + + void mk_ineq(lpvar j, const rational& b, lpvar k, llc cmp); + + void mk_ineq(const rational& a, lpvar j, const rational& b, lpvar k, llc cmp); + + void mk_ineq(const rational& a ,lpvar j, lpvar k, llc cmp, const rational& rs); + + void mk_ineq(lpvar j, lpvar k, llc cmp, const rational& rs); + + void mk_ineq(lpvar j, llc cmp, const rational& rs); + + void mk_ineq(const rational& a, lpvar j, llc cmp, const rational& rs); + void mk_ineq(const rational& a, lpvar j, llc cmp); + + void mk_ineq(lpvar j, llc cmp); + + std::ostream& print_lemma(std::ostream&) const; + + template + std::ostream& print_product(const T & m, std::ostream& out) const; + + std::ostream& print_factor(const factor &, std::ostream& out) const; + std::ostream& print_var(lpvar, std::ostream& out) const; + + std::ostream& print_monomial(const monomial & m, std::ostream& out) const; + std::ostream& print_rooted_monomial(const rooted_mon &, std::ostream& out) const; + std::ostream& print_rooted_monomial_with_vars(const rooted_mon&, std::ostream& out) const; + bool check_monomial(const monomial&) const; + +}; +} diff --git a/src/util/lp/nla_core.cpp b/src/util/lp/nla_core.cpp index 2718139f0..1ab3e41ab 100644 --- a/src/util/lp/nla_core.cpp +++ b/src/util/lp/nla_core.cpp @@ -25,7 +25,8 @@ core::core(lp::lar_solver& s) : m_evars(), m_lar_solver(s), m_tangents(this), - m_basics(this) { + m_basics(this), + m_order(this) { } bool core::compare_holds(const rational& ls, llc cmp, const rational& rs) const { @@ -190,6 +191,8 @@ std::ostream& core::print_product(const T & m, std::ostream& out) const { } return out; } +template std::ostream& core::print_product(const monomial & m, std::ostream& out) const; +template std::ostream& core::print_product(const rooted_mon & m, std::ostream& out) const; std::ostream & core::print_factor(const factor& f, std::ostream& out) const { if (f.is_var()) { @@ -478,19 +481,6 @@ void core:: mk_ineq(const rational& a, lpvar j, llc cmp, const rational& rs) { mk_ineq(t, cmp, rs); } -llc negate(llc cmp) { - switch(cmp) { - case llc::LE: return llc::GT; - case llc::LT: return llc::GE; - case llc::GE: return llc::LT; - case llc::GT: return llc::LE; - case llc::EQ: return llc::NE; - case llc::NE: return llc::EQ; - default: SASSERT(false); - }; - return cmp; // not reachable -} - llc apply_minus(llc cmp) { switch(cmp) { case llc::LE: return llc::GE; @@ -632,84 +622,6 @@ bool core::zero_is_an_inner_point_of_bounds(lpvar j) const { return true; } - -bool core:: try_get_non_strict_sign_from_bounds(lpvar j, int& sign) const { - SASSERT(sign); - if (has_lower_bound(j) && get_lower_bound(j) >= rational(0)) - return true; - if (has_upper_bound(j) && get_upper_bound(j) <= rational(0)) { - sign = -sign; - return true; - } - sign = 0; - return false; -} - -void core:: get_non_strict_sign(lpvar j, int& sign) const { - const rational & v = vvr(j); - if (v.is_zero()) { - try_get_non_strict_sign_from_bounds(j, sign); - } else { - sign *= nla::rat_sign(v); - } -} - - -void core:: add_trival_zero_lemma(lpvar zero_j, const monomial& m) { - add_empty_lemma(); - mk_ineq(zero_j, llc::NE); - mk_ineq(m.var(), llc::EQ); - TRACE("nla_solver", print_lemma(tout);); -} - -void core:: generate_zero_lemmas(const monomial& m) { - SASSERT(!vvr(m).is_zero() && product_value(m).is_zero()); - int sign = nla::rat_sign(vvr(m)); - unsigned_vector fixed_zeros; - lpvar zero_j = find_best_zero(m, fixed_zeros); - SASSERT(is_set(zero_j)); - unsigned zero_power = 0; - for (unsigned j : m){ - if (j == zero_j) { - zero_power++; - continue; - } - get_non_strict_sign(j, sign); - if(sign == 0) - break; - } - if (sign && is_even(zero_power)) - sign = 0; - TRACE("nla_solver_details", tout << "zero_j = " << zero_j << ", sign = " << sign << "\n";); - if (sign == 0) { // have to generate a non-convex lemma - add_trival_zero_lemma(zero_j, m); - } else { - generate_strict_case_zero_lemma(m, zero_j, sign); - } - for (lpvar j : fixed_zeros) - add_fixed_zero_lemma(m, j); -} - -void core:: add_fixed_zero_lemma(const monomial& m, lpvar j) { - add_empty_lemma(); - explain_fixed_var(j); - mk_ineq(m.var(), llc::EQ); - TRACE("nla_solver", print_lemma(tout);); -} - -llc core::negate(llc cmp) { - switch(cmp) { - case llc::LE: return llc::GT; - case llc::LT: return llc::GE; - case llc::GE: return llc::LT; - case llc::GT: return llc::LE; - case llc::EQ: return llc::NE; - case llc::NE: return llc::EQ; - default: SASSERT(false); - }; - return cmp; // not reachable -} - int core::rat_sign(const monomial& m) const { int sign = 1; for (lpvar j : m) { @@ -1890,54 +1802,6 @@ void core::negate_factor_relation(const rational& a_sign, const factor& a, const mk_ineq(a_fs*a_sign, var(a), - b_fs*b_sign, var(b), cmp); } -void core::negate_var_factor_relation(const rational& a_sign, lpvar a, const rational& b_sign, const factor& b) { - rational b_fs = canonize_sign(b); - llc cmp = a_sign*vvr(a) < b_sign*vvr(b)? llc::GE : llc::LE; - mk_ineq(a_sign, a, - b_fs*b_sign, var(b), cmp); -} - -// |c_sign| = 1, and c*c_sign > 0 -// ac > bc => ac/|c| > bc/|c| => a*c_sign > b*c_sign -void core::generate_ol(const rooted_mon& ac, - const factor& a, - int c_sign, - const factor& c, - const rooted_mon& bc, - const factor& b, - llc ab_cmp) { - add_empty_lemma(); - rational rc_sign = rational(c_sign); - mk_ineq(rc_sign * canonize_sign(c), var(c), llc::LE); - mk_ineq(canonize_sign(ac), var(ac), -canonize_sign(bc), var(bc), ab_cmp); - mk_ineq(canonize_sign(a)*rc_sign, var(a), - canonize_sign(b)*rc_sign, var(b), negate(ab_cmp)); - explain(ac, current_expl()); - explain(a, current_expl()); - explain(bc, current_expl()); - explain(b, current_expl()); - explain(c, current_expl()); - TRACE("nla_solver", print_lemma(tout);); -} - -void core::generate_mon_ol(const monomial& ac, - lpvar a, - const rational& c_sign, - lpvar c, - const rooted_mon& bd, - const factor& b, - const rational& d_sign, - lpvar d, - llc ab_cmp) { - add_empty_lemma(); - mk_ineq(c_sign, c, llc::LE); - explain(c, current_expl()); // this explains c == +- d - negate_var_factor_relation(c_sign, a, d_sign, b); - mk_ineq(ac.var(), -canonize_sign(bd), var(bd), ab_cmp); - explain(bd, current_expl()); - explain(b, current_expl()); - explain(d, current_expl()); - TRACE("nla_solver", print_lemma(tout);); -} - std::unordered_set core::collect_vars(const lemma& l) const { std::unordered_set vars; for (const auto& i : current_lemma().ineqs()) { @@ -1966,8 +1830,9 @@ std::unordered_set core::collect_vars(const lemma& l) const { return vars; } -void core::print_lemma(std::ostream& out) { +std::ostream& core::print_lemma(std::ostream& out) const { print_specific_lemma(current_lemma(), out); + return out; } void core::print_specific_lemma(const lemma& l, std::ostream& out) const { @@ -2001,51 +1866,6 @@ void core::trace_print_ol(const rooted_mon& ac, print_factor_with_vars(c, out); } -bool core:: order_lemma_on_ac_and_bc_and_factors(const rooted_mon& ac, - const factor& a, - const factor& c, - const rooted_mon& bc, - const factor& b) { - auto cv = vvr(c); - int c_sign = nla::rat_sign(cv); - SASSERT(c_sign != 0); - auto av_c_s = vvr(a)*rational(c_sign); - auto bv_c_s = vvr(b)*rational(c_sign); - auto acv = vvr(ac); - auto bcv = vvr(bc); - TRACE("nla_solver", trace_print_ol(ac, a, c, bc, b, tout);); - // Suppose ac >= bc, then ac/|c| >= bc/|c|. - // Notice that ac/|c| = a*c_sign , and bc/|c| = b*c_sign, which are correspondingly av_c_s and bv_c_s - if (acv >= bcv && av_c_s < bv_c_s) { - generate_ol(ac, a, c_sign, c, bc, b, llc::LT); - return true; - } else if (acv <= bcv && av_c_s > bv_c_s) { - generate_ol(ac, a, c_sign, c, bc, b, llc::GT); - return true; - } - return false; -} - -// a >< b && c > 0 => ac >< bc -// a >< b && c < 0 => ac <> bc -// ac[k] plays the role of c -bool core:: order_lemma_on_ac_and_bc(const rooted_mon& rm_ac, - const factorization& ac_f, - unsigned k, - const rooted_mon& rm_bd) { - TRACE("nla_solver", tout << "rm_ac = "; - print_rooted_monomial(rm_ac, tout); - tout << "\nrm_bd = "; - print_rooted_monomial(rm_bd, tout); - tout << "\nac_f[k] = "; - print_factor_with_vars(ac_f[k], tout);); - factor b; - if (!divide(rm_bd, ac_f[k], b)){ - return false; - } - - return order_lemma_on_ac_and_bc_and_factors(rm_ac, ac_f[(k + 1) % 2], ac_f[k], rm_bd, b); -} void core::maybe_add_a_factor(lpvar i, const factor& c, std::unordered_set& found_vars, @@ -2070,242 +1890,6 @@ void core::maybe_add_a_factor(lpvar i, } } -bool core:: order_lemma_on_ac_explore(const rooted_mon& rm, const factorization& ac, unsigned k) { - const factor c = ac[k]; - TRACE("nla_solver", tout << "c = "; print_factor_with_vars(c, tout); ); - if (c.is_var()) { - TRACE("nla_solver", tout << "var(c) = " << var(c);); - for (unsigned rm_bc : m_rm_table.var_map()[c.index()]) { - TRACE("nla_solver", ); - if (order_lemma_on_ac_and_bc(rm ,ac, k, m_rm_table.rms()[rm_bc])) { - return true; - } - } - } else { - for (unsigned rm_bc : m_rm_table.proper_multiples()[c.index()]) { - if (order_lemma_on_ac_and_bc(rm , ac, k, m_rm_table.rms()[rm_bc])) { - return true; - } - } - } - return false; -} - -void core::order_lemma_on_factorization(const rooted_mon& rm, const factorization& ab) { - const monomial& m = m_monomials[rm.orig_index()]; - TRACE("nla_solver", tout << "orig_sign = " << rm.orig_sign() << "\n";); - rational sign = rm.orig_sign(); - for(factor f: ab) - sign *= canonize_sign(f); - const rational & fv = vvr(ab[0]) * vvr(ab[1]); - const rational mv = sign * vvr(m); - TRACE("nla_solver", - tout << "ab.size()=" << ab.size() << "\n"; - tout << "we should have sign*vvr(m):" << mv << "=(" << sign << ")*(" << vvr(m) <<") to be equal to " << " vvr(ab[0])*vvr(ab[1]):" << fv << "\n";); - if (mv == fv) - return; - bool gt = mv > fv; - TRACE("nla_solver_f", tout << "m="; print_monomial_with_vars(m, tout); tout << "\nfactorization="; print_factorization(ab, tout);); - for (unsigned j = 0, k = 1; j < 2; j++, k--) { - order_lemma_on_ab(m, sign, var(ab[k]), var(ab[j]), gt); - explain(ab, current_expl()); explain(m, current_expl()); - explain(rm, current_expl()); - TRACE("nla_solver", print_lemma(tout);); - order_lemma_on_ac_explore(rm, ab, j); - } -} - -/** - \brief Add lemma: - a > 0 & b <= value(b) => sign*ab <= value(b)*a if value(a) > 0 - a < 0 & b >= value(b) => sign*ab <= value(b)*a if value(a) < 0 -*/ -void core::order_lemma_on_ab_gt(const monomial& m, const rational& sign, lpvar a, lpvar b) { - SASSERT(sign * vvr(m) > vvr(a) * vvr(b)); - add_empty_lemma(); - if (vvr(a).is_pos()) { - TRACE("nla_solver", tout << "a is pos\n";); - //negate a > 0 - mk_ineq(a, llc::LE); - // negate b <= vvr(b) - mk_ineq(b, llc::GT, vvr(b)); - // ab <= vvr(b)a - mk_ineq(sign, m.var(), -vvr(b), a, llc::LE); - } else { - TRACE("nla_solver", tout << "a is neg\n";); - SASSERT(vvr(a).is_neg()); - //negate a < 0 - mk_ineq(a, llc::GE); - // negate b >= vvr(b) - mk_ineq(b, llc::LT, vvr(b)); - // ab <= vvr(b)a - mk_ineq(sign, m.var(), -vvr(b), a, llc::LE); - } -} -// we need to deduce ab >= vvr(b)*a -/** - \brief Add lemma: - a > 0 & b >= value(b) => sign*ab >= value(b)*a if value(a) > 0 - a < 0 & b <= value(b) => sign*ab >= value(b)*a if value(a) < 0 -*/ -void core::order_lemma_on_ab_lt(const monomial& m, const rational& sign, lpvar a, lpvar b) { - SASSERT(sign * vvr(m) < vvr(a) * vvr(b)); - add_empty_lemma(); - if (vvr(a).is_pos()) { - //negate a > 0 - mk_ineq(a, llc::LE); - // negate b >= vvr(b) - mk_ineq(b, llc::LT, vvr(b)); - // ab <= vvr(b)a - mk_ineq(sign, m.var(), -vvr(b), a, llc::GE); - } else { - SASSERT(vvr(a).is_neg()); - //negate a < 0 - mk_ineq(a, llc::GE); - // negate b <= vvr(b) - mk_ineq(b, llc::GT, vvr(b)); - // ab >= vvr(b)a - mk_ineq(sign, m.var(), -vvr(b), a, llc::GE); - } -} - - -void core::order_lemma_on_ab(const monomial& m, const rational& sign, lpvar a, lpvar b, bool gt) { - if (gt) - order_lemma_on_ab_gt(m, sign, a, b); - else - order_lemma_on_ab_lt(m, sign, a, b); -} - -// void core::order_lemma_on_ab(const monomial& m, const rational& sign, lpvar a, lpvar b, bool gt) { -// add_empty_lemma(); -// if (gt) { -// if (vvr(a).is_pos()) { -// //negate a > 0 -// mk_ineq(a, llc::LE); -// // negate b >= vvr(b) -// mk_ineq(b, llc::LT, vvr(b)); -// // ab <= vvr(b)a -// mk_ineq(sign, m.var(), -vvr(b), a, llc::LE); -// } else { -// SASSERT(vvr(a).is_neg()); -// //negate a < 0 -// mk_ineq(a, llc::GE); -// // negate b <= vvr(b) -// mk_ineq(b, llc::GT, vvr(b)); -// // ab < vvr(b)a -// mk_ineq(sign, m.var(), -vvr(b), a, llc::LE); } -// } -// } - -void core::order_lemma_on_factor_binomial_explore(const monomial& m, unsigned k) { - SASSERT(m.size() == 2); - lpvar c = m[k]; - lpvar d = m_evars.find(c).var(); - auto it = m_rm_table.var_map().find(d); - SASSERT(it != m_rm_table.var_map().end()); - for (unsigned bd_i : it->second) { - order_lemma_on_factor_binomial_rm(m, k, m_rm_table.rms()[bd_i]); - if (done()) - break; - } -} - -void core::order_lemma_on_factor_binomial_rm(const monomial& ac, unsigned k, const rooted_mon& rm_bd) { - factor d(m_evars.find(ac[k]).var(), factor_type::VAR); - factor b; - if (!divide(rm_bd, d, b)) - return; - order_lemma_on_binomial_ac_bd(ac, k, rm_bd, b, d.index()); -} - -void core::order_lemma_on_binomial_ac_bd(const monomial& ac, unsigned k, const rooted_mon& bd, const factor& b, lpvar d) { - TRACE("nla_solver", print_monomial(ac, tout << "ac="); - print_rooted_monomial(bd, tout << "\nrm="); - print_factor(b, tout << ", b="); print_var(d, tout << ", d=") << "\n";); - int p = (k + 1) % 2; - lpvar a = ac[p]; - lpvar c = ac[k]; - SASSERT(m_evars.find(c).var() == d); - rational acv = vvr(ac); - rational av = vvr(a); - rational c_sign = rrat_sign(vvr(c)); - rational d_sign = rrat_sign(vvr(d)); - rational bdv = vvr(bd); - rational bv = vvr(b); - auto av_c_s = av*c_sign; auto bv_d_s = bv*d_sign; - - // suppose ac >= bd, then ac/|c| >= bd/|d|. - // Notice that ac/|c| = a*c_sign , and bd/|d| = b*d_sign - if (acv >= bdv && av_c_s < bv_d_s) - generate_mon_ol(ac, a, c_sign, c, bd, b, d_sign, d, llc::LT); - else if (acv <= bdv && av_c_s > bv_d_s) - generate_mon_ol(ac, a, c_sign, c, bd, b, d_sign, d, llc::GT); - -} - -void core::order_lemma_on_binomial_k(const monomial& m, lpvar k, bool gt) { - SASSERT(gt == (vvr(m) > vvr(m[0]) * vvr(m[1]))); - unsigned p = (k + 1) % 2; - order_lemma_on_binomial_sign(m, m[k], m[p], gt? 1: -1); -} -// sign it the sign of vvr(m) - vvr(m[0]) * vvr(m[1]) -// m = xy -// and val(m) != val(x)*val(y) -// y > 0 and x = a, then xy >= ay -void core::order_lemma_on_binomial_sign(const monomial& ac, lpvar x, lpvar y, int sign) { - SASSERT(!mon_has_zero(ac)); - int sy = rat_sign(y); - add_empty_lemma(); - mk_ineq(y, sy == 1? llc::LE : llc::GE); // negate sy - mk_ineq(x, sy*sign == 1? llc::GT:llc::LT , vvr(x)); // assert x <= vvr(x) if x > 0 - mk_ineq(ac.var(), - vvr(x), y, sign == 1?llc::LE:llc::GE); - TRACE("nla_solver", print_lemma(tout);); -} - -void core::order_lemma_on_binomial(const monomial& ac) { - TRACE("nla_solver", print_monomial(ac, tout);); - SASSERT(!check_monomial(ac) && ac.size() == 2); - const rational & mult_val = vvr(ac[0]) * vvr(ac[1]); - const rational acv = vvr(ac); - bool gt = acv > mult_val; - for (unsigned k = 0; k < 2; k++) { - order_lemma_on_binomial_k(ac, k, gt); - order_lemma_on_factor_binomial_explore(ac, k); - } -} -// a > b && c > 0 => ac > bc -void core::order_lemma_on_rmonomial(const rooted_mon& rm) { - TRACE("nla_solver_details", - tout << "rm = "; print_product(rm, tout); - tout << ", orig = "; print_monomial(m_monomials[rm.orig_index()], tout); - ); - for (auto ac : factorization_factory_imp(rm, *this)) { - if (ac.size() != 2) - continue; - if (ac.is_mon()) - order_lemma_on_binomial(*ac.mon()); - else - order_lemma_on_factorization(rm, ac); - if (done()) - break; - } -} - -void core::order_lemma() { - TRACE("nla_solver", ); - init_rm_to_refine(); - const auto& rm_ref = m_rm_table.to_refine(); - unsigned start = random() % rm_ref.size(); - unsigned i = start; - do { - const rooted_mon& rm = m_rm_table.rms()[rm_ref[i]]; - order_lemma_on_rmonomial(rm); - if (++i == rm_ref.size()) { - i = 0; - } - } while(i != start && !done()); -} std::vector core::get_sorted_key(const rooted_mon& rm) const { std::vector r; @@ -2964,7 +2548,7 @@ lbool core:: inner_check(bool derived) { if (derived) continue; TRACE("nla_solver", tout << "passed derived and basic lemmas\n";); if (search_level == 1) { - order_lemma(); + m_order.order_lemma(); } else { // search_level == 2 monotonicity_lemma(); tangent_lemma(); diff --git a/src/util/lp/nla_core.h b/src/util/lp/nla_core.h index 5a06f16d7..a76506863 100644 --- a/src/util/lp/nla_core.h +++ b/src/util/lp/nla_core.h @@ -1,624 +1,564 @@ -/*++ - Copyright (c) 2017 Microsoft Corporation - - Module Name: - - - - Abstract: - - - - Author: - Nikolaj Bjorner (nbjorner) - Lev Nachmanson (levnach) - - Revision History: - - - --*/ -#pragma once -#include "util/lp/factorization.h" -#include "util/lp/lp_types.h" -#include "util/lp/var_eqs.h" -#include "util/lp/rooted_mons.h" -#include "util/lp/nla_tangent_lemmas.h" -#include "util/lp/nla_basics_lemmas.h" -namespace nla { - -template -bool try_insert(const A& elem, B& collection) { - auto it = collection.find(elem); - if (it != collection.end()) - return false; - collection.insert(elem); - return true; -} - -typedef lp::constraint_index lpci; -typedef lp::lconstraint_kind llc; - -inline int rat_sign(const rational& r) { return r.is_pos()? 1 : ( r.is_neg()? -1 : 0); } -inline rational rrat_sign(const rational& r) { return rational(rat_sign(r)); } -inline bool is_set(unsigned j) { return static_cast(j) != -1; } -inline bool is_even(unsigned k) { return (k >> 1) << 1 == k; } -struct ineq { - lp::lconstraint_kind m_cmp; - lp::lar_term m_term; - rational m_rs; - ineq(lp::lconstraint_kind cmp, const lp::lar_term& term, const rational& rs) : m_cmp(cmp), m_term(term), m_rs(rs) {} - bool operator==(const ineq& a) const { - return m_cmp == a.m_cmp && m_term == a.m_term && m_rs == a.m_rs; - } - const lp::lar_term& term() const { return m_term; }; - lp::lconstraint_kind cmp() const { return m_cmp; }; - const rational& rs() const { return m_rs; }; -}; - -class lemma { - vector m_ineqs; - lp::explanation m_expl; -public: - void push_back(const ineq& i) { m_ineqs.push_back(i);} - size_t size() const { return m_ineqs.size() + m_expl.size(); } - const vector& ineqs() const { return m_ineqs; } - vector& ineqs() { return m_ineqs; } - lp::explanation& expl() { return m_expl; } - const lp::explanation& expl() const { return m_expl; } - bool is_conflict() const { return m_ineqs.empty() && !m_expl.empty(); } -}; - -struct point { - rational x; - rational y; - point(const rational& a, const rational& b): x(a), y(b) {} - point() {} - point& operator*=(rational a) { - x *= a; - y *= a; - return *this; - } -}; - -struct core { - struct bfc { - factor m_x, m_y; - bfc() {} - bfc(const factor& x, const factor& y): m_x(x), m_y(y) {}; - }; - - //fields - var_eqs m_evars; - vector m_monomials; - - rooted_mon_table m_rm_table; - // this field is used for the push/pop operations - unsigned_vector m_monomials_counts; - lp::lar_solver& m_lar_solver; - std::unordered_map> m_monomials_containing_var; - - // m_var_to_its_monomial[m_monomials[i].var()] = i - std::unordered_map m_var_to_its_monomial; - vector * m_lemma_vec; - unsigned_vector m_to_refine; - std::unordered_map m_mkeys; // the key is the sorted vars of a monomial - tangents m_tangents; - basics m_basics; - core(lp::lar_solver& s); - - bool compare_holds(const rational& ls, llc cmp, const rational& rs) const; - - rational value(const lp::lar_term& r) const; - - 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; - - rational vvr(lpvar j) const { return m_lar_solver.get_column_value_rational(j); } - - rational vvr(const monomial& m) const { return m_lar_solver.get_column_value_rational(m.var()); } - - lp::impq vv(lpvar j) const { return m_lar_solver.get_column_value(j); } - - lpvar var(const rooted_mon& rm) const {return m_monomials[rm.orig_index()].var(); } - - rational vvr(const rooted_mon& rm) const { return vvr(m_monomials[rm.orig_index()].var()) * rm.orig_sign(); } - - rational vvr(const factor& f) const { return f.is_var()? vvr(f.index()) : vvr(m_rm_table.rms()[f.index()]); } - - lpvar var(const factor& f) const { return f.is_var()? f.index() : var(m_rm_table.rms()[f.index()]); } - - svector sorted_vars(const factor& f) const; - bool done() const; - - void add_empty_lemma(); - // the value of the factor is equal to the value of the variable multiplied - // by the canonize_sign - rational canonize_sign(const factor& f) const; - - rational canonize_sign_of_var(lpvar j) const; - - // the value of the rooted monomias is equal to the value of the variable multiplied - // by the canonize_sign - rational canonize_sign(const rooted_mon& m) const; - - // returns the monomial index - unsigned add(lpvar v, unsigned sz, lpvar const* vs); - - void push(); - - void deregister_monomial_from_rooted_monomials (const monomial & m, unsigned i); - - void deregister_monomial_from_tables(const monomial & m, unsigned i); - - void pop(unsigned n); - - rational mon_value_by_vars(unsigned i) const; - template - rational product_value(const T & m) const; - - // return true iff the monomial value is equal to the product of the values of the factors - bool check_monomial(const monomial& m) const; - - void explain(const monomial& m, lp::explanation& exp) const; - - void explain(const rooted_mon& rm, lp::explanation& exp) const; - - void explain(const factor& f, lp::explanation& exp) const; - - void explain(lpvar j, lp::explanation& exp) const; - - template - std::ostream& print_product(const T & m, std::ostream& out) const; - - std::ostream & print_factor(const factor& f, std::ostream& out) const; - - std::ostream & print_factor_with_vars(const factor& f, std::ostream& out) const; - - std::ostream& print_monomial(const monomial& m, std::ostream& out) const; - - std::ostream& print_point(const point &a, std::ostream& out) const; - - std::ostream& print_tangent_domain(const point &a, const point &b, std::ostream& out) const; - - std::ostream& print_bfc(const bfc& m, std::ostream& out) const; - - std::ostream& print_monomial(unsigned i, std::ostream& out) const; - - std::ostream& print_monomial_with_vars(unsigned i, std::ostream& out) const; - - template - std::ostream& print_product_with_vars(const T& m, std::ostream& out) const; - - std::ostream& print_monomial_with_vars(const monomial& m, std::ostream& out) const; - - std::ostream& print_rooted_monomial(const rooted_mon& rm, std::ostream& out) const; - - std::ostream& print_rooted_monomial_with_vars(const rooted_mon& rm, std::ostream& out) const; - - std::ostream& print_explanation(const lp::explanation& exp, std::ostream& out) const; - - bool explain_upper_bound(const lp::lar_term& t, const rational& rs, lp::explanation& e) const; - bool explain_lower_bound(const lp::lar_term& t, const rational& rs, lp::explanation& e) const; - - bool explain_coeff_lower_bound(const lp::lar_term::ival& p, rational& bound, lp::explanation& e) const; - - bool explain_coeff_upper_bound(const lp::lar_term::ival& p, rational& bound, lp::explanation& e) const; - - // return true iff the negation of the ineq can be derived from the constraints - bool explain_ineq(const lp::lar_term& t, llc cmp, const rational& rs); - - /** - * \brief - if t is an octagon term -+x -+ y try to explain why the term always - equal zero - */ - bool explain_by_equiv(const lp::lar_term& t, lp::explanation& e); - - void mk_ineq(lp::lar_term& t, llc cmp, const rational& rs); - void mk_ineq(const rational& a, lpvar j, const rational& b, lpvar k, llc cmp, const rational& rs); - - void mk_ineq(lpvar j, const rational& b, lpvar k, llc cmp, const rational& rs); - - void mk_ineq(lpvar j, const rational& b, lpvar k, llc cmp); - - void mk_ineq(const rational& a, lpvar j, const rational& b, lpvar k, llc cmp); - - void mk_ineq(const rational& a ,lpvar j, lpvar k, llc cmp, const rational& rs); - - void mk_ineq(lpvar j, lpvar k, llc cmp, const rational& rs); - - void mk_ineq(lpvar j, llc cmp, const rational& rs); - - void mk_ineq(const rational& a, lpvar j, llc cmp, const rational& rs); - - llc negate(llc cmp); - - llc apply_minus(llc cmp); - - void mk_ineq(const rational& a, lpvar j, llc cmp); - - void mk_ineq(lpvar j, lpvar k, llc cmp, lemma& l); - - void mk_ineq(lpvar j, llc cmp); - - // the monomials should be equal by modulo sign but this is not so in the model - void fill_explanation_and_lemma_sign(const monomial& a, const monomial & b, rational const& sign); - - // Replaces each variable index by the root in the tree and flips the sign if the var comes with a minus. - // Also sorts the result. - // - svector reduce_monomial_to_rooted(const svector & vars, rational & sign) const; - - - // Replaces definition m_v = v1* .. * vn by - // m_v = coeff * w1 * ... * wn, where w1, .., wn are canonical - // representatives, which are the roots of the equivalence tree, under current equations. - // - monomial_coeff canonize_monomial(monomial const& m) const; - - lemma& current_lemma(); - const lemma& current_lemma() const; - vector& current_ineqs(); - lp::explanation& current_expl(); - const lp::explanation& current_expl() const; - - int vars_sign(const svector& v); - - bool has_upper_bound(lpvar j) const; - - bool has_lower_bound(lpvar j) const; - const rational& get_upper_bound(unsigned j) const; - - const rational& get_lower_bound(unsigned j) const; - - - bool zero_is_an_inner_point_of_bounds(lpvar j) const; - - int rat_sign(const monomial& m) const; - inline int rat_sign(lpvar j) const { return nla::rat_sign(vvr(j)); } - - // Returns true if the monomial sign is incorrect - bool sign_contradiction(const monomial& m) const; - - /* - unsigned_vector eq_vars(lpvar j) const; - */ - // Monomials m and n vars have the same values, up to "sign" - // Generate a lemma if values of m.var() and n.var() are not the same up to sign - - bool var_is_fixed_to_zero(lpvar j) const; - bool var_is_fixed_to_val(lpvar j, const rational& v) const; - - bool var_is_fixed(lpvar j) const; - - std::ostream & print_ineq(const ineq & in, std::ostream & out) const; - - std::ostream & print_var(lpvar j, std::ostream & out) const; - - std::ostream & print_monomials(std::ostream & out) const; - - std::ostream & print_ineqs(const lemma& l, std::ostream & out) const; - - std::ostream & print_factorization(const factorization& f, std::ostream& out) const; - - bool find_rm_monomial_of_vars(const svector& vars, unsigned & i) const; - - const monomial* find_monomial_of_vars(const svector& vars) const; - - void explain_existing_lower_bound(lpvar j); - - void explain_existing_upper_bound(lpvar j); - - void explain_separation_from_zero(lpvar j); - - int get_derived_sign(const rooted_mon& rm, const factorization& f) const; - // here we use the fact xy = 0 -> x = 0 or y = 0 - void trace_print_monomial_and_factorization(const rooted_mon& rm, const factorization& f, std::ostream& out) const; - - void explain_var_separated_from_zero(lpvar j); - - void explain_fixed_var(lpvar j); - // x = 0 or y = 0 -> xy = 0 - - bool var_has_positive_lower_bound(lpvar j) const; - - bool var_has_negative_upper_bound(lpvar j) const; - - bool var_is_separated_from_zero(lpvar j) const; - - - bool vars_are_equiv(lpvar a, lpvar b) const; - - - void explain_equiv_vars(lpvar a, lpvar b); - // use the fact that - // |xabc| = |x| and x != 0 -> |a| = |b| = |c| = 1 - void explain(const factorization& f, lp::explanation& exp); - - bool has_zero_factor(const factorization& factorization) const; - - - template - bool has_zero(const T& product) const; - - template - bool mon_has_zero(const T& product) const; - void init_rm_to_refine(); - lp::lp_settings& settings(); - unsigned random(); - void map_monomial_vars_to_monomial_indices(unsigned i); - - void map_vars_to_monomials(); - - // we look for octagon constraints here, with a left part +-x +- y - void collect_equivs(); - - void collect_equivs_of_fixed_vars(); - - // returns true iff the term is in a form +-x-+y. - // the sign is true iff the term is x+y, -x-y. - bool is_octagon_term(const lp::lar_term& t, bool & sign, lpvar& i, lpvar &j) const; - - void add_equivalence_maybe(const lp::lar_term *t, lpci c0, lpci c1); - - // x is equivalent to y if x = +- y - void init_vars_equivalence(); - - bool vars_table_is_ok() const; - - bool rm_table_is_ok() const; - - bool tables_are_ok() const; - - bool var_is_a_root(lpvar j) const; - - template - bool vars_are_roots(const T& v) const; - - void register_monomial_in_tables(unsigned i_mon); - - template - void trace_print_rms(const T& p, std::ostream& out); - - void print_monomial_stats(const monomial& m, std::ostream& out); - - void print_stats(std::ostream& out); - - void register_monomials_in_tables(); - - void clear(); - - void init_search(); - -#if 0 - void init_to_refine(); -#endif - void init_to_refine(); - - bool divide(const rooted_mon& bc, const factor& c, factor & b) const; - - void negate_factor_equality(const factor& c, - const factor& d); - - void negate_factor_relation(const rational& a_sign, const factor& a, const rational& b_sign, const factor& b); - - void negate_var_factor_relation(const rational& a_sign, lpvar a, const rational& b_sign, const factor& b); - - // |c_sign| = 1, and c*c_sign > 0 - // ac > bc => ac/|c| > bc/|c| => a*c_sign > b*c_sign - void generate_ol(const rooted_mon& ac, - const factor& a, - int c_sign, - const factor& c, - const rooted_mon& bc, - const factor& b, - llc ab_cmp); - - void generate_mon_ol(const monomial& ac, - lpvar a, - const rational& c_sign, - lpvar c, - const rooted_mon& bd, - const factor& b, - const rational& d_sign, - lpvar d, - llc ab_cmp); - - std::unordered_set collect_vars(const lemma& l) const; - - void print_lemma(std::ostream& out); - - void print_specific_lemma(const lemma& l, std::ostream& out) const; - - - void trace_print_ol(const rooted_mon& ac, - const factor& a, - const factor& c, - const rooted_mon& bc, - const factor& b, - std::ostream& out); - - bool order_lemma_on_ac_and_bc_and_factors(const rooted_mon& ac, - const factor& a, - const factor& c, - const rooted_mon& bc, - const factor& b); - - // a >< b && c > 0 => ac >< bc - // a >< b && c < 0 => ac <> bc - // ac[k] plays the role of c - bool order_lemma_on_ac_and_bc(const rooted_mon& rm_ac, - const factorization& ac_f, - unsigned k, - const rooted_mon& rm_bd); - void maybe_add_a_factor(lpvar i, - const factor& c, - std::unordered_set& found_vars, - std::unordered_set& found_rm, - vector & r) const; - - bool order_lemma_on_ac_explore(const rooted_mon& rm, const factorization& ac, unsigned k); - - void order_lemma_on_factorization(const rooted_mon& rm, const factorization& ab); - - /** - \brief Add lemma: - a > 0 & b <= value(b) => sign*ab <= value(b)*a if value(a) > 0 - a < 0 & b >= value(b) => sign*ab <= value(b)*a if value(a) < 0 - */ - void order_lemma_on_ab_gt(const monomial& m, const rational& sign, lpvar a, lpvar b); - // we need to deduce ab >= vvr(b)*a - /** - \brief Add lemma: - a > 0 & b >= value(b) => sign*ab >= value(b)*a if value(a) > 0 - a < 0 & b <= value(b) => sign*ab >= value(b)*a if value(a) < 0 - */ - void order_lemma_on_ab_lt(const monomial& m, const rational& sign, lpvar a, lpvar b); - void order_lemma_on_ab(const monomial& m, const rational& sign, lpvar a, lpvar b, bool gt); - bool rm_check(const rooted_mon&) const; - void order_lemma_on_factor_binomial_explore(const monomial& m, unsigned k); - void order_lemma_on_factor_binomial_rm(const monomial& ac, unsigned k, const rooted_mon& rm_bd); - void order_lemma_on_binomial_ac_bd(const monomial& ac, unsigned k, const rooted_mon& bd, const factor& b, lpvar d); - void order_lemma_on_binomial_k(const monomial& m, lpvar k, bool gt); - void order_lemma_on_binomial_sign(const monomial& ac, lpvar x, lpvar y, int sign); - void order_lemma_on_binomial(const monomial& ac); - void order_lemma_on_rmonomial(const rooted_mon& rm); - void order_lemma(); - void print_monotone_array(const vector, unsigned>>& lex_sorted, - std::ostream& out) const; - std::vector get_sorted_key(const rooted_mon& rm) const; - std::unordered_map get_rm_by_arity(); - bool uniform_le(const std::vector& a, - const std::vector& b, - unsigned & strict_i) const; - vector> get_sorted_key_with_vars(const rooted_mon& a) const; - void negate_abs_a_le_abs_b(lpvar a, lpvar b, bool strict); - - void negate_abs_a_lt_abs_b(lpvar a, lpvar b); - void assert_abs_val_a_le_abs_var_b( const rooted_mon& a, const rooted_mon& b, bool strict); - - // strict version - void generate_monl_strict(const rooted_mon& a, - const rooted_mon& b, - unsigned strict); - - // not a strict version - void generate_monl(const rooted_mon& a, - const rooted_mon& b); - - bool monotonicity_lemma_on_lex_sorted_rm_upper(const vector, unsigned>>& lex_sorted, unsigned i, const rooted_mon& rm); - - bool monotonicity_lemma_on_lex_sorted_rm_lower(const vector, unsigned>>& lex_sorted, unsigned i, const rooted_mon& rm); - - bool monotonicity_lemma_on_lex_sorted_rm(const vector, unsigned>>& lex_sorted, unsigned i, const rooted_mon& rm); - bool monotonicity_lemma_on_lex_sorted(const vector, unsigned>>& lex_sorted); - - bool monotonicity_lemma_on_rms_of_same_arity(const unsigned_vector& rms); - - void monotonicity_lemma(); - -#if 0 - void monotonicity_lemma(); - -#endif - - void monotonicity_lemma(unsigned i_mon); - - // add |j| <= |vvr(j)| - void var_abs_val_le(lpvar j); - - // assert that |j| >= |vvr(j)| - void var_abs_val_ge(lpvar j); - - - /** - \brief Add |v| ~ |bound| - where ~ is <, <=, >, >=, - and bound = vvr(v) - - |v| > |bound| - <=> - (v < 0 or v > |bound|) & (v > 0 or -v > |bound|) - => Let s be the sign of vvr(v) - (s*v < 0 or s*v > |bound|) - - |v| < |bound| - <=> - v < |bound| & -v < |bound| - => Let s be the sign of vvr(v) - s*v < |bound| - - */ - - void add_abs_bound(lpvar v, llc cmp); - - void add_abs_bound(lpvar v, llc cmp, rational const& bound); - - /** \brief enforce the inequality |m| <= product |m[i]| . - by enforcing lemma: - /\_i |m[i]| <= |vvr(m[i])| => |m| <= |product_i vvr(m[i])| - <=> - \/_i |m[i]| > |vvr(m[i])} or |m| <= |product_i vvr(m[i])| - */ - - void monotonicity_lemma_gt(const monomial& m, const rational& prod_val); - - /** \brief enforce the inequality |m| >= product |m[i]| . - - /\_i |m[i]| >= |vvr(m[i])| => |m| >= |product_i vvr(m[i])| - <=> - \/_i |m[i]| < |vvr(m[i])} or |m| >= |product_i vvr(m[i])| - */ - void monotonicity_lemma_lt(const monomial& m, const rational& prod_val); - - bool find_bfc_to_refine_on_rmonomial(const rooted_mon& rm, bfc & bf); - - bool find_bfc_to_refine(bfc& bf, lpvar &j, rational& sign, const rooted_mon*& rm_found); - void generate_simple_sign_lemma(const rational& sign, const monomial& m); - - void generate_simple_tangent_lemma(const rooted_mon* rm); - - void tangent_lemma(); - - void generate_explanations_of_tang_lemma(const rooted_mon& rm, const bfc& bf, lp::explanation& exp); - - void tangent_lemma_bf(const bfc& bf, lpvar j, const rational& sign, const rooted_mon* rm); - void generate_tang_plane(const rational & a, const rational& b, const factor& x, const factor& y, bool below, lpvar j, const rational& j_sign); - - void negate_relation(unsigned j, const rational& a); - - void generate_two_tang_lines(const bfc & bf, const point& xy, const rational& sign, lpvar j); - // Get two planes tangent to surface z = xy, one at point a, and another at point b. - // One can show that these planes still create a cut. - void get_initial_tang_points(point &a, point &b, const point& xy, - bool below) const; - - void push_tang_point(point &a, const point& xy, bool below, const rational& correct_val, const rational& val) const; - - void push_tang_points(point &a, point &b, const point& xy, bool below, const rational& correct_val, const rational& val) const; - - rational tang_plane(const point& a, const point& x) const; - - bool plane_is_correct_cut(const point& plane, - const point& xy, - const rational & correct_val, - const rational & val, - bool below) const; - - // "below" means that the val is below the surface xy - void get_tang_points(point &a, point &b, bool below, const rational& val, - const point& xy) const; - - bool conflict_found() const; - lbool inner_check(bool derived); - - lbool check(vector& l_vec); - - bool no_lemmas_hold() const; - - void test_factorization(unsigned /*mon_index*/, unsigned /*number_of_factorizations*/); - - lbool test_check(vector& l); -}; // end of core -} // end of namespace nla +/*++ + Copyright (c) 2017 Microsoft Corporation + + Module Name: + + + + Abstract: + + + + Author: + Nikolaj Bjorner (nbjorner) + Lev Nachmanson (levnach) + + Revision History: + + + --*/ +#pragma once +#include "util/lp/factorization.h" +#include "util/lp/lp_types.h" +#include "util/lp/var_eqs.h" +#include "util/lp/rooted_mons.h" +#include "util/lp/nla_tangent_lemmas.h" +#include "util/lp/nla_basics_lemmas.h" +#include "util/lp/nla_order_lemmas.h" +namespace nla { + +template +bool try_insert(const A& elem, B& collection) { + auto it = collection.find(elem); + if (it != collection.end()) + return false; + collection.insert(elem); + return true; +} + +typedef lp::constraint_index lpci; +typedef lp::lconstraint_kind llc; + +inline int rat_sign(const rational& r) { return r.is_pos()? 1 : ( r.is_neg()? -1 : 0); } +inline rational rrat_sign(const rational& r) { return rational(rat_sign(r)); } +inline bool is_set(unsigned j) { return static_cast(j) != -1; } +inline bool is_even(unsigned k) { return (k >> 1) << 1 == k; } +struct ineq { + lp::lconstraint_kind m_cmp; + lp::lar_term m_term; + rational m_rs; + ineq(lp::lconstraint_kind cmp, const lp::lar_term& term, const rational& rs) : m_cmp(cmp), m_term(term), m_rs(rs) {} + bool operator==(const ineq& a) const { + return m_cmp == a.m_cmp && m_term == a.m_term && m_rs == a.m_rs; + } + const lp::lar_term& term() const { return m_term; }; + lp::lconstraint_kind cmp() const { return m_cmp; }; + const rational& rs() const { return m_rs; }; +}; + +class lemma { + vector m_ineqs; + lp::explanation m_expl; +public: + void push_back(const ineq& i) { m_ineqs.push_back(i);} + size_t size() const { return m_ineqs.size() + m_expl.size(); } + const vector& ineqs() const { return m_ineqs; } + vector& ineqs() { return m_ineqs; } + lp::explanation& expl() { return m_expl; } + const lp::explanation& expl() const { return m_expl; } + bool is_conflict() const { return m_ineqs.empty() && !m_expl.empty(); } +}; + +struct point { + rational x; + rational y; + point(const rational& a, const rational& b): x(a), y(b) {} + point() {} + point& operator*=(rational a) { + x *= a; + y *= a; + return *this; + } +}; + +struct core { + struct bfc { + factor m_x, m_y; + bfc() {} + bfc(const factor& x, const factor& y): m_x(x), m_y(y) {}; + }; + + //fields + var_eqs m_evars; + vector m_monomials; + + rooted_mon_table m_rm_table; + // this field is used for the push/pop operations + unsigned_vector m_monomials_counts; + lp::lar_solver& m_lar_solver; + std::unordered_map> m_monomials_containing_var; + + // m_var_to_its_monomial[m_monomials[i].var()] = i + std::unordered_map m_var_to_its_monomial; + vector * m_lemma_vec; + unsigned_vector m_to_refine; + std::unordered_map m_mkeys; // the key is the sorted vars of a monomial + tangents m_tangents; + basics m_basics; + order m_order; + // methods + unsigned find_monomial(const unsigned_vector& k) const; + core(lp::lar_solver& s); + + bool compare_holds(const rational& ls, llc cmp, const rational& rs) const; + + rational value(const lp::lar_term& r) const; + + 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; + + rational vvr(lpvar j) const { return m_lar_solver.get_column_value_rational(j); } + + rational vvr(const monomial& m) const { return m_lar_solver.get_column_value_rational(m.var()); } + + lp::impq vv(lpvar j) const { return m_lar_solver.get_column_value(j); } + + lpvar var(const rooted_mon& rm) const {return m_monomials[rm.orig_index()].var(); } + + rational vvr(const rooted_mon& rm) const { return vvr(m_monomials[rm.orig_index()].var()) * rm.orig_sign(); } + + rational vvr(const factor& f) const { return f.is_var()? vvr(f.index()) : vvr(m_rm_table.rms()[f.index()]); } + + lpvar var(const factor& f) const { return f.is_var()? f.index() : var(m_rm_table.rms()[f.index()]); } + + svector sorted_vars(const factor& f) const; + bool done() const; + + void add_empty_lemma(); + // the value of the factor is equal to the value of the variable multiplied + // by the canonize_sign + rational canonize_sign(const factor& f) const; + + rational canonize_sign_of_var(lpvar j) const; + + // the value of the rooted monomias is equal to the value of the variable multiplied + // by the canonize_sign + rational canonize_sign(const rooted_mon& m) const; + + // returns the monomial index + unsigned add(lpvar v, unsigned sz, lpvar const* vs); + + void push(); + + void deregister_monomial_from_rooted_monomials (const monomial & m, unsigned i); + + void deregister_monomial_from_tables(const monomial & m, unsigned i); + + void pop(unsigned n); + + rational mon_value_by_vars(unsigned i) const; + template + rational product_value(const T & m) const; + + // return true iff the monomial value is equal to the product of the values of the factors + bool check_monomial(const monomial& m) const; + + void explain(const monomial& m, lp::explanation& exp) const; + + void explain(const rooted_mon& rm, lp::explanation& exp) const; + + void explain(const factor& f, lp::explanation& exp) const; + + void explain(lpvar j, lp::explanation& exp) const; + + template + std::ostream& print_product(const T & m, std::ostream& out) const; + + std::ostream & print_factor(const factor& f, std::ostream& out) const; + + std::ostream & print_factor_with_vars(const factor& f, std::ostream& out) const; + + std::ostream& print_monomial(const monomial& m, std::ostream& out) const; + + std::ostream& print_point(const point &a, std::ostream& out) const; + + std::ostream& print_tangent_domain(const point &a, const point &b, std::ostream& out) const; + + std::ostream& print_bfc(const bfc& m, std::ostream& out) const; + + std::ostream& print_monomial(unsigned i, std::ostream& out) const; + + std::ostream& print_monomial_with_vars(unsigned i, std::ostream& out) const; + + template + std::ostream& print_product_with_vars(const T& m, std::ostream& out) const; + + std::ostream& print_monomial_with_vars(const monomial& m, std::ostream& out) const; + + std::ostream& print_rooted_monomial(const rooted_mon& rm, std::ostream& out) const; + + std::ostream& print_rooted_monomial_with_vars(const rooted_mon& rm, std::ostream& out) const; + + std::ostream& print_explanation(const lp::explanation& exp, std::ostream& out) const; + + bool explain_upper_bound(const lp::lar_term& t, const rational& rs, lp::explanation& e) const; + bool explain_lower_bound(const lp::lar_term& t, const rational& rs, lp::explanation& e) const; + + bool explain_coeff_lower_bound(const lp::lar_term::ival& p, rational& bound, lp::explanation& e) const; + + bool explain_coeff_upper_bound(const lp::lar_term::ival& p, rational& bound, lp::explanation& e) const; + + // return true iff the negation of the ineq can be derived from the constraints + bool explain_ineq(const lp::lar_term& t, llc cmp, const rational& rs); + + /** + * \brief + if t is an octagon term -+x -+ y try to explain why the term always + equal zero + */ + bool explain_by_equiv(const lp::lar_term& t, lp::explanation& e); + + void mk_ineq(lp::lar_term& t, llc cmp, const rational& rs); + void mk_ineq(const rational& a, lpvar j, const rational& b, lpvar k, llc cmp, const rational& rs); + + void mk_ineq(lpvar j, const rational& b, lpvar k, llc cmp, const rational& rs); + + void mk_ineq(lpvar j, const rational& b, lpvar k, llc cmp); + + void mk_ineq(const rational& a, lpvar j, const rational& b, lpvar k, llc cmp); + + void mk_ineq(const rational& a ,lpvar j, lpvar k, llc cmp, const rational& rs); + + void mk_ineq(lpvar j, lpvar k, llc cmp, const rational& rs); + + void mk_ineq(lpvar j, llc cmp, const rational& rs); + + void mk_ineq(const rational& a, lpvar j, llc cmp, const rational& rs); + void maybe_add_a_factor(lpvar i, + const factor& c, + std::unordered_set& found_vars, + std::unordered_set& found_rm, + vector & r) const; + + llc apply_minus(llc cmp); + + void mk_ineq(const rational& a, lpvar j, llc cmp); + + void mk_ineq(lpvar j, lpvar k, llc cmp, lemma& l); + + void mk_ineq(lpvar j, llc cmp); + + // the monomials should be equal by modulo sign but this is not so in the model + void fill_explanation_and_lemma_sign(const monomial& a, const monomial & b, rational const& sign); + + // Replaces each variable index by the root in the tree and flips the sign if the var comes with a minus. + // Also sorts the result. + // + svector reduce_monomial_to_rooted(const svector & vars, rational & sign) const; + + + // Replaces definition m_v = v1* .. * vn by + // m_v = coeff * w1 * ... * wn, where w1, .., wn are canonical + // representatives, which are the roots of the equivalence tree, under current equations. + // + monomial_coeff canonize_monomial(monomial const& m) const; + + lemma& current_lemma(); + const lemma& current_lemma() const; + vector& current_ineqs(); + lp::explanation& current_expl(); + const lp::explanation& current_expl() const; + + int vars_sign(const svector& v); + + bool has_upper_bound(lpvar j) const; + + bool has_lower_bound(lpvar j) const; + const rational& get_upper_bound(unsigned j) const; + + const rational& get_lower_bound(unsigned j) const; + + + bool zero_is_an_inner_point_of_bounds(lpvar j) const; + + int rat_sign(const monomial& m) const; + inline int rat_sign(lpvar j) const { return nla::rat_sign(vvr(j)); } + + // Returns true if the monomial sign is incorrect + bool sign_contradiction(const monomial& m) const; + + /* + unsigned_vector eq_vars(lpvar j) const; + */ + // Monomials m and n vars have the same values, up to "sign" + // Generate a lemma if values of m.var() and n.var() are not the same up to sign + + bool var_is_fixed_to_zero(lpvar j) const; + bool var_is_fixed_to_val(lpvar j, const rational& v) const; + + bool var_is_fixed(lpvar j) const; + + std::ostream & print_ineq(const ineq & in, std::ostream & out) const; + + std::ostream & print_var(lpvar j, std::ostream & out) const; + + std::ostream & print_monomials(std::ostream & out) const; + + std::ostream & print_ineqs(const lemma& l, std::ostream & out) const; + + std::ostream & print_factorization(const factorization& f, std::ostream& out) const; + + bool find_rm_monomial_of_vars(const svector& vars, unsigned & i) const; + + const monomial* find_monomial_of_vars(const svector& vars) const; + + void explain_existing_lower_bound(lpvar j); + + void explain_existing_upper_bound(lpvar j); + + void explain_separation_from_zero(lpvar j); + + int get_derived_sign(const rooted_mon& rm, const factorization& f) const; + // here we use the fact xy = 0 -> x = 0 or y = 0 + void trace_print_monomial_and_factorization(const rooted_mon& rm, const factorization& f, std::ostream& out) const; + + void explain_var_separated_from_zero(lpvar j); + + void explain_fixed_var(lpvar j); + // x = 0 or y = 0 -> xy = 0 + + bool var_has_positive_lower_bound(lpvar j) const; + + bool var_has_negative_upper_bound(lpvar j) const; + + bool var_is_separated_from_zero(lpvar j) const; + + + bool vars_are_equiv(lpvar a, lpvar b) const; + + + void explain_equiv_vars(lpvar a, lpvar b); + // use the fact that + // |xabc| = |x| and x != 0 -> |a| = |b| = |c| = 1 + void explain(const factorization& f, lp::explanation& exp); + + bool has_zero_factor(const factorization& factorization) const; + + + template + bool has_zero(const T& product) const; + + template + bool mon_has_zero(const T& product) const; + void init_rm_to_refine(); + lp::lp_settings& settings(); + unsigned random(); + void map_monomial_vars_to_monomial_indices(unsigned i); + + void map_vars_to_monomials(); + + // we look for octagon constraints here, with a left part +-x +- y + void collect_equivs(); + + void collect_equivs_of_fixed_vars(); + + // returns true iff the term is in a form +-x-+y. + // the sign is true iff the term is x+y, -x-y. + bool is_octagon_term(const lp::lar_term& t, bool & sign, lpvar& i, lpvar &j) const; + + void add_equivalence_maybe(const lp::lar_term *t, lpci c0, lpci c1); + + // x is equivalent to y if x = +- y + void init_vars_equivalence(); + + bool vars_table_is_ok() const; + + bool rm_table_is_ok() const; + + bool tables_are_ok() const; + + bool var_is_a_root(lpvar j) const; + + template + bool vars_are_roots(const T& v) const; + + void register_monomial_in_tables(unsigned i_mon); + + template + void trace_print_rms(const T& p, std::ostream& out); + + void print_monomial_stats(const monomial& m, std::ostream& out); + + void print_stats(std::ostream& out); + + void register_monomials_in_tables(); + + void clear(); + + void init_search(); + +#if 0 + void init_to_refine(); +#endif + void init_to_refine(); + + bool divide(const rooted_mon& bc, const factor& c, factor & b) const; + + void negate_factor_equality(const factor& c, + const factor& d); + + void negate_factor_relation(const rational& a_sign, const factor& a, const rational& b_sign, const factor& b); + + std::unordered_set collect_vars(const lemma& l) const; + + std::ostream& print_lemma(std::ostream& out) const; + + void print_specific_lemma(const lemma& l, std::ostream& out) const; + + + void trace_print_ol(const rooted_mon& ac, + const factor& a, + const factor& c, + const rooted_mon& bc, + const factor& b, + std::ostream& out); + + bool rm_check(const rooted_mon&) const; + void print_monotone_array(const vector, unsigned>>& lex_sorted, + std::ostream& out) const; + std::vector get_sorted_key(const rooted_mon& rm) const; + std::unordered_map get_rm_by_arity(); + bool uniform_le(const std::vector& a, + const std::vector& b, + unsigned & strict_i) const; + vector> get_sorted_key_with_vars(const rooted_mon& a) const; + void negate_abs_a_le_abs_b(lpvar a, lpvar b, bool strict); + + void negate_abs_a_lt_abs_b(lpvar a, lpvar b); + void assert_abs_val_a_le_abs_var_b( const rooted_mon& a, const rooted_mon& b, bool strict); + + // strict version + void generate_monl_strict(const rooted_mon& a, + const rooted_mon& b, + unsigned strict); + + // not a strict version + void generate_monl(const rooted_mon& a, + const rooted_mon& b); + + bool monotonicity_lemma_on_lex_sorted_rm_upper(const vector, unsigned>>& lex_sorted, unsigned i, const rooted_mon& rm); + + bool monotonicity_lemma_on_lex_sorted_rm_lower(const vector, unsigned>>& lex_sorted, unsigned i, const rooted_mon& rm); + + bool monotonicity_lemma_on_lex_sorted_rm(const vector, unsigned>>& lex_sorted, unsigned i, const rooted_mon& rm); + bool monotonicity_lemma_on_lex_sorted(const vector, unsigned>>& lex_sorted); + + bool monotonicity_lemma_on_rms_of_same_arity(const unsigned_vector& rms); + + void monotonicity_lemma(); + +#if 0 + void monotonicity_lemma(); + +#endif + + void monotonicity_lemma(unsigned i_mon); + + // add |j| <= |vvr(j)| + void var_abs_val_le(lpvar j); + + // assert that |j| >= |vvr(j)| + void var_abs_val_ge(lpvar j); + + + /** + \brief Add |v| ~ |bound| + where ~ is <, <=, >, >=, + and bound = vvr(v) + + |v| > |bound| + <=> + (v < 0 or v > |bound|) & (v > 0 or -v > |bound|) + => Let s be the sign of vvr(v) + (s*v < 0 or s*v > |bound|) + + |v| < |bound| + <=> + v < |bound| & -v < |bound| + => Let s be the sign of vvr(v) + s*v < |bound| + + */ + + void add_abs_bound(lpvar v, llc cmp); + + void add_abs_bound(lpvar v, llc cmp, rational const& bound); + + /** \brief enforce the inequality |m| <= product |m[i]| . + by enforcing lemma: + /\_i |m[i]| <= |vvr(m[i])| => |m| <= |product_i vvr(m[i])| + <=> + \/_i |m[i]| > |vvr(m[i])} or |m| <= |product_i vvr(m[i])| + */ + + void monotonicity_lemma_gt(const monomial& m, const rational& prod_val); + + /** \brief enforce the inequality |m| >= product |m[i]| . + + /\_i |m[i]| >= |vvr(m[i])| => |m| >= |product_i vvr(m[i])| + <=> + \/_i |m[i]| < |vvr(m[i])} or |m| >= |product_i vvr(m[i])| + */ + void monotonicity_lemma_lt(const monomial& m, const rational& prod_val); + + bool find_bfc_to_refine_on_rmonomial(const rooted_mon& rm, bfc & bf); + + bool find_bfc_to_refine(bfc& bf, lpvar &j, rational& sign, const rooted_mon*& rm_found); + void generate_simple_sign_lemma(const rational& sign, const monomial& m); + + void generate_simple_tangent_lemma(const rooted_mon* rm); + + void tangent_lemma(); + + void generate_explanations_of_tang_lemma(const rooted_mon& rm, const bfc& bf, lp::explanation& exp); + + void tangent_lemma_bf(const bfc& bf, lpvar j, const rational& sign, const rooted_mon* rm); + void generate_tang_plane(const rational & a, const rational& b, const factor& x, const factor& y, bool below, lpvar j, const rational& j_sign); + + void negate_relation(unsigned j, const rational& a); + + void generate_two_tang_lines(const bfc & bf, const point& xy, const rational& sign, lpvar j); + // Get two planes tangent to surface z = xy, one at point a, and another at point b. + // One can show that these planes still create a cut. + void get_initial_tang_points(point &a, point &b, const point& xy, + bool below) const; + + void push_tang_point(point &a, const point& xy, bool below, const rational& correct_val, const rational& val) const; + + void push_tang_points(point &a, point &b, const point& xy, bool below, const rational& correct_val, const rational& val) const; + + rational tang_plane(const point& a, const point& x) const; + + bool plane_is_correct_cut(const point& plane, + const point& xy, + const rational & correct_val, + const rational & val, + bool below) const; + + // "below" means that the val is below the surface xy + void get_tang_points(point &a, point &b, bool below, const rational& val, + const point& xy) const; + + bool conflict_found() const; + lbool inner_check(bool derived); + + lbool check(vector& l_vec); + + bool no_lemmas_hold() const; + + void test_factorization(unsigned /*mon_index*/, unsigned /*number_of_factorizations*/); + + lbool test_check(vector& l); +}; // end of core +} // end of namespace nla diff --git a/src/util/lp/nla_order_lemmas.cpp b/src/util/lp/nla_order_lemmas.cpp new file mode 100644 index 000000000..6a30ce5b1 --- /dev/null +++ b/src/util/lp/nla_order_lemmas.cpp @@ -0,0 +1,333 @@ +/*++ + Copyright (c) 2017 Microsoft Corporation + + Module Name: + + + + Abstract: + + + + Author: + Nikolaj Bjorner (nbjorner) + Lev Nachmanson (levnach) + + Revision History: + + + --*/ + +#include "util/lp/nla_order_lemmas.h" +#include "util/lp/nla_core.h" +#include "util/lp/nla_common.h" +#include "util/lp/factorization_factory_imp.h" + +namespace nla { + +// a >< b && c > 0 => ac >< bc +// a >< b && c < 0 => ac <> bc +// ac[k] plays the role of c + +bool order::order_lemma_on_ac_and_bc(const rooted_mon& rm_ac, + const factorization& ac_f, + unsigned k, + const rooted_mon& rm_bd) { + TRACE("nla_solver", tout << "rm_ac = "; + c().print_rooted_monomial(rm_ac, tout); + tout << "\nrm_bd = "; + c().print_rooted_monomial(rm_bd, tout); + tout << "\nac_f[k] = "; + c().print_factor_with_vars(ac_f[k], tout);); + factor b; + if (!c().divide(rm_bd, ac_f[k], b)){ + return false; + } + + return order_lemma_on_ac_and_bc_and_factors(rm_ac, ac_f[(k + 1) % 2], ac_f[k], rm_bd, b); +} + +bool order::order_lemma_on_ac_explore(const rooted_mon& rm, const factorization& ac, unsigned k) { + const factor c = ac[k]; + TRACE("nla_solver", tout << "c = "; _().print_factor_with_vars(c, tout); ); + if (c.is_var()) { + TRACE("nla_solver", tout << "var(c) = " << var(c);); + for (unsigned rm_bc : _().m_rm_table.var_map()[c.index()]) { + TRACE("nla_solver", ); + if (order_lemma_on_ac_and_bc(rm ,ac, k, _().m_rm_table.rms()[rm_bc])) { + return true; + } + } + } else { + for (unsigned rm_bc : _().m_rm_table.proper_multiples()[c.index()]) { + if (order_lemma_on_ac_and_bc(rm , ac, k, _().m_rm_table.rms()[rm_bc])) { + return true; + } + } + } + return false; +} + +void order::order_lemma_on_factorization(const rooted_mon& rm, const factorization& ab) { + const monomial& m = _().m_monomials[rm.orig_index()]; + TRACE("nla_solver", tout << "orig_sign = " << rm.orig_sign() << "\n";); + rational sign = rm.orig_sign(); + for(factor f: ab) + sign *= _().canonize_sign(f); + const rational & fv = vvr(ab[0]) * vvr(ab[1]); + const rational mv = sign * vvr(m); + TRACE("nla_solver", + tout << "ab.size()=" << ab.size() << "\n"; + tout << "we should have sign*vvr(m):" << mv << "=(" << sign << ")*(" << vvr(m) <<") to be equal to " << " vvr(ab[0])*vvr(ab[1]):" << fv << "\n";); + if (mv == fv) + return; + bool gt = mv > fv; + TRACE("nla_solver_f", tout << "m="; _().print_monomial_with_vars(m, tout); tout << "\nfactorization="; _().print_factorization(ab, tout);); + for (unsigned j = 0, k = 1; j < 2; j++, k--) { + order_lemma_on_ab(m, sign, var(ab[k]), var(ab[j]), gt); + explain(ab); explain(m); + explain(rm); + TRACE("nla_solver", _().print_lemma(tout);); + order_lemma_on_ac_explore(rm, ab, j); + } +} +// |c_sign| = 1, and c*c_sign > 0 +// ac > bc => ac/|c| > bc/|c| => a*c_sign > b*c_sign +void order::generate_ol(const rooted_mon& ac, + const factor& a, + int c_sign, + const factor& c, + const rooted_mon& bc, + const factor& b, + llc ab_cmp) { + add_empty_lemma(); + rational rc_sign = rational(c_sign); + mk_ineq(rc_sign * canonize_sign(c), var(c), llc::LE); + mk_ineq(canonize_sign(ac), var(ac), -canonize_sign(bc), var(bc), ab_cmp); + mk_ineq(canonize_sign(a)*rc_sign, var(a), - canonize_sign(b)*rc_sign, var(b), negate(ab_cmp)); + explain(ac); + explain(a); + explain(bc); + explain(b); + explain(c); + TRACE("nla_solver", _().print_lemma(tout);); +} + +void order::generate_mon_ol(const monomial& ac, + lpvar a, + const rational& c_sign, + lpvar c, + const rooted_mon& bd, + const factor& b, + const rational& d_sign, + lpvar d, + llc ab_cmp) { + add_empty_lemma(); + mk_ineq(c_sign, c, llc::LE); + explain(c); // this explains c == +- d + negate_var_factor_relation(c_sign, a, d_sign, b); + mk_ineq(ac.var(), -canonize_sign(bd), var(bd), ab_cmp); + explain(bd); + explain(b); + explain(d); + TRACE("nla_solver", print_lemma(tout);); +} + +void order::negate_var_factor_relation(const rational& a_sign, lpvar a, const rational& b_sign, const factor& b) { + rational b_fs = canonize_sign(b); + llc cmp = a_sign*vvr(a) < b_sign*vvr(b)? llc::GE : llc::LE; + mk_ineq(a_sign, a, - b_fs*b_sign, var(b), cmp); +} + +void order::order_lemma() { + TRACE("nla_solver", ); + c().init_rm_to_refine(); + const auto& rm_ref = c().m_rm_table.to_refine(); + unsigned start = random() % rm_ref.size(); + unsigned i = start; + do { + const rooted_mon& rm = c().m_rm_table.rms()[rm_ref[i]]; + order_lemma_on_rmonomial(rm); + if (++i == rm_ref.size()) { + i = 0; + } + } while(i != start && !done()); +} + +bool order::order_lemma_on_ac_and_bc_and_factors(const rooted_mon& ac, + const factor& a, + const factor& c, + const rooted_mon& bc, + const factor& b) { + auto cv = vvr(c); + int c_sign = nla::rat_sign(cv); + SASSERT(c_sign != 0); + auto av_c_s = vvr(a)*rational(c_sign); + auto bv_c_s = vvr(b)*rational(c_sign); + auto acv = vvr(ac); + auto bcv = vvr(bc); + TRACE("nla_solver", _().trace_print_ol(ac, a, c, bc, b, tout);); + // Suppose ac >= bc, then ac/|c| >= bc/|c|. + // Notice that ac/|c| = a*c_sign , and bc/|c| = b*c_sign, which are correspondingly av_c_s and bv_c_s + if (acv >= bcv && av_c_s < bv_c_s) { + generate_ol(ac, a, c_sign, c, bc, b, llc::LT); + return true; + } else if (acv <= bcv && av_c_s > bv_c_s) { + generate_ol(ac, a, c_sign, c, bc, b, llc::GT); + return true; + } + return false; +} +/** + \brief Add lemma: + a > 0 & b <= value(b) => sign*ab <= value(b)*a if value(a) > 0 + a < 0 & b >= value(b) => sign*ab <= value(b)*a if value(a) < 0 +*/ +void order::order_lemma_on_ab_gt(const monomial& m, const rational& sign, lpvar a, lpvar b) { + SASSERT(sign * vvr(m) > vvr(a) * vvr(b)); + add_empty_lemma(); + if (vvr(a).is_pos()) { + TRACE("nla_solver", tout << "a is pos\n";); + //negate a > 0 + mk_ineq(a, llc::LE); + // negate b <= vvr(b) + mk_ineq(b, llc::GT, vvr(b)); + // ab <= vvr(b)a + mk_ineq(sign, m.var(), -vvr(b), a, llc::LE); + } else { + TRACE("nla_solver", tout << "a is neg\n";); + SASSERT(vvr(a).is_neg()); + //negate a < 0 + mk_ineq(a, llc::GE); + // negate b >= vvr(b) + mk_ineq(b, llc::LT, vvr(b)); + // ab <= vvr(b)a + mk_ineq(sign, m.var(), -vvr(b), a, llc::LE); + } +} +// we need to deduce ab >= vvr(b)*a +/** + \brief Add lemma: + a > 0 & b >= value(b) => sign*ab >= value(b)*a if value(a) > 0 + a < 0 & b <= value(b) => sign*ab >= value(b)*a if value(a) < 0 +*/ +void order::order_lemma_on_ab_lt(const monomial& m, const rational& sign, lpvar a, lpvar b) { + SASSERT(sign * vvr(m) < vvr(a) * vvr(b)); + add_empty_lemma(); + if (vvr(a).is_pos()) { + //negate a > 0 + mk_ineq(a, llc::LE); + // negate b >= vvr(b) + mk_ineq(b, llc::LT, vvr(b)); + // ab <= vvr(b)a + mk_ineq(sign, m.var(), -vvr(b), a, llc::GE); + } else { + SASSERT(vvr(a).is_neg()); + //negate a < 0 + mk_ineq(a, llc::GE); + // negate b <= vvr(b) + mk_ineq(b, llc::GT, vvr(b)); + // ab >= vvr(b)a + mk_ineq(sign, m.var(), -vvr(b), a, llc::GE); + } +} + +void order::order_lemma_on_ab(const monomial& m, const rational& sign, lpvar a, lpvar b, bool gt) { + if (gt) + order_lemma_on_ab_gt(m, sign, a, b); + else + order_lemma_on_ab_lt(m, sign, a, b); +} +// a > b && c > 0 => ac > bc +void order::order_lemma_on_rmonomial(const rooted_mon& rm) { + TRACE("nla_solver_details", + tout << "rm = "; print_product(rm, tout); + tout << ", orig = "; print_monomial(c().m_monomials[rm.orig_index()], tout); + ); + for (auto ac : factorization_factory_imp(rm, c())) { + if (ac.size() != 2) + continue; + if (ac.is_mon()) + order_lemma_on_binomial(*ac.mon()); + else + order_lemma_on_factorization(rm, ac); + if (done()) + break; + } +} +void order::order_lemma_on_binomial_k(const monomial& m, lpvar k, bool gt) { + SASSERT(gt == (vvr(m) > vvr(m[0]) * vvr(m[1]))); + unsigned p = (k + 1) % 2; + order_lemma_on_binomial_sign(m, m[k], m[p], gt? 1: -1); +} +// sign it the sign of vvr(m) - vvr(m[0]) * vvr(m[1]) +// m = xy +// and val(m) != val(x)*val(y) +// y > 0 and x = a, then xy >= ay +void order::order_lemma_on_binomial_sign(const monomial& ac, lpvar x, lpvar y, int sign) { + SASSERT(!_().mon_has_zero(ac)); + int sy = rat_sign(vvr(y)); + add_empty_lemma(); + mk_ineq(y, sy == 1? llc::LE : llc::GE); // negate sy + mk_ineq(x, sy*sign == 1? llc::GT:llc::LT , vvr(x)); // assert x <= vvr(x) if x > 0 + mk_ineq(ac.var(), - vvr(x), y, sign == 1?llc::LE:llc::GE); + TRACE("nla_solver", print_lemma(tout);); +} +void order::order_lemma_on_factor_binomial_rm(const monomial& ac, unsigned k, const rooted_mon& rm_bd) { + factor d(_().m_evars.find(ac[k]).var(), factor_type::VAR); + factor b; + if (!_().divide(rm_bd, d, b)) + return; + order_lemma_on_binomial_ac_bd(ac, k, rm_bd, b, d.index()); +} + +void order::order_lemma_on_binomial_ac_bd(const monomial& ac, unsigned k, const rooted_mon& bd, const factor& b, lpvar d) { + TRACE("nla_solver", print_monomial(ac, tout << "ac="); + print_rooted_monomial(bd, tout << "\nrm="); + print_factor(b, tout << ", b="); print_var(d, tout << ", d=") << "\n";); + int p = (k + 1) % 2; + lpvar a = ac[p]; + lpvar c = ac[k]; + SASSERT(_().m_evars.find(c).var() == d); + rational acv = vvr(ac); + rational av = vvr(a); + rational c_sign = rrat_sign(vvr(c)); + rational d_sign = rrat_sign(vvr(d)); + rational bdv = vvr(bd); + rational bv = vvr(b); + auto av_c_s = av*c_sign; auto bv_d_s = bv*d_sign; + + // suppose ac >= bd, then ac/|c| >= bd/|d|. + // Notice that ac/|c| = a*c_sign , and bd/|d| = b*d_sign + if (acv >= bdv && av_c_s < bv_d_s) + generate_mon_ol(ac, a, c_sign, c, bd, b, d_sign, d, llc::LT); + else if (acv <= bdv && av_c_s > bv_d_s) + generate_mon_ol(ac, a, c_sign, c, bd, b, d_sign, d, llc::GT); + +} + +void order::order_lemma_on_factor_binomial_explore(const monomial& m, unsigned k) { + SASSERT(m.size() == 2); + lpvar c = m[k]; + lpvar d = _().m_evars.find(c).var(); + auto it = _().m_rm_table.var_map().find(d); + SASSERT(it != _().m_rm_table.var_map().end()); + for (unsigned bd_i : it->second) { + order_lemma_on_factor_binomial_rm(m, k, _().m_rm_table.rms()[bd_i]); + if (done()) + break; + } +} + +void order::order_lemma_on_binomial(const monomial& ac) { + TRACE("nla_solver", print_monomial(ac, tout);); + SASSERT(!check_monomial(ac) && ac.size() == 2); + const rational & mult_val = vvr(ac[0]) * vvr(ac[1]); + const rational acv = vvr(ac); + bool gt = acv > mult_val; + for (unsigned k = 0; k < 2; k++) { + order_lemma_on_binomial_k(ac, k, gt); + order_lemma_on_factor_binomial_explore(ac, k); + } +} +} // end of namespace nla diff --git a/src/util/lp/nla_order_lemmas.h b/src/util/lp/nla_order_lemmas.h new file mode 100644 index 000000000..0b13f392c --- /dev/null +++ b/src/util/lp/nla_order_lemmas.h @@ -0,0 +1,98 @@ +/*++ + Copyright (c) 2017 Microsoft Corporation + + Module Name: + + + + Abstract: + + + + Author: + Nikolaj Bjorner (nbjorner) + Lev Nachmanson (levnach) + + Revision History: + + + --*/ +#pragma once +#include "util/lp/rooted_mons.h" +#include "util/lp/factorization.h" +#include "util/lp/nla_common.h" + +namespace nla { +struct core; +struct order: common { + + // fields + core * m_core; + core& _() { return *m_core; } + const core& _() const { return *m_core; } + core& c() { return *m_core; } + const core& c() const { return *m_core; } + //constructor + order(core *c) : common(c) {} + bool order_lemma_on_ac_and_bc_and_factors(const rooted_mon& ac, + const factor& a, + const factor& c, + const rooted_mon& bc, + const factor& b); + + // a >< b && c > 0 => ac >< bc + // a >< b && c < 0 => ac <> bc + // ac[k] plays the role of c + bool order_lemma_on_ac_and_bc(const rooted_mon& rm_ac, + const factorization& ac_f, + unsigned k, + const rooted_mon& rm_bd); + + bool order_lemma_on_ac_explore(const rooted_mon& rm, const factorization& ac, unsigned k); + + void order_lemma_on_factorization(const rooted_mon& rm, const factorization& ab); + + /** + \brief Add lemma: + a > 0 & b <= value(b) => sign*ab <= value(b)*a if value(a) > 0 + a < 0 & b >= value(b) => sign*ab <= value(b)*a if value(a) < 0 + */ + void order_lemma_on_ab_gt(const monomial& m, const rational& sign, lpvar a, lpvar b); + // we need to deduce ab >= vvr(b)*a + /** + \brief Add lemma: + a > 0 & b >= value(b) => sign*ab >= value(b)*a if value(a) > 0 + a < 0 & b <= value(b) => sign*ab >= value(b)*a if value(a) < 0 + */ + void order_lemma_on_ab_lt(const monomial& m, const rational& sign, lpvar a, lpvar b); + void order_lemma_on_ab(const monomial& m, const rational& sign, lpvar a, lpvar b, bool gt); + void order_lemma_on_factor_binomial_explore(const monomial& m, unsigned k); + void order_lemma_on_factor_binomial_rm(const monomial& ac, unsigned k, const rooted_mon& rm_bd); + void order_lemma_on_binomial_ac_bd(const monomial& ac, unsigned k, const rooted_mon& bd, const factor& b, lpvar d); + void order_lemma_on_binomial_k(const monomial& m, lpvar k, bool gt); + void order_lemma_on_binomial_sign(const monomial& ac, lpvar x, lpvar y, int sign); + void order_lemma_on_binomial(const monomial& ac); + void order_lemma_on_rmonomial(const rooted_mon& rm); + void order_lemma(); + // |c_sign| = 1, and c*c_sign > 0 + // ac > bc => ac/|c| > bc/|c| => a*c_sign > b*c_sign + void generate_ol(const rooted_mon& ac, + const factor& a, + int c_sign, + const factor& c, + const rooted_mon& bc, + const factor& b, + llc ab_cmp); + + void generate_mon_ol(const monomial& ac, + lpvar a, + const rational& c_sign, + lpvar c, + const rooted_mon& bd, + const factor& b, + const rational& d_sign, + lpvar d, + llc ab_cmp); + void negate_var_factor_relation(const rational& a_sign, lpvar a, const rational& b_sign, const factor& b); +}; +} diff --git a/src/util/lp/nla_tangent_lemmas.cpp b/src/util/lp/nla_tangent_lemmas.cpp index 481d9e48e..3dd6e8c65 100644 --- a/src/util/lp/nla_tangent_lemmas.cpp +++ b/src/util/lp/nla_tangent_lemmas.cpp @@ -23,9 +23,8 @@ namespace nla { template rational tangents::vvr(T const& t) const { return m_core->vvr(t); } template lpvar tangents::var(T const& t) const { return m_core->var(t); } -void tangents::add_empty_lemma() { c().add_empty_lemma(); } -tangents::tangents(core * c) : m_core(c) {} +tangents::tangents(core * c) : common(c) {} std::ostream& tangents::print_point(const point &a, std::ostream& out) const { out << "(" << a.x << ", " << a.y << ")"; return out; diff --git a/src/util/lp/nla_tangent_lemmas.h b/src/util/lp/nla_tangent_lemmas.h new file mode 100644 index 000000000..ce3ede154 --- /dev/null +++ b/src/util/lp/nla_tangent_lemmas.h @@ -0,0 +1,84 @@ +/*++ + Copyright (c) 2017 Microsoft Corporation + + Module Name: + + + + Abstract: + + + + Author: + Nikolaj Bjorner (nbjorner) + Lev Nachmanson (levnach) + + Revision History: + + + --*/ +#pragma once +#include "util/rational.h" +#include "util/lp/rooted_mons.h" +#include "util/lp/factorization.h" +#include "util/lp/nla_common.h" + +namespace nla { +struct core; +struct tangents: common { + struct point { + rational x; + rational y; + point(const rational& a, const rational& b): x(a), y(b) {} + point() {} + inline point& operator*=(rational a) { + x *= a; + y *= a; + return *this; + } + inline point operator+(const point& b) const { + return point(x + b.x, y + b.y); + } + + inline point operator-(const point& b) const { + return point(x - b.x, y - b.y); + } + }; + + core* m_core; + core& c() { return *m_core; } + const core& c() const { return *m_core; } + tangents(core *core); + + void generate_simple_tangent_lemma(const rooted_mon* rm); + + void tangent_lemma(); + + void generate_explanations_of_tang_lemma(const rooted_mon& rm, const bfc& bf, lp::explanation& exp); + + void tangent_lemma_bf(const bfc& bf, lpvar j, const rational& sign, const rooted_mon* rm); + void generate_tang_plane(const rational & a, const rational& b, const factor& x, const factor& y, bool below, lpvar j, const rational& j_sign); + + void generate_two_tang_lines(const bfc & bf, const point& xy, const rational& sign, lpvar j); + // Get two planes tangent to surface z = xy, one at point a, and another at point b. + // One can show that these planes still create a cut. + void get_initial_tang_points(point &a, point &b, const point& xy, bool below) const; + + void push_tang_point(point &a, const point& xy, bool below, const rational& correct_val, const rational& val) const; + + void push_tang_points(point &a, point &b, const point& xy, bool below, const rational& correct_val, const rational& val) const; + + rational tang_plane(const point& a, const point& x) const; + void get_tang_points(point &a, point &b, bool below, const rational& val, const point& xy) const; + std::ostream& print_point(const point &a, std::ostream& out) const; + std::ostream& print_tangent_domain(const point &a, const point &b, std::ostream& out) const; + // "below" means that the val is below the surface xy + bool plane_is_correct_cut(const point& plane, + const point& xy, + const rational & correct_val, + const rational & val, + bool below) const; + template rational vvr(T const& t) const; + template lpvar var(T const& t) const; +}; // end of tangents +}