diff --git a/src/util/lp/int_solver.cpp b/src/util/lp/int_solver.cpp index 938a13fac..f7844e874 100644 --- a/src/util/lp/int_solver.cpp +++ b/src/util/lp/int_solver.cpp @@ -8,14 +8,21 @@ namespace lean { void int_solver::fix_non_base_columns() { + lean_assert(is_feasible() && inf_int_set_is_correct()); auto & lcs = m_lar_solver->m_mpq_lar_core_solver; + bool change = false; for (unsigned j : lcs.m_r_nbasis) { if (column_is_int_inf(j)) { - set_value(j, floor(lcs.m_r_x[j].x)); + change = true; + set_value_for_nbasic_column(j, floor(lcs.m_r_x[j].x)); } } + if (!change) + return; if (m_lar_solver->find_feasible_solution() == INFEASIBLE) failed(); + init_inf_int_set(); + lean_assert(is_feasible() && inf_int_set_is_correct()); } void int_solver::failed() { @@ -392,7 +399,13 @@ void int_solver::set_value(unsigned j, const impq & new_val) { auto delta = new_val - x; x = new_val; m_lar_solver->change_basic_x_by_delta_on_column(j, delta); - update_column_in_inf_set_set(j); + + auto * it = get_column_iterator(j); + update_column_in_int_inf_set(j); + unsigned i; + while (it->next(i)) + update_column_in_int_inf_set(m_lar_solver->m_mpq_lar_core_solver.m_r_basis[i]); + delete it; } void int_solver::patch_int_infeasible_columns() { @@ -707,7 +720,7 @@ bool int_solver::get_freedom_interval_for_column(unsigned x_j, bool & inf_l, imp else { if (x_i_upper) { impq new_l = x_j_val + ((x_i_val - lcs.m_r_upper_bounds()[x_i]) / a_ij); - set_lower(l, inf_u, new_l); + set_lower(l, inf_l, new_l); if (!inf_l && !inf_u && l == u) break;; } if (x_i_lower) { @@ -726,7 +739,11 @@ bool int_solver::get_freedom_interval_for_column(unsigned x_j, bool & inf_l, imp if (inf_l) tout << "-oo"; else tout << l; tout << "; "; if (inf_u) tout << "oo"; else tout << u; - tout << "]\n";); + tout << "]\n"; + tout << "val = " << get_value(x_j) << "\n"; + ); + lean_assert(inf_l || l <= get_value(x_j)); + lean_assert(inf_u || u >= get_value(x_j)); return true; } diff --git a/src/util/lp/lar_core_solver.hpp b/src/util/lp/lar_core_solver.hpp index 62a5c7887..1e7ad6e63 100644 --- a/src/util/lp/lar_core_solver.hpp +++ b/src/util/lp/lar_core_solver.hpp @@ -175,23 +175,7 @@ int lar_core_solver::column_is_out_of_bounds(unsigned j) { void lar_core_solver::calculate_pivot_row(unsigned i) { - SASSERT(!m_r_solver.use_tableau()); - SASSERT(m_r_solver.m_pivot_row.is_OK()); - m_r_solver.m_pivot_row_of_B_1.clear(); - m_r_solver.m_pivot_row_of_B_1.resize(m_r_solver.m_m()); - m_r_solver.m_pivot_row.clear(); - m_r_solver.m_pivot_row.resize(m_r_solver.m_n()); - if (m_r_solver.m_settings.use_tableau()) { - unsigned basis_j = m_r_solver.m_basis[i]; - for (auto & c : m_r_solver.m_A.m_rows[i]) { - if (c.m_j != basis_j) - m_r_solver.m_pivot_row.set_value(c.get_val(), c.m_j); - } - return; - } - - m_r_solver.calculate_pivot_row_of_B_1(i); - m_r_solver.calculate_pivot_row_when_pivot_row_of_B1_is_ready(i); + m_r_solver.calculate_pivot_row(i); } diff --git a/src/util/lp/lar_solver.cpp b/src/util/lp/lar_solver.cpp index def72734c..eff44d1cd 100644 --- a/src/util/lp/lar_solver.cpp +++ b/src/util/lp/lar_solver.cpp @@ -648,6 +648,10 @@ void lar_solver::change_basic_x_by_delta_on_column(unsigned j, const numeric_pai m_basic_columns_with_changed_cost.insert(bj); } m_mpq_lar_core_solver.m_r_solver.update_column_in_inf_set(bj); + TRACE("change_x_del", + tout << "changed basis column " << bj << ", it is " << + ( m_mpq_lar_core_solver.m_r_solver.column_is_feasible(bj)? "feas":"inf") << std::endl;); + } } else { m_column_buffer.clear(); diff --git a/src/util/lp/lp_core_solver_base.h b/src/util/lp/lp_core_solver_base.h index 81f93941c..3b914600f 100644 --- a/src/util/lp/lp_core_solver_base.h +++ b/src/util/lp/lp_core_solver_base.h @@ -703,5 +703,14 @@ public: const unsigned & iters_with_no_cost_growing() const { return m_iters_with_no_cost_growing; } + + linear_combination_iterator * get_iterator_on_row(unsigned i) { + if (m_settings.use_tableau()) + return new iterator_on_row(m_A.m_rows[i]); + calculate_pivot_row(i); + return new iterator_on_pivot_row(m_pivot_row, m_basis[i]); + } + + void calculate_pivot_row(unsigned i); }; } diff --git a/src/util/lp/lp_core_solver_base.hpp b/src/util/lp/lp_core_solver_base.hpp index d870b0c33..090bc62ef 100644 --- a/src/util/lp/lp_core_solver_base.hpp +++ b/src/util/lp/lp_core_solver_base.hpp @@ -971,25 +971,23 @@ template void lp_core_solver_base::pivot_fixed_v // run over basis and non-basis at the same time indexed_vector w(m_basis.size()); // the buffer unsigned i = 0; // points to basis - unsigned j = 0; // points to nonbasis - - for (; i < m_basis.size() && j < m_nbasis.size(); i++) { - unsigned ii = m_basis[i]; - unsigned jj; + for (; i < m_basis.size(); i++) { + unsigned basic_j = m_basis[i]; - if (get_column_type(ii) != column_type::fixed) continue; - while (j < m_nbasis.size()) { - for (; j < m_nbasis.size(); j++) { - jj = m_nbasis[j]; - if (get_column_type(jj) != column_type::fixed) + if (get_column_type(basic_j) != column_type::fixed) continue; + //todo run over the row here!!!!! call get_iterator_on_row(); + T a; + unsigned j; + auto * it = get_iterator_on_row(i); + while (it->next(a, j)) { + if (j == basic_j) + continue; + if (get_column_type(j) != column_type::fixed) { + if (pivot_column_general(j, basic_j, w)) break; } - if (j >= m_nbasis.size()) - break; - j++; - if (pivot_column_general(jj, ii, w)) - break; - } + } + delete it; } } @@ -1047,4 +1045,26 @@ lp_core_solver_base::infeasibility_cost_is_correct_for_column(unsigned j) } } +template +void lp_core_solver_base::calculate_pivot_row(unsigned i) { + lean_assert(!use_tableau()); + lean_assert(m_pivot_row.is_OK()); + m_pivot_row_of_B_1.clear(); + m_pivot_row_of_B_1.resize(m_m()); + m_pivot_row.clear(); + m_pivot_row.resize(m_n()); + if (m_settings.use_tableau()) { + unsigned basis_j = m_basis[i]; + for (auto & c : m_A.m_rows[i]) { + if (c.m_j != basis_j) + m_pivot_row.set_value(c.get_val(), c.m_j); + } + return; + } + + calculate_pivot_row_of_B_1(i); + calculate_pivot_row_when_pivot_row_of_B1_is_ready(i); +} + + }