Skip to content

Commit

Permalink
randomise branches in int_solver
Browse files Browse the repository at this point in the history
Signed-off-by: Lev Nachmanson <[email protected]>
  • Loading branch information
levnach committed Jan 31, 2020
1 parent d2fc88d commit abb13f8
Show file tree
Hide file tree
Showing 3 changed files with 109 additions and 90 deletions.
48 changes: 36 additions & 12 deletions src/math/lp/gomory.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -21,6 +21,8 @@
#include "math/lp/int_solver.h"
#include "math/lp/lar_solver.h"
#include "math/lp/lp_utils.h"

#define SMALL_CUTS 0
namespace lp {

class gomory::imp {
Expand All @@ -35,7 +37,11 @@ class gomory::imp {
mpq m_one_minus_f;
mpq m_fj;
mpq m_one_minus_fj;

#if SMALL_CUTS
mpq m_abs_max, m_big_number;
#endif
struct found_big {};

const impq & get_value(unsigned j) const { return m_int_solver.get_value(j); }
bool is_real(unsigned j) const { return m_int_solver.is_real(j); }
bool at_lower(unsigned j) const { return m_int_solver.at_lower(j); }
Expand Down Expand Up @@ -71,6 +77,9 @@ class gomory::imp {
m_t.add_monomial(new_a, j);
m_lcm_den = lcm(m_lcm_den, denominator(new_a));
TRACE("gomory_cut_detail", tout << "new_a = " << new_a << ", k = " << m_k << ", lcm_den = " << m_lcm_den << "\n";);
#if SMALL_CUTS
if (numerator(new_a) > m_big_number) throw found_big();
#endif
}

void real_case_in_gomory_cut(const mpq & a, unsigned j) {
Expand Down Expand Up @@ -100,6 +109,9 @@ class gomory::imp {
}
TRACE("gomory_cut_detail_real", tout << a << "*v" << j << " k: " << m_k << "\n";);
m_t.add_monomial(new_a, j);
#if SMALL_CUTS
if (numerator(new_a) > m_big_number) throw found_big();
#endif
}

lia_move report_conflict_from_gomory_cut() {
Expand Down Expand Up @@ -280,6 +292,14 @@ class gomory::imp {
tout << "1 - m_f: " << 1 - m_f << ", get_value(m_inf_col).x - m_f = " << get_value(m_inf_col).x - m_f << "\n";);
lp_assert(m_f.is_pos() && (get_value(m_inf_col).x - m_f).is_int());

#if SMALL_CUTS
m_abs_max = 0;
for (const auto & p : m_row) {
mpq t = abs(ceil(p.coeff()));
if (t > m_abs_max) m_abs_max = t;
}
m_big_number = m_abs_max; // .expt(2);
#endif
mpq one_min_m_f = 1 - m_f;
for (const auto & p : m_row) {
unsigned j = p.var();
Expand All @@ -290,20 +310,24 @@ class gomory::imp {
}

// use -p.coeff() to make the format compatible with the format used in: Integrating Simplex with DPLL(T)
if (is_real(j)) {
real_case_in_gomory_cut(- p.coeff(), j);
} else {
if (p.coeff().is_int()) {
// m_fj will be zero and no monomial will be added
continue;
try {
if (is_real(j)) {
real_case_in_gomory_cut(- p.coeff(), j);
}
else if (!p.coeff().is_int()) {
some_int_columns = true;
m_fj = fractional_part(-p.coeff());
m_one_minus_fj = 1 - m_fj;
int_case_in_gomory_cut(j);
}
some_int_columns = true;
m_fj = fractional_part(-p.coeff());
m_one_minus_fj = 1 - m_fj;
int_case_in_gomory_cut(j);
}
catch (found_big) {
m_ex->clear();
m_t.clear();
m_k = 1;
return lia_move::undef;
}
}

if (m_t.is_empty())
return report_conflict_from_gomory_cut();
if (some_int_columns)
Expand Down
141 changes: 65 additions & 76 deletions src/math/lp/int_solver.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -12,26 +12,25 @@
namespace lp {


void int_solver::trace_inf_rows() const {
TRACE("arith_int_rows",
unsigned num = m_lar_solver->A_r().column_count();
for (unsigned v = 0; v < num; v++) {
if (column_is_int(v) && !get_value(v).is_int()) {
display_column(tout, v);
}
}
std::ostream& int_solver::display_inf_rows(std::ostream& out) const {
unsigned num = m_lar_solver->A_r().column_count();
for (unsigned v = 0; v < num; v++) {
if (column_is_int(v) && !get_value(v).is_int()) {
display_column(out, v);
}
}

num = 0;
for (unsigned i = 0; i < m_lar_solver->A_r().row_count(); i++) {
unsigned j = m_lar_solver->m_mpq_lar_core_solver.m_r_basis[i];
if (column_is_int_inf(j)) {
num++;
m_lar_solver->print_row(m_lar_solver->A_r().m_rows[i], tout);
tout << "\n";
}
}
tout << "num of int infeasible: " << num << "\n";
);
num = 0;
for (unsigned i = 0; i < m_lar_solver->A_r().row_count(); i++) {
unsigned j = m_lar_solver->m_mpq_lar_core_solver.m_r_basis[i];
if (column_is_int_inf(j)) {
num++;
m_lar_solver->print_row(m_lar_solver->A_r().m_rows[i], out);
out << "\n";
}
}
out << "num of int infeasible: " << num << "\n";
return out;
}

bool int_solver::has_inf_int() const {
Expand Down Expand Up @@ -134,7 +133,6 @@ constraint_index int_solver::column_lower_bound_constraint(unsigned j) const {
}

lia_move int_solver::mk_gomory_cut( unsigned inf_col, const row_strip<mpq> & row) {

lp_assert(column_is_int_inf(inf_col));
gomory gc(m_t, m_k, m_ex, inf_col, row, *this);
return gc.create_cut();
Expand All @@ -157,18 +155,6 @@ unsigned int_solver::row_of_basic_column(unsigned j) const {
return m_lar_solver->row_of_basic_column(j);
}

// template <typename T>
// void int_solver::fill_cut_solver_for_constraint(constraint_index ci, cut_solver<T> & cs) {
// const lar_base_constraint* c = m_lar_solver->constraints()[ci];
// vector<std::pair<T, var_index>> coeffs;
// T rs;
// get_int_coeffs_from_constraint(c, coeffs, rs);
// vector<constraint_index> explanation;
// explanation.push_back(ci);
// cs.add_ineq(coeffs, -rs, explanation);
// }



// this will allow to enable and disable tracking of the pivot rows
struct check_return_helper {
Expand Down Expand Up @@ -240,8 +226,12 @@ bool int_solver::tighten_terms_for_cube() {
return true;
}

bool int_solver::should_find_cube() {
return m_number_of_calls % settings().m_int_find_cube_period == 0;
}

lia_move int_solver::find_cube() {
if (m_number_of_calls % settings().m_int_find_cube_period != 0)
if (!should_find_cube())
return lia_move::undef;

settings().stats().m_cube_calls++;
Expand Down Expand Up @@ -280,8 +270,12 @@ void int_solver::find_feasible_solution() {
lp_assert(lp_status::OPTIMAL == m_lar_solver->get_status() || lp_status::FEASIBLE == m_lar_solver->get_status());
}

bool int_solver::should_run_gcd_test() {
return settings().m_int_run_gcd_test;
}

lia_move int_solver::run_gcd_test() {
if (settings().m_int_run_gcd_test) {
if (should_run_gcd_test()) {
settings().stats().m_gcd_calls++;
TRACE("int_solver", tout << "gcd-test " << settings().stats().m_gcd_calls << "\n";);
if (!gcd_test()) {
Expand All @@ -293,19 +287,19 @@ lia_move int_solver::run_gcd_test() {
return lia_move::undef;
}

bool int_solver::should_gomory_cut() {
return m_number_of_calls % settings().m_int_gomory_cut_period == 0;
}

lia_move int_solver::gomory_cut() {
TRACE("int_solver", tout << "gomory " << m_number_of_calls << "\n";);
if ((m_number_of_calls) % settings().m_int_gomory_cut_period != 0)
if (!should_gomory_cut())
return lia_move::undef;

if (m_lar_solver->move_non_basic_columns_to_bounds()) {
#if Z3DEBUG
lp_status st =
#endif
m_lar_solver->find_feasible_solution();
#if Z3DEBUG
lp_status st = m_lar_solver->find_feasible_solution();
(void)st;
lp_assert(st == lp_status::FEASIBLE || st == lp_status::OPTIMAL);
#endif
}

int j = find_inf_int_base_column();
Expand All @@ -316,7 +310,6 @@ lia_move int_solver::gomory_cut() {
return proceed_with_gomory_cut(j);
}


void int_solver::try_add_term_to_A_for_hnf(unsigned i) {
mpq rs;
const lar_term* t = m_lar_solver->terms()[i];
Expand Down Expand Up @@ -344,7 +337,7 @@ const lp_settings& int_solver::settings() const {

bool int_solver::hnf_has_var_with_non_integral_value() const {
for (unsigned j : m_hnf_cutter.vars())
if (get_value(j).is_int() == false)
if (!get_value(j).is_int())
return true;
return false;
}
Expand Down Expand Up @@ -393,14 +386,22 @@ lia_move int_solver::make_hnf_cut() {
return r;
}

bool int_solver::should_hnf_cut() {
return settings().m_enable_hnf && m_number_of_calls % m_hnf_cut_period == 0;
}

lia_move int_solver::hnf_cut() {
if (!settings().m_enable_hnf) {
return lia_move::undef;
}
if ((m_number_of_calls) % settings().hnf_cut_period() == 0) {
return make_hnf_cut();
lia_move r = lia_move::undef;
if (should_hnf_cut()) {
r = make_hnf_cut();
if (r == lia_move::undef) {
m_hnf_cut_period *= 2;
}
else {
m_hnf_cut_period = settings().hnf_cut_period();
}
}
return lia_move::undef;
return r;
}

lia_move int_solver::check(lp::explanation * e) {
Expand Down Expand Up @@ -441,12 +442,9 @@ lia_move int_solver::branch_or_sat() {

int int_solver::find_any_inf_int_column_basis_first() {
int j = find_inf_int_base_column();
if (j != -1)
return j;
return find_inf_int_nbasis_column();
return j != -1 ? j : find_inf_int_nbasis_column();
}


void int_solver::set_value_for_nbasic_column_ignore_old_values(unsigned j, const impq & new_val) {
lp_assert(!is_base(j));
auto & x = m_lar_solver->m_mpq_lar_core_solver.m_r_x[j];
Expand All @@ -455,8 +453,6 @@ void int_solver::set_value_for_nbasic_column_ignore_old_values(unsigned j, const
m_lar_solver->change_basic_columns_dependend_on_a_given_nb_column(j, delta);
}



void int_solver::patch_nbasic_column(unsigned j, bool patch_only_int_vals) {
auto & lcs = m_lar_solver->m_mpq_lar_core_solver;
impq & val = lcs.m_r_x[j];
Expand Down Expand Up @@ -694,12 +690,13 @@ linear_combination_iterator<mpq> * int_solver::get_column_iterator(unsigned j) {
int_solver::int_solver(lar_solver* lar_slv) :
m_lar_solver(lar_slv),
m_number_of_calls(0),
m_hnf_cutter(settings()) {
m_hnf_cutter(settings()),
m_hnf_cut_period(settings().hnf_cut_period()) {
m_lar_solver->set_int_solver(this);
}


bool int_solver::has_low(unsigned j) const {
bool int_solver::has_lower(unsigned j) const {
switch (m_lar_solver->m_mpq_lar_core_solver.m_column_types()[j]) {
case column_type::fixed:
case column_type::boxed:
Expand All @@ -722,26 +719,22 @@ bool int_solver::has_upper(unsigned j) const {
}


void set_lower(impq & l,
bool & inf_l,
impq const & v ) {
static void set_lower(impq & l, bool & inf_l, impq const & v ) {
if (inf_l || v > l) {
l = v;
inf_l = false;
}
}


void set_upper(impq & u,
bool & inf_u,
impq const & v) {
static void set_upper(impq & u, bool & inf_u, impq const & v) {
if (inf_u || v < u) {
u = v;
inf_u = false;
}
}

bool int_solver::get_freedom_interval_for_column(unsigned j, bool val_is_int, bool & inf_l, impq & l, bool & inf_u, impq & u, mpq & m) {
bool int_solver::get_freedom_interval_for_column(unsigned j, bool val_is_int, bool & inf_l, impq & l, bool & inf_u, impq & u, mpq & m) {
auto & lcs = m_lar_solver->m_mpq_lar_core_solver;
if (lcs.m_r_heading[j] >= 0) // the basic var
return false;
Expand All @@ -753,7 +746,7 @@ void set_upper(impq & u,
l = u = zero_of_type<impq>();
m = mpq(1);

if (has_low(j)) {
if (has_lower(j)) {
set_lower(l, inf_l, lower_bound(j) - xj);
}
if (has_upper(j)) {
Expand Down Expand Up @@ -794,7 +787,7 @@ void set_upper(impq & u,


if (a.is_neg()) {
if (has_low(i)) {
if (has_lower(i)) {
SET_BOUND(set_lower, l, inf_l, a, xi, lcs.m_r_lower_bounds()[i]);
}
if (has_upper(i)) {
Expand All @@ -805,7 +798,7 @@ void set_upper(impq & u,
if (has_upper(i)) {
SET_BOUND(set_lower, l, inf_l, a, xi, lcs.m_r_upper_bounds()[i]);
}
if (has_low(i)) {
if (has_lower(i)) {
SET_BOUND(set_upper, u, inf_u, a, xi, lcs.m_r_lower_bounds()[i]);
}
}
Expand Down Expand Up @@ -839,9 +832,7 @@ bool int_solver::is_real(unsigned j) const {

bool int_solver::value_is_int(unsigned j) const {
return m_lar_solver->column_value_is_int(j);
}


}

bool int_solver::is_feasible() const {
auto & lcs = m_lar_solver->m_mpq_lar_core_solver;
Expand Down Expand Up @@ -1031,18 +1022,16 @@ lia_move int_solver::create_branch_on_column(int j) {
lp_assert(j != -1);
m_t.add_monomial(mpq(1), m_lar_solver->adjust_column_index_to_term_index(j));
if (is_free(j)) {
m_upper = true;
m_upper = random() % 2;
m_k = mpq(0);
} else {
m_upper = left_branch_is_more_narrow_than_right(j);
m_k = m_upper? floor(get_value(j)) : ceil(get_value(j));
}

TRACE("int_solver", tout << "branching v" << j << " = " << get_value(j) << "\n";
display_column(tout, j);
tout << "k = " << m_k << std::endl;

);
TRACE("int_solver",
display_column(tout << "branching v" << j << " = " << get_value(j) << "\n", j);
tout << "k = " << m_k << std::endl;);
return lia_move::branch;

}
Expand Down
Loading

0 comments on commit abb13f8

Please sign in to comment.