diff --git a/src/math/polysat/forbidden_intervals.cpp b/src/math/polysat/forbidden_intervals.cpp index a06089b07..893cfba33 100644 --- a/src/math/polysat/forbidden_intervals.cpp +++ b/src/math/polysat/forbidden_intervals.cpp @@ -29,9 +29,11 @@ namespace polysat { * \returns True iff a forbidden interval exists and the output parameters were set. */ - bool forbidden_intervals::get_interval(signed_constraint const& c, pvar v, eval_interval& out_interval, vector& out_side_cond) { + bool forbidden_intervals::get_interval(signed_constraint const& c, pvar v, rational & coeff, eval_interval& out_interval, vector& out_side_cond) { if (!c->is_ule()) return false; + + coeff = 1; struct backtrack { bool released = false; @@ -61,21 +63,26 @@ namespace polysat { SASSERT(b1.is_val()); SASSERT(b2.is_val()); + coeff = a1; + _backtrack.released = true; // LOG("add " << c << " " << a1 << " " << b1 << " " << a2 << " " << b2); - if (match_linear1(c, a1, b1, e1, a2, b2, e2, out_interval, out_side_cond)) + if (match_linear1(c, coeff, b1, e1, a2, b2, e2, out_interval, out_side_cond)) return true; - if (match_linear2(c, a1, b1, e1, a2, b2, e2, out_interval, out_side_cond)) + if (match_linear2(c, coeff, b1, e1, a2, b2, e2, out_interval, out_side_cond)) return true; - if (match_linear3(c, a1, b1, e1, a2, b2, e2, out_interval, out_side_cond)) + if (match_linear3(c, coeff, b1, e1, a2, b2, e2, out_interval, out_side_cond)) return true; + + +#if 0 if (match_linear4(c, a1, b1, e1, a2, b2, e2, out_interval, out_side_cond)) return true; if (match_linear5(c, a1, b1, e1, a2, b2, e2, out_interval, out_side_cond)) return true; - +#endif _backtrack.released = false; return false; } @@ -114,7 +121,7 @@ namespace polysat { }; eval_interval forbidden_intervals::to_interval( - signed_constraint const& c, bool is_trivial, rational const& coeff, + signed_constraint const& c, bool is_trivial, rational & coeff, rational & lo_val, pdd & lo, rational & hi_val, pdd & hi) { @@ -131,9 +138,12 @@ namespace polysat { return eval_interval::full(); } - if (!coeff.is_one()) { - rational pow2 = m.max_value() + 1; - SASSERT(coeff == m.max_value()); + rational pow2 = m.max_value() + 1; + + if (coeff > pow2/2) { + + coeff = pow2 - coeff; + SASSERT(coeff > 0); // Transform according to: y \in [l;u[ <=> -y \in [1-u;1-l[ // -y \in [1-u;1-l[ // <=> -y - (1 - u) < (1 - l) - (1 - u) { by: y \in [l;u[ <=> y - l < u - l } @@ -156,14 +166,14 @@ namespace polysat { } /** - * Match e1 + t <= e2, with t = 2^j1*y + * Match e1 + t <= e2, with t = a1*y * condition for empty/full: e2 == -1 */ bool forbidden_intervals::match_linear1(signed_constraint const& c, - rational const& a1, pdd const& b1, pdd const& e1, - rational const& a2, pdd const& b2, pdd const& e2, + rational & a1, pdd const& b1, pdd const& e1, + rational const & a2, pdd const& b2, pdd const& e2, eval_interval& interval, vector& side_cond) { - if (a2.is_zero() && coefficient_is_01(e1.manager(), a1)) { + if (a2.is_zero() && !a1.is_zero()) { SASSERT(!a1.is_zero()); bool is_trivial = (b2 + 1).is_zero(); push_eq(is_trivial, e2 + 1, side_cond); @@ -178,36 +188,37 @@ namespace polysat { } /** - * e1 <= e2 + t, with t = 2^j2*y + * e1 <= e2 + t, with t = a2*y * condition for empty/full: e1 == 0 */ bool forbidden_intervals::match_linear2(signed_constraint const& c, - rational const& a1, pdd const& b1, pdd const& e1, - rational const& a2, pdd const& b2, pdd const& e2, + rational & a1, pdd const& b1, pdd const& e1, + rational const & a2, pdd const& b2, pdd const& e2, eval_interval& interval, vector& side_cond) { - if (a1.is_zero() && coefficient_is_01(e1.manager(), a2)) { + if (a1.is_zero() && !a2.is_zero()) { SASSERT(!a2.is_zero()); + a1 = a2; bool is_trivial = b1.is_zero(); push_eq(is_trivial, e1, side_cond); auto lo = -e2; rational lo_val = (-b2).val(); auto hi = e1 - e2; rational hi_val = (b1 - b2).val(); - interval = to_interval(c, is_trivial, a2, lo_val, lo, hi_val, hi); + interval = to_interval(c, is_trivial, a1, lo_val, lo, hi_val, hi); return true; } return false; } /** - * e1 + t <= e2 + t, with t = 2^j1*y = 2^j2*y - * condition for empty/full: e1 == e2/ + * e1 + t <= e2 + t, with t = a1*y = a2*y + * condition for empty/full: e1 == e2 */ bool forbidden_intervals::match_linear3(signed_constraint const& c, - rational const& a1, pdd const& b1, pdd const& e1, - rational const& a2, pdd const& b2, pdd const& e2, + rational & a1, pdd const& b1, pdd const& e1, + rational const & a2, pdd const& b2, pdd const& e2, eval_interval& interval, vector& side_cond) { - if (coefficient_is_01(e1.manager(), a1) && coefficient_is_01(e1.manager(), a2) && a1 == a2 && !a1.is_zero()) { + if (a1 == a2 && !a1.is_zero()) { bool is_trivial = b1.val() == b2.val(); push_eq(is_trivial, e1 - e2, side_cond); auto lo = -e2; @@ -220,12 +231,13 @@ namespace polysat { return false; } +#if 0 /** * a1*y + e1 = 0, with a1 odd */ bool forbidden_intervals::match_linear4(signed_constraint const& c, - rational const& a1, pdd const& b1, pdd const& e1, - rational const& a2, pdd const& b2, pdd const& e2, + rational & a1, pdd const& b1, pdd const& e1, + rational & a2, pdd const& b2, pdd const& e2, eval_interval& interval, vector& side_cond) { if (a1.is_odd() && a2.is_zero() && b2.val().is_zero()) { push_eq(true, e2, side_cond); @@ -257,8 +269,8 @@ namespace polysat { * - c < ax + b */ bool forbidden_intervals::match_linear5(signed_constraint const& c, - rational const& a1, pdd const& b1, pdd const& e1, - rational const& a2, pdd const& b2, pdd const& e2, + rational & a1, pdd const& b1, pdd const& e1, + rational & a2, pdd const& b2, pdd const& e2, eval_interval& interval, vector& side_cond) { auto& m = e1.manager(); @@ -323,4 +335,5 @@ namespace polysat { } return false; } +#endif } diff --git a/src/math/polysat/forbidden_intervals.h b/src/math/polysat/forbidden_intervals.h index 31ba21bb4..5246ba6c4 100644 --- a/src/math/polysat/forbidden_intervals.h +++ b/src/math/polysat/forbidden_intervals.h @@ -24,40 +24,42 @@ namespace polysat { solver& s; void push_eq(bool is_trivial, pdd const& p, vector& side_cond); - eval_interval to_interval(signed_constraint const& c, bool is_trivial, rational const& coeff, + eval_interval to_interval(signed_constraint const& c, bool is_trivial, rational& coeff, rational & lo_val, pdd & lo, rational & hi_val, pdd & hi); std::tuple linear_decompose(pvar v, pdd const& p, vector& out_side_cond); bool match_linear1(signed_constraint const& c, - rational const& a1, pdd const& b1, pdd const& e1, - rational const& a2, pdd const& b2, pdd const& e2, + rational & a1, pdd const& b1, pdd const& e1, + rational const & a2, pdd const& b2, pdd const& e2, eval_interval& interval, vector& side_cond); bool match_linear2(signed_constraint const& c, - rational const& a1, pdd const& b1, pdd const& e1, - rational const& a2, pdd const& b2, pdd const& e2, + rational & a1, pdd const& b1, pdd const& e1, + rational const & a2, pdd const& b2, pdd const& e2, eval_interval& interval, vector& side_cond); bool match_linear3(signed_constraint const& c, - rational const& a1, pdd const& b1, pdd const& e1, - rational const& a2, pdd const& b2, pdd const& e2, + rational & a1, pdd const& b1, pdd const& e1, + rational const & a2, pdd const& b2, pdd const& e2, eval_interval& interval, vector& side_cond); +#if 0 bool match_linear4(signed_constraint const& c, - rational const& a1, pdd const& b1, pdd const& e1, - rational const& a2, pdd const& b2, pdd const& e2, + rational & a1, pdd const& b1, pdd const& e1, + rational & a2, pdd const& b2, pdd const& e2, eval_interval& interval, vector& side_cond); bool match_linear5(signed_constraint const& c, - rational const& a1, pdd const& b1, pdd const& e1, - rational const& a2, pdd const& b2, pdd const& e2, + rational & a1, pdd const& b1, pdd const& e1, + rational & a2, pdd const& b2, pdd const& e2, eval_interval& interval, vector& side_cond); +#endif - bool coefficient_is_01(dd::pdd_manager& m, rational const& r) { return r.is_zero() || r.is_one() || r == m.max_value(); }; + // bool coefficient_is_01(dd::pdd_manager& m, rational const& r) { return r.is_zero() || r.is_one() || r == m.max_value(); }; public: forbidden_intervals(solver& s) :s(s) {} - bool get_interval(signed_constraint const& c, pvar v, eval_interval& out_interval, vector& side_cond); + bool get_interval(signed_constraint const& c, pvar v, rational & coeff, eval_interval& out_interval, vector& side_cond); }; } diff --git a/src/math/polysat/viable.cpp b/src/math/polysat/viable.cpp index 7d30f8116..362054e59 100644 --- a/src/math/polysat/viable.cpp +++ b/src/math/polysat/viable.cpp @@ -35,50 +35,66 @@ namespace polysat { } viable::entry* viable::alloc_entry() { + rational coeff(1); if (m_alloc.empty()) - return alloc(entry); + return alloc(entry, coeff); auto* e = m_alloc.back(); e->side_cond.reset(); + e->coeff = coeff; m_alloc.pop_back(); return e; } void viable::pop_viable() { - auto& [v, e] = m_trail.back(); - e->remove_from(m_viable[v], e); + auto& [v, is_unit, e] = m_trail.back(); + auto& vec = is_unit ? m_units[v] : m_non_units[v]; + e->remove_from(vec, e); m_alloc.push_back(e); m_trail.pop_back(); } void viable::push_viable() { - auto& [v, e] = m_trail.back(); - SASSERT(e->prev() != e || !m_viable[v]); + auto& [v, is_unit, e] = m_trail.back(); + SASSERT(e->prev() != e || !m_units[v]); SASSERT(e->prev() != e || e->next() == e); + SASSERT(is_unit); + (void)is_unit; if (e->prev() != e) { e->prev()->insert_after(e); if (e->interval.lo_val() < e->next()->interval.lo_val()) - m_viable[v] = e; + m_units[v] = e; } else - m_viable[v] = e; + m_units[v] = e; m_trail.pop_back(); } bool viable::intersect(pvar v, signed_constraint const& c) { auto& fi = s.m_forbidden_intervals; entry* ne = alloc_entry(); - if (!fi.get_interval(c, v, ne->interval, ne->side_cond) || ne->interval.is_currently_empty()) { + if (!fi.get_interval(c, v, ne->coeff, ne->interval, ne->side_cond) || ne->interval.is_currently_empty()) { m_alloc.push_back(ne); return false; } - else { + else if (ne->coeff == 1) { ne->src = c; return intersect(v, ne); } + else { + ne->src = c; + m_trail.push_back({ v, false, ne }); + s.m_trail.push_back(trail_instr_t::viable_add_i); + ne->init(ne); + if (!m_non_units[v]) + m_non_units[v] = ne; + else + ne->insert_after(m_non_units[v]); + return true; + } } bool viable::intersect(pvar v, entry* ne) { - entry* e = m_viable[v]; + entry* e = m_units[v]; if (e && e->interval.is_full()) { m_alloc.push_back(ne); return false; @@ -90,20 +106,20 @@ namespace polysat { } auto create_entry = [&]() { - m_trail.push_back({ v, ne }); + m_trail.push_back({ v, true, ne }); s.m_trail.push_back(trail_instr_t::viable_add_i); ne->init(ne); return ne; }; auto remove_entry = [&](entry* e) { - m_trail.push_back({ v, e }); + m_trail.push_back({ v, true, e }); s.m_trail.push_back(trail_instr_t::viable_rem_i); - e->remove_from(m_viable[v], e); + e->remove_from(m_units[v], e); }; if (!e) - m_viable[v] = create_entry(); + m_units[v] = create_entry(); else { entry* first = e; do { @@ -114,8 +130,8 @@ namespace polysat { while (ne->interval.contains(e->interval)) { entry* n = e->next(); remove_entry(e); - if (!m_viable[v]) { - m_viable[v] = create_entry(); + if (!m_units[v]) { + m_units[v] = create_entry(); return true; } if (e == first) @@ -130,8 +146,8 @@ namespace polysat { } e->insert_before(create_entry()); if (e == first) - m_viable[v] = e->prev(); - SASSERT(well_formed(m_viable[v])); + m_units[v] = e->prev(); + SASSERT(well_formed(m_units[v])); return true; } e = e->next(); @@ -140,32 +156,86 @@ namespace polysat { // otherwise, append to end of list first->insert_before(create_entry()); } - SASSERT(well_formed(m_viable[v])); + SASSERT(well_formed(m_units[v])); return true; } - bool viable::has_viable(pvar v) { - auto* e = m_viable[v]; + /** + * Traverse all interval constraints with coefficients to check whether current value 'val' for + * 'v' is feasible. If not, extract a (maximal) interval to block 'v' from being assigned val. + */ + bool viable::refine_viable(pvar v, rational const& val) { + auto* e = m_non_units[v]; if (!e) return true; entry* first = e; + do { + rational coeff_val = mod(e->coeff * val, s.var2pdd(v).max_value() + 1); + if (e->interval.currently_contains(coeff_val)) { + rational delta_l = floor((coeff_val - e->interval.lo_val()) / e->coeff); + rational delta_u = floor((e->interval.hi_val() - coeff_val - 1) / e->coeff); + rational lo = val - delta_l; + rational hi = val + delta_u + 1; + + if (e->interval.lo_val() < e->interval.hi_val()) { + // pass + } + else if (e->interval.lo_val() <= coeff_val) { + hi = val + 1; + if (hi > s.var2pdd(v).max_value()) + hi = 0; + } + else { + SASSERT(coeff_val < e->interval.hi_val()); + lo = val; + } + SASSERT(hi <= s.var2pdd(v).max_value()); + LOG("forbidden interval [" << lo << ", " << hi << "[\n"); + entry* ne = alloc_entry(); + ne->src = e->src; + ne->side_cond = e->side_cond; + ne->coeff = 1; + pdd lop = s.var2pdd(v).mk_val(lo); // TODO? + pdd hip = s.var2pdd(v).mk_val(hi); + ne->interval = eval_interval::proper(lop, lo, hip, hi); + intersect(v, ne); + return false; + } + e = e->next(); + } + while (e != first); + return true; + } + + bool viable::has_viable(pvar v) { + refined: + auto* e = m_units[v]; + +#define CHECK_RETURN(val) { if (refine_viable(v, val)) return true; else goto refined; } + + if (!e) + CHECK_RETURN(rational::zero()); + entry* first = e; entry* last = e->prev(); // quick check: last interval doesn't wrap around, so hi_val // has not been covered if (last->interval.lo_val() < last->interval.hi_val()) - return true; + CHECK_RETURN(last->interval.hi_val()); do { if (e->interval.is_full()) return false; entry* n = e->next(); if (n == e) - return true; + CHECK_RETURN(e->interval.hi_val()); if (!n->interval.currently_contains(e->interval.hi_val())) - return true; - if (n == first) - return e->interval.lo_val() <= e->interval.hi_val(); + CHECK_RETURN(e->interval.hi_val()); + if (n == first) { + if (e->interval.lo_val() > e->interval.hi_val()) + return false; + CHECK_RETURN(e->interval.hi_val()); + } e = n; } while (e != first); @@ -173,9 +243,9 @@ namespace polysat { } bool viable::is_viable(pvar v, rational const& val) { - auto* e = m_viable[v]; + auto* e = m_units[v]; if (!e) - return true; + return refine_viable(v, val); entry* first = e; entry* last = first->prev(); if (last->interval.currently_contains(val)) @@ -183,15 +253,18 @@ namespace polysat { for (; e != last; e = e->next()) { if (e->interval.currently_contains(val)) return false; - if (val < e->interval.lo_val()) - return true; + if (val < e->interval.lo_val()) + return refine_viable(v, val); } - return true; + return refine_viable(v, val); } rational viable::min_viable(pvar v) { + refined: rational lo(0); - auto* e = m_viable[v]; + auto* e = m_units[v]; + if (!e && !refine_viable(v, lo)) + goto refined; if (!e) return lo; entry* first = e; @@ -205,13 +278,18 @@ namespace polysat { e = e->next(); } while (e != first); - SASSERT(is_viable(v, lo)); + if (!refine_viable(v, lo)) + goto refined; + SASSERT(is_viable(v, lo)); return lo; } rational viable::max_viable(pvar v) { + refined: rational hi = s.var2pdd(v).max_value(); - auto* e = m_viable[v]; + auto* e = m_units[v]; + if (!e && !refine_viable(v, hi)) + goto refined; if (!e) return hi; entry* last = e->prev(); @@ -223,13 +301,20 @@ namespace polysat { e = e->prev(); } while (e != last); + if (!refine_viable(v, hi)) + goto refined; SASSERT(is_viable(v, hi)); return hi; } dd::find_t viable::find_viable(pvar v, rational& lo) { + refined: lo = 0; - auto* e = m_viable[v]; + auto* e = m_units[v]; + if (!e && !refine_viable(v, lo)) + goto refined; + if (!e && !refine_viable(v, rational::one())) + goto refined; if (!e) return dd::find_t::multiple; if (e->interval.is_full()) @@ -244,6 +329,10 @@ namespace polysat { if (last->interval.lo_val() < last->interval.hi_val() && last->interval.hi_val() < max_value) { lo = last->interval.hi_val(); + if (!refine_viable(v, lo)) + goto refined; + if (!refine_viable(v, max_value)) + goto refined; return dd::find_t::multiple; } @@ -271,6 +360,10 @@ namespace polysat { e = e->prev(); } while (e != last); + if (!refine_viable(v, lo)) + goto refined; + if (!refine_viable(v, hi)) + goto refined; if (lo == hi) return dd::find_t::singleton; else @@ -280,7 +373,7 @@ namespace polysat { bool viable::resolve(pvar v, conflict& core) { if (has_viable(v)) return false; - auto* e = m_viable[v]; + auto* e = m_units[v]; entry* first = e; SASSERT(e); core.reset(); @@ -317,9 +410,9 @@ namespace polysat { } void viable::log(pvar v) { - if (!well_formed(m_viable[v])) + if (!well_formed(m_units[v])) LOG("v" << v << " not well formed"); - auto* e = m_viable[v]; + auto* e = m_units[v]; if (!e) return; entry* first = e; @@ -331,16 +424,17 @@ namespace polysat { } void viable::log() { - for (pvar v = 0; v < std::min(10u, m_viable.size()); ++v) + for (pvar v = 0; v < std::min(10u, m_units.size()); ++v) log(v); } - std::ostream& viable::display(std::ostream& out, pvar v) const { - auto* e = m_viable[v]; + std::ostream& viable::display(std::ostream& out, pvar v, entry* e) const { if (!e) return out; entry* first = e; do { + if (e->coeff != 1) + out << e->coeff << " * v" << v << " "; out << e->interval << " " << e->side_cond << " " << e->src << " "; e = e->next(); } @@ -348,8 +442,14 @@ namespace polysat { return out; } + std::ostream& viable::display(std::ostream& out, pvar v) const { + display(out, v, m_units[v]); + display(out, v, m_non_units[v]); + return out; + } + std::ostream& viable::display(std::ostream& out) const { - for (pvar v = 0; v < m_viable.size(); ++v) + for (pvar v = 0; v < m_units.size(); ++v) display(out << "v" << v << ": ", v); return out; } diff --git a/src/math/polysat/viable.h b/src/math/polysat/viable.h index 18b259f14..2822f5806 100644 --- a/src/math/polysat/viable.h +++ b/src/math/polysat/viable.h @@ -32,13 +32,14 @@ namespace polysat { solver& s; struct entry : public dll_base, public fi_record { - public: - entry() : fi_record({ eval_interval::full(), {}, {} }) {} + rational coeff; + entry(rational const& m) : fi_record({ eval_interval::full(), {}, {} }), coeff(m) {} }; ptr_vector m_alloc; - ptr_vector m_viable; // set of viable values. - svector> m_trail; // undo stack + ptr_vector m_units; // set of viable values based on unit multipliers + ptr_vector m_non_units; // entries that have non-unit multipliers + svector> m_trail; // undo stack bool well_formed(entry* e); @@ -46,15 +47,19 @@ namespace polysat { bool intersect(pvar v, entry* e); + bool refine_viable(pvar v, rational const& val); + + std::ostream& display(std::ostream& out, pvar v, entry* e) const; + public: viable(solver& s); ~viable(); // declare and remove var - void push(unsigned) { m_viable.push_back(nullptr); } + void push(unsigned) { m_units.push_back(nullptr); m_non_units.push_back(nullptr); } - void pop() { m_viable.pop_back(); } + void pop() { m_units.pop_back(); m_non_units.pop_back(); } void pop_viable(); @@ -135,8 +140,8 @@ namespace polysat { pvar var; public: constraints(viable& v, pvar var) : v(v), var(var) {} - iterator begin() const { return iterator(v.m_viable[var], false); } - iterator end() const { return iterator(v.m_viable[var], true); } + iterator begin() const { return iterator(v.m_units[var], false); } + iterator end() const { return iterator(v.m_units[var], true); } }; constraints get_constraints(pvar v) {