diff --git a/src/api/z3_api.h b/src/api/z3_api.h index 114490015..5d1d45f37 100644 --- a/src/api/z3_api.h +++ b/src/api/z3_api.h @@ -3068,7 +3068,8 @@ extern "C" { \brief Create a numeral of a given sort. \param c logical context. - \param numeral A string representing the numeral value in decimal notation. If the given sort is a real, then the numeral can be a rational, that is, a string of the form \ccode{[num]* / [num]*}. + \param numeral A string representing the numeral value in decimal notation. The string may be of the form \code{[num]*[.[num]*][E[+|-][num]+]}. + If the given sort is a real, then the numeral can be a rational, that is, a string of the form \ccode{[num]* / [num]*}. \param ty The sort of the numeral. In the current implementation, the given sort can be an int, real, finite-domain, or bit-vectors of arbitrary size. \sa Z3_mk_int diff --git a/src/ast/arith_decl_plugin.h b/src/ast/arith_decl_plugin.h index b72a55e34..7a6281b06 100644 --- a/src/ast/arith_decl_plugin.h +++ b/src/ast/arith_decl_plugin.h @@ -422,5 +422,103 @@ public: expr_ref mk_add_simplify(unsigned sz, expr* const* args); }; + +inline app_ref mk_numeral(rational const& r, app_ref const& x) { + arith_util a(x.get_manager()); + return app_ref(a.mk_numeral(r, r.is_int() && a.is_int(x)), x.get_manager()); +} + +inline app_ref operator+(app_ref const& x, app_ref const& y) { + arith_util a(x.get_manager()); + return app_ref(a.mk_add(x, y), x.get_manager()); +} + +inline app_ref operator+(app_ref const& x, rational const& y) { + return x + mk_numeral(y, x); +} + +inline app_ref operator+(app_ref const& x, int y) { + return x + rational(y); +} + +inline app_ref operator+(rational const& x, app_ref const& y) { + return mk_numeral(x, y) + y; +} + +inline app_ref operator+(int x, app_ref const& y) { + return rational(x) + y; +} + +inline app_ref operator-(app_ref const& x, app_ref const& y) { + arith_util a(x.get_manager()); + return app_ref(a.mk_sub(x, y), x.get_manager()); +} + +inline app_ref operator-(app_ref const& x, rational const& y) { + return x - mk_numeral(y, x); +} + +inline app_ref operator-(app_ref const& x, int y) { + return x - rational(y); +} + +inline app_ref operator-(rational const& x, app_ref const& y) { + return mk_numeral(x, y) - y; +} + +inline app_ref operator-(int x, app_ref const& y) { + return rational(x) - y; +} + + +inline app_ref operator*(app_ref const& x, app_ref const& y) { + arith_util a(x.get_manager()); + return app_ref(a.mk_mul(x, y), x.get_manager()); +} + +inline app_ref operator*(app_ref const& x, rational const& y) { + return x * mk_numeral(y, x); +} + +inline app_ref operator*(rational const& x, app_ref const& y) { + return mk_numeral(x, y) * y; +} + +inline app_ref operator*(app_ref const& x, int y) { + return x * rational(y); +} + +inline app_ref operator*(int x, app_ref const& y) { + return rational(x) * y; +} + +inline app_ref operator<=(app_ref const& x, app_ref const& y) { + arith_util a(x.get_manager()); + return app_ref(a.mk_le(x, y), x.get_manager()); +} + +inline app_ref operator<=(app_ref const& x, rational const& y) { + return x <= mk_numeral(y, x); +} + +inline app_ref operator<=(app_ref const& x, int y) { + return x <= rational(y); +} + +inline app_ref operator>=(app_ref const& x, app_ref const& y) { + arith_util a(x.get_manager()); + return app_ref(a.mk_ge(x, y), x.get_manager()); +} + +inline app_ref operator<(app_ref const& x, app_ref const& y) { + arith_util a(x.get_manager()); + return app_ref(a.mk_lt(x, y), x.get_manager()); +} + +inline app_ref operator>(app_ref const& x, app_ref const& y) { + arith_util a(x.get_manager()); + return app_ref(a.mk_gt(x, y), x.get_manager()); +} + #endif /* ARITH_DECL_PLUGIN_H_ */ diff --git a/src/math/simplex/model_based_opt.cpp b/src/math/simplex/model_based_opt.cpp index c4e4b110e..760c9a941 100644 --- a/src/math/simplex/model_based_opt.cpp +++ b/src/math/simplex/model_based_opt.cpp @@ -7,7 +7,7 @@ Module Name: Abstract: - Model-based optimization for linear real arithmetic. + Model-based optimization and projection for linear real, integer arithmetic. Author: @@ -26,6 +26,7 @@ std::ostream& operator<<(std::ostream& out, opt::ineq_type ie) { case opt::t_eq: return out << " = "; case opt::t_lt: return out << " < "; case opt::t_le: return out << " <= "; + case opt::t_mod: return out << " mod "; } return out; } @@ -55,8 +56,8 @@ namespace opt { vector const& vars = r.m_vars; for (unsigned i = 0; i < vars.size(); ++i) { // variables in each row are sorted and have non-zero coefficients - SASSERT(i + 1 == vars.size() || vars[i].m_id < vars[i+1].m_id); - SASSERT(!vars[i].m_coeff.is_zero()); + PASSERT(i + 1 == vars.size() || vars[i].m_id < vars[i+1].m_id); + PASSERT(!vars[i].m_coeff.is_zero()); } PASSERT(r.m_value == get_row_value(r)); @@ -64,6 +65,7 @@ namespace opt { // values satisfy constraints PASSERT(index == 0 || r.m_type != t_lt || r.m_value.is_neg()); PASSERT(index == 0 || r.m_type != t_le || !r.m_value.is_pos()); + PASSERT(index == 0 || r.m_type != t_mod || (mod(r.m_value, r.m_mod).is_zero())); return true; } @@ -117,7 +119,7 @@ namespace opt { // objective + t2*coeff/a2 <= ub mul_add(false, m_objective_id, - coeff/bound_coeff, bound_row_index); - m_rows[bound_row_index].m_alive = false; + retire_row(bound_row_index); bound_trail.push_back(bound_row_index); bound_vars.push_back(x); } @@ -281,6 +283,12 @@ namespace opt { } return bound_row_index != UINT_MAX; } + + void model_based_opt::retire_row(unsigned row_id) { + m_rows[row_id].m_alive = false; + m_retired_rows.push_back(row_id); + } + rational model_based_opt::get_row_value(row const& r) const { vector const& vars = r.m_vars; @@ -308,7 +316,7 @@ namespace opt { } if (id < var_id) { lo = mid + 1; - } + } else { hi = mid; } @@ -355,10 +363,161 @@ namespace opt { if (m_rows[row_dst].m_alive) { rational a2 = get_coefficient(row_dst, x); - mul_add(row_dst != m_objective_id && a1.is_pos() == a2.is_pos(), row_dst, -a2/a1, row_src); + if (is_int(x)) { + TRACE("opt", + tout << a1 << " " << a2 << ": "; + display(tout, m_rows[row_dst]); + display(tout, m_rows[row_src]);); + if (a1.is_pos() != a2.is_pos()) { + mul_add(x, a1, row_src, a2, row_dst); + } + else { + mul(row_dst, abs(a1)); + mul_add(false, row_dst, -abs(a2), row_src); + } + TRACE("opt", display(tout, m_rows[row_dst]);); + normalize(row_dst); + } + else { + mul_add(row_dst != m_objective_id && a1.is_pos() == a2.is_pos(), row_dst, -a2/a1, row_src); + } } } - + + // resolution for integer rows. + void model_based_opt::mul_add( + unsigned x, rational const& src_c, unsigned row_src, rational const& dst_c, unsigned row_dst) { + row& dst = m_rows[row_dst]; + row const& src = m_rows[row_src]; + SASSERT(is_int(x)); + SASSERT(t_le == dst.m_type && t_le == src.m_type); + SASSERT(src_c.is_int()); + SASSERT(dst_c.is_int()); + + rational abs_src_c = abs(src_c); + rational abs_dst_c = abs(dst_c); + rational x_val = m_var2value[x]; + rational slack = (abs_src_c - rational::one()) * (abs_dst_c - rational::one()); + rational dst_val = dst.m_value - x_val*dst_c; + rational src_val = src.m_value - x_val*src_c; + bool use_case1 = + (src_c * dst_val + dst_c * src_val + slack).is_nonpos() + || abs_src_c.is_one() + || abs_dst_c.is_one(); + + if (use_case1) { + // dst <- abs_src_c*dst + abs_dst_c*src - slack + mul(row_dst, abs_src_c); + sub(row_dst, slack); + mul_add(false, row_dst, abs_dst_c, row_src); + return; + } + + // + // create finite disjunction for |b|. + // exists x, z in [0 .. |b|-2] . b*x + s + z = 0 && ax + t <= 0 && bx + s <= 0 + // <=> + // exists x, z in [0 .. |b|-2] . b*x = -z - s && ax + t <= 0 && bx + s <= 0 + // <=> + // exists x, z in [0 .. |b|-2] . b*x = -z - s && a|b|x + |b|t <= 0 && bx + s <= 0 + // <=> + // exists x, z in [0 .. |b|-2] . b*x = -z - s && a|b|x + |b|t <= 0 && -z - s + s <= 0 + // <=> + // exists x, z in [0 .. |b|-2] . b*x = -z - s && a|b|x + |b|t <= 0 && -z <= 0 + // <=> + // exists x, z in [0 .. |b|-2] . b*x = -z - s && a|b|x + |b|t <= 0 + // <=> + // exists x, z in [0 .. |b|-2] . b*x = -z - s && a*n_sign(b)(s + z) + |b|t <= 0 + // <=> + // exists z in [0 .. |b|-2] . |b| | (z + s) && a*n_sign(b)(s + z) + |b|t <= 0 + // + + vector coeffs; + if (abs_dst_c <= abs_src_c) { + rational z = mod(dst_val, abs_dst_c); + if (!z.is_zero()) z = abs_dst_c - z; + mk_coeffs_without(coeffs, dst.m_vars, x); + add_divides(coeffs, dst.m_coeff + z, abs_dst_c); + add(row_dst, z); + mul(row_dst, src_c * n_sign(dst_c)); + mul_add(false, row_dst, abs_dst_c, row_src); + } + else { + // z := b - (s + bx) mod b + // := b - s mod b + // b | s + z <=> b | s + b - s mod b <=> b | s - s mod b + rational z = mod(src_val, abs_src_c); + if (!z.is_zero()) z = abs_src_c - z; + mk_coeffs_without(coeffs, src.m_vars, x); + add_divides(coeffs, src.m_coeff + z, abs_src_c); + mul(row_dst, abs_src_c); + add(row_dst, z * dst_c * n_sign(src_c)); + mul_add(false, row_dst, dst_c * n_sign(src_c), row_src); + } + } + + void model_based_opt::mk_coeffs_without(vector& dst, vector const src, unsigned x) { + for (unsigned i = 0; i < src.size(); ++i) { + if (src[i].m_id != x) dst.push_back(src[i]); + } + } + + + rational model_based_opt::n_sign(rational const& b) const { + return rational(b.is_pos()?-1:1); + } + + void model_based_opt::mul(unsigned dst, rational const& c) { + if (c.is_one()) return; + row& r = m_rows[dst]; + for (unsigned i = 0; i < r.m_vars.size(); ++i) { + r.m_vars[i].m_coeff *= c; + } + r.m_coeff *= c; + r.m_value *= c; + } + + void model_based_opt::add(unsigned dst, rational const& c) { + row& r = m_rows[dst]; + r.m_coeff += c; + r.m_value += c; + } + + void model_based_opt::sub(unsigned dst, rational const& c) { + row& r = m_rows[dst]; + r.m_coeff -= c; + r.m_value -= c; + } + + void model_based_opt::normalize(unsigned row_id) { + row& r = m_rows[row_id]; + if (r.m_vars.empty()) return; + if (r.m_type == t_mod) return; + rational g(abs(r.m_vars[0].m_coeff)); + bool all_int = g.is_int(); + for (unsigned i = 1; all_int && !g.is_one() && i < r.m_vars.size(); ++i) { + rational const& coeff = r.m_vars[i].m_coeff; + if (coeff.is_int()) { + g = gcd(g, abs(coeff)); + } + else { + all_int = false; + } + } + if (all_int && !r.m_coeff.is_zero()) { + if (r.m_coeff.is_int()) { + g = gcd(g, abs(r.m_coeff)); + } + else { + all_int = false; + } + } + if (all_int && !g.is_one()) { + SASSERT(!g.is_zero()); + mul(row_id, rational::one()/g); + } + } + // // set row1 <- row1 + c*row2 // @@ -452,12 +611,18 @@ namespace opt { else if (r.m_coeff.is_neg()) { out << r.m_coeff << " "; } - out << r.m_type << " 0; value: " << r.m_value << "\n"; + if (r.m_type == opt::t_mod) { + out << r.m_type << " " << r.m_mod << " = 0; value: " << r.m_value << "\n"; + } + else { + out << r.m_type << " 0; value: " << r.m_value << "\n"; + } } - unsigned model_based_opt::add_var(rational const& value) { + unsigned model_based_opt::add_var(rational const& value, bool is_int) { unsigned v = m_var2value.size(); m_var2value.push_back(value); + m_var2is_int.push_back(is_int); m_var2row_ids.push_back(unsigned_vector()); return v; } @@ -466,34 +631,65 @@ namespace opt { return m_var2value[var]; } - void model_based_opt::set_row(unsigned row_id, vector const& coeffs, rational const& c, ineq_type rel) { + void model_based_opt::set_row(unsigned row_id, vector const& coeffs, rational const& c, rational const& m, ineq_type rel) { row& r = m_rows[row_id]; rational val(c); SASSERT(r.m_vars.empty()); r.m_vars.append(coeffs.size(), coeffs.c_ptr()); + bool is_int_row = true; std::sort(r.m_vars.begin(), r.m_vars.end(), var::compare()); for (unsigned i = 0; i < coeffs.size(); ++i) { val += m_var2value[coeffs[i].m_id] * coeffs[i].m_coeff; + SASSERT(!is_int(coeffs[i].m_id) || coeffs[i].m_coeff.is_int()); + is_int_row &= is_int(coeffs[i].m_id); } r.m_alive = true; r.m_coeff = c; r.m_value = val; r.m_type = rel; + r.m_mod = m; + if (is_int_row && rel == t_lt) { + r.m_type = t_le; + r.m_coeff += rational::one(); + r.m_value += rational::one(); + } SASSERT(invariant(row_id, r)); } + unsigned model_based_opt::new_row() { + unsigned row_id = 0; + if (m_retired_rows.empty()) { + row_id = m_rows.size(); + m_rows.push_back(row()); + } + else { + row_id = m_retired_rows.back(); + m_retired_rows.pop_back(); + m_rows[row_id].reset(); + m_rows[row_id].m_alive = true; + } + return row_id; + } + void model_based_opt::add_constraint(vector const& coeffs, rational const& c, ineq_type rel) { - rational val(c); - unsigned row_id = m_rows.size(); - m_rows.push_back(row()); - set_row(row_id, coeffs, c, rel); + add_constraint(coeffs, c, rational::zero(), rel); + } + + void model_based_opt::add_divides(vector const& coeffs, rational const& c, rational const& m) { + add_constraint(coeffs, c, m, t_mod); + } + + void model_based_opt::add_constraint(vector const& coeffs, rational const& c, rational const& m, ineq_type rel) { + unsigned row_id = new_row(); + set_row(row_id, coeffs, c, m, rel); for (unsigned i = 0; i < coeffs.size(); ++i) { m_var2row_ids[coeffs[i].m_id].push_back(row_id); } + SASSERT(invariant(row_id, m_rows[row_id])); } void model_based_opt::set_objective(vector const& coeffs, rational const& c) { - set_row(m_objective_id, coeffs, c, t_le); + set_row(m_objective_id, coeffs, c, rational::zero(), t_le); } void model_based_opt::get_live_rows(vector& rows) { @@ -524,7 +720,8 @@ namespace opt { // void model_based_opt::project(unsigned x) { unsigned_vector& lub_rows = m_lub; - unsigned_vector& glb_rows = m_glb; + unsigned_vector& glb_rows = m_glb; + unsigned_vector& mod_rows = m_mod; unsigned lub_index = UINT_MAX, glb_index = UINT_MAX; bool lub_strict = false, glb_strict = false; rational lub_val, glb_val; @@ -533,6 +730,8 @@ namespace opt { uint_set visited; lub_rows.reset(); glb_rows.reset(); + mod_rows.reset(); + bool lub_is_unit = false, glb_is_unit = false; // select the lub and glb. for (unsigned i = 0; i < row_ids.size(); ++i) { unsigned row_id = row_ids[i]; @@ -552,7 +751,10 @@ namespace opt { solve_for(row_id, x); return; } - if (a.is_pos()) { + if (r.m_type == t_mod) { + mod_rows.push_back(row_id); + } + else if (a.is_pos()) { rational lub_value = x_val - (r.m_value/a); if (lub_rows.empty() || lub_value < lub_val || @@ -562,6 +764,7 @@ namespace opt { lub_strict = r.m_type == t_lt; } lub_rows.push_back(row_id); + lub_is_unit &= a.is_one(); } else { SASSERT(a.is_neg()); @@ -574,36 +777,187 @@ namespace opt { glb_strict = r.m_type == t_lt; } glb_rows.push_back(row_id); + glb_is_unit &= a.is_minus_one(); } } - unsigned row_index = (lub_rows.size() <= glb_rows.size())? lub_index : glb_index; - glb_rows.append(lub_rows); + + if (!mod_rows.empty()) { + solve_mod(x, mod_rows); + return; + } + + unsigned lub_size = lub_rows.size(); + unsigned glb_size = glb_rows.size(); + unsigned row_index = (lub_size <= glb_size) ? lub_index : glb_index; + glb_rows.append(lub_rows); + + // There are only upper or only lower bounds. if (row_index == UINT_MAX) { for (unsigned i = 0; i < glb_rows.size(); ++i) { unsigned row_id = glb_rows[i]; SASSERT(m_rows[row_id].m_alive); SASSERT(!get_coefficient(row_id, x).is_zero()); - m_rows[row_id].m_alive = false; + retire_row(row_id); } + return; } - else { - rational coeff = get_coefficient(row_index, x); - for (unsigned i = 0; i < glb_rows.size(); ++i) { - unsigned row_id = glb_rows[i]; - if (row_id != row_index) { - resolve(row_index, coeff, row_id, x); + + // The number of matching lower and upper bounds is small. + if ((lub_size <= 2 || glb_size <= 2) && + (lub_size <= 3 && glb_size <= 3) && + (!is_int(x) || lub_is_unit || glb_is_unit)) { + for (unsigned i = 0; i < lub_size; ++i) { + unsigned row_id1 = lub_rows[i]; + bool last = i + 1 == lub_rows.size(); + rational coeff = get_coefficient(row_id1, x); + for (unsigned j = 0; j < glb_size; ++j) { + unsigned row_id2 = glb_rows[j]; + if (last) { + resolve(row_id1, coeff, row_id2, x); + } + else { + row r(m_rows[row_id2]); + m_rows.push_back(r); + resolve(row_id1, coeff, m_rows.size()-1, x); + } } } - m_rows[row_index].m_alive = false; + for (unsigned i = 0; i < lub_size; ++i) { + retire_row(lub_rows[i]); + } + return; } + + // General case. + rational coeff = get_coefficient(row_index, x); + for (unsigned i = 0; i < glb_rows.size(); ++i) { + unsigned row_id = glb_rows[i]; + if (row_id != row_index) { + resolve(row_index, coeff, row_id, x); + } + } + retire_row(row_index); } + // + // compute D and u. + // + // D = lcm(d1, d2) + // u = eval(x) mod D + // + // d1 | (a1x + t1) & d2 | (a2x + t2) + // = + // d1 | (a1(D*x' + u) + t1) & d2 | (a2(D*x' + u) + t2) + // = + // d1 | (a1*u + t1) & d2 | (a2*u + t2) + // + // x := D*x' + u + // + + void model_based_opt::solve_mod(unsigned x, unsigned_vector const& mod_rows) { + SASSERT(!mod_rows.empty()); + rational D(1); + for (unsigned i = 0; i < mod_rows.size(); ++i) { + D = lcm(D, m_rows[mod_rows[i]].m_mod); + } + TRACE("opt", display(tout << "lcm: " << D << " tableau\n");); + rational val_x = m_var2value[x]; + rational u = mod(val_x, D); + SASSERT(u.is_nonneg() && u < D); + for (unsigned i = 0; i < mod_rows.size(); ++i) { + replace_var(mod_rows[i], x, u); + SASSERT(invariant(mod_rows[i], m_rows[mod_rows[i]])); + } + // + // update inequalities such that u is added to t and + // D is multiplied to coefficient of x. + // the interpretation of the new version of x is (x-u)/D + // + // a*x + t <= 0 + // a*(D*x' + u) + t <= 0 + // a*D*x' + a*u + t <= 0 + // + rational new_val = (val_x - u) / D; + SASSERT(new_val.is_int()); + unsigned y = add_var(new_val, true); + unsigned_vector const& row_ids = m_var2row_ids[x]; + uint_set visited; + for (unsigned i = 0; i < row_ids.size(); ++i) { + unsigned row_id = row_ids[i]; + if (!visited.contains(row_id)) { + // x |-> D*y + u + replace_var(row_id, x, D, y, u); + visited.insert(row_id); + } + } + project(y); + } + + // update row with: x |-> C + void model_based_opt::replace_var(unsigned row_id, unsigned x, rational const& C) { + row& r = m_rows[row_id]; + SASSERT(!get_coefficient(row_id, x).is_zero()); + unsigned sz = r.m_vars.size(); + unsigned i = 0, j = 0; + rational coeff(0); + for (; i < sz; ++i) { + if (r.m_vars[i].m_id == x) { + coeff = r.m_vars[i].m_coeff; + } + else { + if (i != j) { + r.m_vars[j] = r.m_vars[i]; + } + ++j; + } + } + if (j != sz) { + r.m_vars.shrink(j); + } + r.m_coeff += coeff*C; + r.m_value += coeff*(C - m_var2value[x]); + } + + // update row with: x |-> A*y + B + void model_based_opt::replace_var(unsigned row_id, unsigned x, rational const& A, unsigned y, rational const& B) { + row& r = m_rows[row_id]; + rational coeff = get_coefficient(row_id, x); + if (coeff.is_zero()) return; + if (!r.m_alive) return; + replace_var(row_id, x, B); + r.m_vars.push_back(var(y, coeff*A)); + r.m_value += coeff*A*m_var2value[y]; + if (!r.m_vars.empty() && r.m_vars.back().m_id > y) { + std::sort(r.m_vars.begin(), r.m_vars.end(), var::compare()); + } + m_var2row_ids[y].push_back(row_id); + SASSERT(invariant(row_id, r)); + } + + // 3x + t = 0 & 7 | (c*x + s) & ax <= u + // 3 | -t & 21 | (-ct + 3s) & a-t <= 3u + +#if 0 + void solve_for_int(unsigned row_id1, unsigned x) { + + for (unsigned i = 0; i < num_divs(); ++i) { + add_lit(model, lits, mk_divides(c*div_divisor(i), + mk_sub(mk_mul(c, div_term(i)), mk_mul(div_coeff(i), t)))); + } + } +#endif + void model_based_opt::solve_for(unsigned row_id1, unsigned x) { - rational a = get_coefficient(row_id1, x); - row& r1 = m_rows[row_id1]; + rational a = get_coefficient(row_id1, x), b; SASSERT(!a.is_zero()); - SASSERT(r1.m_type == t_eq); - SASSERT(r1.m_alive); + SASSERT(m_rows[row_id1].m_type == t_eq); + SASSERT(m_rows[row_id1].m_alive); + if (m_var2is_int[x] && !abs(a).is_one()) { + row& r1 = m_rows[row_id1]; + vector coeffs; + mk_coeffs_without(coeffs, r1.m_vars, x); + add_divides(coeffs, r1.m_coeff, abs(a)); + } unsigned_vector const& row_ids = m_var2row_ids[x]; uint_set visited; visited.insert(row_id1); @@ -611,15 +965,20 @@ namespace opt { unsigned row_id2 = row_ids[i]; if (!visited.contains(row_id2)) { visited.insert(row_id2); - resolve(row_id1, a, row_id2, x); + row const& r2 = m_rows[row_id2]; + b = get_coefficient(row_id2, x); + if (!b.is_zero()) { + resolve(row_id1, a, row_id2, x); + } } } - r1.m_alive = false; + retire_row(row_id1); } void model_based_opt::project(unsigned num_vars, unsigned const* vars) { for (unsigned i = 0; i < num_vars; ++i) { project(vars[i]); + TRACE("opt", display(tout << "After projecting: v" << vars[i] << "\n");); } } diff --git a/src/math/simplex/model_based_opt.h b/src/math/simplex/model_based_opt.h index 3bbcc7bf2..92f600522 100644 --- a/src/math/simplex/model_based_opt.h +++ b/src/math/simplex/model_based_opt.h @@ -30,7 +30,8 @@ namespace opt { enum ineq_type { t_eq, t_lt, - t_le + t_le, + t_mod }; @@ -52,9 +53,11 @@ namespace opt { row(): m_type(t_le), m_value(0), m_alive(false) {} vector m_vars; // variables with coefficients rational m_coeff; // constant in inequality + rational m_mod; // value the term divide ineq_type m_type; // inequality type rational m_value; // value of m_vars + m_coeff under interpretation of m_var2value. bool m_alive; // rows can be marked dead if they have been processed. + void reset() { m_vars.reset(); m_coeff.reset(); m_value.reset(); } }; private: @@ -63,9 +66,11 @@ namespace opt { static const unsigned m_objective_id = 0; vector m_var2row_ids; vector m_var2value; + svector m_var2is_int; vector m_new_vars; - unsigned_vector m_lub, m_glb; + unsigned_vector m_lub, m_glb, m_mod; unsigned_vector m_above, m_below; + unsigned_vector m_retired_rows; bool invariant(); bool invariant(unsigned index, row const& r); @@ -82,7 +87,29 @@ namespace opt { void mul_add(bool same_sign, unsigned row_id1, rational const& c, unsigned row_id2); - void set_row(unsigned row_id, vector const& coeffs, rational const& c, ineq_type rel); + void mul_add(unsigned x, rational const& a1, unsigned row_src, rational const& a2, unsigned row_dst); + + void mul(unsigned dst, rational const& c); + + void add(unsigned dst, rational const& c); + + void sub(unsigned dst, rational const& c); + + void set_row(unsigned row_id, vector const& coeffs, rational const& c, rational const& m, ineq_type rel); + + void add_constraint(vector const& coeffs, rational const& c, rational const& m, ineq_type r); + + void replace_var(unsigned row_id, unsigned x, rational const& A, unsigned y, rational const& B); + + void replace_var(unsigned row_id, unsigned x, rational const& C); + + void normalize(unsigned row_id); + + void mk_coeffs_without(vector& dst, vector const src, unsigned x); + + unsigned new_row(); + + rational n_sign(rational const& b) const; void update_values(unsigned_vector const& bound_vars, unsigned_vector const& bound_trail); @@ -92,12 +119,18 @@ namespace opt { void solve_for(unsigned row_id, unsigned x); + void solve_mod(unsigned x, unsigned_vector const& mod_rows); + + bool is_int(unsigned x) const { return m_var2is_int[x]; } + + void retire_row(unsigned row_id); + public: model_based_opt(); // add a fresh variable with value 'value'. - unsigned add_var(rational const& value); + unsigned add_var(rational const& value, bool is_int = false); // retrieve updated value of variable. rational get_value(unsigned var_id); @@ -106,6 +139,9 @@ namespace opt { // satisfied under the values provided to the variables. void add_constraint(vector const& coeffs, rational const& c, ineq_type r); + // add a divisibility constraint. The row should divide m. + void add_divides(vector const& coeffs, rational const& c, rational const& m); + // Set the objective function (linear). void set_objective(vector const& coeffs, rational const& c); diff --git a/src/qe/qe_arith.cpp b/src/qe/qe_arith.cpp index c2d0016c8..dab6325bc 100644 --- a/src/qe/qe_arith.cpp +++ b/src/qe/qe_arith.cpp @@ -15,7 +15,7 @@ Author: Revision History: - + Moved projection functionality to model_based_opt module. 2016-06-26 --*/ @@ -58,38 +58,15 @@ namespace qe { ast_manager& m; arith_util a; th_rewriter m_rw; - expr_ref_vector m_ineq_terms; - vector m_ineq_coeffs; - svector m_ineq_types; - expr_ref_vector m_div_terms; - vector m_div_divisors, m_div_coeffs; - expr_ref_vector m_new_lits; expr_ref_vector m_trail; - rational m_delta, m_u; - scoped_ptr m_var; - unsigned m_num_pos, m_num_neg; - bool m_pos_is_unit, m_neg_is_unit; - - sort* var_sort() const { return m.get_sort(m_var->x()); } - - bool is_int() const { return a.is_int(m_var->x()); } - - void display(std::ostream& out) const { - for (unsigned i = 0; i < num_ineqs(); ++i) { - display_ineq(out, i); - } - for (unsigned i = 0; i < num_divs(); ++i) { - display_div(out, i); - } - } void insert_mul(expr* x, rational const& v, obj_map& ts) { + TRACE("qe", tout << "Adding variable " << mk_pp(x, m) << " " << v << "\n";); rational w; if (ts.find(x, w)) { ts.insert(x, w + v); } else { - TRACE("qe", tout << "Adding variable " << mk_pp(x, m) << "\n";); ts.insert(x, v); } } @@ -99,14 +76,14 @@ namespace qe { // It uses the current model to choose values for conditionals and it primes mbo with the current // interpretation of sub-expressions that are treated as variables for mbo. // - bool linearize(opt::model_based_opt& mbo, model& model, expr* lit, expr_ref_vector& fmls, obj_map& tids) { + bool linearize(opt::model_based_opt& mbo, model_evaluator& eval, expr* lit, expr_ref_vector& fmls, obj_map& tids) { obj_map ts; rational c(0), mul(1); expr_ref t(m); opt::ineq_type ty = opt::t_le; expr* e1, *e2; DEBUG_CODE(expr_ref val(m); - VERIFY(model.eval(lit, val)); + eval(lit, val); CTRACE("qe", !m.is_true(val), tout << mk_pp(lit, m) << " := " << val << "\n";); SASSERT(m.is_true(val));); @@ -115,55 +92,84 @@ namespace qe { mul.neg(); } SASSERT(!m.is_not(lit)); - if ((a.is_le(lit, e1, e2) || a.is_ge(lit, e2, e1)) && a.is_real(e1)) { - linearize(mbo, model, mul, e1, c, fmls, ts, tids); - linearize(mbo, model, -mul, e2, c, fmls, ts, tids); + if ((a.is_le(lit, e1, e2) || a.is_ge(lit, e2, e1))) { + linearize(mbo, eval, mul, e1, c, fmls, ts, tids); + linearize(mbo, eval, -mul, e2, c, fmls, ts, tids); ty = is_not ? opt::t_lt : opt::t_le; } - else if ((a.is_lt(lit, e1, e2) || a.is_gt(lit, e2, e1)) && a.is_real(e1)) { - linearize(mbo, model, mul, e1, c, fmls, ts, tids); - linearize(mbo, model, -mul, e2, c, fmls, ts, tids); + else if ((a.is_lt(lit, e1, e2) || a.is_gt(lit, e2, e1))) { + linearize(mbo, eval, mul, e1, c, fmls, ts, tids); + linearize(mbo, eval, -mul, e2, c, fmls, ts, tids); ty = is_not ? opt::t_le: opt::t_lt; } - else if (m.is_eq(lit, e1, e2) && !is_not && a.is_real(e1)) { - linearize(mbo, model, mul, e1, c, fmls, ts, tids); - linearize(mbo, model, -mul, e2, c, fmls, ts, tids); + else if (m.is_eq(lit, e1, e2) && !is_not && is_arith(e1)) { + linearize(mbo, eval, mul, e1, c, fmls, ts, tids); + linearize(mbo, eval, -mul, e2, c, fmls, ts, tids); ty = opt::t_eq; } - else if (m.is_eq(lit, e1, e2) && is_not && a.is_real(e1)) { + else if (m.is_eq(lit, e1, e2) && is_not && is_arith(e1)) { expr_ref val1(m), val2(m); rational r1, r2; - VERIFY(model.eval(e1, val1) && a.is_numeral(val1, r1)); - VERIFY(model.eval(e2, val2) && a.is_numeral(val2, r2)); + eval(e1, val1); eval(e2, val2); + VERIFY(a.is_numeral(val1, r1)); + VERIFY(a.is_numeral(val2, r2)); SASSERT(r1 != r2); - if (r2 < r1) { + if (r1 < r2) { std::swap(e1, e2); } ty = opt::t_lt; - linearize(mbo, model, mul, e1, c, fmls, ts, tids); - linearize(mbo, model, -mul, e2, c, fmls, ts, tids); + linearize(mbo, eval, mul, e1, c, fmls, ts, tids); + linearize(mbo, eval, -mul, e2, c, fmls, ts, tids); } - else if (m.is_distinct(lit) && !is_not && a.is_real(to_app(lit)->get_arg(0))) { - TRACE("qe", tout << "TBD: handle distinc\n";); - return false; + else if (m.is_distinct(lit) && !is_not && is_arith(to_app(lit)->get_arg(0))) { + expr_ref val(m); + rational r; + app* alit = to_app(lit); + vector > nums; + for (unsigned i = 0; i < alit->get_num_args(); ++i) { + eval(alit->get_arg(i), val); + VERIFY(a.is_numeral(val, r)); + nums.push_back(std::make_pair(alit->get_arg(i), r)); + } + std::sort(nums.begin(), nums.end(), compare_second()); + for (unsigned i = 0; i + 1 < nums.size(); ++i) { + SASSERT(nums[i].second < nums[i+1].second); + expr_ref fml(a.mk_lt(nums[i].first, nums[i+1].first), m); + if (!linearize(mbo, eval, fml, fmls, tids)) { + return false; + } + } + return true; } - else if (m.is_distinct(lit) && is_not && a.is_real(to_app(lit)->get_arg(0))) { - TRACE("qe", tout << "TBD: handle negation of distinc\n";); - return false; + else if (m.is_distinct(lit) && is_not && is_arith(to_app(lit)->get_arg(0))) { + // find the two arguments that are equal. + // linearize these. + map values; + bool found_eq = false; + for (unsigned i = 0; !found_eq && i < to_app(lit)->get_num_args(); ++i) { + expr* arg1 = to_app(lit)->get_arg(i), *arg2 = 0; + expr_ref val(m); + rational r; + eval(arg1, val); + VERIFY(a.is_numeral(val, r)); + if (values.find(r, arg2)) { + ty = opt::t_eq; + linearize(mbo, eval, mul, arg1, c, fmls, ts, tids); + linearize(mbo, eval, -mul, arg2, c, fmls, ts, tids); + found_eq = true; + } + else { + values.insert(r, arg1); + } + } + SASSERT(found_eq); } else { TRACE("qe", tout << "Skipping " << mk_pp(lit, m) << "\n";); return false; } -#if 0 - TBD for integers - if (ty == opt::t_lt && false) { - c += rational(1); - ty = opt::t_le; - } -#endif vars coeffs; - extract_coefficients(mbo, model, ts, tids, coeffs); + extract_coefficients(mbo, eval, ts, tids, coeffs); mbo.add_constraint(coeffs, c, ty); return true; } @@ -171,212 +177,76 @@ namespace qe { // // convert linear arithmetic term into an inequality for mbo. // - void linearize(opt::model_based_opt& mbo, model& model, rational const& mul, expr* t, rational& c, + void linearize(opt::model_based_opt& mbo, model_evaluator& eval, rational const& mul, expr* t, rational& c, expr_ref_vector& fmls, obj_map& ts, obj_map& tids) { expr* t1, *t2, *t3; rational mul1; expr_ref val(m); - if (a.is_mul(t, t1, t2) && is_numeral(model, t1, mul1)) { - linearize(mbo, model, mul* mul1, t2, c, fmls, ts, tids); + if (a.is_mul(t, t1, t2) && is_numeral(eval, t1, mul1)) { + linearize(mbo, eval, mul* mul1, t2, c, fmls, ts, tids); } - else if (a.is_mul(t, t1, t2) && is_numeral(model, t2, mul1)) { - linearize(mbo, model, mul* mul1, t1, c, fmls, ts, tids); + else if (a.is_mul(t, t1, t2) && is_numeral(eval, t2, mul1)) { + linearize(mbo, eval, mul* mul1, t1, c, fmls, ts, tids); } else if (a.is_add(t)) { app* ap = to_app(t); for (unsigned i = 0; i < ap->get_num_args(); ++i) { - linearize(mbo, model, mul, ap->get_arg(i), c, fmls, ts, tids); + linearize(mbo, eval, mul, ap->get_arg(i), c, fmls, ts, tids); } } else if (a.is_sub(t, t1, t2)) { - linearize(mbo, model, mul, t1, c, fmls, ts, tids); - linearize(mbo, model, -mul, t2, c, fmls, ts, tids); + linearize(mbo, eval, mul, t1, c, fmls, ts, tids); + linearize(mbo, eval, -mul, t2, c, fmls, ts, tids); } else if (a.is_uminus(t, t1)) { - linearize(mbo, model, -mul, t1, c, fmls, ts, tids); + linearize(mbo, eval, -mul, t1, c, fmls, ts, tids); } else if (a.is_numeral(t, mul1)) { c += mul*mul1; } else if (m.is_ite(t, t1, t2, t3)) { - VERIFY(model.eval(t1, val)); + eval(t1, val); SASSERT(m.is_true(val) || m.is_false(val)); TRACE("qe", tout << mk_pp(t1, m) << " := " << val << "\n";); if (m.is_true(val)) { - linearize(mbo, model, mul, t2, c, fmls, ts, tids); + linearize(mbo, eval, mul, t2, c, fmls, ts, tids); fmls.push_back(t1); } else { expr_ref not_t1(mk_not(m, t1), m); fmls.push_back(not_t1); - linearize(mbo, model, mul, t3, c, fmls, ts, tids); + linearize(mbo, eval, mul, t3, c, fmls, ts, tids); } } + else if (a.is_mod(t, t1, t2) && is_numeral(eval, t2, mul1)) { + rational r; + eval(t, val); + VERIFY(a.is_numeral(val, r)); + c += mul*r; + // t1 mod mul1 == r + rational c0(-r), mul0(1); + obj_map ts0; + linearize(mbo, eval, mul0, t1, c0, fmls, ts0, tids); + vars coeffs; + extract_coefficients(mbo, eval, ts0, tids, coeffs); + mbo.add_divides(coeffs, c0, mul1); + } else { insert_mul(t, mul, ts); } } - // - // extract linear terms from t into c and ts. - // - void is_linear(model& model, rational const& mul, expr* t, rational& c, expr_ref_vector& ts) { - expr* t1, *t2, *t3; - rational mul1; - expr_ref val(m); - if (t == m_var->x()) { - c += mul; - } - else if (a.is_mul(t, t1, t2) && is_numeral(model, t1, mul1)) { - is_linear(model, mul* mul1, t2, c, ts); - } - else if (a.is_mul(t, t1, t2) && is_numeral(model, t2, mul1)) { - is_linear(model, mul* mul1, t1, c, ts); - } - else if (a.is_add(t)) { - app* ap = to_app(t); - for (unsigned i = 0; i < ap->get_num_args(); ++i) { - is_linear(model, mul, ap->get_arg(i), c, ts); - } - } - else if (a.is_sub(t, t1, t2)) { - is_linear(model, mul, t1, c, ts); - is_linear(model, -mul, t2, c, ts); - } - else if (a.is_uminus(t, t1)) { - is_linear(model, -mul, t1, c, ts); - } - else if (a.is_numeral(t, mul1)) { - ts.push_back(mk_num(mul*mul1)); - } - else if (extract_mod(model, t, val)) { - ts.push_back(mk_mul(mul, val)); - } - else if (m.is_ite(t, t1, t2, t3)) { - VERIFY(model.eval(t1, val)); - SASSERT(m.is_true(val) || m.is_false(val)); - TRACE("qe", tout << mk_pp(t1, m) << " := " << val << "\n";); - if (m.is_true(val)) { - is_linear(model, mul, t2, c, ts); - } - else { - is_linear(model, mul, t3, c, ts); - } - } - else if ((*m_var)(t)) { - TRACE("qe", tout << "can't project:" << mk_pp(t, m) << "\n";); - throw cant_project(); - } - else { - ts.push_back(mk_mul(mul, t)); - } - } - - // - // extract linear inequalities from literal lit. - // - bool is_linear(model& model, expr* lit, bool& found_eq) { - rational c(0), mul(1); - expr_ref t(m); - opt::ineq_type ty = opt::t_le; - expr* e1, *e2; - expr_ref_vector ts(m); - bool is_not = m.is_not(lit, lit); - if (is_not) { - mul.neg(); - } - SASSERT(!m.is_not(lit)); - if (a.is_le(lit, e1, e2) || a.is_ge(lit, e2, e1)) { - is_linear(model, mul, e1, c, ts); - is_linear(model, -mul, e2, c, ts); - ty = is_not? opt::t_lt : opt::t_le; - } - else if (a.is_lt(lit, e1, e2) || a.is_gt(lit, e2, e1)) { - is_linear(model, mul, e1, c, ts); - is_linear(model, -mul, e2, c, ts); - ty = is_not? opt::t_le: opt::t_lt; - } - else if (m.is_eq(lit, e1, e2) && !is_not && is_arith(e1)) { - is_linear(model, mul, e1, c, ts); - is_linear(model, -mul, e2, c, ts); - ty = opt::t_eq; - } - else if (m.is_distinct(lit) && !is_not && is_arith(to_app(lit)->get_arg(0))) { - expr_ref val(m); - rational r; - app* alit = to_app(lit); - vector > nums; - for (unsigned i = 0; i < alit->get_num_args(); ++i) { - VERIFY(model.eval(alit->get_arg(i), val) && a.is_numeral(val, r)); - nums.push_back(std::make_pair(alit->get_arg(i), r)); - } - std::sort(nums.begin(), nums.end(), compare_second()); - for (unsigned i = 0; i + 1 < nums.size(); ++i) { - SASSERT(nums[i].second < nums[i+1].second); - c.reset(); - ts.reset(); - is_linear(model, mul, nums[i+1].first, c, ts); - is_linear(model, -mul, nums[i].first, c, ts); - t = add(ts); - accumulate_linear(model, c, t, opt::t_lt); - } - t = mk_num(0); - c.reset(); - return true; - } - else if (m.is_distinct(lit) && is_not && is_arith(to_app(lit)->get_arg(0))) { - expr_ref eq = project_plugin::pick_equality(m, model, to_app(lit)->get_arg(0)); - return is_linear(model, eq, found_eq); - } - else if (m.is_eq(lit, e1, e2) && is_not && is_arith(e1)) { - expr_ref val1(m), val2(m); - rational r1, r2; - VERIFY(model.eval(e1, val1) && a.is_numeral(val1, r1)); - VERIFY(model.eval(e2, val2) && a.is_numeral(val2, r2)); - SASSERT(r1 != r2); - if (r1 < r2) { - std::swap(e1, e2); - } - ty = opt::t_lt; - is_linear(model, mul, e1, c, ts); - is_linear(model, -mul, e2, c, ts); - } - else { - TRACE("qe", tout << "can't project:" << mk_pp(lit, m) << "\n";); - throw cant_project(); - } - if (ty == opt::t_lt && is_int()) { - ts.push_back(mk_num(1)); - ty = opt::t_le; - } - t = add(ts); - if (ty == opt::t_eq && c.is_neg()) { - t = mk_uminus(t); - c.neg(); - } - if (ty == opt::t_eq && c > rational(1)) { - m_ineq_coeffs.push_back(-c); - m_ineq_terms.push_back(mk_uminus(t)); - m_ineq_types.push_back(opt::t_le); - m_num_neg++; - ty = opt::t_le; - } - accumulate_linear(model, c, t, ty); - found_eq = !c.is_zero() && ty == opt::t_eq; - return true; - } - - bool is_numeral(model& model, expr* t, rational& r) { + bool is_numeral(model_evaluator& eval, expr* t, rational& r) { expr* t1, *t2, *t3; rational r1, r2; expr_ref val(m); if (a.is_numeral(t, r)) return true; - if (a.is_uminus(t, t1) && is_numeral(model, t1, r)) { + if (a.is_uminus(t, t1) && is_numeral(eval, t1, r)) { r.neg(); return true; } - else if (a.is_mul(t, t1, t2) && is_numeral(model, t1, r1) && is_numeral(model, t2, r2)) { + else if (a.is_mul(t, t1, t2) && is_numeral(eval, t1, r1) && is_numeral(eval, t2, r2)) { r = r1*r2; return true; } @@ -384,21 +254,21 @@ namespace qe { app* ap = to_app(t); r = rational(1); for (unsigned i = 0; i < ap->get_num_args(); ++i) { - if (!is_numeral(model, ap->get_arg(i), r1)) return false; + if (!is_numeral(eval, ap->get_arg(i), r1)) return false; r *= r1; } return true; } else if (m.is_ite(t, t1, t2, t3)) { - VERIFY (model.eval(t1, val)); + eval(t1, val); if (m.is_true(val)) { - return is_numeral(model, t1, r); + return is_numeral(eval, t1, r); } else { - return is_numeral(model, t2, r); + return is_numeral(eval, t2, r); } } - else if (a.is_sub(t, t1, t2) && is_numeral(model, t1, r1) && is_numeral(model, t2, r2)) { + else if (a.is_sub(t, t1, t2) && is_numeral(eval, t1, r1) && is_numeral(eval, t2, r2)) { r = r1 - r2; return true; } @@ -413,564 +283,16 @@ namespace qe { } }; - void accumulate_linear(model& model, rational const& c, expr_ref& t, opt::ineq_type ty) { - if (c.is_zero()) { - switch (ty) { - case opt::t_eq: - t = a.mk_eq(t, mk_num(0)); - break; - case opt::t_lt: - t = a.mk_lt(t, mk_num(0)); - break; - case opt::t_le: - t = a.mk_le(t, mk_num(0)); - break; - } - add_lit(model, m_new_lits, t); - } - else { - m_ineq_coeffs.push_back(c); - m_ineq_terms.push_back(t); - m_ineq_types.push_back(ty); - if (ty == opt::t_eq) { - // skip - } - else if (c.is_pos()) { - ++m_num_pos; - m_pos_is_unit &= c.is_one(); - } - else { - ++m_num_neg; - m_neg_is_unit &= c.is_minus_one(); - } - } - } - bool is_arith(expr* e) { return a.is_int(e) || a.is_real(e); } - - expr_ref add(expr_ref_vector const& ts) { - switch (ts.size()) { - case 0: - return mk_num(0); - case 1: - return expr_ref(ts[0], m); - default: - return expr_ref(a.mk_add(ts.size(), ts.c_ptr()), m); - } - } - - - // e is of the form (ax + t) mod k - bool is_mod(model& model, expr* e, rational& k, expr_ref& t, rational& c) { - expr* t1, *t2; - expr_ref_vector ts(m); - if (a.is_mod(e, t1, t2) && - a.is_numeral(t2, k) && - (*m_var)(t1)) { - c.reset(); - rational mul(1); - is_linear(model, mul, t1, c, ts); - t = add(ts); - return true; - } - return false; - } - - bool extract_mod(model& model, expr* e, expr_ref& val) { - rational c, k; - expr_ref t(m); - if (is_mod(model, e, k, t, c)) { - VERIFY (model.eval(e, val)); - SASSERT (a.is_numeral(val)); - TRACE("qe", tout << "extract: " << mk_pp(e, m) << " evals: " << val << " c: " << c << " t: " << t << "\n";); - if (!c.is_zero()) { - t = mk_sub(t, val); - m_div_terms.push_back(t); - m_div_divisors.push_back(k); - m_div_coeffs.push_back(c); - } - else { - t = m.mk_eq(a.mk_mod(t, mk_num(k)), val); - add_lit(model, m_new_lits, t); - } - return true; - } - return false; - } - - bool lit_is_true(model& model, expr* e) { - expr_ref val(m); - VERIFY(model.eval(e, val)); - CTRACE("qe", !m.is_true(val), tout << "eval: " << mk_pp(e, m) << " " << val << "\n";); - return m.is_true(val); - } - - expr_ref mk_num(unsigned n) { - rational r(n); - return mk_num(r); - } - - expr_ref mk_num(rational const& r) const { - return expr_ref(a.mk_numeral(r, var_sort()), m); - } - - expr_ref mk_divides(rational const& k, expr* t) { - return expr_ref(m.mk_eq(a.mk_mod(t, mk_num(abs(k))), mk_num(0)), m); - } - - void reset() { - reset_ineqs(); - reset_divs(); - m_delta = rational(1); - m_u = rational(0); - m_new_lits.reset(); - } - - void reset_divs() { - m_div_terms.reset(); - m_div_coeffs.reset(); - m_div_divisors.reset(); - } - - void reset_ineqs() { - m_ineq_terms.reset(); - m_ineq_coeffs.reset(); - m_ineq_types.reset(); - } - - expr* ineq_term(unsigned i) const { return m_ineq_terms[i]; } - rational const& ineq_coeff(unsigned i) const { return m_ineq_coeffs[i]; } - opt::ineq_type ineq_ty(unsigned i) const { return m_ineq_types[i]; } - app_ref mk_ineq_pred(unsigned i) { - app_ref result(m); - result = to_app(mk_add(mk_mul(ineq_coeff(i), m_var->x()), ineq_term(i))); - switch (ineq_ty(i)) { - case opt::t_lt: - result = a.mk_lt(result, mk_num(0)); - break; - case opt::t_le: - result = a.mk_le(result, mk_num(0)); - break; - case opt::t_eq: - result = m.mk_eq(result, mk_num(0)); - break; - } - return result; - } - void display_ineq(std::ostream& out, unsigned i) const { - out << mk_pp(ineq_term(i), m) << " " << ineq_coeff(i) << "*" << mk_pp(m_var->x(), m); - switch (ineq_ty(i)) { - case opt::t_eq: out << " = 0\n"; break; - case opt::t_le: out << " <= 0\n"; break; - case opt::t_lt: out << " < 0\n"; break; - } - } - unsigned num_ineqs() const { return m_ineq_terms.size(); } - expr* div_term(unsigned i) const { return m_div_terms[i]; } - rational const& div_coeff(unsigned i) const { return m_div_coeffs[i]; } - rational const& div_divisor(unsigned i) const { return m_div_divisors[i]; } - void display_div(std::ostream& out, unsigned i) const { - out << div_divisor(i) << " | ( " << mk_pp(div_term(i), m) << " " << div_coeff(i) << "*" - << mk_pp(m_var->x(), m) << ")\n"; - } - unsigned num_divs() const { return m_div_terms.size(); } - - void project(model& model, expr_ref_vector& lits) { - TRACE("qe", - tout << "project: " << mk_pp(m_var->x(), m) << "\n"; - tout << lits; - model_v2_pp(tout, model); ); - - m_num_pos = 0; m_num_neg = 0; - m_pos_is_unit = true; m_neg_is_unit = true; - unsigned eq_index = 0; - reset(); - bool found_eq = false; - for (unsigned i = 0; i < lits.size(); ++i) { - bool found_eq0 = false; - expr* e = lits[i].get(); - if (!(*m_var)(e)) { - m_new_lits.push_back(e); - } - else if (!is_linear(model, e, found_eq0)) { - TRACE("qe", tout << "can't project:" << mk_pp(e, m) << "\n";); - throw cant_project(); - } - if (found_eq0 && !found_eq) { - found_eq = true; - eq_index = num_ineqs()-1; - } - } - TRACE("qe", display(tout << mk_pp(m_var->x(), m) << ":\n"); - tout << "found eq: " << found_eq << " @ " << eq_index << "\n"; - tout << "num pos: " << m_num_pos << " num neg: " << m_num_neg << " num divs " << num_divs() << "\n"; - ); - lits.reset(); - lits.append(m_new_lits); - if (found_eq) { - apply_equality(model, eq_index, lits); - return; - } - if (num_divs() == 0 && (m_num_pos == 0 || m_num_neg == 0)) { - return; - } - if (num_divs() > 0) { - apply_divides(model, lits); - TRACE("qe", display(tout << "after division " << mk_pp(m_var->x(), m) << "\n");); - } - if (m_num_pos == 0 || m_num_neg == 0) { - return; - } - if ((m_num_pos <= 2 || m_num_neg <= 2) && - (m_num_pos == 1 || m_num_neg == 1 || (m_num_pos <= 3 && m_num_neg <= 3)) && - (!is_int() || m_pos_is_unit || m_neg_is_unit)) { - - unsigned index1 = num_ineqs(); - unsigned index2 = num_ineqs(); - bool is_pos = m_num_pos <= m_num_neg; - for (unsigned i = 0; i < num_ineqs(); ++i) { - if (ineq_coeff(i).is_pos() == is_pos) { - if (index1 == num_ineqs()) { - index1 = i; - } - else { - SASSERT(index2 == num_ineqs()); - index2 = i; - } - } - } - for (unsigned i = 0; i < num_ineqs(); ++i) { - if (ineq_coeff(i).is_pos() != is_pos) { - SASSERT(index1 != num_ineqs()); - mk_lt(model, lits, i, index1); - if (index2 != num_ineqs()) { - mk_lt(model, lits, i, index2); - } - } - } - } - else { - expr_ref t(m); - bool use_pos = m_num_pos < m_num_neg; - unsigned max_t = find_max(model, use_pos); - - for (unsigned i = 0; i < num_ineqs(); ++i) { - if (i != max_t) { - if (ineq_coeff(i).is_pos() == use_pos) { - t = mk_le(i, max_t); - add_lit(model, lits, t); - } - else { - mk_lt(model, lits, i, max_t); - } - } - } - } - TRACE("qe", tout << lits;); - } - - unsigned find_max(model& mdl, bool do_pos) { - unsigned result = 0; - bool new_max = true; - rational max_r, r; - expr_ref val(m); - model_evaluator eval(mdl); - eval.set_model_completion(true); - - bool is_int = a.is_int(m_var->x()); - for (unsigned i = 0; i < num_ineqs(); ++i) { - rational const& ac = m_ineq_coeffs[i]; - SASSERT(!is_int || opt::t_le == ineq_ty(i)); - - // - // ac*x + t < 0 - // ac > 0: x + max { t/ac | ac > 0 } < 0 <=> x < - max { t/ac | ac > 0 } - // ac < 0: x + t/ac > 0 <=> x > max { - t/ac | ac < 0 } = max { t/|ac| | ac < 0 } - // - if (ac.is_pos() == do_pos) { - eval(ineq_term(i), val); - VERIFY(a.is_numeral(val, r)); - r /= abs(ac); - new_max = - new_max || - (r > max_r) || - (r == max_r && opt::t_lt == ineq_ty(i)) || - (r == max_r && is_int && ac.is_one()); - TRACE("qe", tout << "max: " << mk_pp(ineq_term(i), m) << "/" << abs(ac) << " := " << r << " " - << (new_max?"":"not ") << "new max\n";); - if (new_max) { - result = i; - max_r = r; - } - new_max = false; - } - } - SASSERT(!new_max); - return result; - } - - // ax + t <= 0 - // bx + s <= 0 - // a and b have different signs. - // Infer: a|b|x + |b|t + |a|bx + |a|s <= 0 - // e.g. |b|t + |a|s <= 0 - void mk_lt(model& model, expr_ref_vector& lits, unsigned i, unsigned j) { - rational const& ac = ineq_coeff(i); - rational const& bc = ineq_coeff(j); - SASSERT(ac.is_pos() != bc.is_pos()); - SASSERT(ac.is_neg() != bc.is_neg()); - - TRACE("qe", display_ineq(tout, i); display_ineq(tout, j);); - - if (is_int() && !abs(ac).is_one() && !abs(bc).is_one()) { - return mk_int_lt(model, lits, i, j); - } - expr* t = ineq_term(i); - expr* s = ineq_term(j); - expr_ref bt = mk_mul(abs(bc), t); - expr_ref as = mk_mul(abs(ac), s); - expr_ref ts = mk_add(bt, as); - expr_ref z = mk_num(0); - expr_ref fml(m); - if (opt::t_lt == ineq_ty(i) || opt::t_lt == ineq_ty(j)) { - fml = a.mk_lt(ts, z); - } - else { - fml = a.mk_le(ts, z); - } - add_lit(model, lits, fml); - } - - void mk_int_lt(model& model, expr_ref_vector& lits, unsigned i, unsigned j) { - TRACE("qe", display_ineq(tout, i); display_ineq(tout, j);); - - expr* t = ineq_term(i); - expr* s = ineq_term(j); - rational ac = ineq_coeff(i); - rational bc = ineq_coeff(j); - expr_ref tmp(m); - SASSERT(opt::t_le == ineq_ty(i) && opt::t_le == ineq_ty(j)); - SASSERT(ac.is_pos() == bc.is_neg()); - rational abs_a = abs(ac); - rational abs_b = abs(bc); - expr_ref as(mk_mul(abs_a, s), m); - expr_ref bt(mk_mul(abs_b, t), m); - - rational slack = (abs_a - rational(1))*(abs_b-rational(1)); - rational sval, tval; - VERIFY (model.eval(ineq_term(i), tmp) && a.is_numeral(tmp, tval)); - VERIFY (model.eval(ineq_term(j), tmp) && a.is_numeral(tmp, sval)); - bool use_case1 = ac*sval + bc*tval + slack <= rational(0); - if (use_case1) { - expr_ref_vector ts(m); - ts.push_back(as); - ts.push_back(bt); - ts.push_back(mk_num(-slack)); - tmp = a.mk_le(add(ts), mk_num(0)); - add_lit(model, lits, tmp); - return; - } - - if (abs_a < abs_b) { - std::swap(abs_a, abs_b); - std::swap(ac, bc); - std::swap(s, t); - std::swap(as, bt); - std::swap(sval, tval); - } - SASSERT(abs_a >= abs_b); - - // create finite disjunction for |b|. - // exists x, z in [0 .. |b|-2] . b*x + s + z = 0 && ax + t <= 0 && bx + s <= 0 - // <=> - // exists x, z in [0 .. |b|-2] . b*x = -z - s && ax + t <= 0 && bx + s <= 0 - // <=> - // exists x, z in [0 .. |b|-2] . b*x = -z - s && a|b|x + |b|t <= 0 && bx + s <= 0 - // <=> - // exists x, z in [0 .. |b|-2] . b*x = -z - s && a|b|x + |b|t <= 0 && -z - s + s <= 0 - // <=> - // exists x, z in [0 .. |b|-2] . b*x = -z - s && a|b|x + |b|t <= 0 && -z <= 0 - // <=> - // exists x, z in [0 .. |b|-2] . b*x = -z - s && a|b|x + |b|t <= 0 - // <=> - // exists x, z in [0 .. |b|-2] . b*x = -z - s && a*n_sign(b)(s + z) + |b|t <= 0 - // <=> - // exists z in [0 .. |b|-2] . |b| | (z + s) && a*n_sign(b)(s + z) + |b|t <= 0 - // - - rational z = mod(sval, abs_b); - if (!z.is_zero()) z = abs_b - z; - expr_ref s_plus_z(mk_add(z, s), m); - - tmp = mk_divides(abs_b, s_plus_z); - add_lit(model, lits, tmp); - tmp = a.mk_le(mk_add(mk_mul(ac*n_sign(bc), s_plus_z), - mk_mul(abs_b, t)), mk_num(0)); - add_lit(model, lits, tmp); - } - rational n_sign(rational const& b) { return rational(b.is_pos()?-1:1); } - // ax + t <= 0 - // bx + s <= 0 - // a and b have same signs. - // encode: - // t/|a| <= s/|b| - // e.g. |b|t <= |a|s - expr_ref mk_le(unsigned i, unsigned j) { - rational const& ac = ineq_coeff(i); - rational const& bc = ineq_coeff(j); - SASSERT(ac.is_pos() == bc.is_pos()); - SASSERT(ac.is_neg() == bc.is_neg()); - expr* t = ineq_term(i); - expr* s = ineq_term(j); - expr_ref bt = mk_mul(abs(bc), t); - expr_ref as = mk_mul(abs(ac), s); - if (opt::t_lt == ineq_ty(i) && opt::t_le == ineq_ty(j)) { - return expr_ref(a.mk_lt(bt, as), m); - } - else { - return expr_ref(a.mk_le(bt, as), m); - } - } - - expr_ref mk_add(expr* t1, expr* t2) { - rational r; - if (a.is_numeral(t1, r) && r.is_zero()) return expr_ref(t2, m); - if (a.is_numeral(t2, r) && r.is_zero()) return expr_ref(t1, m); - return expr_ref(a.mk_add(t1, t2), m); - } - expr_ref mk_add(rational const& r, expr* e) { - if (r.is_zero()) return expr_ref(e, m); - return mk_add(mk_num(r), e); - } - - expr_ref mk_mul(rational const& r, expr* t) { - if (r.is_one()) return expr_ref(t, m); - return expr_ref(a.mk_mul(mk_num(r), t), m); - } - - expr_ref mk_sub(expr* t1, expr* t2) { - rational r1, r2; - if (a.is_numeral(t2, r2) && r2.is_zero()) return expr_ref(t1, m); - if (a.is_numeral(t1, r1) && a.is_numeral(t2, r2)) return mk_num(r1 - r2); - return expr_ref(a.mk_sub(t1, t2), m); - } - - expr_ref mk_uminus(expr* t) { - rational r; - if (a.is_numeral(t, r)) { - return mk_num(-r); - } - return expr_ref(a.mk_uminus(t), m); - } - - void add_lit(model& model, expr_ref_vector& lits, expr* e) { - expr_ref orig(e, m), result(m); - m_rw(orig, result); - TRACE("qe", tout << mk_pp(orig, m) << " -> " << result << "\n";); - SASSERT(lit_is_true(model, orig)); - SASSERT(lit_is_true(model, result)); - if (!m.is_true(result)) { - lits.push_back(result); - } - } - - - // 3x + t = 0 & 7 | (c*x + s) & ax <= u - // 3 | -t & 21 | (-ct + 3s) & a-t <= 3u - - void apply_equality(model& model, unsigned eq_index, expr_ref_vector& lits) { - rational c = ineq_coeff(eq_index); - expr* t = ineq_term(eq_index); - SASSERT(c.is_pos()); - if (is_int()) { - add_lit(model, lits, mk_divides(c, ineq_term(eq_index))); - } - - for (unsigned i = 0; i < num_divs(); ++i) { - add_lit(model, lits, mk_divides(c*div_divisor(i), - mk_sub(mk_mul(c, div_term(i)), mk_mul(div_coeff(i), t)))); - } - for (unsigned i = 0; i < num_ineqs(); ++i) { - if (eq_index != i) { - expr_ref lhs(m), val(m); - lhs = mk_sub(mk_mul(c, ineq_term(i)), mk_mul(ineq_coeff(i), t)); - switch (ineq_ty(i)) { - case opt::t_lt: lhs = a.mk_lt(lhs, mk_num(0)); break; - case opt::t_le: lhs = a.mk_le(lhs, mk_num(0)); break; - case opt::t_eq: lhs = m.mk_eq(lhs, mk_num(0)); break; - } - TRACE("qe", tout << lhs << "\n";); - SASSERT (model.eval(lhs, val) && m.is_true(val)); - add_lit(model, lits, lhs); - } - } - } - - // - // compute D and u. - // - // D = lcm(d1, d2) - // u = eval(x) mod D - // - // d1 | (a1x + t1) & d2 | (a2x + t2) - // = - // D | (D/d1)(a1x + t1) & D | (D/d2)(a2x + t2) - // = - // D | D1(a1*u + t1) & D | D2(a2*u + t2) & x = D*x' + u & 0 <= u < D - // = - // D | D1(a1*u + t1) & D | D2(a2*u + t2) & x = D*x' + u & 0 <= u < D - // - // x := D*x' + u - // - void apply_divides(model& model, expr_ref_vector& lits) { - SASSERT(m_delta.is_one()); - unsigned n = num_divs(); - if (n == 0) { - return; - } - for (unsigned i = 0; i < n; ++i) { - m_delta = lcm(m_delta, div_divisor(i)); - } - expr_ref val(m); - rational r; - VERIFY (model.eval(m_var->x(), val) && a.is_numeral(val, r)); - m_u = mod(r, m_delta); - SASSERT(m_u < m_delta && rational(0) <= m_u); - for (unsigned i = 0; i < n; ++i) { - add_lit(model, lits, mk_divides(div_divisor(i), - mk_add(mk_num(div_coeff(i) * m_u), div_term(i)))); - } - reset_divs(); - // - // update inequalities such that u is added to t and - // D is multiplied to coefficient of x. - // the interpretation of the new version of x is (x-u)/D - // - // a*x + t <= 0 - // a*(D*x' + u) + t <= 0 - // a*D*x' + a*u + t <= 0 - for (unsigned i = 0; i < num_ineqs(); ++i) { - if (!m_u.is_zero()) { - m_ineq_terms[i] = a.mk_add(ineq_term(i), mk_num(m_ineq_coeffs[i]*m_u)); - } - m_ineq_coeffs[i] *= m_delta; - } - r = (r - m_u) / m_delta; - SASSERT(r.is_int()); - val = a.mk_numeral(r, true); - model.register_decl(m_var->x()->get_decl(), val); - TRACE("qe", model_v2_pp(tout, model);); - } - imp(ast_manager& m): - m(m), a(m), m_rw(m), m_ineq_terms(m), m_div_terms(m), m_new_lits(m), m_trail(m) { + m(m), a(m), m_rw(m), m_trail(m) { params_ref params; params.set_bool("gcd_rouding", true); m_rw.updt_params(params); @@ -984,43 +306,43 @@ namespace qe { } bool operator()(model& model, app* v, app_ref_vector& vars, expr_ref_vector& lits) { - SASSERT(a.is_real(v) || a.is_int(v)); - m_var = alloc(contains_app, m, v); - try { - project(model, lits); - } - catch (cant_project) { - TRACE("qe", tout << "can't project:" << mk_pp(v, m) << "\n";); - return false; - } - return true; + app_ref_vector vs(m); + vs.push_back(v); + (*this)(model, vs, lits); + return vs.empty(); } typedef opt::model_based_opt::var var; typedef opt::model_based_opt::row row; typedef vector vars; - void operator()(model& model, app_ref_vector& vars, expr_ref_vector& fmls) { - bool has_real = false; - for (unsigned i = 0; !has_real && i < vars.size(); ++i) { - has_real = a.is_real(vars[i].get()); + bool has_arith = false; + for (unsigned i = 0; !has_arith && i < vars.size(); ++i) { + expr* v = vars[i].get(); + has_arith |= is_arith(v); } - if (!has_real) { + if (!has_arith) { return; } + model_evaluator eval(model); + // eval.set_model_completion(true); opt::model_based_opt mbo; obj_map tids; m_trail.reset(); unsigned j = 0; for (unsigned i = 0; i < fmls.size(); ++i) { - if (!linearize(mbo, model, fmls[i].get(), fmls, tids)) { + expr* fml = fmls[i].get(); + if (!linearize(mbo, eval, fml, fmls, tids)) { if (i != j) { fmls[j] = fmls[i].get(); } ++j; } + else { + TRACE("qe", tout << mk_pp(fml, m) << "\n";); + } } fmls.resize(j); @@ -1036,8 +358,13 @@ namespace qe { for (unsigned i = 0; i < vars.size(); ++i) { app* v = vars[i].get(); var_mark.mark(v); - if (a.is_real(v) && !tids.contains(v)) { - tids.insert(v, tids.size()); + if (is_arith(v) && !tids.contains(v)) { + expr_ref val(m); + rational r; + eval(v, val); + a.is_numeral(val, r); + TRACE("qe", tout << mk_pp(v, m) << " " << val << "\n";); + tids.insert(v, mbo.add_var(r, a.is_int(v))); } } for (unsigned i = 0; i < fmls.size(); ++i) { @@ -1056,7 +383,7 @@ namespace qe { unsigned_vector real_vars; for (unsigned i = 0; i < vars.size(); ++i) { app* v = vars[i].get(); - if (a.is_real(v) && !fmls_mark.is_marked(v)) { + if (is_arith(v) && !fmls_mark.is_marked(v)) { real_vars.push_back(tids.find(v)); } else { @@ -1067,7 +394,14 @@ namespace qe { } } vars.resize(j); + TRACE("qe", tout << "remaining vars: " << vars << "\n"; + for (unsigned i = 0; i < real_vars.size(); ++i) { + unsigned v = real_vars[i]; + tout << "v" << v << " " << mk_pp(index2expr[v], m) << "\n"; + } + mbo.display(tout);); mbo.project(real_vars.size(), real_vars.c_ptr()); + TRACE("qe", mbo.display(tout);); vector rows; mbo.get_live_rows(rows); @@ -1078,17 +412,18 @@ namespace qe { if (r.m_vars.size() == 0) { continue; } - if (r.m_vars.size() == 1 && r.m_vars[0].m_coeff.is_neg()) { + if (r.m_vars.size() == 1 && r.m_vars[0].m_coeff.is_neg() && r.m_type != opt::t_mod) { var const& v = r.m_vars[0]; t = index2expr[v.m_id]; if (!v.m_coeff.is_minus_one()) { - t = a.mk_mul(t, a.mk_numeral(-v.m_coeff, a.is_int(t))); + t = a.mk_mul(a.mk_numeral(-v.m_coeff, a.is_int(t)), t); } s = a.mk_numeral(r.m_coeff, a.is_int(t)); switch (r.m_type) { case opt::t_lt: t = a.mk_gt(t, s); break; case opt::t_le: t = a.mk_ge(t, s); break; case opt::t_eq: t = a.mk_eq(t, s); break; + default: UNREACHABLE(); } fmls.push_back(t); VERIFY(model.eval(t, val)); @@ -1099,7 +434,7 @@ namespace qe { var const& v = r.m_vars[j]; t = index2expr[v.m_id]; if (!v.m_coeff.is_one()) { - t = a.mk_mul(t, a.mk_numeral(v.m_coeff, a.is_int(t))); + t = a.mk_mul(a.mk_numeral(v.m_coeff, a.is_int(t)), t); } ts.push_back(t); } @@ -1114,10 +449,17 @@ namespace qe { case opt::t_lt: t = a.mk_lt(t, s); break; case opt::t_le: t = a.mk_le(t, s); break; case opt::t_eq: t = a.mk_eq(t, s); break; + case opt::t_mod: { + rational sval; + if (!a.is_numeral(s, sval) || !sval.is_zero()) { + t = a.mk_sub(t, s); + } + t = a.mk_eq(a.mk_mod(t, a.mk_numeral(r.m_mod, true)), a.mk_int(0)); + break; + } } fmls.push_back(t); - - + VERIFY(model.eval(t, val)); CTRACE("qe", !m.is_true(val), tout << "Evaluated " << t << " to " << val << "\n";); @@ -1133,18 +475,18 @@ namespace qe { opt::inf_eps value; obj_map ts; obj_map tids; - + model_evaluator eval(mdl); // extract objective function. vars coeffs; rational c(0), mul(1); - linearize(mbo, mdl, mul, t, c, fmls, ts, tids); - extract_coefficients(mbo, mdl, ts, tids, coeffs); + linearize(mbo, eval, mul, t, c, fmls, ts, tids); + extract_coefficients(mbo, eval, ts, tids, coeffs); mbo.set_objective(coeffs, c); // extract linear constraints for (unsigned i = 0; i < fmls.size(); ++i) { - linearize(mbo, mdl, fmls[i].get(), fmls, tids); + linearize(mbo, eval, fmls[i].get(), fmls, tids); } // find optimal value @@ -1168,7 +510,7 @@ namespace qe { } expr_ref val(a.mk_numeral(value.get_rational(), false), m); expr_ref tval(m); - VERIFY (mdl.eval(t, tval)); + eval(t, tval); // update the predicate 'bound' which forces larger values when 'strict' is true. // strict: bound := valuue < t @@ -1202,25 +544,31 @@ namespace qe { return valid; } - void extract_coefficients(opt::model_based_opt& mbo, model& model, obj_map const& ts, obj_map& tids, vars& coeffs) { + void extract_coefficients(opt::model_based_opt& mbo, model_evaluator& eval, obj_map const& ts, obj_map& tids, vars& coeffs) { coeffs.reset(); + eval.set_model_completion(true); obj_map::iterator it = ts.begin(), end = ts.end(); for (; it != end; ++it) { unsigned id; - if (!tids.find(it->m_key, id)) { + expr* v = it->m_key; + if (!tids.find(v, id)) { rational r; expr_ref val(m); - if (model.eval(it->m_key, val) && a.is_numeral(val, r)) { - id = mbo.add_var(r); + eval(v, val); + if (a.is_numeral(val, r)) { + id = mbo.add_var(r, a.is_int(v)); } else { TRACE("qe", tout << "extraction of coefficients cancelled\n";); return; } - tids.insert(it->m_key, id); - m_trail.push_back(it->m_key); + tids.insert(v, id); + m_trail.push_back(v); + } + CTRACE("qe", it->m_value.is_zero(), tout << mk_pp(v, m) << " has coefficeint 0\n";); + if (!it->m_value.is_zero()) { + coeffs.push_back(var(id, it->m_value)); } - coeffs.push_back(var(id, it->m_value)); } } diff --git a/src/qe/qe_mbp.cpp b/src/qe/qe_mbp.cpp index 85371d8ef..eaedc470b 100644 --- a/src/qe/qe_mbp.cpp +++ b/src/qe/qe_mbp.cpp @@ -270,6 +270,8 @@ class mbp::impl { sub(fmls[i].get(), val); m_rw(val); if (!m.is_true(val)) { + TRACE("qe", tout << mk_pp(fmls[i].get(), m) << " -> " << val << "\n";); + fmls[i] = val; if (j != i) { fmls[j] = fmls[i].get(); } diff --git a/src/test/model_based_opt.cpp b/src/test/model_based_opt.cpp index b5b8e2953..64b2df285 100644 --- a/src/test/model_based_opt.cpp +++ b/src/test/model_based_opt.cpp @@ -287,9 +287,76 @@ static void test7() { mbo.display(std::cout); } +static void test8() { + opt::model_based_opt mbo; + unsigned x0 = mbo.add_var(rational(2)); + unsigned x = mbo.add_var(rational(1)); + unsigned y = mbo.add_var(rational(3)); + unsigned z = mbo.add_var(rational(4)); + unsigned u = mbo.add_var(rational(5)); + unsigned v = mbo.add_var(rational(6)); + unsigned w = mbo.add_var(rational(6)); + + add_ineq(mbo, x0, 1, y, -1, 0, opt::t_le); + add_ineq(mbo, x, 1, y, -1, 0, opt::t_lt); + add_ineq(mbo, y, 1, u, -1, 0, opt::t_le); + add_ineq(mbo, y, 1, z, -1, 1, opt::t_le); + add_ineq(mbo, y, 1, v, -1, 1, opt::t_le); + + mbo.display(std::cout); + mbo.project(1, &y); + mbo.display(std::cout); +} + + +static void test9() { + opt::model_based_opt mbo; + unsigned x0 = mbo.add_var(rational(2), true); + unsigned x = mbo.add_var(rational(1), true); + unsigned y = mbo.add_var(rational(3), true); + unsigned z = mbo.add_var(rational(4), true); + unsigned u = mbo.add_var(rational(5), true); + unsigned v = mbo.add_var(rational(6), true); + + add_ineq(mbo, x0, 1, y, -1, 0, opt::t_le); + add_ineq(mbo, x, 1, y, -1, 0, opt::t_lt); + add_ineq(mbo, y, 1, u, -1, 0, opt::t_le); + add_ineq(mbo, y, 1, z, -1, 1, opt::t_le); + add_ineq(mbo, y, 1, v, -1, 1, opt::t_le); + + mbo.display(std::cout); + mbo.project(1, &y); + mbo.display(std::cout); +} + + +static void test10() { + opt::model_based_opt mbo; + unsigned x0 = mbo.add_var(rational(2), true); + unsigned x = mbo.add_var(rational(1), true); + unsigned y = mbo.add_var(rational(3), true); + unsigned z = mbo.add_var(rational(4), true); + unsigned u = mbo.add_var(rational(5), true); + unsigned v = mbo.add_var(rational(6), true); + + add_ineq(mbo, x0, 1, y, -2, 0, opt::t_le); + add_ineq(mbo, x, 1, y, -2, 0, opt::t_lt); + add_ineq(mbo, y, 3, u, -4, 0, opt::t_le); + add_ineq(mbo, y, 3, z, -5, 1, opt::t_le); + add_ineq(mbo, y, 3, v, -6, 1, opt::t_le); + + mbo.display(std::cout); + mbo.project(1, &y); + mbo.display(std::cout); + mbo.project(1, &x0); + mbo.display(std::cout); +} + // test with mix of upper and lower bounds void tst_model_based_opt() { + test10(); + return; check_random_ineqs(); test1(); test2(); @@ -298,5 +365,6 @@ void tst_model_based_opt() { test5(); test6(); test7(); - + test8(); + test9(); } diff --git a/src/test/qe_arith.cpp b/src/test/qe_arith.cpp index f9cd9aa23..79495e24f 100644 --- a/src/test/qe_arith.cpp +++ b/src/test/qe_arith.cpp @@ -438,10 +438,158 @@ static void check_random_ineqs() { } } +static void test_project() { + ast_manager m; + reg_decl_plugins(m); + qe::arith_project_plugin plugin(m); + arith_util a(m); + app_ref_vector vars(m); + expr_ref_vector lits(m), ds(m); + model mdl(m); + app_ref x(m), y(m), z(m), u(m); + x = m.mk_const(symbol("x"), a.mk_int()); + y = m.mk_const(symbol("y"), a.mk_int()); + z = m.mk_const(symbol("z"), a.mk_int()); + u = m.mk_const(symbol("u"), a.mk_int()); + func_decl_ref f(m); + sort* int_sort = a.mk_int(); + f = m.mk_func_decl(symbol("f"), 1, &int_sort, int_sort); + // test non-projection + mdl.register_decl(x->get_decl(), a.mk_int(0)); + mdl.register_decl(y->get_decl(), a.mk_int(1)); + mdl.register_decl(z->get_decl(), a.mk_int(2)); + mdl.register_decl(u->get_decl(), a.mk_int(3)); + func_interp* fi = alloc(func_interp, m, 1); + expr_ref_vector nums(m); + nums.push_back(a.mk_int(0)); + nums.push_back(a.mk_int(1)); + nums.push_back(a.mk_int(2)); + fi->insert_new_entry(nums.c_ptr(), a.mk_int(1)); + fi->insert_new_entry(nums.c_ptr()+1, a.mk_int(2)); + fi->insert_new_entry(nums.c_ptr()+2, a.mk_int(3)); + fi->set_else(a.mk_int(10)); + mdl.register_decl(f, fi); + vars.reset(); + lits.reset(); + vars.push_back(x); + lits.push_back(x <= app_ref(m.mk_app(f, (expr*)x), m)); + lits.push_back(x < y); + plugin(mdl, vars, lits); + std::cout << lits << "\n"; + + // test not-equals + vars.reset(); + lits.reset(); + vars.push_back(x); + lits.push_back(m.mk_not(m.mk_eq(x, y))); + plugin(mdl, vars, lits); + std::cout << lits << "\n"; + + // test negation of distinct using bound variables + mdl.register_decl(x->get_decl(), a.mk_int(0)); + mdl.register_decl(y->get_decl(), a.mk_int(1)); + mdl.register_decl(z->get_decl(), a.mk_int(0)); + mdl.register_decl(u->get_decl(), a.mk_int(6)); + vars.reset(); + lits.reset(); + ds.reset(); + vars.push_back(x); + vars.push_back(y); + ds.push_back(x); + ds.push_back(y); + ds.push_back(z + 2); + ds.push_back(u); + ds.push_back(z); + lits.push_back(m.mk_not(m.mk_distinct(ds.size(), ds.c_ptr()))); + plugin(mdl, vars, lits); + std::cout << lits << "\n"; + + // test negation of distinct, not using bound variables + mdl.register_decl(x->get_decl(), a.mk_int(0)); + mdl.register_decl(y->get_decl(), a.mk_int(1)); + mdl.register_decl(z->get_decl(), a.mk_int(0)); + mdl.register_decl(u->get_decl(), a.mk_int(6)); + vars.reset(); + lits.reset(); + ds.reset(); + vars.push_back(x); + vars.push_back(y); + ds.push_back(x); + ds.push_back(y); + ds.push_back(z + 2); + ds.push_back(u); + ds.push_back(z + 10); + ds.push_back(u + 4); + lits.push_back(m.mk_not(m.mk_distinct(ds.size(), ds.c_ptr()))); + plugin(mdl, vars, lits); + std::cout << lits << "\n"; + + + // test distinct + mdl.register_decl(x->get_decl(), a.mk_int(0)); + mdl.register_decl(y->get_decl(), a.mk_int(1)); + mdl.register_decl(z->get_decl(), a.mk_int(0)); + mdl.register_decl(u->get_decl(), a.mk_int(6)); + vars.reset(); + lits.reset(); + ds.reset(); + vars.push_back(x); + vars.push_back(y); + ds.push_back(x); + ds.push_back(y); + ds.push_back(z + 2); + ds.push_back(u); + lits.push_back(m.mk_distinct(ds.size(), ds.c_ptr())); + plugin(mdl, vars, lits); + std::cout << lits << "\n"; + + // equality over modulus + mdl.register_decl(y->get_decl(), a.mk_int(4)); + mdl.register_decl(z->get_decl(), a.mk_int(8)); + lits.reset(); + vars.reset(); + vars.push_back(y); + lits.push_back(m.mk_eq(a.mk_mod(y, a.mk_int(3)), a.mk_int(1))); + lits.push_back(m.mk_eq(2*y, z)); + plugin(mdl, vars, lits); + std::cout << lits << "\n"; + + // inequality test + mdl.register_decl(x->get_decl(), a.mk_int(0)); + mdl.register_decl(y->get_decl(), a.mk_int(1)); + mdl.register_decl(z->get_decl(), a.mk_int(0)); + mdl.register_decl(u->get_decl(), a.mk_int(6)); + vars.reset(); + lits.reset(); + vars.push_back(x); + vars.push_back(y); + lits.push_back(z <= (x + (2*y))); + lits.push_back(2*x < u + 3); + lits.push_back(2*y <= u); + plugin(mdl, vars, lits); + std::cout << lits << "\n"; + + // non-unit equalities + mdl.register_decl(x->get_decl(), a.mk_int(1)); + mdl.register_decl(z->get_decl(), a.mk_int(2)); + mdl.register_decl(u->get_decl(), a.mk_int(3)); + mdl.register_decl(y->get_decl(), a.mk_int(4)); + lits.reset(); + vars.reset(); + vars.push_back(x); + lits.push_back(m.mk_eq(2*x, z)); + lits.push_back(m.mk_eq(3*x, u)); + plugin(mdl, vars, lits); + std::cout << lits << "\n"; + + +} void tst_qe_arith() { + test_project(); + return; check_random_ineqs(); return; // enable_trace("qe");