diff --git a/src/math/lp/lar_constraints.h b/src/math/lp/lar_constraints.h index 19a389f05..158c5a780 100644 --- a/src/math/lp/lar_constraints.h +++ b/src/math/lp/lar_constraints.h @@ -52,14 +52,16 @@ class lar_base_constraint { lconstraint_kind m_kind; mpq m_right_side; bool m_active; + unsigned m_j; public: virtual vector> coeffs() const = 0; - lar_base_constraint(lconstraint_kind kind, const mpq& right_side) :m_kind(kind), m_right_side(right_side), m_active(false) {} + lar_base_constraint(unsigned j, lconstraint_kind kind, const mpq& right_side) :m_kind(kind), m_right_side(right_side), m_active(false), m_j(j) {} virtual ~lar_base_constraint() {} lconstraint_kind kind() const { return m_kind; } mpq const& rhs() const { return m_right_side; } + unsigned column() const { return m_j; } void activate() { m_active = true; } void deactivate() { m_active = false; } @@ -70,14 +72,13 @@ public: }; class lar_var_constraint: public lar_base_constraint { - unsigned m_j; public: lar_var_constraint(unsigned j, lconstraint_kind kind, const mpq& right_side) : - lar_base_constraint(kind, right_side), m_j(j) {} + lar_base_constraint(j, kind, right_side) {} vector> coeffs() const override { vector> ret; - ret.push_back(std::make_pair(one_of_type(), m_j)); + ret.push_back(std::make_pair(one_of_type(), column())); return ret; } unsigned size() const override { return 1;} @@ -87,8 +88,8 @@ public: class lar_term_constraint: public lar_base_constraint { const lar_term * m_term; public: - lar_term_constraint(const lar_term *t, lconstraint_kind kind, const mpq& right_side) : - lar_base_constraint(kind, right_side), m_term(t) {} + lar_term_constraint(unsigned j, const lar_term *t, lconstraint_kind kind, const mpq& right_side) : + lar_base_constraint(j, kind, right_side), m_term(t) {} vector> coeffs() const override { return m_term->coeffs_as_vector(); } unsigned size() const override { return m_term->size();} @@ -150,7 +151,7 @@ public: m_constraint_count = m_constraints.size(); m_constraint_count.push(); m_region.push_scope(); -#if 0 +#if 1 m_active_lim = m_active.size(); m_active_lim.push(); #endif @@ -161,7 +162,7 @@ public: for (unsigned i = m_constraints.size(); i-- > m_constraint_count; ) m_constraints[i]->~lar_base_constraint(); m_constraints.shrink(m_constraint_count); -#if 0 +#if 1 m_active_lim.pop(k); for (unsigned i = m_active.size(); i-- > m_active_lim; ) { m_constraints[m_active[i]]->deactivate(); @@ -175,11 +176,11 @@ public: return add(new (m_region) lar_var_constraint(j, k, rhs)); } - constraint_index add_term_constraint(const lar_term* t, lconstraint_kind k, mpq const& rhs) { - return add(new (m_region) lar_term_constraint(t, k, rhs)); + constraint_index add_term_constraint(unsigned j, const lar_term* t, lconstraint_kind k, mpq const& rhs) { + return add(new (m_region) lar_term_constraint(j, t, k, rhs)); } -#if 1 +#if 0 bool is_active(constraint_index ci) const { return true; } void activate(constraint_index ci) {} @@ -188,7 +189,7 @@ public: // future behavior uses activation bit. bool is_active(constraint_index ci) const { return m_constraints[ci]->is_active(); } - void activate(constraint_index ci) { m_constraints[ci]->activate(); } + void activate(constraint_index ci) { auto& c = *m_constraints[ci]; if (!c.is_active()) { c.activate(); m_active.push_back(ci); } } #endif lar_base_constraint const& operator[](constraint_index ci) const { return *m_constraints[ci]; } diff --git a/src/math/lp/lar_solver.cpp b/src/math/lp/lar_solver.cpp index 9505ff2a0..f7e98028d 100644 --- a/src/math/lp/lar_solver.cpp +++ b/src/math/lp/lar_solver.cpp @@ -1764,12 +1764,22 @@ bool lar_solver::bound_is_integer_for_integer_column(unsigned j, const mpq & rig } constraint_index lar_solver::add_var_bound(var_index j, lconstraint_kind kind, const mpq & right_side) { + constraint_index ci = mk_var_bound(j, kind, right_side); + activate(ci); + return ci; +} + +void lar_solver::activate(constraint_index ci) { + auto const& c = m_constraints[ci]; + update_column_type_and_bound(c.column(), c.kind(), c.rhs(), ci); +} + +constraint_index lar_solver::mk_var_bound(var_index j, lconstraint_kind kind, const mpq & right_side) { TRACE("lar_solver", tout << "j = " << j << " " << lconstraint_kind_string(kind) << " " << right_side<< std::endl;); constraint_index ci; if (!is_term(j)) { // j is a var lp_assert(bound_is_integer_for_integer_column(j, right_side)); ci = m_constraints.add_var_constraint(j, kind, right_side); - update_column_type_and_bound(j, kind, right_side, ci); } else { ci = add_var_bound_on_constraint_for_term(j, kind, right_side); @@ -1811,26 +1821,20 @@ constraint_index lar_solver::add_var_bound_on_constraint_for_term(var_index j, l unsigned adjusted_term_index = adjust_term_index(j); // lp_assert(!term_is_int(m_terms[adjusted_term_index]) || right_side.is_int()); unsigned term_j; - constraint_index ci; + lar_term const* term = m_terms[adjusted_term_index]; if (m_var_register.external_is_used(j, term_j)) { - ci = m_constraints.add_term_constraint(m_terms[adjusted_term_index], kind, right_side); - update_column_type_and_bound(term_j, kind, right_side, ci); + return m_constraints.add_term_constraint(term_j, term, kind, right_side); } else { - ci = add_constraint_from_term_and_create_new_column_row(j, m_terms[adjusted_term_index], kind, right_side); + return add_constraint_from_term_and_create_new_column_row(j, term, kind, right_side); } - return ci; } -constraint_index lar_solver::add_constraint_from_term_and_create_new_column_row(unsigned term_j, const lar_term* term, - lconstraint_kind kind, const mpq & right_side) { - +constraint_index lar_solver::add_constraint_from_term_and_create_new_column_row( + unsigned term_j, const lar_term* term, lconstraint_kind kind, const mpq & right_side) { add_row_from_term_no_constraint(term, term_j); unsigned j = A_r().column_count() - 1; - constraint_index ci = m_constraints.add_term_constraint(term, kind, right_side); - update_column_type_and_bound(j, kind, right_side, ci); - lp_assert(A_r().column_count() == m_mpq_lar_core_solver.m_r_solver.m_costs.size()); - return ci; + return m_constraints.add_term_constraint(j, term, kind, right_side); } void lar_solver::decide_on_strategy_and_adjust_initial_state() { diff --git a/src/math/lp/lar_solver.h b/src/math/lp/lar_solver.h index b2ea8288e..707225384 100644 --- a/src/math/lp/lar_solver.h +++ b/src/math/lp/lar_solver.h @@ -183,6 +183,10 @@ public: void add_basic_var_to_core_fields(); + constraint_index mk_var_bound(var_index j, lconstraint_kind kind, const mpq & right_side); + + void activate(constraint_index ci); + constraint_index add_var_bound(var_index j, lconstraint_kind kind, const mpq & right_side) ; bool compare_values(var_index j, lconstraint_kind kind, const mpq & right_side); diff --git a/src/math/lp/stacked_vector.h b/src/math/lp/stacked_vector.h index 90131fbeb..c38433016 100644 --- a/src/math/lp/stacked_vector.h +++ b/src/math/lp/stacked_vector.h @@ -19,22 +19,25 @@ Revision History: --*/ #pragma once -#include -#include -#include #include "util/vector.h" namespace lp { template < typename B> class stacked_vector { + struct log_entry { + unsigned m_i; unsigned m_ts; B b; + log_entry(unsigned i, unsigned t, B const& b): m_i(i), m_ts(t), b(b) {} + log_entry():m_i(UINT_MAX), m_ts(0) {} + }; svector m_stack_of_vector_sizes; svector m_stack_of_change_sizes; - vector> m_changes; - vector m_vector; + vector m_log; + vector m_vector; + svector m_last_update; public: class ref { stacked_vector & m_vec; unsigned m_i; public: - ref(stacked_vector &m, unsigned key) :m_vec(m), m_i(key) { + ref(stacked_vector &m, unsigned key): m_vec(m), m_i(key) { lp_assert(key < m.size()); } ref & operator=(const B & b) { @@ -77,14 +80,24 @@ public: }; private: + + unsigned now() const { return m_stack_of_change_sizes.size(); } + void emplace_replace(unsigned i,const B & b) { - if (m_vector[i] != b) { - m_changes.push_back(std::make_pair(i, m_vector[i])); + unsigned n = now(); + if (m_last_update[i] == n) { m_vector[i] = b; } + else if (m_vector[i] != b) { + m_log.push_back(log_entry(i, m_last_update[i], m_vector[i])); + m_vector[i] = b; + m_last_update[i] = n; + } } public: + stacked_vector() { } + ref operator[] (unsigned a) { return ref(*this, a); } @@ -102,11 +115,10 @@ public: unsigned size() const { return m_vector.size(); } - void push() { - m_stack_of_change_sizes.push_back(m_changes.size()); - m_stack_of_vector_sizes.push_back(m_vector.size()); + m_stack_of_change_sizes.push_back(m_log.size()); + m_stack_of_vector_sizes.push_back(m_vector.size()); } void pop() { @@ -133,62 +145,31 @@ public: void pop(unsigned k) { lp_assert(m_stack_of_vector_sizes.size() >= k); lp_assert(k > 0); - resize(m_vector, m_stack_of_vector_sizes[m_stack_of_vector_sizes.size() - k]); + m_vector.resize(m_stack_of_vector_sizes[m_stack_of_vector_sizes.size() - k]); + m_last_update.resize(m_stack_of_vector_sizes[m_stack_of_vector_sizes.size() - k]); pop_tail(m_stack_of_vector_sizes, k); unsigned first_change = m_stack_of_change_sizes[m_stack_of_change_sizes.size() - k]; pop_tail(m_stack_of_change_sizes, k); - for (unsigned j = m_changes.size(); j-- > first_change; ) { - const auto & p = m_changes[j]; - unsigned jc = p.first; + for (unsigned j = m_log.size(); j-- > first_change; ) { + const auto & p = m_log[j]; + unsigned jc = p.m_i; if (jc < m_vector.size()) { - m_vector[jc] = p.second; // restore the old value - m_vector[jc] = p.second; // restore the old value - } + m_vector[jc] = p.b; // restore the old value + m_last_update[jc] = p.m_ts; + } } - resize(m_changes, first_change); + resize(m_log, first_change); - /* - while (k-- > 0) { - - if (m_stack.empty()) - return; - - delta & d = m_stack.back(); - lp_assert(m_vector.size() >= d.m_size); - while (m_vector.size() > d.m_size) - m_vector.pop_back(); - - for (auto & t : d.m_original_changed) { - lp_assert(t.first < m_vector.size()); - m_vector[t.first] = t.second; - } - // lp_assert(d.m_deb_copy == m_vector); - m_stack.pop_back();*/ } - - - // void clear() { - // if (m_stack.empty()) { - // m_vector.clear(); - // return; - // } - - // delta & d = m_stack.top(); - // auto & oc = d.m_original_changed; - // for (auto & p : m_vector) { - // const auto & it = oc.find(p.first); - // if (it == oc.end() && d.m_new.find(p.first) == d.m_new.end()) - // oc.emplace(p.first, p.second); - // } - // m_vector.clear(); - // } - + void push_back(const B & b) { m_vector.push_back(b); + m_last_update.push_back(now()); } void increase_size_by_one() { m_vector.resize(m_vector.size() + 1); + m_last_update.push_back(now()); } unsigned peek_size(unsigned k) const { diff --git a/src/smt/theory_lra.cpp b/src/smt/theory_lra.cpp index fba955473..04afd7152 100644 --- a/src/smt/theory_lra.cpp +++ b/src/smt/theory_lra.cpp @@ -65,24 +65,31 @@ std::ostream& operator<<(std::ostream& out, bound_kind const& k) { class bound { smt::bool_var m_bv; smt::theory_var m_var; + lpvar m_vi; bool m_is_int; rational m_value; bound_kind m_bound_kind; + lp::constraint_index m_constraints[2]; public: - bound(smt::bool_var bv, smt::theory_var v, bool is_int, rational const & val, bound_kind k): + bound(smt::bool_var bv, smt::theory_var v, lpvar vi, bool is_int, rational const & val, bound_kind k, lp::constraint_index ct, lp::constraint_index cf): m_bv(bv), m_var(v), + m_vi(vi), m_is_int(is_int), m_value(val), m_bound_kind(k) { + m_constraints[0] = cf; + m_constraints[1] = ct; } virtual ~bound() {} smt::theory_var get_var() const { return m_var; } + lpvar lp_solver_var() const { return m_vi; } smt::bool_var get_bv() const { return m_bv; } bound_kind get_bound_kind() const { return m_bound_kind; } bool is_int() const { return m_is_int; } rational const& get_value() const { return m_value; } + lp::constraint_index get_constraint(bool b) const { return m_constraints[b]; } inf_rational get_value(bool is_true) const { if (is_true) return inf_rational(m_value); // v >= value or v <= value if (m_is_int) { @@ -1050,7 +1057,7 @@ public: if (is_int(v) && !r.is_int()) { r = (k == lp_api::upper_t) ? floor(r) : ceil(r); } - lp_api::bound* b = alloc(lp_api::bound, bv, v, is_int(v), r, k); + lp_api::bound* b = mk_var_bound(bv, v, k, r); m_bounds[v].push_back(b); updt_unassigned_bounds(v, +1); m_bounds_trail.push_back(v); @@ -2587,7 +2594,7 @@ public: rational const& k2 = b2.get_value(); lp_api::bound_kind kind1 = b1.get_bound_kind(); lp_api::bound_kind kind2 = b2.get_bound_kind(); - bool v_is_int = is_int(v); + bool v_is_int = b1.is_int(); SASSERT(v == b2.get_var()); if (k1 == k2 && kind1 == kind2) return; SASSERT(k1 != k2 || kind1 != kind2); @@ -2794,7 +2801,7 @@ public: SASSERT(!bounds.empty()); if (bounds.size() == 1) return; if (m_unassigned_bounds[v] == 0) return; - bool v_is_int = is_int(v); + bool v_is_int = b.is_int(); literal lit1(bv, !is_true); literal lit2 = null_literal; bool find_glb = (is_true == (k == lp_api::lower_t)); @@ -3007,47 +3014,59 @@ public: return true; } - void assert_bound(bool_var bv, bool is_true, lp_api::bound& b) { - if (is_infeasible()) return; - scoped_internalize_state st(*this); - st.vars().push_back(b.get_var()); - st.coeffs().push_back(rational::one()); - init_left_side(st); - lp::lconstraint_kind k = lp::EQ; - bool is_int = b.is_int(); - switch (b.get_bound_kind()) { + lp::lconstraint_kind bound2constraint_kind(bool is_int, lp_api::bound_kind bk, bool is_true) { + switch (bk) { case lp_api::lower_t: - k = is_true ? lp::GE : (is_int ? lp::LE : lp::LT); - break; + return is_true ? lp::GE : (is_int ? lp::LE : lp::LT); case lp_api::upper_t: - k = is_true ? lp::LE : (is_int ? lp::GE : lp::GT); - break; - } + return is_true ? lp::LE : (is_int ? lp::GE : lp::GT); + } + UNREACHABLE(); + return lp::EQ; + } + + void assert_bound(bool_var bv, bool is_true, lp_api::bound& b) { + lp::constraint_index ci = b.get_constraint(is_true); + m_solver->activate(ci); + if (is_infeasible()) { + return; + } + lp::lconstraint_kind k = bound2constraint_kind(b.is_int(), b.get_bound_kind(), is_true); if (k == lp::LT || k == lp::LE) { ++m_stats.m_assert_lower; } else { ++m_stats.m_assert_upper; } - auto vi = register_theory_var_in_lar_solver(b.get_var()); - rational bound = b.get_value(); - lp::constraint_index ci; - TRACE("arith", tout << "v" << b.get_var() << ", vi = " << vi;); - if (is_int && !is_true) { - rational bound = b.get_value(false).get_rational(); - ci = m_solver->add_var_bound(vi, k, bound); - TRACE("arith", tout << "\bbound = " << bound << ", ci = " << ci << "\n";); + propagate_eqs(b.lp_solver_var(), ci, k, b); + } + + lp_api::bound* mk_var_bound(bool_var bv, theory_var v, lp_api::bound_kind bk, rational const& bound) { + scoped_internalize_state st(*this); + st.vars().push_back(v); + st.coeffs().push_back(rational::one()); + init_left_side(st); + lp::constraint_index cT, cF; + bool v_is_int = is_int(v); + auto vi = register_theory_var_in_lar_solver(v); + + lp::lconstraint_kind kT = bound2constraint_kind(v_is_int, bk, true); + lp::lconstraint_kind kF = bound2constraint_kind(v_is_int, bk, false); + + cT = m_solver->mk_var_bound(vi, kT, bound); + if (v_is_int) { + rational boundF = (bk == lp_api::lower_t) ? bound - 1 : bound + 1; + cF = m_solver->mk_var_bound(vi, kF, boundF); } else { - ci = m_solver->add_var_bound(vi, k, b.get_value()); - TRACE("arith", tout << "\nbound = " << bound << ", ci = " << ci << "\n";); + cF = m_solver->mk_var_bound(vi, kF, bound); } - add_ineq_constraint(ci, literal(bv, !is_true)); - if (is_infeasible()) { - return; - } - propagate_eqs(vi, ci, k, b); + add_ineq_constraint(cT, literal(bv, false)); + add_ineq_constraint(cF, literal(bv, true)); + + return alloc(lp_api::bound, bv, v, vi, v_is_int, bound, bk, cT, cF); } + // // fixed equalities. // A fixed equality is inferred if there are two variables v1, v2 whose @@ -3761,7 +3780,7 @@ public: // ctx().set_enode_flag(bv, true); lp_api::bound_kind bkind = lp_api::bound_kind::lower_t; if (is_strict) bkind = lp_api::bound_kind::upper_t; - lp_api::bound* a = alloc(lp_api::bound, bv, v, is_int, r, bkind); + lp_api::bound* a = mk_var_bound(bv, v, bkind, r); mk_bound_axioms(*a); updt_unassigned_bounds(v, +1); m_bounds[v].push_back(a);