z3-z3-4.13.0.src.smt.theory_arith_aux.h Maven / Gradle / Ivy
The newest version!
/*++
Copyright (c) 2006 Microsoft Corporation
Module Name:
theory_arith_aux.h
Abstract:
Author:
Leonardo de Moura (leonardo) 2008-04-29.
Revision History:
--*/
#pragma once
#include "util/inf_eps_rational.h"
#include "smt/theory_arith.h"
#include "smt/smt_farkas_util.h"
#include "ast/rewriter/th_rewriter.h"
#include "ast/converters/generic_model_converter.h"
namespace smt {
// -----------------------------------
//
// Rows
//
// -----------------------------------
template
theory_arith::row::row():
m_size(0),
m_base_var(null_theory_var),
m_first_free_idx(-1) {
}
template
void theory_arith::row::reset() {
m_entries.reset();
m_size = 0;
m_base_var = -1;
m_first_free_idx = -1;
}
/**
\brief Add a new row_entry. The result is a reference to the new row_entry. The position of the new row_entry in the
row is stored in pos_idx.
*/
template
typename theory_arith::row_entry & theory_arith::row::add_row_entry(int & pos_idx) {
m_size++;
if (m_first_free_idx == -1) {
pos_idx = m_entries.size();
m_entries.push_back(row_entry());
return m_entries.back();
}
else {
pos_idx = m_first_free_idx;
row_entry & result = m_entries[pos_idx];
SASSERT(result.is_dead());
m_first_free_idx = result.m_next_free_row_entry_idx;
return result;
}
}
/**
\brief Delete row_entry at position idx.
*/
template
void theory_arith::row::del_row_entry(unsigned idx) {
row_entry & t = m_entries[idx];
SASSERT(!t.is_dead());
t.m_next_free_row_entry_idx = m_first_free_idx;
t.m_var = null_theory_var;
m_size--;
SASSERT(t.is_dead());
}
/**
\brief Remove holes (i.e., dead entries) from the row.
*/
template
void theory_arith::row::compress(vector & cols) {
unsigned i = 0;
unsigned j = 0;
unsigned sz = m_entries.size();
for (; i < sz; i++) {
row_entry & t1 = m_entries[i];
if (!t1.is_dead()) {
if (i != j) {
row_entry & t2 = m_entries[j];
t2.m_coeff.swap(t1.m_coeff);
t2.m_var = t1.m_var;
t2.m_col_idx = t1.m_col_idx;
SASSERT(!t2.is_dead());
column & col = cols[t2.m_var];
col[t2.m_col_idx].m_row_idx = j;
}
j++;
}
}
SASSERT(j == m_size);
m_entries.shrink(m_size);
m_first_free_idx = -1;
}
/**
\brief Invoke compress if the row contains too many holes (i.e., dead entries).
*/
template
inline void theory_arith::row::compress_if_needed(vector & cols) {
if (size() * 2 < num_entries()) {
compress(cols);
}
}
/**
\brief Fill the map var -> pos/idx
*/
template
inline void theory_arith::row::save_var_pos(svector & result_map) const {
typename vector::const_iterator it = m_entries.begin();
typename vector::const_iterator end = m_entries.end();
unsigned idx = 0;
for (; it != end; ++it, ++idx) {
if (!it->is_dead()) {
result_map[it->m_var] = idx;
}
}
}
/**
\brief Reset the map var -> pos/idx. That is for all variables v in the row, set result[v] = -1
This method can be viewed as the "inverse" of save_var_pos.
*/
template
inline void theory_arith::row::reset_var_pos(svector & result_map) const {
typename vector::const_iterator it = m_entries.begin();
typename vector::const_iterator end = m_entries.end();
for (; it != end; ++it) {
if (!it->is_dead()) {
result_map[it->m_var] = -1;
}
}
};
#ifdef Z3DEBUG
/**
\brief Return true if the coefficient of v in the row is equals to 'expected'.
*/
template
bool theory_arith::row::is_coeff_of(theory_var v, numeral const & expected) const {
typename vector::const_iterator it = m_entries.begin();
typename vector::const_iterator end = m_entries.end();
for (; it != end; ++it) {
if (!it->is_dead() && it->m_var == v) {
return it->m_coeff == expected;
}
}
return false;
}
#endif
template
void theory_arith::row::display(std::ostream & out) const {
out << "v" << m_base_var << ", ";
typename vector::const_iterator it = m_entries.begin();
typename vector::const_iterator end = m_entries.end();
for (; it != end; ++it) {
if (!it->is_dead()) {
out << it->m_coeff << "*v" << it->m_var << " ";
}
}
out << "\n";
}
template
typename theory_arith::numeral theory_arith::row::get_denominators_lcm() const {
numeral r(1);
TRACE("lcm_bug", tout << "starting get_denominators_lcm...\n";);
typename vector::const_iterator it = m_entries.begin();
typename vector::const_iterator end = m_entries.end();
for (; it != end; ++it) {
if (!it->is_dead()) {
r = lcm(r, denominator(it->m_coeff));
TRACE("lcm_bug", tout << "it->m_coeff: " << it->m_coeff << ", denominator(it->m_coeff): " << denominator(it->m_coeff) << ", r: " << r << "\n";);
}
}
return r;
}
template
int theory_arith::row::get_idx_of(theory_var v) const {
typename vector::const_iterator it = m_entries.begin();
typename vector::const_iterator end = m_entries.end();
for (unsigned idx = 0; it != end; ++it, ++idx) {
if (!it->is_dead() && it->m_var == v)
return idx;
}
return -1;
}
// -----------------------------------
//
// Columns
//
// -----------------------------------
template
void theory_arith::column::reset() {
m_entries.reset();
m_size = 0;
m_first_free_idx = -1;
}
/**
\brief Remove holes (i.e., dead entries) from the column.
*/
template
void theory_arith::column::compress(vector & rows) {
unsigned i = 0;
unsigned j = 0;
unsigned sz = m_entries.size();
for (; i < sz; i++) {
col_entry & e1 = m_entries[i];
if (!e1.is_dead()) {
if (i != j) {
m_entries[j] = e1;
row & r = rows[e1.m_row_id];
r[e1.m_row_idx].m_col_idx = j;
}
j++;
}
}
SASSERT(j == m_size);
m_entries.shrink(m_size);
m_first_free_idx = -1;
}
/**
\brief Invoke compress if the column contains too many holes (i.e., dead entries).
*/
template
inline void theory_arith::column::compress_if_needed(vector & rows) {
if (size() * 2 < num_entries()) {
compress(rows);
}
}
/**
\brief Special version of compress, that is used when the column contain
only one entry located at position singleton_pos.
*/
template
void theory_arith::column::compress_singleton(vector & rows, unsigned singleton_pos) {
SASSERT(m_size == 1);
if (singleton_pos != 0) {
col_entry & s = m_entries[singleton_pos];
m_entries[0] = s;
row & r = rows[s.m_row_id];
r[s.m_row_idx].m_col_idx = 0;
}
m_first_free_idx = -1;
m_entries.shrink(1);
}
template
const typename theory_arith::col_entry * theory_arith::column::get_first_col_entry() const {
typename svector::const_iterator it = m_entries.begin();
typename svector::const_iterator end = m_entries.end();
for (; it != end; ++it) {
if (!it->is_dead()) {
return it;
}
}
return nullptr;
}
template
typename theory_arith::col_entry & theory_arith::column::add_col_entry(int & pos_idx) {
m_size++;
if (m_first_free_idx == -1) {
pos_idx = m_entries.size();
m_entries.push_back(col_entry());
return m_entries.back();
}
else {
pos_idx = m_first_free_idx;
col_entry & result = m_entries[pos_idx];
SASSERT(result.is_dead());
m_first_free_idx = result.m_next_free_row_entry_idx;
return result;
}
}
template
void theory_arith::column::del_col_entry(unsigned idx) {
col_entry & c = m_entries[idx];
SASSERT(!c.is_dead());
c.m_row_id = dead_row_id;
c.m_next_free_row_entry_idx = m_first_free_idx;
m_first_free_idx = idx;
m_size--;
}
// -----------------------------------
//
// Antecedents
//
// -----------------------------------
template
void theory_arith::antecedents_t::init() {
if (!m_init && !empty()) {
m_params.push_back(parameter(symbol("unknown-arith")));
for (unsigned i = 0; i < m_lits.size(); i++) {
m_params.push_back(parameter(m_lit_coeffs[i].to_rational()));
}
for (unsigned i = 0; i < m_eqs.size(); i++) {
m_params.push_back(parameter(m_eq_coeffs[i].to_rational()));
}
m_init = true;
}
}
template
void theory_arith::antecedents_t::reset() {
m_init = false;
m_eq_coeffs.reset();
m_lit_coeffs.reset();
m_eqs.reset();
m_lits.reset();
m_params.reset();
}
template
void theory_arith::antecedents_t::push_lit(literal l, numeral const& r, bool proofs_enabled) {
m_lits.push_back(l);
if (proofs_enabled) {
m_lit_coeffs.push_back(r);
}
}
template
void theory_arith::antecedents_t::push_eq(enode_pair const& p, numeral const& r, bool proofs_enabled) {
m_eqs.push_back(p);
if (proofs_enabled) {
m_eq_coeffs.push_back(r);
}
}
template
parameter * theory_arith::antecedents_t::params(char const* name) {
if (empty()) return nullptr;
init();
m_params[0] = parameter(symbol(name));
return m_params.data();
}
// -----------------------------------
//
// Bounds
//
// -----------------------------------
template
std::ostream& theory_arith::bound::display(theory_arith const& th, std::ostream& out) const {
return out << "v" << get_var() << " " << get_bound_kind() << " " << get_value();
}
// -----------------------------------
//
// Atoms
//
// -----------------------------------
template
theory_arith::atom::atom(bool_var bv, theory_var v, inf_numeral const & k, atom_kind kind):
bound(v, inf_numeral::zero(), B_LOWER, true),
m_bvar(bv),
m_k(k),
m_atom_kind(kind),
m_is_true(false) {
}
template
void theory_arith::atom::assign_eh(bool is_true, inf_numeral const & epsilon) {
m_is_true = is_true;
if (is_true) {
this->m_value = m_k;
this->m_bound_kind = static_cast(m_atom_kind);
SASSERT((this->m_bound_kind == B_LOWER) == (m_atom_kind == A_LOWER));
}
else {
if (get_atom_kind() == A_LOWER) {
this->m_value = m_k;
this->m_value -= epsilon;
this->m_bound_kind = B_UPPER;
}
else {
SASSERT(get_atom_kind() == A_UPPER);
this->m_value = m_k;
this->m_value += epsilon;
this->m_bound_kind = B_LOWER;
}
}
}
template
std::ostream& theory_arith::atom::display(theory_arith const& th, std::ostream& out) const {
literal l(get_bool_var(), !m_is_true);
th.ctx.display_detailed_literal(out, l);
return out;
}
// -----------------------------------
//
// eq_bound
//
// -----------------------------------
template
std::ostream& theory_arith::eq_bound::display(theory_arith const& th, std::ostream& out) const {
ast_manager& m = th.get_manager();
return out << "#" << m_lhs->get_owner_id() << " " << mk_pp(m_lhs->get_expr(), m) << " = "
<< "#" << m_rhs->get_owner_id() << " " << mk_pp(m_rhs->get_expr(), m);
}
// -----------------------------------
//
// Auxiliary methods
//
// -----------------------------------
template
bool theory_arith::at_bound(theory_var v) const {
bound * l = lower(v);
if (l != nullptr && get_value(v) == l->get_value())
return true;
bound * u = upper(v);
return u != nullptr && get_value(v) == u->get_value();
}
template
bool theory_arith::is_fixed(theory_var v) const {
bound * l = lower(v);
if (l == nullptr)
return false;
bound * u = upper(v);
if (u == nullptr)
return false;
return l->get_value() == u->get_value();
}
template
void theory_arith::set_bound(bound * new_bound, bool upper) {
SASSERT(new_bound);
SASSERT(!upper || new_bound->get_bound_kind() == B_UPPER);
SASSERT(upper || new_bound->get_bound_kind() == B_LOWER);
theory_var v = new_bound->get_var();
set_bound_core(v, new_bound, upper);
if ((propagate_eqs() || propagate_diseqs()) && is_fixed(v))
fixed_var_eh(v);
}
/**
\brief Return the col_entry that points to a base row that contains the given variable.
Return 0 if no row contains v.
*/
template
typename theory_arith::col_entry const * theory_arith::get_a_base_row_that_contains(theory_var v) {
while (true) {
column const & c = m_columns[v];
if (c.size() == 0)
return nullptr;
int quasi_base_rid = -1;
typename svector::const_iterator it = c.begin_entries();
typename svector::const_iterator end = c.end_entries();
for (; it != end; ++it) {
if (!it->is_dead()) {
unsigned rid = it->m_row_id;
row & r = m_rows[rid];
theory_var v = r.get_base_var();
if (v == null_theory_var) {
// skip
}
else if (is_base(v)) {
return it;
}
else if (quasi_base_rid == -1)
quasi_base_rid = rid;
}
}
if (quasi_base_rid == -1) {
return nullptr;
}
quasi_base_row2base_row(quasi_base_rid);
// There is no guarantee that v is still a variable of row quasi_base_rid.
// However, this loop will always terminate since I'm creating
// a base row that contains v, or decreasing c.size().
}
}
/**
\brief Return true if all coefficients of the given row are int.
*/
template
bool theory_arith::all_coeff_int(row const & r) const {
typename vector::const_iterator it = r.begin_entries();
typename vector::const_iterator end = r.end_entries();
for (; it != end; ++it) {
if (!it->is_dead() && !it->m_coeff.is_int()) {
TRACE("gomory_cut", display_row(tout, r, true););
return false;
}
}
return true;
}
/**
\brief Return the col_entry that points to row that contains the given variable.
This row should not be owned by an unconstrained quasi-base variable.
Return 0 if failed.
This method is used by move_unconstrained_to_base
*/
template
typename theory_arith::col_entry const * theory_arith::get_row_for_eliminating(theory_var v) const {
column const & c = m_columns[v];
if (c.size() == 0)
return nullptr;
typename svector::const_iterator it = c.begin_entries();
typename svector::const_iterator end = c.end_entries();
for (; it != end; ++it) {
if (!it->is_dead()) {
row const & r = m_rows[it->m_row_id];
theory_var s = r.get_base_var();
if (is_quasi_base(s) && m_var_occs[s].empty())
continue;
if (is_int(v)) {
numeral const & c = r[it->m_row_idx].m_coeff;
// If c == 1 or c == -1, and all other coefficients of r are integer,
// then if we pivot v with the base var of r, we will produce a row
// that will guarantee an integer assignment for v, when the
// non-base vars have integer assignment.
if (!c.is_one() && !c.is_minus_one())
continue;
if (!all_coeff_int(r))
continue;
}
return it;
}
}
return nullptr;
}
template
void theory_arith::move_unconstrained_to_base() {
if (lazy_pivoting_lvl() == 0)
return;
TRACE("move_unconstrained_to_base", tout << "before...\n"; display(tout););
int num = get_num_vars();
for (theory_var v = 0; v < num; v++) {
if (m_var_occs[v].empty() && is_free(v)) {
switch (get_var_kind(v)) {
case QUASI_BASE:
break;
case BASE:
if (is_int(v) && !all_coeff_int(m_rows[get_var_row(v)]))
// If the row contains non integer coefficients, then v may be assigned
// to a non-integer value even if all non-base variables are integer.
// So, v should not be "eliminated"
break;
eliminate(v, m_eager_gcd);
break;
case NON_BASE: {
col_entry const * entry = get_row_for_eliminating(v);
if (entry) {
TRACE("move_unconstrained_to_base", tout << "moving v" << v << " to the base\n";);
row & r = m_rows[entry->m_row_id];
SASSERT(r[entry->m_row_idx].m_var == v);
pivot(r.get_base_var(), v, r[entry->m_row_idx].m_coeff, m_eager_gcd);
SASSERT(is_base(v));
set_var_kind(v, QUASI_BASE);
SASSERT(is_quasi_base(v));
}
break;
} }
}
}
TRACE("move_unconstrained_to_base", tout << "after...\n"; display(tout););
CASSERT("arith", wf_rows());
CASSERT("arith", wf_columns());
CASSERT("arith", valid_row_assignment());
}
/**
\brief Force all quasi_base rows to become base rows.
*/
template
void theory_arith::elim_quasi_base_rows() {
int num = get_num_vars();
for (theory_var v = 0; v < num; v++) {
if (is_quasi_base(v)) {
quasi_base_row2base_row(get_var_row(v));
}
}
}
/**
\brief Remove fixed vars from the base.
*/
template
void theory_arith::remove_fixed_vars_from_base() {
int num = get_num_vars();
for (theory_var v = 0; v < num; v++) {
if (is_base(v) && is_fixed(v)) {
row const & r = m_rows[get_var_row(v)];
typename vector::const_iterator it = r.begin_entries();
typename vector::const_iterator end = r.end_entries();
for (; it != end; ++it) {
if (!it->is_dead() && it->m_var != v && !is_fixed(it->m_var)) {
break;
}
}
if (it != end) {
pivot(v, it->m_var, it->m_coeff, false);
}
}
}
}
/**
\brief Try to minimize the number of rational coefficients.
The idea is to pivot x_i and x_j whenever there is a row
x_i + 1/n * x_j + ... = 0
where
- x_i is a base variable
- x_j is a non-base variables
- x_j is not a fixed variable
- The denominator of any other coefficient a_ik divides n (I only consider the coefficient of non-fixed variables)
remark if there are more than one variable with such properties, we give preference to free variables,
then to variables where upper - lower is maximal.
*/
template
void theory_arith::try_to_minimize_rational_coeffs() {
int num = get_num_vars();
for (theory_var v = 0; v < num; v++) {
if (!is_base(v) || !is_int(v))
continue;
numeral max_den;
row const & r = m_rows[get_var_row(v)];
typename vector::const_iterator it = r.begin_entries();
typename vector::const_iterator end = r.end_entries();
for (; it != end; ++it) {
if (it->is_dead())
continue;
if (it->m_var == v)
continue;
if (is_fixed(it->m_var))
continue;
numeral num = numerator(it->m_coeff);
if (!num.is_one() && !num.is_minus_one())
continue;
numeral den = denominator(it->m_coeff);
if (den > max_den)
max_den = den;
}
if (max_den <= numeral(1))
continue;
// check whether all a_ik denominators divide max_den
it = r.begin_entries();
for (; it != end; ++it) {
if (it->is_dead())
continue;
if (is_fixed(it->m_var))
continue;
numeral den = denominator(it->m_coeff);
if (!(max_den / den).is_int())
break;
}
if (it != end)
continue;
// pick best candidate
theory_var x_j = null_theory_var;
numeral a_ij;
it = r.begin_entries();
for (; it != end; ++it) {
if (it->is_dead())
continue;
if (it->m_var == v)
continue;
if (is_fixed(it->m_var))
continue;
numeral num = numerator(it->m_coeff);
if (!num.is_one() && !num.is_minus_one())
continue;
numeral den = denominator(it->m_coeff);
if (den != max_den)
continue;
if (x_j == null_theory_var ||
// TODO: add extra cases...
is_free(it->m_var) ||
(is_bounded(x_j) && !is_bounded(it->m_var)) ||
(is_bounded(x_j) && is_bounded(it->m_var) && (upper_bound(x_j) - lower_bound(x_j) > upper_bound(it->m_var) - lower_bound(it->m_var)))) {
x_j = it->m_var;
a_ij = it->m_coeff;
if (is_free(x_j))
break;
}
}
if (x_j != null_theory_var)
pivot(v, x_j, a_ij, false);
}
}
// -----------------------------------
//
// Derived bounds
//
// -----------------------------------
template
void theory_arith::derived_bound::push_justification(antecedents& a, numeral const& coeff, bool proofs_enabled) {
TRACE("arith", tout << m_lits << " " << m_eqs.size() << "\n";);
if (proofs_enabled) {
for (literal l : m_lits)
a.push_lit(l, coeff, proofs_enabled);
for (auto const& e : m_eqs)
a.push_eq(e, coeff, proofs_enabled);
}
else {
a.append(m_lits.size(), m_lits.data());
a.append(m_eqs.size(), m_eqs.data());
}
}
template
std::ostream& theory_arith::derived_bound::display(theory_arith const& th, std::ostream& out) const {
ast_manager& m = th.get_manager();
out << "v" << bound::get_var() << " " << bound::get_bound_kind() << " " << bound::get_value() << "\n";
out << "expr: " << mk_pp(th.var2expr(bound::get_var()), m) << "\n";
for (auto const& e : m_eqs) {
enode* a = e.first;
enode* b = e.second;
out << " ";
out << "#" << a->get_owner_id() << " " << mk_pp(a->get_expr(), m) << " = "
<< "#" << b->get_owner_id() << " " << mk_pp(b->get_expr(), m) << "\n";
}
for (literal l : m_lits) {
out << l << ":"; th.ctx.display_detailed_literal(out, l) << "\n";
}
return out;
}
template
void theory_arith::justified_derived_bound::push_justification(antecedents& a, numeral const& coeff, bool proofs_enabled) {
for (unsigned i = 0; i < this->m_lits.size(); ++i) {
a.push_lit(this->m_lits[i], coeff*m_lit_coeffs[i], proofs_enabled);
}
for (unsigned i = 0; i < this->m_eqs.size(); ++i) {
a.push_eq(this->m_eqs[i], coeff*m_eq_coeffs[i], proofs_enabled);
}
}
template
void theory_arith::justified_derived_bound::push_lit(literal l, numeral const& coeff) {
for (unsigned i = 0; i < this->m_lits.size(); ++i) {
if (this->m_lits[i] == l) {
m_lit_coeffs[i] += coeff;
return;
}
}
this->m_lits.push_back(l);
m_lit_coeffs.push_back(coeff);
}
template
void theory_arith::justified_derived_bound::push_eq(enode_pair const& p, numeral const& coeff) {
for (unsigned i = 0; i < this->m_eqs.size(); ++i) {
if (this->m_eqs[i] == p) {
m_eq_coeffs[i] += coeff;
return;
}
}
this->m_eqs.push_back(p);
m_eq_coeffs.push_back(coeff);
}
/**
\brief Copy the justification of b to new_bound. Only literals and equalities not in lits and eqs are copied.
The justification of b is also copied to lits and eqs.
*/
template
void theory_arith::accumulate_justification(bound & b, derived_bound& new_bound, numeral const& coeff, literal_idx_set & lits, eq_set & eqs) {
antecedents ante(*this);
b.push_justification(ante, coeff, proofs_enabled());
unsigned num_lits = ante.lits().size();
for (unsigned i = 0; i < num_lits; ++i) {
literal l = ante.lits()[i];
if (lits.contains(l.index()))
continue;
if (proofs_enabled()) {
new_bound.push_lit(l, ante.lit_coeffs()[i]);
}
else {
new_bound.push_lit(l, numeral::zero());
lits.insert(l.index());
}
}
unsigned num_eqs = ante.eqs().size();
for (unsigned i = 0; i < num_eqs; ++i) {
enode_pair const & p = ante.eqs()[i];
if (eqs.contains(p))
continue;
if (proofs_enabled()) {
new_bound.push_eq(p, ante.eq_coeffs()[i]);
}
else {
new_bound.push_eq(p, numeral::zero());
eqs.insert(p);
}
}
}
template
typename theory_arith::inf_numeral theory_arith::normalize_bound(theory_var v, inf_numeral const & k, bound_kind kind) {
if (is_real(v)) {
return k;
}
if (kind == B_LOWER)
return inf_numeral(ceil(k));
SASSERT(kind == B_UPPER);
return inf_numeral(floor(k));
}
/**
\brief Create a derived bound for v using the given row as an explanation.
*/
template
void theory_arith::mk_bound_from_row(theory_var v, inf_numeral const & k, bound_kind kind, row const & r) {
inf_numeral k_norm = normalize_bound(v, k, kind);
derived_bound * new_bound = proofs_enabled()?alloc(justified_derived_bound, v, k_norm, kind):alloc(derived_bound, v, k_norm, kind);
m_bounds_to_delete.push_back(new_bound);
m_asserted_bounds.push_back(new_bound);
m_tmp_lit_set.reset();
m_tmp_eq_set.reset();
#ifdef Z3DEBUG
inf_numeral val;
#endif
typename vector::const_iterator it = r.begin_entries();
typename vector::const_iterator end = r.end_entries();
for (; it != end; ++it) {
if (!it->is_dead()) {
theory_var v = it->m_var;
bool use_upper = (kind == B_UPPER);
if (!it->m_coeff.is_pos())
use_upper = !use_upper;
TRACE("derived_bound", tout << "using " << (use_upper ? "upper" : "lower") << " bound of v" << v << "\n";);
bound * b = get_bound(v, use_upper);
SASSERT(b);
DEBUG_CODE({
inf_numeral tmp = b->get_value();
tmp *= it->m_coeff;
val += tmp;
});
accumulate_justification(*b, *new_bound, it->m_coeff, m_tmp_lit_set, m_tmp_eq_set);
}
}
TRACE("derived_bound",
tout << "explanation:\n";
for (literal l : new_bound->m_lits) tout << l << " ";
tout << " ";
for (auto const& e : new_bound->m_eqs)
tout << "#" << e.first->get_owner_id() << "=#" << e.second->get_owner_id() << " ";
tout << "\n";);
DEBUG_CODE(CTRACE("derived_bound", k != val, tout << "k: " << k << ", k_norm: " << k_norm << ", val: " << val << "\n";););
SASSERT(k == val);
}
// -----------------------------------
//
// Maximization/Minimization
//
// -----------------------------------
/**
\brief Set: row1 <- row1 + coeff * row2,
where row1 is a temporary row.
\remark Columns do not need to be updated when updating a temporary row.
*/
template
void theory_arith::add_tmp_row(row & r1, numeral const & coeff, row const & r2) {
r1.save_var_pos(m_var_pos);
//
// loop over variables in row2,
// add terms in row2 to row1.
//
#define ADD_TMP_ROW(_SET_COEFF_, _ADD_COEFF_) \
typename vector::const_iterator it = r2.begin_entries(); \
typename vector::const_iterator end = r2.end_entries(); \
for (; it != end; ++it) { \
if (!it->is_dead()) { \
theory_var v = it->m_var; \
int pos = m_var_pos[v]; \
if (pos == -1) { \
/* variable v is not in row1 */ \
int row_idx; \
row_entry & r_entry = r1.add_row_entry(row_idx); \
r_entry.m_var = v; \
_SET_COEFF_; \
} \
else { \
/* variable v is in row1 */ \
row_entry & r_entry = r1[pos]; \
SASSERT(r_entry.m_var == v); \
_ADD_COEFF_; \
if (r_entry.m_coeff.is_zero()) { \
r1.del_row_entry(pos); \
} \
m_var_pos[v] = -1; \
} \
} \
} ((void) 0)
if (coeff.is_one()) {
ADD_TMP_ROW(r_entry.m_coeff = it->m_coeff,
r_entry.m_coeff += it->m_coeff);
}
else if (coeff.is_minus_one()) {
ADD_TMP_ROW(r_entry.m_coeff = it->m_coeff; r_entry.m_coeff.neg(),
r_entry.m_coeff -= it->m_coeff);
}
else {
ADD_TMP_ROW(r_entry.m_coeff = it->m_coeff; r_entry.m_coeff *= coeff,
r_entry.m_coeff += it->m_coeff * coeff);
}
r1.reset_var_pos(m_var_pos);
}
template
bool theory_arith::is_safe_to_leave(theory_var x, bool inc, bool& has_int, bool& shared) {
shared |= ctx.is_shared(get_enode(x));
column & c = m_columns[x];
typename svector::iterator it = c.begin_entries();
typename svector::iterator end = c.end_entries();
has_int = false;
bool unbounded = (inc && !upper(x)) || (!inc && !lower(x));
bool was_unsafe = false;
for (; it != end; ++it) {
if (it->is_dead()) continue;
row const & r = m_rows[it->m_row_id];
theory_var s = r.get_base_var();
numeral const & coeff = r[it->m_row_idx].m_coeff;
if (s != null_theory_var && is_int(s)) has_int = true;
bool is_unsafe = (s != null_theory_var && is_int(s) && !coeff.is_int());
shared |= (s != null_theory_var && ctx.is_shared(get_enode(s)));
was_unsafe |= is_unsafe;
bool inc_s = coeff.is_neg() ? inc : !inc;
unbounded &= !get_bound(s, inc_s);
TRACE("opt", tout << "is v" << x << " safe to leave for v" << s
<< "? " << (is_unsafe?"no":"yes") << " " << (has_int?"int":"real") << " " << (unbounded?"unbounded":"bounded") << "\n";
display_row(tout, r, true););
if (was_unsafe && !unbounded) return false;
}
return !was_unsafe || unbounded;
}
template
bool theory_arith::get_theory_vars(expr * n, uint_set & vars) {
rational r;
expr* x, *y;
if (m_util.is_numeral(n, r)) {
return true;
}
else if (m_util.is_add(n)) {
for (unsigned i = 0; i < to_app(n)->get_num_args(); ++i) {
if (!get_theory_vars(to_app(n)->get_arg(i), vars)) {
return false;
}
}
}
else if (m_util.is_to_real(n, x) || m_util.is_to_int(n, x)) {
return get_theory_vars(x, vars);
}
else if (m_util.is_mul(n, x, y) && m_util.is_numeral(x, r)) {
return get_theory_vars(y, vars);
}
else if (m_util.is_mul(n, y, x) && m_util.is_numeral(x, r)) {
return get_theory_vars(y, vars);
}
else if (!is_app(n)) {
return false;
}
else if (to_app(n)->get_family_id() == m_util.get_family_id()) {
return false;
}
else {
SASSERT(ctx.e_internalized(n));
enode * e = ctx.get_enode(n);
if (is_attached_to_var(e)) {
vars.insert(e->get_th_var(get_id()));
}
return true;
}
return true;
}
//
// add_objective(expr* term) internalizes the arithmetic term and creates
// a row for it if it is not already internalized.
// Then return the variable corresponding to the term.
//
template
theory_var theory_arith::add_objective(app* term) {
theory_var v = internalize_term_core(term);
TRACE("opt", tout << mk_pp(term, get_manager()) << " |-> v" << v << "\n";);
SASSERT(!is_quasi_base(v));
if (!is_linear(get_manager(), term)) {
v = null_theory_var;
}
return v;
}
template
inf_eps_rational theory_arith::value(theory_var v) {
return inf_eps_rational(get_value(v));
}
template
inf_eps_rational theory_arith::maximize(theory_var v, expr_ref& blocker, bool& has_shared) {
TRACE("bound_bug", display_var(tout, v); display(tout););
if (ctx.get_fparams().m_threads > 1)
throw default_exception("multi-threaded optimization is not supported");
has_shared = false;
if (!m_nl_monomials.empty()) {
has_shared = true;
blocker = mk_gt(v);
return inf_eps_rational(get_value(v));
}
max_min_t r = max_min(v, true, true, has_shared);
if (r == UNBOUNDED) {
has_shared = false;
blocker = get_manager().mk_false();
return inf_eps_rational::infinity();
}
else {
blocker = mk_gt(v);
return inf_eps_rational(get_value(v));
}
}
/**
\brief: Create an atom that enforces the inequality v > val
The arithmetical expression encoding the inequality suffices
for the theory of arithmetic.
*/
template
expr_ref theory_arith::mk_gt(theory_var v) {
ast_manager& m = get_manager();
inf_numeral const& val = get_value(v);
expr* obj = get_enode(v)->get_expr();
expr_ref e(m);
rational r = val.get_rational();
if (m_util.is_int(obj->get_sort())) {
if (r.is_int()) {
r += rational::one();
}
else {
r = ceil(r);
}
e = m_util.mk_numeral(r, obj->get_sort());
e = m_util.mk_ge(obj, e);
}
else {
// obj is over the reals.
e = m_util.mk_numeral(r, obj->get_sort());
if (val.get_infinitesimal().is_neg()) {
e = m_util.mk_ge(obj, e);
}
else {
e = m_util.mk_gt(obj, e);
}
}
TRACE("opt", tout << e << "\n";);
return e;
}
/**
\brief create atom that enforces: val <= v
The atom that enforces the inequality is created directly
as opposed to using arithmetical terms.
This allows to handle inequalities with non-standard numbers.
*/
template
expr_ref theory_arith::mk_ge(generic_model_converter& fm, theory_var v, inf_numeral const& val) {
ast_manager& m = get_manager();
std::ostringstream strm;
strm << val << " <= " << mk_pp(get_enode(v)->get_expr(), get_manager());
app* b = m.mk_const(symbol(strm.str()), m.mk_bool_sort());
expr_ref result(b, m);
TRACE("opt", tout << result << "\n";);
if (!ctx.b_internalized(b)) {
fm.hide(b->get_decl());
bool_var bv = ctx.mk_bool_var(b);
ctx.set_var_theory(bv, get_id());
// ctx.set_enode_flag(bv, true);
atom* a = alloc(atom, bv, v, val, A_LOWER);
mk_bound_axioms(a);
m_unassigned_atoms[v]++;
m_var_occs[v].push_back(a);
m_atoms.push_back(a);
insert_bv2a(bv, a);
TRACE("arith", tout << mk_pp(b, m) << "\n";
display_atom(tout, a, false););
}
return result;
}
/**
\brief enable watching bound atom.
*/
template
void theory_arith::enable_record_conflict(expr* bound) {
m_params.m_arith_bound_prop = bound_prop_mode::BP_NONE;
SASSERT(propagation_mode() == bound_prop_mode::BP_NONE); // bound propagation rules are not (yet) handled.
if (bound) {
m_bound_watch = ctx.get_bool_var(bound);
}
else {
m_bound_watch = null_bool_var;
}
m_upper_bound = -inf_eps_rational::infinity();
}
/**
\brief
pos < 0
==
r(Ax <= b) + q(v <= val)
==
val' <= q*v & q*v <= q*val
q*v - val' >= 0
=>
(q*v - val' - q*v)/q >= -v
==
val/q <= v
*/
template
void theory_arith::record_conflict(
unsigned num_lits, literal const * lits,
unsigned num_eqs, enode_pair const * eqs,
unsigned num_params, parameter* params) {
ast_manager& m = get_manager();
expr_ref tmp(m), vq(m);
expr* x, *y, *e;
if (null_bool_var == m_bound_watch) {
return;
}
unsigned idx = num_lits;
for (unsigned i = 0; i < num_lits; ++i) {
if (m_bound_watch == lits[i].var()) {
//SASSERT(!lits[i].sign());
idx = i;
break;
}
}
if (idx == num_lits || num_params == 0) {
return;
}
for (unsigned i = 0; i < num_lits; ++i) {
ctx.literal2expr(lits[i], tmp);
}
for (unsigned i = 0; i < num_eqs; ++i) {
enode_pair const& p = eqs[i];
x = p.first->get_expr();
y = p.second->get_expr();
tmp = m.mk_eq(x,y);
}
SASSERT(num_params == 1 + num_lits + num_eqs);
SASSERT(params[0].is_symbol());
SASSERT(params[0].get_symbol() == symbol("farkas")); // for now, just handle this rule.
farkas_util farkas(m);
rational q;
for (unsigned i = 0; i < num_lits; ++i) {
parameter const& pa = params[i+1];
SASSERT(pa.is_rational());
if (idx == i) {
q = abs(pa.get_rational());
continue;
}
ctx.literal2expr(lits[i], tmp);
if (!farkas.add(abs(pa.get_rational()), to_app(tmp)))
return;
}
for (unsigned i = 0; i < num_eqs; ++i) {
enode_pair const& p = eqs[i];
x = p.first->get_expr();
y = p.second->get_expr();
tmp = m.mk_eq(x,y);
parameter const& pa = params[1 + num_lits + i];
SASSERT(pa.is_rational());
if (!farkas.add(abs(pa.get_rational()), to_app(tmp)))
return;
}
tmp = farkas.get();
if (m.has_trace_stream()) {
log_axiom_instantiation(tmp);
m.trace_stream() << "[end-of-instance]\n";
}
// IF_VERBOSE(1, verbose_stream() << "Farkas result: " << tmp << "\n";);
atom* a = get_bv2a(m_bound_watch);
SASSERT(a);
expr_ref_vector terms(m);
vector mults;
bool strict = false;
if (m_util.is_le(tmp, x, y) || m_util.is_ge(tmp, y, x)) {
}
else if (m.is_not(tmp, e) && (m_util.is_le(e, y, x) || m_util.is_ge(e, x, y))) {
strict = true;
}
else if (m.is_eq(tmp, x, y)) {
}
else {
UNREACHABLE();
}
e = var2expr(a->get_var());
q *= farkas.get_normalize_factor();
SASSERT(!m_util.is_int(e) || q.is_int()); // TBD: not fully handled.
if (q.is_one()) {
vq = e;
}
else {
vq = m_util.mk_mul(m_util.mk_numeral(q, q.is_int()), e);
}
vq = m_util.mk_add(m_util.mk_sub(x, y), vq);
if (!q.is_one()) {
vq = m_util.mk_div(vq, m_util.mk_numeral(q, q.is_int()));
}
th_rewriter rw(m);
rw(vq, tmp);
if (m_util.is_numeral(tmp, q) && m_upper_bound < q) {
m_upper_bound = q;
if (strict) {
m_upper_bound -= get_epsilon(a->get_var());
}
IF_VERBOSE(1, verbose_stream() << "new upper bound: " << m_upper_bound << "\n";);
}
}
/**
\brief find the minimal upper bound on the variable that was last enabled
for conflict recording.
*/
template
inf_eps_rational theory_arith::conflict_minimize() {
return m_upper_bound;
}
/**
\brief Select tightest variable x_i to pivot with x_j. The goal
is to select a x_i such that the value of x_j is increased
(decreased) if inc = true (inc = false), and the tableau
remains feasible. Store the gain in x_j of the pivoting
operation in 'gain'. Note the gain can be too much. That is,
it may make x_i infeasible. In this case, instead of pivoting
we move x_j to its upper bound (lower bound) when inc = true (inc = false).
If no x_i imposes a restriction on x_j, then return null_theory_var.
That is, x_j is free to move to its upper bound (lower bound).
Get the equations for x_j:
x_i1 = coeff_1 * x_j + rest_1
...
x_in = coeff_n * x_j + rest_n
gain_k := (upper_bound(x_ik) - value(x_ik))/coeff_k
*/
template
bool theory_arith::pick_var_to_leave(
theory_var x_j, // non-base variable to increment/decrement
bool inc,
numeral & a_ij, // coefficient of x_i
inf_numeral& min_gain, // minimal required gain on x_j (integral value on integers)
inf_numeral& max_gain, // maximal possible gain on x_j
bool& has_shared, // determine if pivot involves shared variable
theory_var& x_i) { // base variable to pivot with x_j
x_i = null_theory_var;
init_gains(x_j, inc, min_gain, max_gain);
has_shared |= ctx.is_shared(get_enode(x_j));
if (is_int(x_j) && !get_value(x_j).is_int()) {
return false;
}
column & c = m_columns[x_j];
typename svector::iterator it = c.begin_entries();
typename svector::iterator end = c.end_entries();
bool empty_column = true;
for (; it != end; ++it) {
if (it->is_dead()) continue;
empty_column = false;
row const & r = m_rows[it->m_row_id];
theory_var s = r.get_base_var();
numeral const & coeff_ij = r[it->m_row_idx].m_coeff;
if (update_gains(inc, s, coeff_ij, min_gain, max_gain) ||
(x_i == null_theory_var && !unbounded_gain(max_gain))) {
x_i = s;
a_ij = coeff_ij;
}
has_shared |= ctx.is_shared(get_enode(s));
}
TRACE("opt",
tout << (safe_gain(min_gain, max_gain)?"safe":"unsafe") << "\n";
tout << "min gain: " << min_gain;
tout << " max gain: " << max_gain << "\n";
tout << "v" << x_i << " ";
tout << (has_shared?"shared":"not shared") << "\n";);
(void) empty_column;
SASSERT(!safe_gain(min_gain, max_gain) ||
empty_column ||
(unbounded_gain(max_gain) == (x_i == null_theory_var)));
return safe_gain(min_gain, max_gain);
}
template
bool theory_arith::unbounded_gain(inf_numeral const & max_gain) const {
return max_gain.is_minus_one();
}
/*
A gain is 'safe' with respect to the tableau if:
- the selected variable is unbounded and every base variable where it occurs is unbounded
in the direction of the gain. max_gain == -1 is used to indicate unbounded variables.
- the selected variable is a rational (min_gain == -1, max_gain >= 0).
-
*/
template
bool theory_arith::safe_gain(inf_numeral const& min_gain, inf_numeral const & max_gain) const {
return
unbounded_gain(max_gain) ||
min_gain <= max_gain;
}
/**
\brief ensure that maximal gain is divisible by divisor.
*/
template
void theory_arith::normalize_gain(numeral const& divisor, inf_numeral & max_gain) const {
SASSERT(divisor.is_int());
if (!divisor.is_minus_one() && !max_gain.is_minus_one()) {
max_gain = floor(max_gain/divisor)*divisor;
}
}
/**
\brief initialize gains for x_j based on the bounds for x_j.
*/
template
void theory_arith::init_gains(
theory_var x, // non-base variable to increment/decrement
bool inc,
inf_numeral& min_gain, // min value to increment, -1 if rational
inf_numeral& max_gain) { // max value to decrement, -1 if unbounded
min_gain = -inf_numeral::one();
max_gain = -inf_numeral::one();
if (inc && upper(x)) {
max_gain = upper_bound(x) - get_value(x);
}
else if (!inc && lower(x)) {
max_gain = get_value(x) - lower_bound(x);
}
if (is_int(x)) {
min_gain = inf_numeral::one();
}
TRACE("opt",
tout << "v" << x << " := " << get_value(x) << " "
<< "min gain: " << min_gain << " "
<< "max gain: " << max_gain << "\n";);
SASSERT(max_gain.is_minus_one() || !max_gain.is_neg());
SASSERT(min_gain.is_minus_one() || min_gain.is_one());
SASSERT(is_int(x) == min_gain.is_one());
}
template
bool theory_arith::update_gains(
bool inc, // increment/decrement x_j
theory_var x_i, // potential base variable to pivot
numeral const& a_ij, // coefficient of x_j in row where x_i is base.
inf_numeral& min_gain, // min value to increment, -1 if rational
inf_numeral& max_gain) { // max value to decrement, -1 if unbounded
// x_i = row + a_ij*x_j
// a_ij > 0, inc -> decrement x_i
// a_ij < 0, !inc -> decrement x_i
// a_ij denominator
SASSERT(!a_ij.is_zero());
if (!safe_gain(min_gain, max_gain)) return false;
inf_numeral max_inc = inf_numeral::minus_one();
bool decrement_x_i = (inc && a_ij.is_pos()) || (!inc && a_ij.is_neg());
if (decrement_x_i && lower(x_i)) {
max_inc = abs((get_value(x_i) - lower_bound(x_i))/a_ij);
}
else if (!decrement_x_i && upper(x_i)) {
max_inc = abs((upper_bound(x_i) - get_value(x_i))/a_ij);
}
numeral den_aij(1);
bool is_tighter = false;
if (is_int(x_i)) den_aij = denominator(a_ij);
SASSERT(den_aij.is_pos() && den_aij.is_int());
if (is_int(x_i) && !den_aij.is_one()) {
if (min_gain.is_neg()) {
min_gain = inf_numeral(den_aij);
}
else {
min_gain = inf_numeral(lcm(min_gain.get_rational(), den_aij));
}
normalize_gain(min_gain.get_rational(), max_gain);
}
if (is_int(x_i) && !max_gain.is_int()) {
max_gain = inf_numeral(floor(max_gain));
normalize_gain(min_gain.get_rational(), max_gain);
}
if (!max_inc.is_minus_one()) {
if (is_int(x_i)) {
TRACE("opt",
tout << "v" << x_i << " a_ij " << a_ij << " "
<< "min gain: " << min_gain << " "
<< "max gain: " << max_gain << "\n";);
max_inc = floor(max_inc);
normalize_gain(min_gain.get_rational(), max_inc);
}
if (unbounded_gain(max_gain)) {
max_gain = max_inc;
is_tighter = true;
}
else if (max_gain > max_inc) {
max_gain = max_inc;
is_tighter = true;
}
}
TRACE("opt",
tout << "v" << x_i << (is_int(x_i)?" int":" real") << " a_ij " << a_ij << " "
<< "min gain: " << min_gain << " "
<< "max gain: " << max_gain << " tighter: "
<< (is_tighter?"true":"false") << "\n";);
SASSERT(max_gain.is_minus_one() || !max_gain.is_neg());
SASSERT(min_gain.is_minus_one() || !min_gain.is_neg());
//SASSERT(!is_int(x_i) || min_gain.is_pos());
//SASSERT(!is_int(x_i) || min_gain.is_int());
//SASSERT(!is_int(x_i) || max_gain.is_int());
return is_tighter;
}
/**
\brief Check if bound change affects interface equality.
*/
template
bool theory_arith::has_interface_equality(theory_var x) {
theory_var num = get_num_vars();
enode* r = get_enode(x)->get_root();
for (theory_var v = 0; v < num; v++) {
if (v == x) continue;
enode* n = get_enode(v);
if (ctx.is_shared(n) && n->get_root() == r) {
return true;
}
}
return false;
}
/**
\brief Maximize (Minimize) the given temporary row.
Return true if succeeded.
*/
template
typename theory_arith::max_min_t theory_arith::max_min(
row & r,
bool max,
bool maintain_integrality,
bool& has_shared) {
m_stats.m_max_min++;
unsigned best_efforts = 0;
bool inc = false;
SASSERT(!maintain_integrality || valid_assignment());
SASSERT(satisfy_bounds());
numeral a_ij, curr_a_ij, coeff, curr_coeff;
inf_numeral min_gain, max_gain, curr_min_gain, curr_max_gain;
unsigned round = 0;
max_min_t result = OPTIMIZED;
has_shared = false;
unsigned max_efforts = 10 + (ctx.get_random_value() % 20);
while (best_efforts < max_efforts && !ctx.get_cancel_flag()) {
theory_var x_j = null_theory_var;
theory_var x_i = null_theory_var;
bool has_bound = false;
max_gain.reset();
min_gain.reset();
++round;
TRACE("opt", tout << "round: " << round << ", max: " << max << "\n"; display_row(tout, r, true); tout << "state:\n"; display(tout););
typename vector::const_iterator it = r.begin_entries();
typename vector::const_iterator end = r.end_entries();
for (; it != end; ++it) {
if (it->is_dead()) continue;
theory_var curr_x_j = it->m_var;
theory_var curr_x_i = null_theory_var;
SASSERT(is_non_base(curr_x_j));
curr_coeff = it->m_coeff;
bool curr_inc = curr_coeff.is_pos() ? max : !max;
if ((curr_inc && upper(curr_x_j)) || (!curr_inc && lower(curr_x_j))) {
has_bound = true;
}
if ((curr_inc && at_upper(curr_x_j)) || (!curr_inc && at_lower(curr_x_j))) {
// variable cannot be used for max/min.
continue;
}
bool safe_to_leave = pick_var_to_leave(curr_x_j, curr_inc, curr_a_ij,
curr_min_gain, curr_max_gain,
has_shared, curr_x_i);
if (!safe_to_leave) {
TRACE("opt", tout << "no variable picked\n";);
has_bound = true;
best_efforts++;
}
else if (curr_x_i == null_theory_var) {
TRACE("opt", tout << "v" << curr_x_j << " is unrestricted by other variables\n";);
// we can increase/decrease curr_x_j as much as we want.
x_i = null_theory_var; // unbounded
x_j = curr_x_j;
inc = curr_inc;
min_gain = curr_min_gain;
max_gain = curr_max_gain;
break;
}
else if (curr_max_gain > max_gain) {
x_i = curr_x_i;
x_j = curr_x_j;
a_ij = curr_a_ij;
coeff = curr_coeff;
max_gain = curr_max_gain;
min_gain = curr_min_gain;
inc = curr_inc;
}
else if (curr_max_gain.is_zero() && (x_i == null_theory_var || curr_x_i < x_i)) {
x_i = curr_x_i;
x_j = curr_x_j;
a_ij = curr_a_ij;
coeff = curr_coeff;
max_gain = curr_max_gain;
min_gain = curr_min_gain;
inc = curr_inc;
// continue
}
}
TRACE("opt", tout << "after traversing row:\nx_i: v" << x_i << ", x_j: v" << x_j << ", gain: " << max_gain << "\n";
tout << "best efforts: " << best_efforts << " has shared: " << has_shared << "\n";);
if (!has_bound && x_i == null_theory_var && x_j == null_theory_var) {
has_shared = false;
best_efforts = 0;
result = UNBOUNDED;
break;
}
if (x_j == null_theory_var) {
TRACE("opt", tout << "row is " << (max ? "maximized" : "minimized") << "\n";
display_row(tout, r, true););
SASSERT(!maintain_integrality || valid_assignment());
SASSERT(satisfy_bounds());
result = OPTIMIZED;
break;
}
if (min_gain.is_pos() && !min_gain.is_one()) {
++best_efforts;
}
if (x_i == null_theory_var) {
// can increase/decrease x_j as much as we want.
if (inc && upper(x_j)) {
if (max_gain.is_zero()) return BEST_EFFORT;
SASSERT(!unbounded_gain(max_gain));
update_value(x_j, max_gain);
TRACE("opt", tout << "moved v" << x_j << " to upper bound\n";);
SASSERT(!maintain_integrality || valid_assignment());
SASSERT(satisfy_bounds());
continue;
}
if (!inc && lower(x_j)) {
if (max_gain.is_zero()) return BEST_EFFORT;
SASSERT(!unbounded_gain(max_gain));
SASSERT(max_gain.is_pos());
max_gain.neg();
update_value(x_j, max_gain);
TRACE("opt", tout << "moved v" << x_j << " to lower bound\n";);
SASSERT(!maintain_integrality || valid_assignment());
SASSERT(satisfy_bounds());
continue;
}
#if 0
if (ctx.is_shared(get_enode(x_j)) && has_interface_equality(x_j)) {
++best_efforts;
}
else {
SASSERT(unbounded_gain(max_gain));
has_shared = false;
best_efforts = 0;
}
#endif
//
// NB. As it stands this is a possibly unsound conclusion for shared theories.
// the tradeoff is non-termination for unbounded objectives in the
// presence of sharing.
//
has_shared = false;
best_efforts = 0;
result = UNBOUNDED;
break;
}
if (!is_fixed(x_j) && is_bounded(x_j) &&
(upper_bound(x_j) - lower_bound(x_j) == max_gain)) {
// can increase/decrease x_j up to upper/lower bound.
if (inc) {
TRACE("opt", tout << "moved v" << x_j << " to upper bound\n";);
}
else {
max_gain.neg();
TRACE("opt", tout << "moved v" << x_j << " to lower bound\n";);
}
update_value(x_j, max_gain);
SASSERT(!maintain_integrality || valid_assignment());
SASSERT(satisfy_bounds());
continue;
}
TRACE("opt", tout << "max: " << max << ", x_i: v" << x_i << ", x_j: v" << x_j << ", a_ij: " << a_ij << ", coeff: " << coeff << "\n";
if (upper(x_i)) tout << "upper x_i: " << upper_bound(x_i) << " ";
if (lower(x_i)) tout << "lower x_i: " << lower_bound(x_i) << " ";
tout << "value x_i: " << get_value(x_i) << "\n";
if (upper(x_j)) tout << "upper x_j: " << upper_bound(x_j) << " ";
if (lower(x_j)) tout << "lower x_j: " << lower_bound(x_j) << " ";
tout << "value x_j: " << get_value(x_j) << "\n";
);
pivot(x_i, x_j, a_ij, false);
SASSERT(is_non_base(x_i));
SASSERT(is_base(x_j));
bool inc_xi = inc?a_ij.is_neg():a_ij.is_pos();
if (!move_to_bound(x_i, inc_xi, best_efforts, has_shared)) {
TRACE("opt", tout << "can't move bound fully\n";);
// break; // break;
}
row & r2 = m_rows[get_var_row(x_j)];
coeff.neg();
add_tmp_row(r, coeff, r2);
SASSERT(r.get_idx_of(x_j) == -1);
SASSERT(!maintain_integrality || valid_assignment());
SASSERT(satisfy_bounds());
}
TRACE("opt_verbose", display(tout););
return (best_efforts>0 || ctx.get_cancel_flag())?BEST_EFFORT:result;
}
/**
Move the variable x_i maximally towards its bound as long as
bounds of other variables are not violated.
Returns false if an integer bound was truncated and no
progress was made.
*/
template
bool theory_arith::move_to_bound(
theory_var x_i, // variable to move
bool inc, // increment variable or decrement
unsigned& best_efforts, // is bound move a best effort?
bool& has_shared) { // does move include shared variables?
inf_numeral min_gain, max_gain;
if (is_int(x_i) && !get_value(x_i).is_int()) {
++best_efforts;
return false;
}
init_gains(x_i, inc, min_gain, max_gain);
column & c = m_columns[x_i];
typename svector::iterator it = c.begin_entries();
typename svector::iterator end = c.end_entries();
for (; it != end; ++it) {
if (it->is_dead()) continue;
row const & r = m_rows[it->m_row_id];
theory_var s = r.get_base_var();
numeral const & coeff = r[it->m_row_idx].m_coeff;
update_gains(inc, s, coeff, min_gain, max_gain);
has_shared |= ctx.is_shared(get_enode(s));
}
bool result = false;
if (safe_gain(min_gain, max_gain)) {
TRACE("opt", tout << "Safe delta: " << max_gain << "\n";);
SASSERT(!unbounded_gain(max_gain));
if (!inc) {
max_gain.neg();
}
update_value(x_i, max_gain);
if (!min_gain.is_pos() || min_gain.is_one()) {
++best_efforts;
}
result = !max_gain.is_zero();
}
if (!result) {
++best_efforts;
}
return result;
}
/**
\brief Add an entry to a temporary row.
\remark Columns do not need to be updated when updating a temporary row.
*/
template
template
void theory_arith::add_tmp_row_entry(row & r, numeral const & coeff, theory_var v) {
int r_idx;
row_entry & r_entry = r.add_row_entry(r_idx);
r_entry.m_var = v;
r_entry.m_coeff = coeff;
if (invert)
r_entry.m_coeff .neg();
}
/**
\brief Maximize/Minimize the given variable. The bounds of v are update if procedure succeeds.
*/
template
typename theory_arith::max_min_t theory_arith::max_min(theory_var v, bool max, bool maintain_integrality, bool& has_shared) {
expr* e = get_enode(v)->get_expr();
(void)e;
SASSERT(!maintain_integrality || valid_assignment());
SASSERT(satisfy_bounds());
SASSERT(!is_quasi_base(v));
if ((max && at_upper(v)) || (!max && at_lower(v))) {
TRACE("opt", display_var(tout << "At " << (max?"max: ":"min: ") << mk_pp(e, get_manager()) << " \n", v););
return AT_BOUND; // nothing to be done...
}
m_tmp_row.reset();
if (is_non_base(v)) {
add_tmp_row_entry(m_tmp_row, numeral(1), v);
}
else {
row & r = m_rows[get_var_row(v)];
typename vector::const_iterator it = r.begin_entries();
typename vector::const_iterator end = r.end_entries();
for (; it != end; ++it) {
if (!it->is_dead() && it->m_var != v)
add_tmp_row_entry(m_tmp_row, it->m_coeff, it->m_var);
}
}
max_min_t r = max_min(m_tmp_row, max, maintain_integrality, has_shared);
if (r == OPTIMIZED) {
TRACE("opt", tout << mk_pp(e, get_manager()) << " " << (max ? "max" : "min") << " value is: " << get_value(v) << "\n";
display_row(tout, m_tmp_row, true); display_row_info(tout, m_tmp_row););
mk_bound_from_row(v, get_value(v), max ? B_UPPER : B_LOWER, m_tmp_row);
}
else if (r == UNBOUNDED) {
TRACE("opt", display_var(tout << "unbounded: " << mk_pp(e, get_manager()) << "\n", v););
}
else {
TRACE("opt", display_var(tout << "not optimized: " << mk_pp(e, get_manager()) << "\n", v););
}
return r;
}
/**
\brief Maximize & Minimize variables in vars.
Return false if an inconsistency was detected.
*/
template
bool theory_arith::max_min(svector const & vars) {
bool succ = false;
bool has_shared = false;
svector::const_iterator it = vars.begin();
svector::const_iterator end = vars.end();
for (; it != end; ++it) {
if (max_min(*it, true, false, has_shared) == OPTIMIZED && !has_shared)
succ = true;
if (max_min(*it, false, false, has_shared) == OPTIMIZED && !has_shared)
succ = true;
}
if (succ) {
// process new bounds
bool r = propagate_core();
TRACE("opt", tout << "after max/min round:\n"; display(tout););
return r;
}
return true;
}
// -----------------------------------
//
// Freedom intervals
//
// -----------------------------------
/**
\brief See Model-based theory combination paper.
Return false if failed to build the freedom interval.
\remark If x_j is an integer variable, then m will contain the lcm of the denominators of a_ij.
We only consider the a_ij coefficients for x_i
*/
template
bool theory_arith::get_freedom_interval(theory_var x_j, bool & inf_l, inf_numeral & l, bool & inf_u, inf_numeral & u, numeral & m) {
if (is_base(x_j))
return false;
inf_numeral const & x_j_val = get_value(x_j);
column & c = m_columns[x_j];
typename svector::iterator it = c.begin_entries();
typename svector::iterator end = c.end_entries();
inf_l = true;
inf_u = true;
l.reset();
u.reset();
m = numeral(1);
#define IS_FIXED() { if (!inf_l && !inf_u && l == u) goto fi_succeeded; }
#define SET_LOWER(VAL) { inf_numeral const & _VAL = VAL; if (inf_l || _VAL > l) { l = _VAL; inf_l = false; } IS_FIXED(); }
#define SET_UPPER(VAL) { inf_numeral const & _VAL = VAL; if (inf_u || _VAL < u) { u = _VAL; inf_u = false; } IS_FIXED(); }
if (lower(x_j)) {
SET_LOWER(lower_bound(x_j));
}
if (upper(x_j)) {
SET_UPPER(upper_bound(x_j));
}
for (; it != end; ++it) {
if (!it->is_dead()) {
row & r = m_rows[it->m_row_id];
theory_var x_i = r.get_base_var();
if (x_i != null_theory_var && !is_quasi_base(x_i)) {
numeral const & a_ij = r[it->m_row_idx].m_coeff;
inf_numeral const & x_i_val = get_value(x_i);
if (is_int(x_i) && is_int(x_j) && !a_ij.is_int())
m = lcm(m, denominator(a_ij));
bound * x_i_lower = lower(x_i);
bound * x_i_upper = upper(x_i);
if (a_ij.is_neg()) {
if (x_i_lower) {
inf_numeral new_l = x_j_val + ((x_i_val - x_i_lower->get_value()) / a_ij);
SET_LOWER(new_l);
}
if (x_i_upper) {
inf_numeral new_u = x_j_val + ((x_i_val - x_i_upper->get_value()) / a_ij);
SET_UPPER(new_u);
}
}
else {
if (x_i_upper) {
inf_numeral new_l = x_j_val + ((x_i_val - x_i_upper->get_value()) / a_ij);
SET_LOWER(new_l);
}
if (x_i_lower) {
inf_numeral new_u = x_j_val + ((x_i_val - x_i_lower->get_value()) / a_ij);
SET_UPPER(new_u);
}
}
}
}
}
fi_succeeded:
TRACE("freedom_interval",
tout << "freedom variable for:\n";
display_var(tout, x_j);
tout << "[";
if (inf_l) tout << "-oo"; else tout << l;
tout << "; ";
if (inf_u) tout << "oo"; else tout << u;
tout << "]\n";);
return true;
}
// -----------------------------------
//
// Implied eqs
//
// -----------------------------------
/**
\brief Try to check whether v1 == v2 is implied by the current state.
If it is return true.
*/
template
bool theory_arith::try_to_imply_eq(theory_var v1, theory_var v2) {
SASSERT(v1 != v2);
SASSERT(get_value(v1) == get_value(v2));
SASSERT(valid_assignment());
if (is_quasi_base(v1) || is_quasi_base(v2))
return false;
m_tmp_row.reset();
if (is_non_base(v1)) {
add_tmp_row_entry(m_tmp_row, numeral(1), v1);
}
else {
row & r = m_rows[get_var_row(v1)];
typename vector::const_iterator it = r.begin_entries();
typename vector::const_iterator end = r.end_entries();
for (; it != end; ++it) {
if (!it->is_dead() && it->m_var != v1)
add_tmp_row_entry(m_tmp_row, it->m_coeff, it->m_var);
}
}
m_tmp_row.save_var_pos(m_var_pos);
#define ADD_ENTRY(COEFF, VAR) { \
int pos = m_var_pos[VAR]; \
if (pos == -1) { \
add_tmp_row_entry(m_tmp_row, COEFF, VAR); \
} \
else { \
row_entry & r_entry = m_tmp_row[pos]; \
SASSERT(r_entry.m_var == VAR); \
r_entry.m_coeff += COEFF; \
if (r_entry.m_coeff.is_zero()) \
m_tmp_row.del_row_entry(pos); \
m_var_pos[VAR] = -1; \
} \
}
if (is_non_base(v2)) {
ADD_ENTRY(numeral(-1), v2);
}
else {
row & r = m_rows[get_var_row(v2)];
typename vector::const_iterator it = r.begin_entries();
typename vector::const_iterator end = r.end_entries();
for (; it != end; ++it) {
if (!it->is_dead() && it->m_var != v2) {
numeral c = it->m_coeff;
c.neg();
ADD_ENTRY(c, it->m_var);
}
}
}
m_tmp_row.reset_var_pos(m_var_pos);
SASSERT(m_tmp_row.size() > 0);
#if 0
TRACE("imply_eq", display_row_info(tout, m_tmp_row););
m_tmp_acc_lits.reset();
m_tmp_acc_eqs.reset();
m_tmp_lit_set.reset();
m_tmp_eq_set.reset();
if ((OPTIMIZED == max_min(m_tmp_row, true)) &&
is_zero_row(m_tmp_row, true, m_tmp_acc_lits, m_tmp_acc_eqs, m_tmp_lit_set, m_tmp_eq_set) &&
(OPTIMIZED == max_min(m_tmp_row, false)) &&
is_zero_row(m_tmp_row, false, m_tmp_acc_lits, m_tmp_acc_eqs, m_tmp_lit_set, m_tmp_eq_set)) {
// v1 == v2
TRACE("imply_eq", tout << "found new implied equality:\n";
display_var(tout, v1); display_var(tout, v2););
// TODO: assert implied equality
// return true;
}
#endif
return false;
}
// -----------------------------------
//
// Assume eqs
//
// The revamped assume eqs try to perturbate the
// current assignment using pivoting operations.
//
// -----------------------------------
#define RANGE 10000
/**
\brief Performs a random update on v using its freedom interval.
Return true if it was possible to change.
*/
template
bool theory_arith::random_update(theory_var v) {
if (is_fixed(v) || !is_non_base(v))
return false;
bool inf_l, inf_u;
inf_numeral l, u;
numeral m;
get_freedom_interval(v, inf_l, l, inf_u, u, m);
if (inf_l && inf_u) {
inf_numeral new_val = inf_numeral(m_random() % (RANGE + 1));
set_value(v, new_val);
return true;
}
if (is_int(v)) {
if (!inf_l) {
l = ceil(l);
if (!m.is_one())
l = m*ceil(l/m);
}
if (!inf_u) {
u = floor(u);
if (!m.is_one())
u = m*floor(u/m);
}
}
if (!inf_l && !inf_u && l >= u)
return false;
if (inf_u) {
SASSERT(!inf_l);
inf_numeral delta = inf_numeral(m_random() % (RANGE + 1));
inf_numeral new_val = l + m*delta;
set_value(v, new_val);
return true;
}
if (inf_l) {
SASSERT(!inf_u);
inf_numeral delta = inf_numeral(m_random() % (RANGE + 1));
inf_numeral new_val = u - m*delta;
set_value(v, new_val);
return true;
}
if (!is_int(v)) {
SASSERT(!inf_l && !inf_u);
numeral delta = numeral(m_random() % (RANGE + 1));
inf_numeral new_val = l + ((delta * (u - l)) / numeral(RANGE));
set_value(v, new_val);
return true;
}
else {
unsigned range = RANGE;
numeral r = (u.get_rational() - l.get_rational()) / m;
if (r < numeral(RANGE))
range = static_cast(r.get_uint64());
inf_numeral new_val = l + m * (inf_numeral(m_random() % (range + 1)));
set_value(v, new_val);
return true;
}
}
template
void theory_arith::mutate_assignment() {
SASSERT(m_to_patch.empty());
remove_fixed_vars_from_base();
int num_vars = get_num_vars();
m_var_value_table.reset();
m_tmp_var_set.reset();
sbuffer candidates;
for (theory_var v = 0; v < num_vars; v++) {
enode * n1 = get_enode(v);
if (!is_relevant_and_shared(n1))
continue;
theory_var other = m_var_value_table.insert_if_not_there(v);
if (other == v)
continue; // first node with the given value...
enode * n2 = get_enode(other);
if (n1->get_root() == n2->get_root())
continue;
if (!is_fixed(v)) {
candidates.push_back(v);
}
else if (!is_fixed(other) && !m_tmp_var_set.contains(other)) {
m_tmp_var_set.insert(other);
candidates.push_back(other);
}
}
TRACE("arith_rand", tout << "candidates.size() == " << candidates.size() << "\n";);
if (candidates.empty())
return;
m_tmp_var_set.reset();
m_tmp_var_set2.reset();
for (theory_var v : candidates) {
SASSERT(!is_fixed(v));
if (is_base(v)) {
row & r = m_rows[get_var_row(v)];
typename vector::const_iterator it2 = r.begin_entries();
typename vector::const_iterator end2 = r.end_entries();
for (; it2 != end2; ++it2) {
if (!it2->is_dead() && it2->m_var != v && !is_fixed(it2->m_var) && random_update(it2->m_var))
break;
}
}
else {
random_update(v);
}
}
SASSERT(m_to_patch.empty());
}
/**
\brief We must redefine this method, because theory of arithmetic contains
underspecified operators such as division by 0.
(/ a b) is essentially an uninterpreted function when b = 0.
Thus, 'a' must be considered a shared var if it is the child of an underspecified operator.
*/
template
bool theory_arith::is_shared(theory_var v) const {
if (m_underspecified_ops.empty())
return false;
enode * n = get_enode(v);
enode * r = n->get_root();
enode_vector::const_iterator it = r->begin_parents();
enode_vector::const_iterator end = r->end_parents();
TRACE("shared", tout << ctx.get_scope_level() << " " << v << " " << r->get_num_parents() << "\n";);
for (; it != end; ++it) {
enode * parent = *it;
app * o = parent->get_expr();
if (o->get_family_id() == get_id()) {
switch (o->get_decl_kind()) {
case OP_DIV:
case OP_IDIV:
case OP_REM:
case OP_MOD:
return true;
default:
break;
}
}
}
return false;
}
template
bool theory_arith::assume_eqs() {
// See comment in m_liberal_final_check declaration
if (m_liberal_final_check)
mutate_assignment();
TRACE("assume_eq_int", display(tout););
unsigned old_sz = m_assume_eq_candidates.size();
m_var_value_table.reset();
bool result = false;
int num = get_num_vars();
for (theory_var v = 0; v < num; v++) {
enode * n = get_enode(v);
if (!is_relevant_and_shared(n))
continue;
theory_var other = null_theory_var;
other = m_var_value_table.insert_if_not_there(v);
if (other == v)
continue;
enode * n2 = get_enode(other);
if (n->get_root() == n2->get_root())
continue;
m_assume_eq_candidates.push_back({ other , v });
result = true;
}
if (result)
ctx.push_trail(restore_vector(m_assume_eq_candidates, old_sz));
return delayed_assume_eqs();
}
template
bool theory_arith::delayed_assume_eqs() {
if (m_assume_eq_head == m_assume_eq_candidates.size())
return false;
ctx.push_trail(value_trail(m_assume_eq_head));
while (m_assume_eq_head < m_assume_eq_candidates.size()) {
auto const& [v1, v2] = m_assume_eq_candidates[m_assume_eq_head];
enode* n1 = get_enode(v1);
enode* n2 = get_enode(v2);
m_assume_eq_head++;
CTRACE("arith",
get_value(v1) == get_value(v2) && n1->get_root() != n2->get_root(),
tout << "assuming eq: " << ctx.pp(n1) << " = #" << ctx.pp(n2) << "\n";);
if (get_value(v1) == get_value(v2) &&
n1->get_root() != n2->get_root() &&
assume_eq(n1, n2)) {
++m_stats.m_assume_eqs;
return true;
}
}
return false;
}
#if 0
/**
\brief Check if the given row is implying a zero upper/lower bound.
Accumulate the justification in the given vectors.
Return true if it is.
*/
template
bool theory_arith::is_zero_row(row const & r, bool upper, literal_vector & lit_vect, eq_vector & eq_vect, literal_idx_set & lits, eq_set & eqs) {
inf_numeral val;
typename vector