mirrors/z3
3
0
Fork 0
mirror of https://github.com/Z3Prover/z3 synced 2025-04-08 02:15:19 +00:00
Code Activity

Arith min max (#6864)

Browse source 
* prepare for dependencies

Signed-off-by: Nikolaj Bjorner <nbjorner@microsoft.com>

* snapshot

Signed-off-by: Nikolaj Bjorner <nbjorner@microsoft.com>

* more refactoring

Signed-off-by: Nikolaj Bjorner <nbjorner@microsoft.com>

* more refactoring

Signed-off-by: Nikolaj Bjorner <nbjorner@microsoft.com>

* build

Signed-off-by: Nikolaj Bjorner <nbjorner@microsoft.com>

* pass in u_dependency_manager

Signed-off-by: Nikolaj Bjorner <nbjorner@microsoft.com>

* address NYIs

Signed-off-by: Nikolaj Bjorner <nbjorner@microsoft.com>

* more refactoring names

Signed-off-by: Nikolaj Bjorner <nbjorner@microsoft.com>

* eq_explanation update

Signed-off-by: Nikolaj Bjorner <nbjorner@microsoft.com>

* add outline of bounds improvement functionality

Signed-off-by: Nikolaj Bjorner <nbjorner@microsoft.com>

* fix unit tests

Signed-off-by: Nikolaj Bjorner <nbjorner@microsoft.com>

* remove unused structs

Signed-off-by: Nikolaj Bjorner <nbjorner@microsoft.com>

* more bounds

Signed-off-by: Nikolaj Bjorner <nbjorner@microsoft.com>

* more bounds

Signed-off-by: Nikolaj Bjorner <nbjorner@microsoft.com>

* convert more internals to use u_dependency instead of constraint_index

Signed-off-by: Nikolaj Bjorner <nbjorner@microsoft.com>

* convert more internals to use u_dependency instead of constraint_index

Signed-off-by: Nikolaj Bjorner <nbjorner@microsoft.com>

* remember to push/pop scopes

Signed-off-by: Nikolaj Bjorner <nbjorner@microsoft.com>

* use the main function for updating bounds

Signed-off-by: Nikolaj Bjorner <nbjorner@microsoft.com>

* na

Signed-off-by: Nikolaj Bjorner <nbjorner@microsoft.com>

* na

Signed-off-by: Nikolaj Bjorner <nbjorner@microsoft.com>

* remove reset of shared dep manager

Signed-off-by: Nikolaj Bjorner <nbjorner@microsoft.com>

* disable improve-bounds, add statistics

Signed-off-by: Nikolaj Bjorner <nbjorner@microsoft.com>

---------

Signed-off-by: Nikolaj Bjorner <nbjorner@microsoft.com>
This commit is contained in:
Nikolaj Bjorner 2023-08-19 17:44:09 -07:00 committed by GitHub
parent c3b344ec47
commit 5e3df9ee77
No known key found for this signature in database
GPG key ID: 4AEE18F83AFDEB23
40 changed files with 630 additions and 529 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

3
src/ast/simplifiers/bound_simplifier.h
View file

@ -37,6 +37,7 @@ class bound_simplifier : public dependent_expr_simplifier {
unsynch_mpq_manager nm;
small_object_allocator m_alloc;
bound_propagator bp;
u_dependency_manager m_dep_manager;
dep_intervals m_interval;
ptr_vector<expr> m_var2expr;
unsigned_vector m_expr2var;
@ -105,7 +106,7 @@ public:
a(m),
m_rewriter(m),
bp(nm, m_alloc, p),
m_interval(m.limit()),
m_interval(m_dep_manager, m.limit()),
m_trail(m),
m_num_buffer(nm) {
updt_params(p);

5
src/math/grobner/pdd_solver.cpp
View file

@ -59,9 +59,10 @@ namespace dd {
*/
solver::solver(reslimit& lim, pdd_manager& m) :
solver::solver(reslimit& lim, u_dependency_manager& dm, pdd_manager& m) :
m(m),
m_limit(lim)
m_limit(lim),
m_dep_manager(dm)
{}
solver::~solver() {

5
src/math/grobner/pdd_solver.h
View file

@ -112,6 +112,7 @@ private:
pdd_manager& m;
reslimit& m_limit;
u_dependency_manager& m_dep_manager;
stats m_stats;
config m_config;
print_dep_t m_print_dep;
@ -119,12 +120,11 @@ private:
equation_vector m_processed;
equation_vector m_to_simplify;
vector<std::tuple<unsigned, pdd, u_dependency*>> m_subst;
mutable u_dependency_manager m_dep_manager;
equation_vector m_all_eqs;
equation* m_conflict = nullptr;
bool m_too_complex;
public:
solver(reslimit& lim, pdd_manager& m);
solver(reslimit& lim, u_dependency_manager& dm, pdd_manager& m);
~solver();
pdd_manager& get_manager() { return m; }
@ -144,7 +144,6 @@ public:
void saturate();
equation_vector const& equations();
u_dependency_manager& dep() const { return m_dep_manager; }
void collect_statistics(statistics & st) const;
std::ostream& display(std::ostream& out, const equation& eq) const;

17
src/math/interval/dep_intervals.h
View file

@ -27,6 +27,7 @@
#include "math/interval/interval.h"
class dep_intervals {
public:
enum with_deps_t { with_deps, without_deps };
@ -142,8 +143,9 @@ private:
public:
typedef interval_manager<im_config>::interval interval;
u_dependency_manager& m_dep_manager;
mutable unsynch_mpq_manager m_num_manager;
mutable u_dependency_manager m_dep_manager;
im_config m_config;
mutable interval_manager<im_config> m_imanager;
@ -158,9 +160,10 @@ public:
public:
u_dependency_manager& dep_manager() { return m_dep_manager; }
dep_intervals(reslimit& lim) :
m_config(m_num_manager, m_dep_manager),
m_imanager(lim, im_config(m_num_manager, m_dep_manager))
dep_intervals(u_dependency_manager& dm, reslimit& lim) :
m_dep_manager(dm),
m_config(m_num_manager, dm),
m_imanager(lim, im_config(m_num_manager, dm))
{}
std::ostream& display(std::ostream& out, const interval& i) const;
@ -335,15 +338,17 @@ public:
bool is_empty(interval const& a) const;
void set_interval_for_scalar(interval&, const rational&);
template <typename T>
void linearize(u_dependency* dep, T& expl) const {
vector<unsigned, false> v;
m_dep_manager.linearize(dep, v);
for (unsigned ci: v)
for (auto ci: v)
expl.push_back(ci);
}
void reset() { m_dep_manager.reset(); }
void reset() { }
void del(interval& i) { m_imanager.del(i); }

7
src/math/lp/emonics.cpp
View file

@ -517,11 +517,10 @@ bool emonics::invariant() const {
TRACE("nla_solver_mons", display(tout););
// the variable index contains exactly the active monomials
unsigned mons = 0;
for (lpvar v = 0; v < m_var2index.size(); v++) {
if (is_monic_var(v)) {
for (lpvar v = 0; v < m_var2index.size(); v++)
if (is_monic_var(v))
mons++;
}
}
if (m_monics.size() != mons) {
TRACE("nla_solver_mons", tout << "missmatch of monic vars\n";);
return false;

20
src/math/lp/gomory.cpp
View file

@ -50,9 +50,13 @@ class create_cut {
bool at_upper(unsigned j) const { return lia.at_upper(j); }
const impq & lower_bound(unsigned j) const { return lia.lower_bound(j); }
const impq & upper_bound(unsigned j) const { return lia.upper_bound(j); }
constraint_index column_lower_bound_constraint(unsigned j) const { return lia.column_lower_bound_constraint(j); }
constraint_index column_upper_bound_constraint(unsigned j) const { return lia.column_upper_bound_constraint(j); }
u_dependency* column_lower_bound_constraint(unsigned j) const { return lia.column_lower_bound_constraint(j); }
u_dependency* column_upper_bound_constraint(unsigned j) const { return lia.column_upper_bound_constraint(j); }
bool column_is_fixed(unsigned j) const { return lia.lra.column_is_fixed(j); }
void push_explanation(u_dependency* d) {
for (auto ci : lia.lra.flatten(d))
m_ex->push_back(ci);
}
void int_case_in_gomory_cut(unsigned j) {
lp_assert(is_int(j) && m_fj.is_pos());
@ -67,7 +71,7 @@ class create_cut {
new_a = m_fj <= m_one_minus_f ? m_fj / m_one_minus_f : ((1 - m_fj) / m_f);
lp_assert(new_a.is_pos());
m_k.addmul(new_a, lower_bound(j).x);
m_ex->push_back(column_lower_bound_constraint(j));
push_explanation(column_lower_bound_constraint(j));
}
else {
lp_assert(at_upper(j));
@ -75,7 +79,7 @@ class create_cut {
new_a = - (m_fj <= m_f ? m_fj / m_f : ((1 - m_fj) / m_one_minus_f));
lp_assert(new_a.is_neg());
m_k.addmul(new_a, upper_bound(j).x);
m_ex->push_back(column_upper_bound_constraint(j));
push_explanation(column_upper_bound_constraint(j));
}
m_t.add_monomial(new_a, j);
m_lcm_den = lcm(m_lcm_den, denominator(new_a));
@ -99,7 +103,7 @@ class create_cut {
new_a = - a / m_f;
m_k.addmul(new_a, lower_bound(j).x); // is it a faster operation than
// k += lower_bound(j).x * new_a;
m_ex->push_back(column_lower_bound_constraint(j));
push_explanation(column_lower_bound_constraint(j));
}
else {
lp_assert(at_upper(j));
@ -110,7 +114,7 @@ class create_cut {
// the delta is positive works again m_one_minus_f
new_a = a / m_one_minus_f;
m_k.addmul(new_a, upper_bound(j).x); // k += upper_bound(j).x * new_a;
m_ex->push_back(column_upper_bound_constraint(j));
push_explanation(column_upper_bound_constraint(j));
}
m_t.add_monomial(new_a, j);
TRACE("gomory_cut_detail_real", tout << "add " << new_a << "*v" << j << ", k: " << m_k << "\n";
@ -313,8 +317,8 @@ public:
// use -p.coeff() to make the format compatible with the format used in: Integrating Simplex with DPLL(T)
try {
if (lia.is_fixed(j)) {
m_ex->push_back(column_lower_bound_constraint(j));
m_ex->push_back(column_upper_bound_constraint(j));
push_explanation(column_lower_bound_constraint(j));
push_explanation(column_upper_bound_constraint(j));
continue;
}
if (is_real(j))

27
src/math/lp/hnf_cutter.cpp
View file

@ -39,7 +39,7 @@ namespace lp {
m_overflow = false;
}
void hnf_cutter::add_term(const lar_term* t, const mpq &rs, constraint_index ci, bool upper_bound) {
void hnf_cutter::add_term(const lar_term* t, const mpq &rs, u_dependency* ci, bool upper_bound) {
m_terms.push_back(t);
m_terms_upper.push_back(upper_bound);
if (upper_bound)
@ -146,7 +146,7 @@ namespace lp {
}
#endif
void hnf_cutter::shrink_explanation(const svector<unsigned>& basis_rows) {
svector<unsigned> new_expl;
ptr_vector<u_dependency> new_expl;
for (unsigned i : basis_rows) {
new_expl.push_back(m_constraints_for_explanation[i]);
}
@ -232,10 +232,10 @@ branch y_i >= ceil(y0_i) is impossible.
void hnf_cutter::try_add_term_to_A_for_hnf(tv const &i) {
mpq rs;
const lar_term& t = lra.get_term(i);
constraint_index ci;
u_dependency* dep;
bool upper_bound;
if (!is_full() && lra.get_equality_and_right_side_for_term_on_current_x(i, rs, ci, upper_bound)) {
add_term(&t, rs, ci, upper_bound);
if (!is_full() && lra.get_equality_and_right_side_for_term_on_current_x(i, rs, dep, upper_bound)) {
add_term(&t, rs, dep, upper_bound);
}
}
@ -259,8 +259,9 @@ branch y_i >= ceil(y0_i) is impossible.
}
lia.settings().stats().m_hnf_cutter_calls++;
TRACE("hnf_cut", tout << "settings().stats().m_hnf_cutter_calls = " << lia.settings().stats().m_hnf_cutter_calls << "\n";
for (unsigned i : constraints_for_explanation()) {
lra.constraints().display(tout, i);
for (u_dependency* d : constraints_for_explanation()) {
for (auto ci : lra.flatten(d))
lra.constraints().display(tout, ci);
}
tout << lra.constraints();
);
@ -277,16 +278,16 @@ branch y_i >= ceil(y0_i) is impossible.
TRACE("hnf_cut",
lra.print_term(lia.m_t, tout << "cut:");
tout << " <= " << lia.m_k << std::endl;
for (unsigned i : constraints_for_explanation()) {
lra.constraints().display(tout, i);
}
for (auto* dep : constraints_for_explanation())
for (auto ci : lra.flatten(dep))
lra.constraints().display(tout, ci);
);
lp_assert(lia.current_solution_is_inf_on_cut());
lia.settings().stats().m_hnf_cuts++;
lia.m_ex->clear();
for (unsigned i : constraints_for_explanation()) {
lia.m_ex->push_back(i);
}
for (u_dependency* dep : constraints_for_explanation())
for (auto ci : lia.lra.flatten(dep))
lia.m_ex->push_back(ci);
}
return r;
}

6
src/math/lp/hnf_cutter.h
View file

@ -34,7 +34,7 @@ class hnf_cutter {
general_matrix m_A;
vector<const lar_term*> m_terms;
vector<bool> m_terms_upper;
svector<constraint_index> m_constraints_for_explanation;
ptr_vector<u_dependency> m_constraints_for_explanation;
vector<mpq> m_right_sides;
mpq m_abs_max;
bool m_overflow;
@ -55,13 +55,13 @@ private:
unsigned terms_count() const { return m_terms.size(); }
const mpq & abs_max() const { return m_abs_max; }
const vector<const lar_term*>& terms() const { return m_terms; }
const svector<unsigned>& constraints_for_explanation() const { return m_constraints_for_explanation; }
const ptr_vector<u_dependency>& constraints_for_explanation() const { return m_constraints_for_explanation; }
const vector<mpq> & right_sides() const { return m_right_sides; }
bool is_full() const;
void clear();
void add_term(const lar_term* t, const mpq &rs, constraint_index ci, bool upper_bound);
void add_term(const lar_term* t, const mpq &rs, u_dependency* ci, bool upper_bound);
void print(std::ostream & out);

2
src/math/lp/horner.cpp
View file

@ -103,7 +103,7 @@ bool horner::horner_lemmas() {
return false;
}
c().lp_settings().stats().m_horner_calls++;
const auto& matrix = c().m_lar_solver.A_r();
const auto& matrix = c().lra.A_r();
// choose only rows that depend on m_to_refine variables
std::set<unsigned> rows_to_check;
for (lpvar j : c().m_to_refine) {

7
src/math/lp/int_gcd_test.cpp
View file

@ -250,10 +250,9 @@ namespace lp {
}
void int_gcd_test::add_to_explanation_from_fixed_or_boxed_column(unsigned j) {
constraint_index lc, uc;
lra.get_bound_constraint_witnesses_for_column(j, lc, uc);
lia.m_ex->push_back(lc);
lia.m_ex->push_back(uc);
auto* deps = lra.get_bound_constraint_witnesses_for_column(j);
for (auto d : lra.flatten(deps))
lia.m_ex->push_back(d);
}
bool int_gcd_test::accumulate_parity(const row_strip<mpq> & row, unsigned least_idx) {

4
src/math/lp/int_solver.cpp
View file

@ -242,11 +242,11 @@ namespace lp {
return lra.has_inf_int();
}
constraint_index int_solver::column_upper_bound_constraint(unsigned j) const {
u_dependency* int_solver::column_upper_bound_constraint(unsigned j) const {
return lra.get_column_upper_bound_witness(j);
}
constraint_index int_solver::column_lower_bound_constraint(unsigned j) const {
u_dependency* int_solver::column_lower_bound_constraint(unsigned j) const {
return lra.get_column_lower_bound_witness(j);
}

4
src/math/lp/int_solver.h
View file

@ -113,8 +113,8 @@ private:
public:
std::ostream& display_column(std::ostream & out, unsigned j) const;
constraint_index column_upper_bound_constraint(unsigned j) const;
constraint_index column_lower_bound_constraint(unsigned j) const;
u_dependency* column_upper_bound_constraint(unsigned j) const;
u_dependency* column_lower_bound_constraint(unsigned j) const;
bool current_solution_is_inf_on_cut() const;
bool shift_var(unsigned j, unsigned range);

42
src/math/lp/lar_constraints.h
View file

@ -1,21 +1,10 @@
/*++
Copyright (c) 2017 Microsoft Corporation
Module Name:
<name>
Abstract:
<abstract>
Author:
Lev Nachmanson (levnach)
Revision History:
--*/
#pragma once
@ -53,15 +42,19 @@ class lar_base_constraint {
mpq m_right_side;
bool m_active;
unsigned m_j;
public:
u_dependency* m_dep;
public:
virtual vector<std::pair<mpq, var_index>> coeffs() const = 0;
lar_base_constraint(unsigned j, lconstraint_kind kind, const mpq& right_side) :m_kind(kind), m_right_side(right_side), m_active(false), m_j(j) {}
lar_base_constraint(unsigned j, lconstraint_kind kind, u_dependency* dep, const mpq& right_side) :
m_kind(kind), m_right_side(right_side), m_active(false), m_j(j), m_dep(dep) {}
virtual ~lar_base_constraint() = default;
lconstraint_kind kind() const { return m_kind; }
mpq const& rhs() const { return m_right_side; }
unsigned column() const { return m_j; }
u_dependency* dep() const { return m_dep; }
void activate() { m_active = true; }
void deactivate() { m_active = false; }
@ -73,8 +66,8 @@ public:
class lar_var_constraint: public lar_base_constraint {
public:
lar_var_constraint(unsigned j, lconstraint_kind kind, const mpq& right_side) :
lar_base_constraint(j, kind, right_side) {}
lar_var_constraint(unsigned j, lconstraint_kind kind, u_dependency* dep, const mpq& right_side) :
lar_base_constraint(j, kind, dep, right_side) {}
vector<std::pair<mpq, var_index>> coeffs() const override {
vector<std::pair<mpq, var_index>> ret;
@ -88,8 +81,8 @@ public:
class lar_term_constraint: public lar_base_constraint {
const lar_term * m_term;
public:
lar_term_constraint(unsigned j, const lar_term *t, lconstraint_kind kind, const mpq& right_side) :
lar_base_constraint(j, kind, right_side), m_term(t) {}
lar_term_constraint(unsigned j, const lar_term* t, lconstraint_kind kind, u_dependency* dep, const mpq& right_side) :
lar_base_constraint(j, kind, dep, right_side), m_term(t) {}
vector<std::pair<mpq, var_index>> coeffs() const override { return m_term->coeffs_as_vector(); }
unsigned size() const override { return m_term->size();}
@ -98,10 +91,11 @@ public:
class constraint_set {
region m_region;
column_namer& m_namer;
u_dependency_manager& m_dep_manager;
vector<lar_base_constraint*> m_constraints;
stacked_value<unsigned> m_constraint_count;
unsigned_vector m_active;
stacked_value<unsigned> m_active_lim;
stacked_value<unsigned> m_active_lim;
constraint_index add(lar_base_constraint* c) {
constraint_index ci = m_constraints.size();
@ -137,8 +131,13 @@ class constraint_set {
return out << "constraint " << T_to_string(ci) << " is not found" << std::endl;
}
u_dependency* mk_dep() {
return m_dep_manager.mk_leaf(m_constraints.size());
}
public:
constraint_set(column_namer& cn):
constraint_set(u_dependency_manager& d, column_namer& cn):
m_dep_manager(d),
m_namer(cn) {}
~constraint_set() {
@ -169,11 +168,12 @@ public:
}
constraint_index add_var_constraint(var_index j, lconstraint_kind k, mpq const& rhs) {
return add(new (m_region) lar_var_constraint(j, k, rhs));
return add(new (m_region) lar_var_constraint(j, k, mk_dep(), rhs));
}
constraint_index add_term_constraint(unsigned j, const lar_term* t, lconstraint_kind k, mpq const& rhs) {
return add(new (m_region) lar_term_constraint(j, t, k, rhs));
auto* dep = mk_dep();
return add(new (m_region) lar_term_constraint(j, t, k, dep, rhs));
}
// future behavior uses activation bit.

23
src/math/lp/lar_core_solver.h
View file

@ -21,9 +21,10 @@ class lar_core_solver {
int m_infeasible_sum_sign; // todo: get rid of this field
vector<numeric_pair<mpq>> m_right_sides_dummy;
vector<mpq> m_costs_dummy;
public:
stacked_value<simplex_strategy_enum> m_stacked_simplex_strategy;
public:
stacked_vector<column_type> m_column_types;
// r - solver fields, for rational numbers
vector<numeric_pair<mpq>> m_r_x; // the solution
@ -43,8 +44,6 @@ public:
const column_namer & column_names
);
lp_settings & settings() { return m_r_solver.m_settings;}
const lp_settings & settings() const { return m_r_solver.m_settings;}
int get_infeasible_sum_sign() const { return m_infeasible_sum_sign; }
@ -57,8 +56,7 @@ public:
void fill_not_improvable_zero_sum_from_inf_row();
column_type get_column_type(unsigned j) { return m_column_types[j];}
void print_pivot_row(std::ostream & out, unsigned row_index) const {
for (unsigned j : m_r_solver.m_pivot_row.m_index) {
if (numeric_traits<mpq>::is_pos(m_r_solver.m_pivot_row.m_data[j]))
@ -68,9 +66,9 @@ public:
out << " +" << m_r_solver.column_name(m_r_solver.m_basis[row_index]) << std::endl;
for (unsigned j : m_r_solver.m_pivot_row.m_index) {
for (unsigned j : m_r_solver.m_pivot_row.m_index)
m_r_solver.print_column_bound_info(j, out);
}
m_r_solver.print_column_bound_info(m_r_solver.m_basis[row_index], out);
}
@ -110,10 +108,7 @@ public:
m_column_types.push();
// rational
m_r_lower_bounds.push();
m_r_upper_bounds.push();
m_r_upper_bounds.push();
}
void pop(unsigned k) {
@ -127,11 +122,9 @@ public:
m_r_solver.m_d.resize(m_r_A.column_count());
m_stacked_simplex_strategy.pop(k);
settings().set_simplex_strategy(m_stacked_simplex_strategy);
m_r_solver.m_settings.set_simplex_strategy(m_stacked_simplex_strategy);
lp_assert(m_r_solver.basis_heading_is_correct());
}
bool r_basis_is_OK() const {
#ifdef Z3DEBUG

15
src/math/lp/lar_core_solver_def.h
View file

@ -88,19 +88,18 @@ void lar_core_solver::solve() {
lp_assert(m_r_solver.inf_heap_is_correct());
TRACE("find_feas_stats", tout << "infeasibles = " << m_r_solver.inf_heap_size() << ", int_infs = " << get_number_of_non_ints() << std::endl;);
if (m_r_solver.current_x_is_feasible() && m_r_solver.m_look_for_feasible_solution_only) {
m_r_solver.set_status(lp_status::OPTIMAL);
TRACE("lar_solver", tout << m_r_solver.get_status() << "\n";);
return;
m_r_solver.set_status(lp_status::OPTIMAL);
TRACE("lar_solver", tout << m_r_solver.get_status() << "\n";);
return;
}
++settings().stats().m_need_to_solve_inf;
++m_r_solver.m_settings.stats().m_need_to_solve_inf;
lp_assert( r_basis_is_OK());
if (m_r_solver.m_look_for_feasible_solution_only) //todo : should it be set?
m_r_solver.find_feasible_solution();
else {
else
m_r_solver.solve();
}
lp_assert(r_basis_is_OK());
switch (m_r_solver.get_status())

163
src/math/lp/lar_solver.cpp
View file

@ -26,7 +26,7 @@ namespace lp {
m_mpq_lar_core_solver(m_settings, *this),
m_var_register(false),
m_term_register(true),
m_constraints(*this) {}
m_constraints(m_dependencies, *this) {}
// start or ends tracking the rows that were changed by solve()
void lar_solver::track_touched_rows(bool v) {
@ -218,8 +218,14 @@ namespace lp {
// this is the case when the lower bound is in conflict with the upper one
const ul_pair& ul = m_columns_to_ul_pairs[m_crossed_bounds_column];
evidence.add_pair(ul.upper_bound_witness(), numeric_traits<mpq>::one());
evidence.add_pair(ul.lower_bound_witness(), -numeric_traits<mpq>::one());
svector<constraint_index> deps;
m_dependencies.linearize(ul.upper_bound_witness(), deps);
for (auto d : deps)
evidence.add_pair(d, numeric_traits<mpq>::one());
deps.reset();
m_dependencies.linearize(ul.lower_bound_witness(), deps);
for (auto d : deps)
evidence.add_pair(d, -numeric_traits<mpq>::one());
}
void lar_solver::push() {
@ -232,6 +238,7 @@ namespace lp {
m_term_count.push();
m_constraints.push();
m_usage_in_terms.push();
m_dependencies.push_scope();
}
void lar_solver::clean_popped_elements(unsigned n, indexed_uint_set& set) {
@ -297,6 +304,7 @@ namespace lp {
lp_assert(sizes_are_correct());
lp_assert(m_mpq_lar_core_solver.m_r_solver.reduced_costs_are_correct_tableau());
m_usage_in_terms.pop(k);
m_dependencies.pop_scope(k);
set_status(lp_status::UNKNOWN);
}
@ -320,6 +328,39 @@ namespace lp {
}
}
bool lar_solver::improve_bound(lpvar j, bool improve_lower_bound) {
lar_term term = get_term_to_maximize(j);
if (improve_lower_bound)
term.negate();
impq bound;
if (!maximize_term_on_tableau(term, bound))
return false;
return false;
// TODO
if (improve_lower_bound) {
if (column_has_lower_bound(j) && bound.x == column_lower_bound(j).x)
return false;
SASSERT(!column_has_lower_bound(j) || column_lower_bound(j).x < bound.x);
// TODO - explain new lower bound.
// Seems the relevant information is in the "costs" that are used when
// setting the optimization objecive. The costs are cleared after a call so
// maybe have some way of extracting bound dependencies from the costs.
u_dependency* dep = nullptr;
update_column_type_and_bound(j, bound.y > 0 ? lconstraint_kind::GT : lconstraint_kind::GE, bound.x, dep);
}
else {
if (column_has_upper_bound(j) && bound.x == column_upper_bound(j).x)
return false;
SASSERT(!column_has_upper_bound(j) || column_upper_bound(j).x > bound.x);
// similar for upper bounds
u_dependency* dep = nullptr;
update_column_type_and_bound(j, bound.y < 0 ? lconstraint_kind::LT : lconstraint_kind::LE, bound.x, dep);
}
return true;
}
bool lar_solver::costs_are_zeros_for_r_solver() const {
for (unsigned j = 0; j < m_mpq_lar_core_solver.m_r_solver.m_costs.size(); j++) {
lp_assert(is_zero(m_mpq_lar_core_solver.m_r_solver.m_costs[j]));
@ -349,7 +390,7 @@ namespace lp {
jset.insert(rc.var());
}
}
for (unsigned j : jset)
rslv.m_d[j] = zero_of_type<mpq>();
@ -577,15 +618,15 @@ namespace lp {
}
void lar_solver::set_upper_bound_witness(var_index j, constraint_index ci) {
void lar_solver::set_upper_bound_witness(var_index j, u_dependency* dep) {
ul_pair ul = m_columns_to_ul_pairs[j];
ul.upper_bound_witness() = ci;
ul.upper_bound_witness() = dep;
m_columns_to_ul_pairs[j] = ul;
}
void lar_solver::set_lower_bound_witness(var_index j, constraint_index ci) {
void lar_solver::set_lower_bound_witness(var_index j, u_dependency* dep) {
ul_pair ul = m_columns_to_ul_pairs[j];
ul.lower_bound_witness() = ci;
ul.lower_bound_witness() = dep;
m_columns_to_ul_pairs[j] = ul;
}
@ -920,7 +961,7 @@ namespace lp {
return ret;
}
bool lar_solver::has_lower_bound(var_index var, constraint_index& ci, mpq& value, bool& is_strict) const {
bool lar_solver::has_lower_bound(var_index var, u_dependency*& ci, mpq& value, bool& is_strict) const {
if (var >= m_columns_to_ul_pairs.size()) {
// TBD: bounds on terms could also be used, caller may have to track these.
@ -928,7 +969,7 @@ namespace lp {
}
const ul_pair& ul = m_columns_to_ul_pairs[var];
ci = ul.lower_bound_witness();
if (ci != null_ci) {
if (ci != nullptr) {
auto& p = m_mpq_lar_core_solver.m_r_lower_bounds()[var];
value = p.x;
is_strict = p.y.is_pos();
@ -939,7 +980,7 @@ namespace lp {
}
}
bool lar_solver::has_upper_bound(var_index var, constraint_index& ci, mpq& value, bool& is_strict) const {
bool lar_solver::has_upper_bound(var_index var, u_dependency*& ci, mpq& value, bool& is_strict) const {
if (var >= m_columns_to_ul_pairs.size()) {
// TBD: bounds on terms could also be used, caller may have to track these.
@ -947,7 +988,7 @@ namespace lp {
}
const ul_pair& ul = m_columns_to_ul_pairs[var];
ci = ul.upper_bound_witness();
if (ci != null_ci) {
if (ci != nullptr) {
auto& p = m_mpq_lar_core_solver.m_r_upper_bounds()[var];
value = p.x;
is_strict = p.y.is_neg();
@ -1005,9 +1046,13 @@ namespace lp {
int adj_sign = coeff.is_pos() ? inf_sign : -inf_sign;
const ul_pair& ul = m_columns_to_ul_pairs[j];
constraint_index bound_constr_i = adj_sign < 0 ? ul.upper_bound_witness() : ul.lower_bound_witness();
lp_assert(m_constraints.valid_index(bound_constr_i));
exp.add_pair(bound_constr_i, coeff);
u_dependency* bound_constr_i = adj_sign < 0 ? ul.upper_bound_witness() : ul.lower_bound_witness();
svector<constraint_index> deps;
m_dependencies.linearize(bound_constr_i, deps);
for (auto d : deps) {
lp_assert(m_constraints.valid_index(d));
exp.add_pair(d, coeff);
}
}
}
@ -1698,12 +1743,12 @@ namespace lp {
void lar_solver::activate_check_on_equal(constraint_index ci, unsigned& equal_column) {
auto const& c = m_constraints[ci];
update_column_type_and_bound_check_on_equal(c.column(), c.kind(), c.rhs(), ci, equal_column);
update_column_type_and_bound_check_on_equal(c.column(), c.rhs(), ci, equal_column);
}
void lar_solver::activate(constraint_index ci) {
auto const& c = m_constraints[ci];
update_column_type_and_bound(c.column(), c.kind(), c.rhs(), ci);
update_column_type_and_bound(c.column(), c.rhs(), ci);
}
mpq lar_solver::adjust_bound_for_int(lpvar j, lconstraint_kind& k, const mpq& bound) {
@ -1765,31 +1810,37 @@ namespace lp {
}
void lar_solver::update_column_type_and_bound(unsigned j,
lconstraint_kind kind,
const mpq& right_side,
constraint_index constr_index) {
TRACE("lar_solver_feas", tout << "j = " << j << " was " << (this->column_is_feasible(j)?"feas":"non-feas") << std::endl;);
m_constraints.activate(constr_index);
if (column_has_upper_bound(j))
update_column_type_and_bound_with_ub(j, kind, right_side, constr_index);
else
update_column_type_and_bound_with_no_ub(j, kind, right_side, constr_index);
if (is_base(j) && column_is_fixed(j))
m_fixed_base_var_set.insert(j);
TRACE("lar_solver_feas", tout << "j = " << j << " became " << (this->column_is_feasible(j)?"feas":"non-feas") << ", and " << (this->column_is_bounded(j)? "bounded":"non-bounded") << std::endl;);
lconstraint_kind kind = m_constraints[constr_index].kind();
u_dependency* dep = m_constraints[constr_index].dep();
update_column_type_and_bound(j, kind, right_side, dep);
}
void lar_solver::update_column_type_and_bound(unsigned j, lconstraint_kind kind, const mpq& right_side, u_dependency* dep) {
if (column_has_upper_bound(j))
update_column_type_and_bound_with_ub(j, kind, right_side, dep);
else
update_column_type_and_bound_with_no_ub(j, kind, right_side, dep);
if (is_base(j) && column_is_fixed(j))
m_fixed_base_var_set.insert(j);
TRACE("lar_solver_feas", tout << "j = " << j << " became " << (this->column_is_feasible(j) ? "feas" : "non-feas") << ", and " << (this->column_is_bounded(j) ? "bounded" : "non-bounded") << std::endl;);
}
void lar_solver::insert_to_columns_with_changed_bounds(unsigned j) {
m_columns_with_changed_bounds.insert(j);
TRACE("lar_solver", tout << "column " << j << (column_is_feasible(j) ? " feas" : " non-feas") << "\n";);
}
void lar_solver::update_column_type_and_bound_check_on_equal(unsigned j,
lconstraint_kind kind,
const mpq& right_side,
constraint_index constr_index,
unsigned& equal_to_j) {
update_column_type_and_bound(j, kind, right_side, constr_index);
update_column_type_and_bound(j, right_side, constr_index);
equal_to_j = null_lpvar;
if (column_is_fixed(j)) {
register_in_fixed_var_table(j, equal_to_j);
@ -1852,28 +1903,28 @@ namespace lp {
}
void lar_solver::update_column_type_and_bound_with_ub(unsigned j, lp::lconstraint_kind kind, const mpq& right_side, unsigned constraint_index) {
void lar_solver::update_column_type_and_bound_with_ub(unsigned j, lp::lconstraint_kind kind, const mpq& right_side, u_dependency* dep) {
SASSERT(column_has_upper_bound(j));
if (column_has_lower_bound(j)) {
update_bound_with_ub_lb(j, kind, right_side, constraint_index);
update_bound_with_ub_lb(j, kind, right_side, dep);
}
else {
update_bound_with_ub_no_lb(j, kind, right_side, constraint_index);
update_bound_with_ub_no_lb(j, kind, right_side, dep);
}
}
void lar_solver::update_column_type_and_bound_with_no_ub(unsigned j, lp::lconstraint_kind kind, const mpq& right_side, unsigned constraint_index) {
void lar_solver::update_column_type_and_bound_with_no_ub(unsigned j, lp::lconstraint_kind kind, const mpq& right_side, u_dependency* dep) {
SASSERT(!column_has_upper_bound(j));
if (column_has_lower_bound(j)) {
update_bound_with_no_ub_lb(j, kind, right_side, constraint_index);
update_bound_with_no_ub_lb(j, kind, right_side, dep);
}
else {
update_bound_with_no_ub_no_lb(j, kind, right_side, constraint_index);
update_bound_with_no_ub_no_lb(j, kind, right_side, dep);
}
}
// clang-format on
void lar_solver::update_bound_with_ub_lb(var_index j, lconstraint_kind kind, const mpq& right_side, constraint_index ci) {
void lar_solver::update_bound_with_ub_lb(var_index j, lconstraint_kind kind, const mpq& right_side, u_dependency* dep) {
lp_assert(column_has_lower_bound(j) && column_has_upper_bound(j));
lp_assert(m_mpq_lar_core_solver.m_column_types[j] == column_type::boxed ||
m_mpq_lar_core_solver.m_column_types[j] == column_type::fixed);
@ -1889,7 +1940,7 @@ namespace lp {
}
if (up >= m_mpq_lar_core_solver.m_r_upper_bounds[j]) return;
m_mpq_lar_core_solver.m_r_upper_bounds[j] = up;
set_upper_bound_witness(j, ci);
set_upper_bound_witness(j, dep);
insert_to_columns_with_changed_bounds(j);
break;
}
@ -1904,7 +1955,7 @@ namespace lp {
return;
}
m_mpq_lar_core_solver.m_r_lower_bounds[j] = low;
set_lower_bound_witness(j, ci);
set_lower_bound_witness(j, dep);
m_mpq_lar_core_solver.m_column_types[j] = (low == m_mpq_lar_core_solver.m_r_upper_bounds[j] ? column_type::fixed : column_type::boxed);
insert_to_columns_with_changed_bounds(j);
break;
@ -1914,8 +1965,8 @@ namespace lp {
if (v > m_mpq_lar_core_solver.m_r_upper_bounds[j] || v < m_mpq_lar_core_solver.m_r_lower_bounds[j]) {
set_infeasible_column(j);
}
set_upper_bound_witness(j, ci);
set_lower_bound_witness(j, ci);
set_upper_bound_witness(j, dep);
set_lower_bound_witness(j, dep);
m_mpq_lar_core_solver.m_r_upper_bounds[j] = m_mpq_lar_core_solver.m_r_lower_bounds[j] = v;
break;
}
@ -1928,7 +1979,7 @@ namespace lp {
}
}
// clang-format off
void lar_solver::update_bound_with_no_ub_lb(var_index j, lconstraint_kind kind, const mpq& right_side, constraint_index ci) {
void lar_solver::update_bound_with_no_ub_lb(var_index j, lconstraint_kind kind, const mpq& right_side, u_dependency* dep) {
lp_assert(column_has_lower_bound(j) && !column_has_upper_bound(j));
lp_assert(m_mpq_lar_core_solver.m_column_types[j] == column_type::lower_bound);
@ -1942,7 +1993,7 @@ namespace lp {
set_infeasible_column(j);
}
m_mpq_lar_core_solver.m_r_upper_bounds[j] = up;
set_upper_bound_witness(j, ci);
set_upper_bound_witness(j, dep);
m_mpq_lar_core_solver.m_column_types[j] = (up == m_mpq_lar_core_solver.m_r_lower_bounds[j] ? column_type::fixed : column_type::boxed);
insert_to_columns_with_changed_bounds(j);
break;
@ -1955,7 +2006,7 @@ namespace lp {
return;
}
m_mpq_lar_core_solver.m_r_lower_bounds[j] = low;
set_lower_bound_witness(j, ci);
set_lower_bound_witness(j, dep);
insert_to_columns_with_changed_bounds(j);
break;
}
@ -1965,8 +2016,8 @@ namespace lp {
set_infeasible_column(j);
}
set_upper_bound_witness(j, ci);
set_lower_bound_witness(j, ci);
set_upper_bound_witness(j, dep);
set_lower_bound_witness(j, dep);
m_mpq_lar_core_solver.m_r_upper_bounds[j] = m_mpq_lar_core_solver.m_r_lower_bounds[j] = v;
m_mpq_lar_core_solver.m_column_types[j] = column_type::fixed;
break;
@ -1977,7 +2028,7 @@ namespace lp {
}
}
// clang-format off
void lar_solver::update_bound_with_ub_no_lb(var_index j, lconstraint_kind kind, const mpq& right_side, constraint_index ci) {
void lar_solver::update_bound_with_ub_no_lb(var_index j, lconstraint_kind kind, const mpq& right_side, u_dependency* dep) {
lp_assert(!column_has_lower_bound(j) && column_has_upper_bound(j));
lp_assert(m_mpq_lar_core_solver.m_column_types[j] == column_type::upper_bound);
mpq y_of_bound(0);
@ -1989,7 +2040,7 @@ namespace lp {
auto up = numeric_pair<mpq>(right_side, y_of_bound);
if (up >= m_mpq_lar_core_solver.m_r_upper_bounds[j]) return;
m_mpq_lar_core_solver.m_r_upper_bounds[j] = up;
set_upper_bound_witness(j, ci);
set_upper_bound_witness(j, dep);
insert_to_columns_with_changed_bounds(j);
}
break;
@ -2002,7 +2053,7 @@ namespace lp {
set_infeasible_column(j);
}
m_mpq_lar_core_solver.m_r_lower_bounds[j] = low;
set_lower_bound_witness(j, ci);
set_lower_bound_witness(j, dep);
m_mpq_lar_core_solver.m_column_types[j] = (low == m_mpq_lar_core_solver.m_r_upper_bounds[j] ? column_type::fixed : column_type::boxed);
insert_to_columns_with_changed_bounds(j);
@ -2015,8 +2066,8 @@ namespace lp {
set_infeasible_column(j);
}
set_upper_bound_witness(j, ci);
set_lower_bound_witness(j, ci);
set_upper_bound_witness(j, dep);
set_lower_bound_witness(j, dep);
m_mpq_lar_core_solver.m_r_upper_bounds[j] = m_mpq_lar_core_solver.m_r_lower_bounds[j] = v;
m_mpq_lar_core_solver.m_column_types[j] = column_type::fixed;
break;
@ -2027,7 +2078,7 @@ namespace lp {
}
}
// clang-format on
void lar_solver::update_bound_with_no_ub_no_lb(var_index j, lconstraint_kind kind, const mpq& right_side, constraint_index ci) {
void lar_solver::update_bound_with_no_ub_no_lb(var_index j, lconstraint_kind kind, const mpq& right_side, u_dependency* dep) {
lp_assert(!column_has_lower_bound(j) && !column_has_upper_bound(j));
mpq y_of_bound(0);
@ -2037,7 +2088,7 @@ namespace lp {
case LE: {
auto up = numeric_pair<mpq>(right_side, y_of_bound);
m_mpq_lar_core_solver.m_r_upper_bounds[j] = up;
set_upper_bound_witness(j, ci);
set_upper_bound_witness(j, dep);
m_mpq_lar_core_solver.m_column_types[j] = column_type::upper_bound;
} break;
case GT:
@ -2045,14 +2096,14 @@ namespace lp {
case GE: {
auto low = numeric_pair<mpq>(right_side, y_of_bound);
m_mpq_lar_core_solver.m_r_lower_bounds[j] = low;
set_lower_bound_witness(j, ci);
set_lower_bound_witness(j, dep);
m_mpq_lar_core_solver.m_column_types[j] = column_type::lower_bound;
} break;
case EQ: {
auto v = numeric_pair<mpq>(right_side, zero_of_type<mpq>());
set_upper_bound_witness(j, ci);
set_lower_bound_witness(j, ci);
set_upper_bound_witness(j, dep);
set_lower_bound_witness(j, dep);
m_mpq_lar_core_solver.m_r_upper_bounds[j] = m_mpq_lar_core_solver.m_r_lower_bounds[j] = v;
m_mpq_lar_core_solver.m_column_types[j] = column_type::fixed;
break;
@ -2165,7 +2216,7 @@ namespace lp {
return true;
}
bool lar_solver::get_equality_and_right_side_for_term_on_current_x(tv const& t, mpq& rs, constraint_index& ci, bool& upper_bound) const {
bool lar_solver::get_equality_and_right_side_for_term_on_current_x(tv const& t, mpq& rs, u_dependency*& ci, bool& upper_bound) const {
lp_assert(t.is_term());
unsigned j;
bool is_int;

63
src/math/lp/lar_solver.h
View file

@ -47,6 +47,8 @@ namespace lp {
class int_branch;
class int_solver;
class lar_solver : public column_namer {
struct term_hasher {
std::size_t operator()(const lar_term& t) const {
@ -88,6 +90,8 @@ class lar_solver : public column_namer {
indexed_uint_set m_columns_with_changed_bounds;
indexed_uint_set m_touched_rows;
unsigned_vector m_row_bounds_to_replay;
u_dependency_manager m_dependencies;
svector<constraint_index> m_tmp_dependencies;
indexed_uint_set m_basic_columns_with_changed_cost;
// these are basic columns with the value changed, so the corresponding row in the tableau
@ -136,14 +140,15 @@ class lar_solver : public column_namer {
inline void clear_columns_with_changed_bounds() { m_columns_with_changed_bounds.reset(); }
void insert_to_columns_with_changed_bounds(unsigned j);
void update_column_type_and_bound_check_on_equal(unsigned j, lconstraint_kind kind, const mpq& right_side, constraint_index constr_index, unsigned&);
void update_column_type_and_bound(unsigned j, lconstraint_kind kind, const mpq& right_side, constraint_index constr_index);
void update_column_type_and_bound_with_ub(var_index j, lconstraint_kind kind, const mpq& right_side, constraint_index constr_index);
void update_column_type_and_bound_with_no_ub(var_index j, lconstraint_kind kind, const mpq& right_side, constraint_index constr_index);
void update_bound_with_ub_lb(var_index j, lconstraint_kind kind, const mpq& right_side, constraint_index constr_index);
void update_bound_with_no_ub_lb(var_index j, lconstraint_kind kind, const mpq& right_side, constraint_index constr_index);
void update_bound_with_ub_no_lb(var_index j, lconstraint_kind kind, const mpq& right_side, constraint_index constr_index);
void update_bound_with_no_ub_no_lb(var_index j, lconstraint_kind kind, const mpq& right_side, constraint_index constr_index);
void update_column_type_and_bound_check_on_equal(unsigned j, const mpq& right_side, constraint_index ci, unsigned&);
void update_column_type_and_bound(unsigned j, const mpq& right_side, constraint_index ci);
void update_column_type_and_bound(unsigned j, lconstraint_kind kind, const mpq& right_side, u_dependency* dep);
void update_column_type_and_bound_with_ub(var_index j, lconstraint_kind kind, const mpq& right_side, u_dependency* dep);
void update_column_type_and_bound_with_no_ub(var_index j, lconstraint_kind kind, const mpq& right_side, u_dependency* dep);
void update_bound_with_ub_lb(var_index j, lconstraint_kind kind, const mpq& right_side, u_dependency* dep);
void update_bound_with_no_ub_lb(var_index j, lconstraint_kind kind, const mpq& right_side, u_dependency* dep);
void update_bound_with_ub_no_lb(var_index j, lconstraint_kind kind, const mpq& right_side, u_dependency* dep);
void update_bound_with_no_ub_no_lb(var_index j, lconstraint_kind kind, const mpq& right_side, u_dependency* dep);
void register_in_fixed_var_table(unsigned, unsigned&);
void remove_non_fixed_from_fixed_var_table();
constraint_index add_var_bound_on_constraint_for_term(var_index j, lconstraint_kind kind, const mpq& right_side);
@ -186,8 +191,8 @@ class lar_solver : public column_namer {
bool maximize_term_on_corrected_r_solver(lar_term& term, impq& term_max);
void pop_core_solver_params();
void pop_core_solver_params(unsigned k);
void set_upper_bound_witness(var_index j, constraint_index ci);
void set_lower_bound_witness(var_index j, constraint_index ci);
void set_upper_bound_witness(var_index j, u_dependency* ci);
void set_lower_bound_witness(var_index j, u_dependency* ci);
void substitute_terms_in_linear_expression(const vector<std::pair<mpq, var_index>>& left_side_with_terms,
vector<std::pair<mpq, var_index>>& left_side) const;
@ -283,6 +288,8 @@ class lar_solver : public column_namer {
lp_status maximize_term(unsigned j_or_term, impq& term_max);
bool improve_bound(lpvar j, bool is_lower);
inline core_solver_pretty_printer<lp::mpq, lp::impq> pp(std::ostream& out) const {
return core_solver_pretty_printer<lp::mpq, lp::impq>(m_mpq_lar_core_solver.m_r_solver, out);
}
@ -310,9 +317,10 @@ class lar_solver : public column_namer {
int a_sign = is_pos(a) ? 1 : -1;
int sign = j_sign * a_sign;
const ul_pair& ul = m_columns_to_ul_pairs[j];
auto witness = sign > 0 ? ul.upper_bound_witness() : ul.lower_bound_witness();
lp_assert(is_valid(witness));
bp.consume(a, witness);
auto* witness = sign > 0 ? ul.upper_bound_witness() : ul.lower_bound_witness();
lp_assert(witness);
for (auto ci : flatten(witness))
bp.consume(a, ci);
}
}
@ -460,7 +468,20 @@ class lar_solver : public column_namer {
return m_mpq_lar_core_solver.m_r_solver.column_has_lower_bound(j);
}
inline constraint_index get_column_upper_bound_witness(unsigned j) const {
svector<constraint_index> const& flatten(u_dependency* d) {
m_tmp_dependencies.reset();
m_dependencies.linearize(d, m_tmp_dependencies);
return m_tmp_dependencies;
}
void push_explanation(u_dependency* d, explanation& ex) {
for (auto ci : flatten(d))
ex.push_back(ci);
}
u_dependency_manager& dep_manager() { return m_dependencies; }
inline u_dependency* get_column_upper_bound_witness(unsigned j) const {
if (tv::is_term(j)) {
j = m_var_register.external_to_local(j);
}
@ -474,8 +495,8 @@ class lar_solver : public column_namer {
inline const impq& get_lower_bound(column_index j) const {
return m_mpq_lar_core_solver.m_r_solver.m_lower_bounds[j];
}
bool has_lower_bound(var_index var, constraint_index& ci, mpq& value, bool& is_strict) const;
bool has_upper_bound(var_index var, constraint_index& ci, mpq& value, bool& is_strict) const;
bool has_lower_bound(var_index var, u_dependency*& ci, mpq& value, bool& is_strict) const;
bool has_upper_bound(var_index var, u_dependency*& ci, mpq& value, bool& is_strict) const;
bool has_value(var_index var, mpq& value) const;
bool fetch_normalized_term_column(const lar_term& t, std::pair<mpq, lpvar>&) const;
unsigned map_term_index_to_column_index(unsigned j) const;
@ -530,15 +551,15 @@ class lar_solver : public column_namer {
std::pair<constraint_index, constraint_index> add_equality(lpvar j, lpvar k);
inline void get_bound_constraint_witnesses_for_column(unsigned j, constraint_index& lc, constraint_index& uc) const {
u_dependency* get_bound_constraint_witnesses_for_column(unsigned j) {
const ul_pair& ul = m_columns_to_ul_pairs[j];
lc = ul.lower_bound_witness();
uc = ul.upper_bound_witness();
return m_dependencies.mk_join(ul.lower_bound_witness(), ul.upper_bound_witness());
}
inline constraint_set const& constraints() const { return m_constraints; }
void push();
void pop();
inline constraint_index get_column_lower_bound_witness(unsigned j) const {
inline u_dependency* get_column_lower_bound_witness(unsigned j) const {
if (tv::is_term(j)) {
j = m_var_register.external_to_local(j);
}
@ -601,7 +622,7 @@ class lar_solver : public column_namer {
inline const column_strip& get_column(unsigned i) const { return A_r().m_columns[i]; }
bool row_is_correct(unsigned i) const;
bool ax_is_correct() const;
bool get_equality_and_right_side_for_term_on_current_x(tv const& t, mpq& rs, constraint_index& ci, bool& upper_bound) const;
bool get_equality_and_right_side_for_term_on_current_x(tv const& t, mpq& rs, u_dependency*& ci, bool& upper_bound) const;
bool var_is_int(var_index v) const;
inline const vector<int>& r_heading() const { return m_mpq_lar_core_solver.m_r_heading; }
inline const vector<unsigned>& r_basis() const { return m_mpq_lar_core_solver.r_basis(); }

17
src/math/lp/lp_bound_propagator.h
View file

@ -228,33 +228,32 @@ class lp_bound_propagator {
return lp().column_is_int(j);
}
void explain_fixed_in_row(unsigned row, explanation& ex) const {
void explain_fixed_in_row(unsigned row, explanation& ex) {
TRACE("eq", tout << lp().get_row(row) << std::endl);
for (const auto& c : lp().get_row(row))
if (lp().is_fixed(c.var()))
explain_fixed_column(c.var(), ex);
}
unsigned explain_fixed_in_row_and_get_base(unsigned row, explanation& ex) const {
unsigned explain_fixed_in_row_and_get_base(unsigned row, explanation& ex) {
unsigned base = UINT_MAX;
TRACE("eq", tout << lp().get_row(row) << std::endl);
for (const auto& c : lp().get_row(row)) {
if (lp().is_fixed(c.var())) {
explain_fixed_column(c.var(), ex);
} else if (lp().is_base(c.var())) {
}
else if (lp().is_base(c.var())) {
base = c.var();
}
}
return base;
}
void explain_fixed_column(unsigned j, explanation& ex) const {
void explain_fixed_column(unsigned j, explanation& ex) {
SASSERT(column_is_fixed(j));
constraint_index lc, uc;
lp().get_bound_constraint_witnesses_for_column(j, lc, uc);
ex.push_back(lc);
if (lc != uc)
ex.push_back(uc);
auto* deps = lp().get_bound_constraint_witnesses_for_column(j);
for (auto ci : lp().flatten(deps))
ex.push_back(ci);
}
#ifdef Z3DEBUG
bool all_fixed_in_row(unsigned row) const {

171
src/math/lp/lp_primal_core_solver.h
View file

@ -138,95 +138,98 @@ namespace lp {
return false;
}
/**
* Return the number of base non-free variables depending on the column j,
* different from bj,
* but take the min with the (bound+1).
* This function is used to select the pivot variable.
*/
unsigned get_num_of_not_free_basic_dependent_vars(
unsigned j, unsigned bound, unsigned bj) const { // consider looking at the signs here: todo
unsigned r = 0;
for (const auto &cc : this->m_A.m_columns[j]) {
unsigned basic_for_row = this->m_basis[cc.var()];
if (basic_for_row == bj)
continue;
// std::cout << this->m_A.m_rows[cc.var()] << std::endl;
if (this->m_column_types[basic_for_row] != column_type::free_column)
if (r++ > bound) return r;
}
return r;
}
int find_beneficial_entering_in_row_tableau_rows_bland_mode(int i, T &a_ent) {
int j = -1;
unsigned bj = this->m_basis[i];
bool bj_needs_to_grow = needs_to_grow(bj);
for (const row_cell<T> &rc : this->m_A.m_rows[i]) {
if (rc.var() == bj)
continue;
if (bj_needs_to_grow) {
if (!monoid_can_decrease(rc))
continue;
} else {
if (!monoid_can_increase(rc))
/**
* Return the number of base non-free variables depending on the column j,
* different from bj,
* but take the min with the (bound+1).
* This function is used to select the pivot variable.
*/
unsigned get_num_of_not_free_basic_dependent_vars(unsigned j, unsigned bound, unsigned bj) const {
// consider looking at the signs here: todo
unsigned r = 0;
for (const auto &cc : this->m_A.m_columns[j]) {
unsigned basic_for_row = this->m_basis[cc.var()];
if (basic_for_row == bj)
continue;
// std::cout << this->m_A.m_rows[cc.var()] << std::endl;
if (this->m_column_types[basic_for_row] != column_type::free_column)
if (r++ > bound) return r;
}
if (rc.var() < static_cast<unsigned>(j)) {
j = rc.var();
a_ent = rc.coeff();
}
}
if (j == -1)
m_inf_row_index_for_tableau = i;
return j;
}
int find_beneficial_entering_tableau_rows(int i, T &a_ent) {
if (m_bland_mode_tableau)
return find_beneficial_entering_in_row_tableau_rows_bland_mode(i, a_ent);
// a short row produces short infeasibility explanation and benefits at
// least one pivot operation
int choice = -1;
int nchoices = 0;
unsigned min_non_free_so_far = -1;
unsigned best_col_sz = -1;
unsigned bj = this->m_basis[i];
bool bj_needs_to_grow = needs_to_grow(bj);
for (unsigned k = 0; k < this->m_A.m_rows[i].size(); k++) {
const row_cell<T> &rc = this->m_A.m_rows[i][k];
unsigned j = rc.var();
if (j == bj)
continue;
if (bj_needs_to_grow) {
if (!monoid_can_decrease(rc))
continue;
} else {
if (!monoid_can_increase(rc))
continue;
}
unsigned not_free = get_num_of_not_free_basic_dependent_vars(j, min_non_free_so_far, bj);
unsigned col_sz = this->m_A.m_columns[j].size();
if (not_free < min_non_free_so_far || (not_free == min_non_free_so_far && col_sz < best_col_sz)) {
min_non_free_so_far = not_free;
best_col_sz = this->m_A.m_columns[j].size();
choice = k;
nchoices = 1;
} else if (not_free == min_non_free_so_far &&
col_sz == best_col_sz) {
if (this->m_settings.random_next(++nchoices) == 0) {
choice = k;
}
}
return r;
}
if (choice == -1) {
m_inf_row_index_for_tableau = i;
return -1;
int find_beneficial_entering_in_row_tableau_rows_bland_mode(int i, T &a_ent) {
int j = -1;
unsigned bj = this->m_basis[i];
bool bj_needs_to_grow = needs_to_grow(bj);
for (const row_cell<T> &rc : this->m_A.m_rows[i]) {
if (rc.var() == bj)
continue;
if (bj_needs_to_grow) {
if (!monoid_can_decrease(rc))
continue;
}
else {
if (!monoid_can_increase(rc))
continue;
}
if (rc.var() < static_cast<unsigned>(j)) {
j = rc.var();
a_ent = rc.coeff();
}
}
if (j == -1)
m_inf_row_index_for_tableau = i;
return j;
}
int find_beneficial_entering_tableau_rows(int i, T &a_ent) {
if (m_bland_mode_tableau)
return find_beneficial_entering_in_row_tableau_rows_bland_mode(i, a_ent);
// a short row produces short infeasibility explanation and benefits at
// least one pivot operation
int choice = -1;
int nchoices = 0;
unsigned min_non_free_so_far = -1;
unsigned best_col_sz = -1;
unsigned bj = this->m_basis[i];
bool bj_needs_to_grow = needs_to_grow(bj);
for (unsigned k = 0; k < this->m_A.m_rows[i].size(); k++) {
const row_cell<T> &rc = this->m_A.m_rows[i][k];
unsigned j = rc.var();
if (j == bj)
continue;
if (bj_needs_to_grow) {
if (!monoid_can_decrease(rc))
continue;
} else {
if (!monoid_can_increase(rc))
continue;
}
unsigned not_free = get_num_of_not_free_basic_dependent_vars(j, min_non_free_so_far, bj);
unsigned col_sz = this->m_A.m_columns[j].size();
if (not_free < min_non_free_so_far || (not_free == min_non_free_so_far && col_sz < best_col_sz)) {
min_non_free_so_far = not_free;
best_col_sz = this->m_A.m_columns[j].size();
choice = k;
nchoices = 1;
}
else if (not_free == min_non_free_so_far &&
col_sz == best_col_sz) {
if (this->m_settings.random_next(++nchoices) == 0)
choice = k;
}
}
if (choice == -1) {
m_inf_row_index_for_tableau = i;
return -1;
}
const row_cell<T> &rc = this->m_A.m_rows[i][choice];
a_ent = rc.coeff();
return rc.var();
}
const row_cell<T> &rc = this->m_A.m_rows[i][choice];
a_ent = rc.coeff();
return rc.var();
}
bool try_jump_to_another_bound_on_entering(unsigned entering, const X &theta, X &t, bool &unlimited);

3
src/math/lp/lp_types.h
View file

@ -22,6 +22,8 @@ Revision History:
#include <sstream>
#include <limits.h>
#include "util/debug.h"
#include "util/dependency.h"
namespace nla {
class core;
}
@ -36,6 +38,7 @@ typedef unsigned lpvar;
const lpvar null_lpvar = UINT_MAX;
const constraint_index null_ci = UINT_MAX;
class column_index {
unsigned m_index;
friend class lar_solver;

10
src/math/lp/monomial_bounds.cpp
View file

@ -149,22 +149,22 @@ namespace nla {
}
void monomial_bounds::var2interval(lpvar v, scoped_dep_interval& i) {
lp::constraint_index ci;
u_dependency* d = nullptr;
rational bound;
bool is_strict;
if (c().has_lower_bound(v, ci, bound, is_strict)) {
if (c().has_lower_bound(v, d, bound, is_strict)) {
dep.set_lower_is_open(i, is_strict);
dep.set_lower(i, bound);
dep.set_lower_dep(i, dep.mk_leaf(ci));
dep.set_lower_dep(i, d);
dep.set_lower_is_inf(i, false);
}
else {
dep.set_lower_is_inf(i, true);
}
if (c().has_upper_bound(v, ci, bound, is_strict)) {
if (c().has_upper_bound(v, d, bound, is_strict)) {
dep.set_upper_is_open(i, is_strict);
dep.set_upper(i, bound);
dep.set_upper_dep(i, dep.mk_leaf(ci));
dep.set_upper_dep(i, d);
dep.set_upper_is_inf(i, false);
}
else {

12
src/math/lp/nla_common.cpp
View file

@ -60,11 +60,9 @@ unsigned common::random() {
}
void common::add_deps_of_fixed(lpvar j, u_dependency*& dep) {
unsigned lc, uc;
auto& dep_manager = c().m_intervals.get_dep_intervals().dep_manager();
c().m_lar_solver.get_bound_constraint_witnesses_for_column(j, lc, uc);
dep = dep_manager.mk_join(dep, dep_manager.mk_leaf(lc));
dep = dep_manager.mk_join(dep, dep_manager.mk_leaf(uc));
auto& dm = c().lra.dep_manager();
auto* deps = c().lra.get_bound_constraint_witnesses_for_column(j);
dep = dm.mk_join(dep, deps);
}
@ -73,7 +71,7 @@ nex * common::nexvar(const rational & coeff, lpvar j, nex_creator& cn, u_depende
SASSERT(!coeff.is_zero());
if (c().params().arith_nl_horner_subs_fixed() == 1 && c().var_is_fixed(j)) {
add_deps_of_fixed(j, dep);
return cn.mk_scalar(coeff * c().m_lar_solver.column_lower_bound(j).x);
return cn.mk_scalar(coeff * c().lra.column_lower_bound(j).x);
}
if (c().params().arith_nl_horner_subs_fixed() == 2 && c().var_is_fixed_to_zero(j)) {
add_deps_of_fixed(j, dep);
@ -91,7 +89,7 @@ nex * common::nexvar(const rational & coeff, lpvar j, nex_creator& cn, u_depende
for (lpvar k : m.vars()) {
if (c().params().arith_nl_horner_subs_fixed() == 1 && c().var_is_fixed(k)) {
add_deps_of_fixed(k, dep);
mf *= c().m_lar_solver.column_lower_bound(k).x;
mf *= c().lra.column_lower_bound(k).x;
} else if (c().params().arith_nl_horner_subs_fixed() == 2 &&
c().var_is_fixed_to_zero(k)) {
dep = initial_dep;

15
src/math/lp/nla_common.h
View file

@ -80,21 +80,6 @@ struct common {
bool check_monic(const monic&) const;
unsigned random();
void add_deps_of_fixed(lpvar j, u_dependency*& dep);
class ci_value_manager {
public:
void inc_ref(lp::constraint_index const & v) {
}
void dec_ref(lp::constraint_index const & v) {
}
};
struct u_dependency_config {
typedef ci_value_manager value_manager;
typedef region allocator;
static const bool ref_count = false;
typedef lp::constraint_index value;
};
nex* nexvar(const rational& coeff, lpvar j, nex_creator&, u_dependency*&);
template <typename T>

180
src/math/lp/nla_core.cpp
View file

@ -23,7 +23,7 @@ typedef lp::lar_term term;
core::core(lp::lar_solver& s, params_ref const& p, reslimit & lim) :
m_evars(),
m_lar_solver(s),
lra(s),
m_reslim(lim),
m_params(p),
m_tangents(this),
@ -69,7 +69,7 @@ lp::lar_term core::subs_terms_to_columns(const lp::lar_term& t) const {
for (lp::lar_term::ival p : t) {
lpvar j = p.column();
if (lp::tv::is_term(j))
j = m_lar_solver.map_term_index_to_column_index(j);
j = lra.map_term_index_to_column_index(j);
r.add_monomial(p.coeff(), j);
}
return r;
@ -133,7 +133,7 @@ void core::add_monic(lpvar v, unsigned sz, lpvar const* vs) {
for (unsigned i = 0; i < sz; i++) {
lpvar j = vs[i];
if (lp::tv::is_term(j))
j = m_lar_solver.map_term_index_to_column_index(j);
j = lra.map_term_index_to_column_index(j);
m_add_buffer[i] = j;
}
m_emons.add(v, m_add_buffer);
@ -154,7 +154,7 @@ void core::pop(unsigned n) {
rational core::product_value(const monic& m) const {
rational r(1);
for (auto j : m.vars()) {
r *= m_lar_solver.get_column_value(j).x;
r *= lra.get_column_value(j).x;
}
return r;
}
@ -167,10 +167,10 @@ bool core::check_monic(const monic& m) const {
if (!is_relevant(m.var()))
return true;
#endif
if (m_lar_solver.column_is_int(m.var()) && !m_lar_solver.get_column_value(m.var()).is_int())
if (lra.column_is_int(m.var()) && !lra.get_column_value(m.var()).is_int())
return true;
bool ret = product_value(m) == m_lar_solver.get_column_value(m.var()).x;
bool ret = product_value(m) == lra.get_column_value(m.var()).x;
CTRACE("nla_solver_check_monic", !ret, print_monic(m, tout) << '\n';);
return ret;
}
@ -182,7 +182,7 @@ std::ostream& core::print_product(const T & m, std::ostream& out) const {
for (lpvar v : m) {
if (!first) out << "*"; else first = false;
if (lp_settings().print_external_var_name())
out << "(" << m_lar_solver.get_variable_name(v) << "=" << val(v) << ")";
out << "(" << lra.get_variable_name(v) << "=" << val(v) << ")";
else
out << "(j" << v << " = " << val(v) << ")";
@ -228,7 +228,7 @@ std::ostream & core::print_factor_with_vars(const factor& f, std::ostream& out)
std::ostream& core::print_monic(const monic& m, std::ostream& out) const {
if (lp_settings().print_external_var_name())
out << "([" << m.var() << "] = " << m_lar_solver.get_variable_name(m.var()) << " = " << val(m.var()) << " = ";
out << "([" << m.var() << "] = " << lra.get_variable_name(m.var()) << " = " << val(m.var()) << " = ";
else
out << "(j" << m.var() << " = " << val(m.var()) << " = ";
print_product(m.vars(), out) << ")\n";
@ -271,7 +271,7 @@ std::ostream& core::print_explanation(const lp::explanation& exp, std::ostream&
unsigned i = 0;
for (auto p : exp) {
out << "(" << p.ci() << ")";
m_lar_solver.constraints().display(out, [this](lpvar j) { return var_str(j);}, p.ci());
lra.constraints().display(out, [this](lpvar j) { return var_str(j);}, p.ci());
if (++i < exp.size())
out << " ";
}
@ -316,21 +316,20 @@ bool core::explain_lower_bound(const lp::lar_term& t, const rational& rs, lp::ex
bool core::explain_coeff_lower_bound(const lp::lar_term::ival& p, rational& bound, lp::explanation& e) const {
const rational& a = p.coeff();
SASSERT(!a.is_zero());
unsigned c; // the index for the lower or the upper bound
if (a.is_pos()) {
unsigned c = m_lar_solver.get_column_lower_bound_witness(p.column());
if (c + 1 == 0)
auto* dep = lra.get_column_lower_bound_witness(p.column());
if (!dep)
return false;
bound = a * m_lar_solver.get_lower_bound(p.column()).x;
e.push_back(c);
bound = a * lra.get_lower_bound(p.column()).x;
lra.push_explanation(dep, e);
return true;
}
// a.is_neg()
c = m_lar_solver.get_column_upper_bound_witness(p.column());
if (c + 1 == 0)
auto* dep = lra.get_column_upper_bound_witness(p.column());
if (!dep)
return false;
bound = a * m_lar_solver.get_upper_bound(p.column()).x;
e.push_back(c);
bound = a * lra.get_upper_bound(p.column()).x;
lra.push_explanation(dep, e);
return true;
}
@ -338,21 +337,20 @@ bool core::explain_coeff_upper_bound(const lp::lar_term::ival& p, rational& boun
const rational& a = p.coeff();
lpvar j = p.column();
SASSERT(!a.is_zero());
unsigned c; // the index for the lower or the upper bound
if (a.is_neg()) {
unsigned c = m_lar_solver.get_column_lower_bound_witness(j);
if (c + 1 == 0)
auto *dep = lra.get_column_lower_bound_witness(j);
if (!dep)
return false;
bound = a * m_lar_solver.get_lower_bound(j).x;
e.push_back(c);
bound = a * lra.get_lower_bound(j).x;
lra.push_explanation(dep, e);
return true;
}
// a.is_pos()
c = m_lar_solver.get_column_upper_bound_witness(j);
if (c + 1 == 0)
auto* dep = lra.get_column_upper_bound_witness(j);
if (!dep)
return false;
bound = a * m_lar_solver.get_upper_bound(j).x;
e.push_back(c);
bound = a * lra.get_upper_bound(j).x;
lra.push_explanation(dep, e);
return true;
}
@ -413,12 +411,12 @@ bool core::explain_by_equiv(const lp::lar_term& t, lp::explanation& e) const {
return false;
m_evars.explain(signed_var(i, false), signed_var(j, sign), e);
TRACE("nla_solver", tout << "explained :"; m_lar_solver.print_term_as_indices(t, tout););
TRACE("nla_solver", tout << "explained :"; lra.print_term_as_indices(t, tout););
return true;
}
void core::mk_ineq_no_expl_check(new_lemma& lemma, lp::lar_term& t, llc cmp, const rational& rs) {
TRACE("nla_solver_details", m_lar_solver.print_term_as_indices(t, tout << "t = "););
TRACE("nla_solver_details", lra.print_term_as_indices(t, tout << "t = "););
lemma |= ineq(cmp, t, rs);
CTRACE("nla_solver", ineq_holds(ineq(cmp, t, rs)), print_ineq(ineq(cmp, t, rs), tout) << "\n";);
SASSERT(!ineq_holds(ineq(cmp, t, rs)));
@ -485,18 +483,18 @@ int core::vars_sign(const svector<lpvar>& v) {
}
bool core::has_upper_bound(lpvar j) const {
return m_lar_solver.column_has_upper_bound(j);
return lra.column_has_upper_bound(j);
}
bool core::has_lower_bound(lpvar j) const {
return m_lar_solver.column_has_lower_bound(j);
return lra.column_has_lower_bound(j);
}
const rational& core::get_upper_bound(unsigned j) const {
return m_lar_solver.get_upper_bound(j).x;
return lra.get_upper_bound(j).x;
}
const rational& core::get_lower_bound(unsigned j) const {
return m_lar_solver.get_lower_bound(j).x;
return lra.get_lower_bound(j).x;
}
bool core::zero_is_an_inner_point_of_bounds(lpvar j) const {
@ -540,26 +538,26 @@ bool core::sign_contradiction(const monic& m) const {
bool core::var_is_fixed_to_zero(lpvar j) const {
return
m_lar_solver.column_is_fixed(j) &&
m_lar_solver.get_lower_bound(j) == lp::zero_of_type<lp::impq>();
lra.column_is_fixed(j) &&
lra.get_lower_bound(j) == lp::zero_of_type<lp::impq>();
}
bool core::var_is_fixed_to_val(lpvar j, const rational& v) const {
return
m_lar_solver.column_is_fixed(j) &&
m_lar_solver.get_lower_bound(j) == lp::impq(v);
lra.column_is_fixed(j) &&
lra.get_lower_bound(j) == lp::impq(v);
}
bool core::var_is_fixed(lpvar j) const {
return m_lar_solver.column_is_fixed(j);
return lra.column_is_fixed(j);
}
bool core::var_is_free(lpvar j) const {
return m_lar_solver.column_is_free(j);
return lra.column_is_free(j);
}
std::ostream & core::print_ineq(const ineq & in, std::ostream & out) const {
m_lar_solver.print_term_as_indices(in.term(), out);
lra.print_term_as_indices(in.term(), out);
out << " " << lconstraint_kind_string(in.cmp()) << " " << in.rs();
return out;
}
@ -569,14 +567,14 @@ std::ostream & core::print_var(lpvar j, std::ostream & out) const {
print_monic(m_emons[j], out);
}
m_lar_solver.print_column_info(j, out);
lra.print_column_info(j, out);
signed_var jr = m_evars.find(j);
out << "root=";
if (jr.sign()) {
out << "-";
}
out << m_lar_solver.get_variable_name(jr.var()) << "\n";
out << lra.get_variable_name(jr.var()) << "\n";
return out;
}
@ -643,11 +641,11 @@ void core::trace_print_monic_and_factorization(const monic& rm, const factorizat
bool core::var_has_positive_lower_bound(lpvar j) const {
return m_lar_solver.column_has_lower_bound(j) && m_lar_solver.get_lower_bound(j) > lp::zero_of_type<lp::impq>();
return lra.column_has_lower_bound(j) && lra.get_lower_bound(j) > lp::zero_of_type<lp::impq>();
}
bool core::var_has_negative_upper_bound(lpvar j) const {
return m_lar_solver.column_has_upper_bound(j) && m_lar_solver.get_upper_bound(j) < lp::zero_of_type<lp::impq>();
return lra.column_has_upper_bound(j) && lra.get_upper_bound(j) < lp::zero_of_type<lp::impq>();
}
bool core::var_is_separated_from_zero(lpvar j) const {
@ -684,10 +682,10 @@ template bool core::mon_has_zero<unsigned_vector>(const unsigned_vector& product
lp::lp_settings& core::lp_settings() {
return m_lar_solver.settings();
return lra.settings();
}
const lp::lp_settings& core::lp_settings() const {
return m_lar_solver.settings();
return lra.settings();
}
unsigned core::random() { return lp_settings().random_next(); }
@ -695,7 +693,7 @@ unsigned core::random() { return lp_settings().random_next(); }
// we look for octagon constraints here, with a left part +-x +- y
void core::collect_equivs() {
const lp::lar_solver& s = m_lar_solver;
const lp::lar_solver& s = lra;
for (unsigned i = 0; i < s.terms().size(); i++) {
if (!s.term_is_used_as_row(i))
@ -737,7 +735,7 @@ bool core::is_octagon_term(const lp::lar_term& t, bool & sign, lpvar& i, lpvar &
return true;
}
void core::add_equivalence_maybe(const lp::lar_term *t, lpci c0, lpci c1) {
void core::add_equivalence_maybe(const lp::lar_term* t, u_dependency* c0, u_dependency* c1) {
bool sign;
lpvar i, j;
if (!is_octagon_term(*t, sign, i, j))
@ -866,7 +864,7 @@ std::unordered_set<lpvar> core::collect_vars(const lemma& l) const {
}
}
for (auto p : l.expl()) {
const auto& c = m_lar_solver.constraints()[p.ci()];
const auto& c = lra.constraints()[p.ci()];
for (const auto& r : c.coeffs()) {
insert_j(r.second);
}
@ -1125,8 +1123,8 @@ new_lemma& new_lemma::explain_equiv(lpvar a, lpvar b) {
new_lemma& new_lemma::explain_var_separated_from_zero(lpvar j) {
SASSERT(c.var_is_separated_from_zero(j));
if (c.m_lar_solver.column_has_upper_bound(j) &&
(c.m_lar_solver.get_upper_bound(j)< lp::zero_of_type<lp::impq>()))
if (c.lra.column_has_upper_bound(j) &&
(c.lra.get_upper_bound(j)< lp::zero_of_type<lp::impq>()))
explain_existing_upper_bound(j);
else
explain_existing_lower_bound(j);
@ -1136,7 +1134,7 @@ new_lemma& new_lemma::explain_var_separated_from_zero(lpvar j) {
new_lemma& new_lemma::explain_existing_lower_bound(lpvar j) {
SASSERT(c.has_lower_bound(j));
lp::explanation ex;
ex.push_back(c.m_lar_solver.get_column_lower_bound_witness(j));
c.lra.push_explanation(c.lra.get_column_lower_bound_witness(j), ex);
*this &= ex;
TRACE("nla_solver", tout << j << ": " << *this << "\n";);
return *this;
@ -1145,7 +1143,7 @@ new_lemma& new_lemma::explain_existing_lower_bound(lpvar j) {
new_lemma& new_lemma::explain_existing_upper_bound(lpvar j) {
SASSERT(c.has_upper_bound(j));
lp::explanation ex;
ex.push_back(c.m_lar_solver.get_column_upper_bound_witness(j));
c.lra.push_explanation(c.lra.get_column_upper_bound_witness(j), ex);
*this &= ex;
return *this;
}
@ -1155,7 +1153,7 @@ std::ostream& new_lemma::display(std::ostream & out) const {
for (auto p : lemma.expl()) {
out << "(" << p.ci() << ") ";
c.m_lar_solver.constraints().display(out, [this](lpvar j) { return c.var_str(j);}, p.ci());
c.lra.constraints().display(out, [this](lpvar j) { return c.var_str(j);}, p.ci());
}
out << " ==> ";
if (lemma.ineqs().empty()) {
@ -1303,7 +1301,7 @@ bool core::has_real(const monic& m) const {
bool core::is_patch_blocked(lpvar u, const lp::impq& ival) const {
TRACE("nla_solver", tout << "u = " << u << '\n';);
if (m_cautious_patching &&
(!m_lar_solver.inside_bounds(u, ival) || (var_is_int(u) && ival.is_int() == false))) {
(!lra.inside_bounds(u, ival) || (var_is_int(u) && ival.is_int() == false))) {
TRACE("nla_solver", tout << "u = " << u << " blocked, for feas or integr\n";);
return true; // block
}
@ -1334,7 +1332,7 @@ bool core::is_patch_blocked(lpvar u, const lp::impq& ival) const {
bool core::try_to_patch(const rational& v) {
auto is_blocked = [this](lpvar u, const lp::impq& iv) { return is_patch_blocked(u, iv); };
auto change_report = [this](lpvar u) { update_to_refine_of_var(u); };
return m_lar_solver.try_to_patch(m_patched_var, v, is_blocked, change_report);
return lra.try_to_patch(m_patched_var, v, is_blocked, change_report);
}
bool in_power(const svector<lpvar>& vs, unsigned l) {
@ -1343,7 +1341,7 @@ bool in_power(const svector<lpvar>& vs, unsigned l) {
}
bool core::to_refine_is_correct() const {
for (unsigned j = 0; j < m_lar_solver.number_of_vars(); j++) {
for (unsigned j = 0; j < lra.number_of_vars(); j++) {
if (!is_monic_var(j)) continue;
bool valid = check_monic(emons()[j]);
if (valid == m_to_refine.contains(j)) {
@ -1435,17 +1433,17 @@ void core::patch_monomials() {
NOT_IMPLEMENTED_YET();
m_cautious_patching = false;
patch_monomials_on_to_refine();
m_lar_solver.push();
lra.push();
save_tableau();
constrain_nl_in_tableau();
if (solve_tableau() && integrality_holds()) {
m_lar_solver.pop(1);
lra.pop(1);
} else {
m_lar_solver.pop();
lra.pop();
restore_tableau();
m_lar_solver.clear_inf_heap();
lra.clear_inf_heap();
}
SASSERT(m_lar_solver.ax_is_correct());
SASSERT(lra.ax_is_correct());
}
void core::constrain_nl_in_tableau() {
@ -1507,11 +1505,11 @@ void core::check_bounded_divisions(vector<lemma>& l_vec) {
lbool core::check(vector<lemma>& l_vec) {
lp_settings().stats().m_nla_calls++;
TRACE("nla_solver", tout << "calls = " << lp_settings().stats().m_nla_calls << "\n";);
m_lar_solver.get_rid_of_inf_eps();
lra.get_rid_of_inf_eps();
m_lemma_vec = &l_vec;
if (!(m_lar_solver.get_status() == lp::lp_status::OPTIMAL ||
m_lar_solver.get_status() == lp::lp_status::FEASIBLE)) {
TRACE("nla_solver", tout << "unknown because of the m_lar_solver.m_status = " << m_lar_solver.get_status() << "\n";);
if (!(lra.get_status() == lp::lp_status::OPTIMAL ||
lra.get_status() == lp::lp_status::FEASIBLE)) {
TRACE("nla_solver", tout << "unknown because of the lra.m_status = " << lra.get_status() << "\n";);
return l_undef;
}
@ -1528,6 +1526,9 @@ lbool core::check(vector<lemma>& l_vec) {
if (l_vec.empty() && !done())
m_monomial_bounds();
if (l_vec.empty() && !done() && improve_bounds())
return l_false;
if (l_vec.empty() && !done() && run_horner)
m_horner.horner_lemmas();
@ -1635,21 +1636,21 @@ bool core::no_lemmas_hold() const {
}
lbool core::test_check(vector<lemma>& l) {
m_lar_solver.set_status(lp::lp_status::OPTIMAL);
lra.set_status(lp::lp_status::OPTIMAL);
return check(l);
}
std::ostream& core::print_terms(std::ostream& out) const {
for (unsigned i = 0; i< m_lar_solver.terms().size(); i++) {
for (unsigned i = 0; i< lra.terms().size(); i++) {
unsigned ext = lp::tv::mask_term(i);
if (!m_lar_solver.var_is_registered(ext)) {
if (!lra.var_is_registered(ext)) {
out << "term is not registered\n";
continue;
}
const lp::lar_term & t = *m_lar_solver.terms()[i];
const lp::lar_term & t = *lra.terms()[i];
out << "term:"; print_term(t, out) << std::endl;
lpvar j = m_lar_solver.external_to_local(ext);
lpvar j = lra.external_to_local(ext);
print_var(j, out);
}
return out;
@ -1677,7 +1678,7 @@ std::ostream& core::print_term( const lp::lar_term& t, std::ostream& out) const
std::unordered_set<lpvar> core::get_vars_of_expr_with_opening_terms(const nex *e ) {
auto ret = get_vars_of_expr(e);
auto & ls = m_lar_solver;
auto & ls = lra;
svector<lpvar> added;
for (auto j : ret) {
added.push_back(j);
@ -1685,7 +1686,7 @@ std::unordered_set<lpvar> core::get_vars_of_expr_with_opening_terms(const nex *e
for (unsigned i = 0; i < added.size(); ++i) {
lpvar j = added[i];
if (ls.column_corresponds_to_term(j)) {
const auto& t = m_lar_solver.get_term(lp::tv::raw(ls.local_to_external(j)));
const auto& t = lra.get_term(lp::tv::raw(ls.local_to_external(j)));
for (auto p : t) {
if (ret.find(p.column()) == ret.end()) {
added.push_back(p.column());
@ -1705,7 +1706,7 @@ bool core::is_nl_var(lpvar j) const {
unsigned core::get_var_weight(lpvar j) const {
unsigned k;
switch (m_lar_solver.get_column_type(j)) {
switch (lra.get_column_type(j)) {
case lp::column_type::fixed:
k = 0;
@ -1734,7 +1735,7 @@ unsigned core::get_var_weight(lpvar j) const {
void core::set_active_vars_weights(nex_creator& nc) {
nc.set_number_of_vars(m_lar_solver.column_count());
nc.set_number_of_vars(lra.column_count());
for (lpvar j : active_var_set())
nc.set_var_weight(j, get_var_weight(j));
}
@ -1744,8 +1745,8 @@ bool core::influences_nl_var(lpvar j) const {
j = lp::tv::unmask_term(j);
if (is_nl_var(j))
return true;
for (const auto & c : m_lar_solver.A_r().m_columns[j]) {
lpvar basic_in_row = m_lar_solver.r_basis()[c.var()];
for (const auto & c : lra.A_r().m_columns[j]) {
lpvar basic_in_row = lra.r_basis()[c.var()];
if (is_nl_var(basic_in_row))
return true;
}
@ -1762,9 +1763,32 @@ void core::set_use_nra_model(bool m) {
void core::collect_statistics(::statistics & st) {
st.update("arith-nla-explanations", m_stats.m_nla_explanations);
st.update("arith-nla-lemmas", m_stats.m_nla_lemmas);
st.update("arith-nra-calls", m_stats.m_nra_calls);
st.update("arith-nra-calls", m_stats.m_nra_calls);
st.update("arith-bounds-improvements", m_stats.m_bounds_improvements);
}
bool core::improve_bounds() {
return false;
uint_set seen;
bool bounds_improved = false;
auto insert = [&](lpvar v) {
if (seen.contains(v))
return;
seen.insert(v);
if (lra.improve_bound(v, false))
bounds_improved = true, m_stats.m_bounds_improvements++;
if (lra.improve_bound(v, true))
bounds_improved = true, m_stats.m_bounds_improvements++;
};
for (auto & m : m_emons) {
insert(m.var());
for (auto v : m.vars())
insert(v);
}
return bounds_improved;
}
} // end of nla

24
src/math/lp/nla_core.h
View file

@ -64,6 +64,7 @@ class core {
unsigned m_nla_explanations;
unsigned m_nla_lemmas;
unsigned m_nra_calls;
unsigned m_bounds_improvements;
stats() { reset(); }
void reset() {
memset(this, 0, sizeof(*this));
@ -79,7 +80,7 @@ class core {
var_eqs<emonics> m_evars;
lp::lar_solver& m_lar_solver;
lp::lar_solver& lra;
reslimit& m_reslim;
smt_params_helper m_params;
std::function<bool(lpvar)> m_relevant;
@ -110,6 +111,9 @@ class core {
void check_weighted(unsigned sz, std::pair<unsigned, std::function<void(void)>>* checks);
// try to improve bounds for variables in monomials.
bool improve_bounds();
public:
// constructor
core(lp::lar_solver& s, params_ref const& p, reslimit&);
@ -142,13 +146,13 @@ public:
bool ineq_holds(const ineq& n) const;
bool lemma_holds(const lemma& l) const;
bool is_monic_var(lpvar j) const { return m_emons.is_monic_var(j); }
const rational& val(lpvar j) const { return m_lar_solver.get_column_value(j).x; }
const rational& val(lpvar j) const { return lra.get_column_value(j).x; }
const rational& var_val(const monic& m) const { return m_lar_solver.get_column_value(m.var()).x; }
const rational& var_val(const monic& m) const { return lra.get_column_value(m.var()).x; }
rational mul_val(const monic& m) const {
rational r(1);
for (lpvar v : m.vars()) r *= m_lar_solver.get_column_value(v).x;
for (lpvar v : m.vars()) r *= lra.get_column_value(v).x;
return r;
}
@ -287,16 +291,16 @@ public:
}
const rational& get_upper_bound(unsigned j) const;
const rational& get_lower_bound(unsigned j) const;
bool has_lower_bound(lp::var_index var, lp::constraint_index& ci, lp::mpq& value, bool& is_strict) const {
return m_lar_solver.has_lower_bound(var, ci, value, is_strict);
bool has_lower_bound(lp::var_index var, u_dependency*& ci, lp::mpq& value, bool& is_strict) const {
return lra.has_lower_bound(var, ci, value, is_strict);
}
bool has_upper_bound(lp::var_index var, lp::constraint_index& ci, lp::mpq& value, bool& is_strict) const {
return m_lar_solver.has_upper_bound(var, ci, value, is_strict);
bool has_upper_bound(lp::var_index var, u_dependency*& ci, lp::mpq& value, bool& is_strict) const {
return lra.has_upper_bound(var, ci, value, is_strict);
}
bool zero_is_an_inner_point_of_bounds(lpvar j) const;
bool var_is_int(lpvar j) const { return m_lar_solver.column_is_int(j); }
bool var_is_int(lpvar j) const { return lra.column_is_int(j); }
int rat_sign(const monic& m) const;
inline int rat_sign(lpvar j) const { return nla::rat_sign(val(j)); }
@ -338,7 +342,7 @@ public:
bool is_octagon_term(const lp::lar_term& t, bool & sign, lpvar& i, lpvar &j) const;
void add_equivalence_maybe(const lp::lar_term *t, lpci c0, lpci c1);
void add_equivalence_maybe(const lp::lar_term* t, u_dependency* c0, u_dependency* c1);
void init_vars_equivalence();

9
src/math/lp/nla_defs.h
View file

@ -14,11 +14,10 @@
#include "math/lp/explanation.h"
namespace nla {
typedef lp::constraint_index lpci;
typedef lp::lconstraint_kind llc;
typedef lp::constraint_index lpci;
typedef lp::explanation expl_set;
typedef lp::var_index lpvar;
typedef lp::constraint_index lpci;
typedef lp::lconstraint_kind llc;
typedef lp::explanation expl_set;
typedef lp::var_index lpvar;
struct from_index_dummy{};
class signed_var {

41
src/math/lp/nla_grobner.cpp
View file

@ -22,9 +22,9 @@ namespace nla {
grobner::grobner(core* c):
common(c),
m_pdd_manager(m_core.m_lar_solver.number_of_vars()),
m_solver(m_core.m_reslim, m_pdd_manager),
m_lar_solver(m_core.m_lar_solver),
m_pdd_manager(m_core.lra.number_of_vars()),
m_solver(m_core.m_reslim, m_core.lra.dep_manager(), m_pdd_manager),
lra(m_core.lra),
m_quota(m_core.params().arith_nl_gr_q())
{}
@ -211,12 +211,12 @@ namespace nla {
TRACE("grobner",
tout << "base vars: ";
for (lpvar j : c().active_var_set())
if (m_lar_solver.is_base(j))
if (lra.is_base(j))
tout << "j" << j << " ";
tout << "\n");
for (lpvar j : c().active_var_set()) {
if (m_lar_solver.is_base(j))
add_row(m_lar_solver.basic2row(j));
if (lra.is_base(j))
add_row(lra.basic2row(j));
if (c().is_monic_var(j) && c().var_is_fixed(j))
add_fixed_monic(j);
@ -267,12 +267,12 @@ namespace nla {
}
}
for (unsigned j = 0; j < m_lar_solver.number_of_vars(); ++j) {
if (m_lar_solver.column_has_lower_bound(j) || m_lar_solver.column_has_upper_bound(j)) {
for (unsigned j = 0; j < lra.number_of_vars(); ++j) {
if (lra.column_has_lower_bound(j) || lra.column_has_upper_bound(j)) {
out << j << ": [";
if (m_lar_solver.column_has_lower_bound(j)) out << m_lar_solver.get_lower_bound(j);
if (lra.column_has_lower_bound(j)) out << lra.get_lower_bound(j);
out << "..";
if (m_lar_solver.column_has_upper_bound(j)) out << m_lar_solver.get_upper_bound(j);
if (lra.column_has_upper_bound(j)) out << lra.get_upper_bound(j);
out << "]\n";
}
}
@ -343,14 +343,14 @@ namespace nla {
if (c().var_is_fixed(j))
return;
const auto& matrix = m_lar_solver.A_r();
const auto& matrix = lra.A_r();
for (auto & s : matrix.m_columns[j]) {
unsigned row = s.var();
if (m_rows.contains(row))
continue;
m_rows.insert(row);
unsigned k = m_lar_solver.get_base_column_in_row(row);
if (m_lar_solver.column_is_free(k) && k != j)
unsigned k = lra.get_base_column_in_row(row);
if (lra.column_is_free(k) && k != j)
continue;
CTRACE("grobner", matrix.m_rows[row].size() > c().params().arith_nl_grobner_row_length_limit(),
tout << "ignore the row " << row << " with the size " << matrix.m_rows[row].size() << "\n";);
@ -362,11 +362,10 @@ namespace nla {
}
const rational& grobner::val_of_fixed_var_with_deps(lpvar j, u_dependency*& dep) {
unsigned lc, uc;
m_lar_solver.get_bound_constraint_witnesses_for_column(j, lc, uc);
dep = c().m_intervals.mk_join(dep, c().m_intervals.mk_leaf(lc));
dep = c().m_intervals.mk_join(dep, c().m_intervals.mk_leaf(uc));
return m_lar_solver.column_lower_bound(j).x;
auto* d = lra.get_bound_constraint_witnesses_for_column(j);
if (d)
dep = c().m_intervals.mk_join(dep, d);
return lra.column_lower_bound(j).x;
}
dd::pdd grobner::pdd_expr(const rational& coeff, lpvar j, u_dependency*& dep) {
@ -415,7 +414,7 @@ namespace nla {
SASSERT(r.hi().is_val());
v = r.var();
rational val = r.hi().val();
switch (m_lar_solver.get_column_type(v)) {
switch (lra.get_column_type(v)) {
case lp::column_type::lower_bound:
if (val > 0) num_lo++, lo = v, lc = val; else num_hi++, hi = v, hc = val;
break;
@ -508,7 +507,7 @@ namespace nla {
void grobner::display_matrix_of_m_rows(std::ostream & out) const {
const auto& matrix = m_lar_solver.A_r();
const auto& matrix = lra.A_r();
out << m_rows.size() << " rows" << "\n";
out << "the matrix\n";
for (const auto & r : matrix.m_rows)
@ -516,7 +515,7 @@ namespace nla {
}
void grobner::set_level2var() {
unsigned n = m_lar_solver.column_count();
unsigned n = lra.column_count();
unsigned_vector sorted_vars(n), weighted_vars(n);
for (unsigned j = 0; j < n; j++) {
sorted_vars[j] = j;

2
src/math/lp/nla_grobner.h
View file

@ -21,7 +21,7 @@ namespace nla {
class grobner : common {
dd::pdd_manager m_pdd_manager;
dd::solver m_solver;
lp::lar_solver& m_lar_solver;
lp::lar_solver& lra;
indexed_uint_set m_rows;
unsigned m_quota = 0;

38
src/math/lp/nla_intervals.cpp
View file

@ -4,6 +4,11 @@
#include "util/mpq.h"
namespace nla {
intervals::intervals(core* c, reslimit& lim):
m_dep_intervals(c->lra.dep_manager(), lim),
m_core(c) {}
typedef enum dep_intervals::with_deps_t e_with_deps;
const nex* intervals::get_inf_interval_child(const nex_sum& e) const {
@ -173,7 +178,7 @@ lp::lar_term intervals::expression_to_normalized_term(const nex_sum* e, rational
// where m_terms[k] corresponds to the returned lpvar
lpvar intervals::find_term_column(const lp::lar_term & norm_t, rational& a) const {
std::pair<rational, lpvar> a_j;
if (m_core->m_lar_solver.fetch_normalized_term_column(norm_t, a_j)) {
if (m_core->lra.fetch_normalized_term_column(norm_t, a_j)) {
a /= a_j.first;
return a_j.second;
}
@ -206,19 +211,10 @@ void intervals::set_zero_interval_deps_for_mult(interval& a) {
a.m_upper_dep = a.m_lower_dep;
}
u_dependency *intervals::mk_dep(lp::constraint_index ci) {
return m_dep_intervals.mk_leaf(ci);
}
u_dependency *intervals::mk_dep(const lp::explanation& expl) {
u_dependency* intervals::mk_dep(const lp::explanation& expl) {
u_dependency * r = nullptr;
for (auto p : expl) {
if (r == nullptr) {
r = m_dep_intervals.mk_leaf(p.ci());
} else {
r = m_dep_intervals.mk_join(r, m_dep_intervals.mk_leaf(p.ci()));
}
}
for (auto p : expl)
r = m_dep_intervals.mk_join(r, m_dep_intervals.mk_leaf(p.ci()));
return r;
}
@ -249,25 +245,25 @@ std::ostream& intervals::display(std::ostream& out, const interval& i) const {
template <e_with_deps wd>
void intervals::set_var_interval(lpvar v, interval& b) {
TRACE("nla_intervals_details", m_core->print_var(v, tout) << "\n";);
lp::constraint_index ci;
u_dependency* dep = nullptr;
rational val;
bool is_strict;
if (ls().has_lower_bound(v, ci, val, is_strict)) {
if (ls().has_lower_bound(v, dep, val, is_strict)) {
m_dep_intervals.set_lower(b, val);
m_dep_intervals.set_lower_is_open(b, is_strict);
m_dep_intervals.set_lower_is_inf(b, false);
if (wd == e_with_deps::with_deps) b.m_lower_dep = mk_dep(ci);
if (wd == e_with_deps::with_deps) b.m_lower_dep = dep;
}
else {
m_dep_intervals.set_lower_is_open(b, true);
m_dep_intervals.set_lower_is_inf(b, true);
if (wd == e_with_deps::with_deps) b.m_lower_dep = nullptr;
}
if (ls().has_upper_bound(v, ci, val, is_strict)) {
if (ls().has_upper_bound(v, dep, val, is_strict)) {
m_dep_intervals.set_upper(b, val);
m_dep_intervals.set_upper_is_open(b, is_strict);
m_dep_intervals.set_upper_is_inf(b, false);
if (wd == e_with_deps::with_deps) b.m_upper_dep = mk_dep(ci);
if (wd == e_with_deps::with_deps) b.m_upper_dep = dep;
}
else {
m_dep_intervals.set_upper_is_open(b, true);
@ -303,7 +299,7 @@ bool intervals::interval_from_term(const nex& e, scoped_dep_interval& i) {
m_dep_intervals.set<wd>(i, bi);
TRACE("nla_intervals",
m_core->m_lar_solver.print_column_info(j, tout) << "\n";
m_core->lra.print_column_info(j, tout) << "\n";
tout << "a=" << a << ", b=" << b << "\n";
tout << e << ", interval = "; display(tout, i););
return true;
@ -476,9 +472,9 @@ bool intervals::interval_of_expr(const nex* e, unsigned p, scoped_dep_interval&
}
lp::lar_solver& intervals::ls() { return m_core->m_lar_solver; }
lp::lar_solver& intervals::ls() { return m_core->lra; }
const lp::lar_solver& intervals::ls() const { return m_core->m_lar_solver; }
const lp::lar_solver& intervals::ls() const { return m_core->lra; }
} // end of nla namespace

7
src/math/lp/nla_intervals.h
View file

@ -25,16 +25,13 @@ class intervals {
public:
typedef dep_intervals::interval interval;
private:
u_dependency* mk_dep(lp::constraint_index ci);
u_dependency* mk_dep(lp::explanation const&);
lp::lar_solver& ls();
const lp::lar_solver& ls() const;
public:
intervals(core* c, reslimit& lim) :
m_dep_intervals(lim),
m_core(c)
{}
intervals(core* c, reslimit& lim);
dep_intervals& get_dep_intervals() { return m_dep_intervals; }
u_dependency* mk_join(u_dependency* a, u_dependency* b) { return m_dep_intervals.mk_join(a, b); }
u_dependency* mk_leaf(lp::constraint_index ci) { return m_dep_intervals.mk_leaf(ci); }

4
src/math/lp/nra_solver.cpp
View file

@ -186,7 +186,7 @@ struct solver::imp {
for (auto const& eq : eqs)
add_eq(eq);
for (auto const& [v, w] : m_lp2nl) {
auto& ls = m_nla_core.m_lar_solver;
auto& ls = m_nla_core.lra;
if (ls.column_has_lower_bound(v))
add_lb(ls.get_lower_bound(v), w);
if (ls.column_has_upper_bound(v))
@ -209,7 +209,7 @@ struct solver::imp {
IF_VERBOSE(0, verbose_stream() << "check-nra " << r << "\n";
m_nlsat->display(verbose_stream());
for (auto const& [v, w] : m_lp2nl) {
auto& ls = m_nla_core.m_lar_solver;
auto& ls = m_nla_core.lra;
if (ls.column_has_lower_bound(v))
verbose_stream() << w << " >= " << ls.get_lower_bound(v) << "\n";
if (ls.column_has_upper_bound(v))

37
src/math/lp/ul_pair.h
View file

@ -1,25 +1,19 @@
/*++
Copyright (c) 2017 Microsoft Corporation
Module Name:
<name>
Abstract:
<abstract>
justifications for upper or lower bounds
Author:
Lev Nachmanson (levnach)
Revision History:
--*/
#pragma once
#include "util/vector.h"
#include "util/dependency.h"
#include <string>
#include <algorithm>
#include <utility>
@ -48,14 +42,20 @@ inline bool compare(const std::pair<mpq, var_index> & a, const std::pair<mpq, va
}
class ul_pair {
constraint_index m_lower_bound_witness;
constraint_index m_upper_bound_witness;
bool m_associated_with_row;
u_dependency* m_lower_bound_witness = nullptr;
u_dependency* m_upper_bound_witness = nullptr;
bool m_associated_with_row = false;
public:
constraint_index& lower_bound_witness() {return m_lower_bound_witness;}
constraint_index lower_bound_witness() const {return m_lower_bound_witness;}
constraint_index& upper_bound_witness() { return m_upper_bound_witness;}
constraint_index upper_bound_witness() const {return m_upper_bound_witness;}
// TODO - seems more straight-forward to just expose ul_pair as a struct with direct access to attributes.
u_dependency*& lower_bound_witness() { return m_lower_bound_witness; }
u_dependency* lower_bound_witness() const { return m_lower_bound_witness; }
u_dependency*& upper_bound_witness() { return m_upper_bound_witness; }
u_dependency* upper_bound_witness() const { return m_upper_bound_witness; }
// equality is used by stackedvector operations.
// this appears to be a low level reason
bool operator!=(const ul_pair & p) const {
return !(*this == p);
}
@ -66,14 +66,9 @@ public:
m_associated_with_row == p.m_associated_with_row;
}
// empty constructor
ul_pair() :
m_lower_bound_witness(UINT_MAX),
m_upper_bound_witness(UINT_MAX),
m_associated_with_row(false) {}
ul_pair() {}
ul_pair(bool associated_with_row) :
m_lower_bound_witness(UINT_MAX),
m_upper_bound_witness(UINT_MAX),
m_associated_with_row(associated_with_row) {}
bool associated_with_row() const { return m_associated_with_row; }

33
src/math/lp/var_eqs.h
View file

@ -28,22 +28,24 @@
namespace nla {
class eq_justification {
lpci m_cs[4];
u_dependency* m_cs[4];
public:
eq_justification(std::initializer_list<lpci> cs) {
eq_justification(std::initializer_list<u_dependency*> cs) {
int i = 0;
for (lpci c: cs) {
for (auto c: cs) {
m_cs[i++] = c;
}
for (; i < 4; i++) {
m_cs[i] = -1;
m_cs[i] = nullptr;
}
}
void explain(lp::explanation& e) const {
for (lpci c : m_cs)
if (c + 1 != 0) // c != -1
e.push_back(c);
u_dependency* const* begin() const { return m_cs; }
u_dependency* const* end() const {
unsigned i = 0;
for (; i < 4 && m_cs[i]; ++i)
;
return m_cs + i;
}
};
@ -202,7 +204,7 @@ public:
}
for (eq_justification const& j : m_justtrail) {
j.explain(e);
explain_eq(j, e);
}
m_stats.m_num_explains += m_justtrail.size();
m_stats.m_num_explain_calls++;
@ -216,6 +218,17 @@ public:
// IF_VERBOSE(2, verbose_stream() << (double)m_stats.m_num_explains / m_stats.m_num_explain_calls << "\n");
}
void explain_eq(eq_justification const& eq, lp::explanation& e) const {
u_dependency_manager dm;
unsigned_vector deps;
for (auto* dep : eq) {
deps.reset();
dm.linearize(dep, deps);
for (auto ci : deps)
e.push_back(ci);
}
}
void explain_bfs(signed_var v1, signed_var v2, lp::explanation& e) const {
SASSERT(find(v1) == find(v2));
if (v1 == v2) {
@ -249,7 +262,7 @@ public:
}
while (head != 0) {
m_justtrail[head].explain(e);
explain_eq(m_justtrail[head], e);
head = m_todo[head].m_index;
++m_stats.m_num_explains;
}

3
src/sat/sat_anf_simplifier.cpp
View file

@ -42,7 +42,8 @@ namespace sat {
void anf_simplifier::operator()() {
dd::pdd_manager m(20, dd::pdd_manager::semantics::mod2_e);
pdd_solver solver(s.rlimit(), m);
u_dependency_manager dm;
pdd_solver solver(s.rlimit(), dm, m);
report _report(*this);
configure_solver(solver);
clauses2anf(solver);

6
src/sat/smt/arith_axioms.cpp
View file

@ -524,7 +524,7 @@ namespace arith {
return all_divs_valid;
}
void solver::fixed_var_eh(theory_var v, lp::constraint_index ci1, lp::constraint_index ci2, rational const& bound) {
void solver::fixed_var_eh(theory_var v, u_dependency* dep, rational const& bound) {
theory_var w = euf::null_theory_var;
enode* x = var2enode(v);
if (bound.is_zero())
@ -540,8 +540,8 @@ namespace arith {
return;
reset_evidence();
m_explanation.clear();
consume(rational::one(), ci1);
consume(rational::one(), ci2);
for (auto ci : lp().flatten(dep))
consume(rational::one(), ci);
++m_stats.m_fixed_eqs;
auto* hint = explain_implied_eq(m_explanation, x, y);
auto* jst = euf::th_explain::propagate(*this, m_core, m_eqs, x, y, hint);

46
src/sat/smt/arith_solver.cpp
View file

@ -402,12 +402,13 @@ namespace arith {
}
void solver::propagate_eqs(lp::tv t, lp::constraint_index ci1, lp::lconstraint_kind k, api_bound& b, rational const& value) {
lp::constraint_index ci2;
if (k == lp::GE && set_lower_bound(t, ci1, value) && has_upper_bound(t.index(), ci2, value)) {
fixed_var_eh(b.get_var(), ci1, ci2, value);
u_dependency* dep;
auto& dm = lp().dep_manager();
if (k == lp::GE && set_lower_bound(t, ci1, value) && has_upper_bound(t.index(), dep, value)) {
fixed_var_eh(b.get_var(), dm.mk_join(dm.mk_leaf(ci1), dep), value);
}
else if (k == lp::LE && set_upper_bound(t, ci1, value) && has_lower_bound(t.index(), ci2, value)) {
fixed_var_eh(b.get_var(), ci1, ci2, value);
else if (k == lp::LE && set_upper_bound(t, ci1, value) && has_lower_bound(t.index(), dep, value)) {
fixed_var_eh(b.get_var(), dm.mk_join(dm.mk_leaf(ci1), dep), value);
}
}
@ -433,11 +434,12 @@ namespace arith {
// m_solver already tracks bounds on proper variables, but not on terms.
bool is_strict = false;
rational b;
u_dependency* dep = nullptr;
if (is_lower) {
return lp().has_lower_bound(tv.id(), ci, b, is_strict) && !is_strict && b == v;
return lp().has_lower_bound(tv.id(), dep, b, is_strict) && !is_strict && b == v;
}
else {
return lp().has_upper_bound(tv.id(), ci, b, is_strict) && !is_strict && b == v;
return lp().has_upper_bound(tv.id(), dep, b, is_strict) && !is_strict && b == v;
}
}
}
@ -692,7 +694,7 @@ namespace arith {
void solver::report_equality_of_fixed_vars(unsigned vi1, unsigned vi2) {
rational bound;
lp::constraint_index ci1, ci2, ci3, ci4;
u_dependency* ci1 = nullptr, *ci2 = nullptr, *ci3 = nullptr, *ci4 = nullptr;
theory_var v1 = lp().local_to_external(vi1);
theory_var v2 = lp().local_to_external(vi2);
TRACE("arith", tout << "fixed: " << mk_pp(var2expr(v1), m) << " " << mk_pp(var2expr(v2), m) << "\n";);
@ -714,10 +716,10 @@ namespace arith {
++m_stats.m_fixed_eqs;
reset_evidence();
m_explanation.clear();
consume(rational::one(), ci1);
consume(rational::one(), ci2);
consume(rational::one(), ci3);
consume(rational::one(), ci4);
auto& dm = lp().dep_manager();
auto* d = dm.mk_join(dm.mk_join(ci1, ci2), dm.mk_join(ci3, ci4));
for (auto ci : lp().flatten(d))
consume(rational::one(), ci);
enode* x = var2enode(v1);
enode* y = var2enode(v2);
auto* ex = explain_implied_eq(m_explanation, x, y);
@ -729,26 +731,28 @@ namespace arith {
return x == y || var2enode(x)->get_root() == var2enode(y)->get_root();
}
bool solver::has_upper_bound(lpvar vi, lp::constraint_index& ci, rational const& bound) { return has_bound(vi, ci, bound, false); }
bool solver::has_upper_bound(lpvar vi, u_dependency*& ci, rational const& bound) { return has_bound(vi, ci, bound, false); }
bool solver::has_lower_bound(lpvar vi, lp::constraint_index& ci, rational const& bound) { return has_bound(vi, ci, bound, true); }
bool solver::has_lower_bound(lpvar vi, u_dependency*& ci, rational const& bound) { return has_bound(vi, ci, bound, true); }
bool solver::has_bound(lpvar vi, lp::constraint_index& ci, rational const& bound, bool is_lower) {
bool solver::has_bound(lpvar vi, u_dependency*& dep, rational const& bound, bool is_lower) {
if (lp::tv::is_term(vi)) {
theory_var v = lp().local_to_external(vi);
rational val;
TRACE("arith", tout << lp().get_variable_name(vi) << " " << v << "\n";);
if (v != euf::null_theory_var && a.is_numeral(var2expr(v), val) && bound == val) {
ci = UINT_MAX;
dep = nullptr;
return bound == val;
}
auto& vec = is_lower ? m_lower_terms : m_upper_terms;
lpvar ti = lp::tv::unmask_term(vi);
if (vec.size() > ti) {
constraint_bound& b = vec[ti];
ci = b.first;
return ci != UINT_MAX && bound == b.second;
auto& [ci, coeff] = vec[ti];
if (ci == UINT_MAX)
return false;
dep = lp().dep_manager().mk_leaf(ci);
return bound == coeff;
}
else {
return false;
@ -758,10 +762,10 @@ namespace arith {
bool is_strict = false;
rational b;
if (is_lower) {
return lp().has_lower_bound(vi, ci, b, is_strict) && b == bound && !is_strict;
return lp().has_lower_bound(vi, dep, b, is_strict) && b == bound && !is_strict;
}
else {
return lp().has_upper_bound(vi, ci, b, is_strict) && b == bound && !is_strict;
return lp().has_upper_bound(vi, dep, b, is_strict) && b == bound && !is_strict;
}
}
}

8
src/sat/smt/arith_solver.h
View file

@ -359,7 +359,7 @@ namespace arith {
void assert_idiv_mod_axioms(theory_var u, theory_var v, theory_var w, rational const& r);
api_bound* mk_var_bound(sat::literal lit, theory_var v, lp_api::bound_kind bk, rational const& bound);
lp::lconstraint_kind bound2constraint_kind(bool is_int, lp_api::bound_kind bk, bool is_true);
void fixed_var_eh(theory_var v1, lp::constraint_index ci1, lp::constraint_index ci2, rational const& bound);
void fixed_var_eh(theory_var v1, u_dependency* dep, rational const& bound);
bool set_upper_bound(lp::tv t, lp::constraint_index ci, rational const& v) { return set_bound(t, ci, v, false); }
bool set_lower_bound(lp::tv t, lp::constraint_index ci, rational const& v) { return set_bound(t, ci, v, true); }
bool set_bound(lp::tv tv, lp::constraint_index ci, rational const& v, bool is_lower);
@ -412,9 +412,9 @@ namespace arith {
nlsat::anum const& nl_value(theory_var v, scoped_anum& r) const;
bool has_bound(lpvar vi, lp::constraint_index& ci, rational const& bound, bool is_lower);
bool has_lower_bound(lpvar vi, lp::constraint_index& ci, rational const& bound);
bool has_upper_bound(lpvar vi, lp::constraint_index& ci, rational const& bound);
bool has_bound(lpvar vi, u_dependency*& ci, rational const& bound, bool is_lower);
bool has_lower_bound(lpvar vi, u_dependency*& ci, rational const& bound);
bool has_upper_bound(lpvar vi, u_dependency*& ci, rational const& bound);
/*
* Facility to put a small box around integer variables used in branch and bounds.

58
src/smt/theory_lra.cpp
View file

@ -2863,7 +2863,7 @@ public:
lp::lar_term const& term = lp().get_term(ti);
for (auto const mono : term) {
auto wi = lp().column2tv(mono.column());
lp::constraint_index ci;
u_dependency* ci = nullptr;
rational value;
bool is_strict;
if (wi.is_term()) {
@ -2966,12 +2966,13 @@ public:
vector<constraint_bound> m_upper_terms;
void propagate_eqs(lp::tv t, lp::constraint_index ci1, lp::lconstraint_kind k, api_bound& b, rational const& value) {
lp::constraint_index ci2;
u_dependency* ci2 = nullptr;
auto pair = [&]() { return lp().dep_manager().mk_join(lp().dep_manager().mk_leaf(ci1), ci2); };
if (k == lp::GE && set_lower_bound(t, ci1, value) && has_upper_bound(t.index(), ci2, value)) {
fixed_var_eh(b.get_var(), t, ci1, ci2, value);
fixed_var_eh(b.get_var(), t, pair(), value);
}
else if (k == lp::LE && set_upper_bound(t, ci1, value) && has_lower_bound(t.index(), ci2, value)) {
fixed_var_eh(b.get_var(), t, ci1, ci2, value);
fixed_var_eh(b.get_var(), t, pair(), value);
}
}
@ -3012,11 +3013,12 @@ public:
// m_solver already tracks bounds on proper variables, but not on terms.
bool is_strict = false;
rational b;
u_dependency* dep = nullptr;
if (is_lower) {
return lp().has_lower_bound(tv.id(), ci, b, is_strict) && !is_strict && b == v;
return lp().has_lower_bound(tv.id(), dep, b, is_strict) && !is_strict && b == v;
}
else {
return lp().has_upper_bound(tv.id(), ci, b, is_strict) && !is_strict && b == v;
return lp().has_upper_bound(tv.id(), dep, b, is_strict) && !is_strict && b == v;
}
}
}
@ -3024,35 +3026,37 @@ public:
bool var_has_bound(lpvar vi, bool is_lower) {
bool is_strict = false;
rational b;
lp::constraint_index ci;
u_dependency* dep;
if (is_lower) {
return lp().has_lower_bound(vi, ci, b, is_strict);
return lp().has_lower_bound(vi, dep, b, is_strict);
}
else {
return lp().has_upper_bound(vi, ci, b, is_strict);
return lp().has_upper_bound(vi, dep, b, is_strict);
}
}
bool has_upper_bound(lpvar vi, lp::constraint_index& ci, rational const& bound) { return has_bound(vi, ci, bound, false); }
bool has_upper_bound(lpvar vi, u_dependency*& ci, rational const& bound) { return has_bound(vi, ci, bound, false); }
bool has_lower_bound(lpvar vi, lp::constraint_index& ci, rational const& bound) { return has_bound(vi, ci, bound, true); }
bool has_lower_bound(lpvar vi, u_dependency*& ci, rational const& bound) { return has_bound(vi, ci, bound, true); }
bool has_bound(lpvar vi, lp::constraint_index& ci, rational const& bound, bool is_lower) {
bool has_bound(lpvar vi, u_dependency*& dep, rational const& bound, bool is_lower) {
if (lp::tv::is_term(vi)) {
theory_var v = lp().local_to_external(vi);
rational val;
TRACE("arith", tout << lp().get_variable_name(vi) << " " << v << "\n";);
if (v != null_theory_var && a.is_numeral(get_owner(v), val) && bound == val) {
ci = UINT_MAX;
dep = nullptr;
return bound == val;
}
auto& vec = is_lower ? m_lower_terms : m_upper_terms;
lpvar ti = lp::tv::unmask_term(vi);
if (vec.size() > ti) {
constraint_bound& b = vec[ti];
ci = b.first;
return ci != UINT_MAX && bound == b.second;
auto const& [ci, coeff] = vec[ti];
if (ci == UINT_MAX)
return false;
dep = lp().dep_manager().mk_leaf(ci);
return bound == coeff;
}
else {
return false;
@ -3062,10 +3066,10 @@ public:
bool is_strict = false;
rational b;
if (is_lower) {
return lp().has_lower_bound(vi, ci, b, is_strict) && b == bound && !is_strict;
return lp().has_lower_bound(vi, dep, b, is_strict) && b == bound && !is_strict;
}
else {
return lp().has_upper_bound(vi, ci, b, is_strict) && b == bound && !is_strict;
return lp().has_upper_bound(vi, dep, b, is_strict) && b == bound && !is_strict;
}
}
}
@ -3078,7 +3082,7 @@ public:
void report_equality_of_fixed_vars(unsigned vi1, unsigned vi2) {
rational bound(0);
lp::constraint_index ci1, ci2, ci3, ci4;
u_dependency* ci1 = nullptr, *ci2 = nullptr, *ci3 = nullptr, *ci4 = nullptr;
theory_var v1 = lp().local_to_external(vi1);
theory_var v2 = lp().local_to_external(vi2);
TRACE("arith", tout << "fixed: " << mk_pp(get_owner(v1), m) << " " << mk_pp(get_owner(v2), m) << "\n";);
@ -3130,7 +3134,7 @@ public:
ctx().assign_eq(x, y, eq_justification(js));
}
void fixed_var_eh(theory_var v, lp::tv t, lp::constraint_index ci1, lp::constraint_index ci2, rational const& bound) {
void fixed_var_eh(theory_var v, lp::tv t, u_dependency* dep, rational const& bound) {
theory_var w = null_theory_var;
enode* x = get_enode(v);
if (m_value2var.find(bound, w))
@ -3147,8 +3151,7 @@ public:
if (x->get_root() == y->get_root())
return;
reset_evidence();
set_evidence(ci1, m_core, m_eqs);
set_evidence(ci2, m_core, m_eqs);
set_evidence(dep, m_core, m_eqs);
++m_stats.m_fixed_eqs;
assign_eq(v, w);
}
@ -3182,6 +3185,11 @@ public:
// lp::constraint_index const null_constraint_index = UINT_MAX; // not sure what a correct fix is
void set_evidence(u_dependency* dep, literal_vector& core, svector<enode_pair>& eqs) {
for (auto ci : lp().flatten(dep))
set_evidence(ci, core, eqs);
}
void set_evidence(lp::constraint_index idx, literal_vector& core, svector<enode_pair>& eqs) {
if (idx == UINT_MAX) {
return;
@ -3402,7 +3410,7 @@ public:
if (!is_registered_var(v))
return false;
lpvar vi = get_lpvar(v);
lp::constraint_index ci;
u_dependency* ci;
return lp().has_lower_bound(vi, ci, val, is_strict);
}
@ -3421,8 +3429,8 @@ public:
if (!is_registered_var(v))
return false;
lpvar vi = get_lpvar(v);
lp::constraint_index ci;
return lp().has_upper_bound(vi, ci, val, is_strict);
u_dependency* dep = nullptr;
return lp().has_upper_bound(vi, dep, val, is_strict);
}

7
src/test/pdd_solver.cpp
View file

@ -20,13 +20,15 @@ namespace dd {
}
void test1() {
pdd_manager m(4);
u_dependency_manager dm;
reslimit lim;
pdd v0 = m.mk_var(0);
pdd v1 = m.mk_var(1);
pdd v2 = m.mk_var(2);
pdd v3 = m.mk_var(3);
solver gb(lim, m);
solver gb(lim, dm, m);
gb.add(v1*v2 + v1*v3);
gb.add(v1 - 1);
gb.display(std::cout);
@ -198,10 +200,11 @@ namespace dd {
void test_simplify(expr_ref_vector& fmls, bool use_mod2) {
ast_manager& m = fmls.get_manager();
unsigned_vector id2var;
u_dependency_manager dm;
collect_id2var(id2var, fmls);
pdd_manager p(id2var.size(), use_mod2 ? pdd_manager::mod2_e : pdd_manager::zero_one_vars_e);
solver g(m.limit(), p);
solver g(m.limit(), dm, p);
for (expr* e : subterms::ground(fmls)) {
add_def(id2var, to_app(e), m, p, g);

15
src/util/dependency.h
View file

@ -69,7 +69,7 @@ private:
value_manager & m_vmanager;
allocator & m_allocator;
ptr_vector<dependency> m_todo;
mutable ptr_vector<dependency> m_todo;
void inc_ref(value const & v) {
if (C::ref_count)
@ -106,12 +106,9 @@ private:
}
}
void unmark_todo() {
typename ptr_vector<dependency>::iterator it = m_todo.begin();
typename ptr_vector<dependency>::iterator end = m_todo.end();
for (; it != end; ++it) {
(*it)->unmark();
}
void unmark_todo() const {
for (auto* d : m_todo)
d->unmark();
m_todo.reset();
}
@ -193,7 +190,7 @@ public:
return false;
}
void linearize(dependency * d, vector<value, false> & vs) {
void linearize(dependency * d, vector<value, false> & vs) const {
if (d) {
m_todo.reset();
d->mark();
@ -300,7 +297,7 @@ public:
return m_dep_manager.contains(d, v);
}
void linearize(dependency * d, vector<value, false> & vs) {
void linearize(dependency * d, vector<value, false> & vs) const {
return m_dep_manager.linearize(d, vs);
}