diff --git a/src/math/lp/lar_core_solver.h b/src/math/lp/lar_core_solver.h index e024b199c1d..52c5b6f1a84 100644 --- a/src/math/lp/lar_core_solver.h +++ b/src/math/lp/lar_core_solver.h @@ -43,12 +43,6 @@ class lar_core_solver { stacked_vector m_r_columns_nz; stacked_vector m_r_rows_nz; - // d - solver fields, for doubles - stacked_vector m_d_pushed_basis; - vector m_d_basis; - vector m_d_nbasis; - vector m_d_heading; - lp_primal_core_solver> m_r_solver; // solver in rational numbers @@ -123,14 +117,6 @@ class lar_core_solver { void fill_not_improvable_zero_sum(); - void pop_basis(unsigned k) { - - m_d_basis = m_r_basis; - m_d_nbasis = m_r_nbasis; - m_d_heading = m_r_heading; - - } - void push() { lp_assert(m_r_solver.basis_heading_is_correct()); lp_assert(m_column_types.size() == m_r_A.column_count()); @@ -180,7 +166,6 @@ class lar_core_solver { m_r_solver.m_costs.resize(m_r_A.column_count()); m_r_solver.m_d.resize(m_r_A.column_count()); - pop_basis(k); m_stacked_simplex_strategy.pop(k); settings().set_simplex_strategy(m_stacked_simplex_strategy); lp_assert(m_r_solver.basis_heading_is_correct()); @@ -356,10 +341,6 @@ class lar_core_solver { return delta; } - void init_column_row_nz_for_r_solver() { - m_r_solver.init_column_row_non_zeroes(); - } - bool column_is_fixed(unsigned j) const { return m_column_types()[j] == column_type::fixed || ( m_column_types()[j] == column_type::boxed && diff --git a/src/math/lp/lar_solver.cpp b/src/math/lp/lar_solver.cpp index 67cf90bd081..500ee23e060 100644 --- a/src/math/lp/lar_solver.cpp +++ b/src/math/lp/lar_solver.cpp @@ -245,11 +245,7 @@ namespace lp { set.erase(j); } - void lar_solver::shrink_inf_set_after_pop(unsigned n, u_set& set) { - clean_popped_elements(n, set); - set.resize(n); - } - + void lar_solver::pop(unsigned k) { TRACE("lar_solver", tout << "k = " << k << std::endl;); @@ -714,11 +710,6 @@ namespace lp { detect_rows_with_changed_bounds_for_column(j); } - void lar_solver::update_x_and_inf_costs_for_columns_with_changed_bounds() { - for (auto j : m_columns_with_changed_bounds) - update_x_and_inf_costs_for_column_with_changed_bounds(j); - } - void lar_solver::update_x_and_inf_costs_for_columns_with_changed_bounds_tableau() { for (auto j : m_columns_with_changed_bounds) update_x_and_inf_costs_for_column_with_changed_bounds(j); @@ -792,31 +783,6 @@ namespace lp { } - void lar_solver::fill_last_row_of_A_r(static_matrix>& A, const lar_term* ls) { - lp_assert(A.row_count() > 0); - lp_assert(A.column_count() > 0); - unsigned last_row = A.row_count() - 1; - lp_assert(A.m_rows[last_row].size() == 0); - for (auto t : *ls) { - lp_assert(!is_zero(t.coeff())); - var_index j = t.column(); - A.set(last_row, j, -t.coeff()); - } - unsigned basis_j = A.column_count() - 1; - A.set(last_row, basis_j, mpq(1)); - } - - template - void lar_solver::copy_from_mpq_matrix(static_matrix& matr) { - matr.m_rows.resize(A_r().row_count()); - matr.m_columns.resize(A_r().column_count()); - for (unsigned i = 0; i < matr.row_count(); i++) { - for (auto& it : A_r().m_rows[i]) { - matr.set(i, it.var(), convert_struct::convert(it.coeff())); - } - } - } - bool lar_solver::all_constrained_variables_are_registered(const vector>& left_side) { for (auto it : left_side) { if (!var_is_registered(it.second)) diff --git a/src/math/lp/lar_solver.h b/src/math/lp/lar_solver.h index 237b9de810a..83d2a25fe29 100644 --- a/src/math/lp/lar_solver.h +++ b/src/math/lp/lar_solver.h @@ -122,7 +122,6 @@ class lar_solver : public column_namer { bool term_is_int(const lar_term * t) const; bool term_is_int(const vector> & coeffs) const; void add_non_basic_var_to_core_fields(unsigned ext_j, bool is_int); - void add_new_var_to_core_fields_for_doubles(bool register_in_basis); void add_new_var_to_core_fields_for_mpq(bool register_in_basis); mpq adjust_bound_for_int(lpvar j, lconstraint_kind&, const mpq&); @@ -131,7 +130,6 @@ class lar_solver : public column_namer { var_index add_term_undecided(const vector> & coeffs); bool term_coeffs_are_ok(const vector> & coeffs); void push_term(lar_term* t); - void add_row_for_term(const lar_term * term, unsigned term_ext_index); void add_row_from_term_no_constraint(const lar_term * term, unsigned term_ext_index); void add_basic_var_to_core_fields(); bool compare_values(impq const& lhs, lconstraint_kind k, const mpq & rhs); @@ -187,7 +185,6 @@ class lar_solver : public column_namer { analyze_new_bounds_on_row_tableau(i, bp); } static void clean_popped_elements(unsigned n, u_set& set); - static void shrink_inf_set_after_pop(unsigned n, u_set & set); bool maximize_term_on_tableau(const lar_term & term, impq &term_max); bool costs_are_zeros_for_r_solver() const; @@ -213,17 +210,10 @@ class lar_solver : public column_namer { void detect_rows_with_changed_bounds_for_column(unsigned j); void detect_rows_with_changed_bounds(); - void update_x_and_inf_costs_for_columns_with_changed_bounds(); void update_x_and_inf_costs_for_columns_with_changed_bounds_tableau(); void solve_with_core_solver(); numeric_pair get_basic_var_value_from_row(unsigned i); bool x_is_correct() const; - void fill_last_row_of_A_r(static_matrix> & A, const lar_term * ls); - template - void create_matrix_A(static_matrix & matr); - template - void copy_from_mpq_matrix(static_matrix & matr); - bool try_to_set_fixed(column_info & ci); bool all_constrained_variables_are_registered(const vector>& left_side); bool all_constraints_hold() const; bool constraint_holds(const lar_base_constraint & constr, std::unordered_map & var_map) const; @@ -231,7 +221,6 @@ class lar_solver : public column_namer { static void register_in_map(std::unordered_map & coeffs, const lar_base_constraint & cn, const mpq & a); static void register_monoid_in_map(std::unordered_map & coeffs, const mpq & a, unsigned j); bool the_left_sides_sum_to_zero(const vector> & evidence) const; - bool the_right_sides_do_not_sum_to_zero(const vector> & evidence); bool explanation_is_correct(explanation&) const; bool inf_explanation_is_correct() const; mpq sum_of_right_sides_of_explanation(explanation &) const; @@ -251,21 +240,16 @@ class lar_solver : public column_namer { void remove_last_column_from_tableau(); void pop_tableau(); void clean_inf_set_of_r_solver_after_pop(); - void shrink_explanation_to_minimum(vector> & explanation) const; inline bool column_value_is_integer(unsigned j) const { return get_column_value(j).is_int(); } bool model_is_int_feasible() const; bool bound_is_integer_for_integer_column(unsigned j, const mpq & right_side) const; inline lar_core_solver & get_core_solver() { return m_mpq_lar_core_solver; } - void catch_up_in_updating_int_solver(); var_index to_column(unsigned ext_j) const; void fix_terms_with_rounded_columns(); - void update_delta_for_terms(const impq & delta, unsigned j, const vector&); - void fill_vars_to_terms(vector> & vars_to_terms); bool remove_from_basis(unsigned); lar_term get_term_to_maximize(unsigned ext_j) const; bool sum_first_coords(const lar_term& t, mpq & val) const; - void collect_rounded_rows_to_fix(); void register_normalized_term(const lar_term&, lpvar); void deregister_normalized_term(const lar_term&); diff --git a/src/math/lp/lp_core_solver_base.h b/src/math/lp/lp_core_solver_base.h index 964e447a125..96f8729a009 100644 --- a/src/math/lp/lp_core_solver_base.h +++ b/src/math/lp/lp_core_solver_base.h @@ -80,7 +80,7 @@ class lp_core_solver_base { vector & m_basis; vector& m_nbasis; vector& m_basis_heading; - vector & m_x; // a feasible solution, the fist time set in the constructor + vector & m_x; // a feasible solution, the first time set in the constructor vector & m_costs; lp_settings & m_settings; diff --git a/src/math/lp/lp_primal_core_solver.h b/src/math/lp/lp_primal_core_solver.h index aacabf94426..fecf66e68cf 100644 --- a/src/math/lp/lp_primal_core_solver.h +++ b/src/math/lp/lp_primal_core_solver.h @@ -49,7 +49,6 @@ class lp_primal_core_solver:public lp_core_solver_base { indexed_vector m_beta; // see Swietanowski working vector beta for column norms T m_epsilon_of_reduced_cost; vector m_costs_backup; - T m_converted_harris_eps; unsigned m_inf_row_index_for_tableau; bool m_bland_mode_tableau; u_set m_left_basis_tableau; @@ -277,13 +276,8 @@ class lp_primal_core_solver:public lp_core_solver_base { return convert_struct::convert(std::numeric_limits::max()); } - bool get_harris_theta(X & theta); - - void zero_harris_eps() { m_converted_harris_eps = zero_of_type(); } - int find_leaving_on_harris_theta(X const & harris_theta, X & t); bool try_jump_to_another_bound_on_entering(unsigned entering, const X & theta, X & t, bool & unlimited); bool try_jump_to_another_bound_on_entering_unlimited(unsigned entering, X & t); - int find_leaving_and_t_precise(unsigned entering, X & t); int find_leaving_and_t_tableau(unsigned entering, X & t); void limit_theta(const X & lim, X & theta, bool & unlimited) { @@ -317,9 +311,6 @@ class lp_primal_core_solver:public lp_core_solver_base { limit_inf_on_bound_m_pos(m, this->m_x[j], this->m_upper_bounds[j], theta, unlimited); }; - X harris_eps_for_bound(const X & bound) const { return ( convert_struct::convert(1) + abs(bound)/10) * m_converted_harris_eps/3; - } - void get_bound_on_variable_and_update_leaving_precisely(unsigned j, vector & leavings, T m, X & t, T & abs_of_d_of_leaving); vector m_lower_bounds_dummy; // needed for the base class only @@ -489,20 +480,9 @@ class lp_primal_core_solver:public lp_core_solver_base { this->set_status(this->current_x_is_feasible()? lp_status::OPTIMAL: lp_status::INFEASIBLE); } - // void limit_theta_on_basis_column_for_feas_case_m_neg(unsigned j, const T & m, X & theta) { - // lp_assert(m < 0); - // lp_assert(this->m_column_type[j] == lower_bound || this->m_column_type[j] == boxed); - // const X & eps = harris_eps_for_bound(this->m_lower_bounds[j]); - // if (this->above_bound(this->m_x[j], this->m_lower_bounds[j])) { - // theta = std::min((this->m_lower_bounds[j] -this->m_x[j] - eps) / m, theta); - // if (theta < zero_of_type()) theta = zero_of_type(); - // } - // } - void limit_theta_on_basis_column_for_feas_case_m_neg_no_check(unsigned j, const T & m, X & theta, bool & unlimited) { lp_assert(m < 0); - const X& eps = harris_eps_for_bound(this->m_lower_bounds[j]); - limit_theta((this->m_lower_bounds[j] - this->m_x[j] - eps) / m, theta, unlimited); + limit_theta((this->m_lower_bounds[j] - this->m_x[j]) / m, theta, unlimited); if (theta < zero_of_type()) theta = zero_of_type(); } @@ -545,25 +525,21 @@ class lp_primal_core_solver:public lp_core_solver_base { void limit_inf_on_upper_bound_m_neg(const T & m, const X & x, const X & bound, X & theta, bool & unlimited) { // x gets smaller lp_assert(m < 0); - const X& eps = harris_eps_for_bound(bound); if (this->above_bound(x, bound)) { - limit_theta((bound - x - eps) / m, theta, unlimited); + limit_theta((bound - x) / m, theta, unlimited); } } void limit_theta_on_basis_column_for_inf_case_m_pos_boxed(unsigned j, const T & m, X & theta, bool & unlimited) { - // lp_assert(m > 0 && this->m_column_type[j] == column_type::boxed); const X & x = this->m_x[j]; const X & lbound = this->m_lower_bounds[j]; if (this->below_bound(x, lbound)) { - const X& eps = harris_eps_for_bound(this->m_upper_bounds[j]); - limit_theta((lbound - x + eps) / m, theta, unlimited); + limit_theta((lbound - x) / m, theta, unlimited); } else { const X & ubound = this->m_upper_bounds[j]; if (this->below_bound(x, ubound)){ - const X& eps = harris_eps_for_bound(ubound); - limit_theta((ubound - x + eps) / m, theta, unlimited); + limit_theta((ubound - x) / m, theta, unlimited); } else if (!this->above_bound(x, ubound)) { theta = zero_of_type(); unlimited = false; @@ -576,13 +552,11 @@ class lp_primal_core_solver:public lp_core_solver_base { const X & x = this->m_x[j]; const X & ubound = this->m_upper_bounds[j]; if (this->above_bound(x, ubound)) { - const X& eps = harris_eps_for_bound(ubound); - limit_theta((ubound - x - eps) / m, theta, unlimited); + limit_theta((ubound - x) / m, theta, unlimited); } else { const X & lbound = this->m_lower_bounds[j]; if (this->above_bound(x, lbound)){ - const X& eps = harris_eps_for_bound(lbound); - limit_theta((lbound - x - eps) / m, theta, unlimited); + limit_theta((lbound - x) / m, theta, unlimited); } else if (!this->below_bound(x, lbound)) { theta = zero_of_type(); unlimited = false; @@ -591,9 +565,8 @@ class lp_primal_core_solver:public lp_core_solver_base { } void limit_theta_on_basis_column_for_feas_case_m_pos(unsigned j, const T & m, X & theta, bool & unlimited) { lp_assert(m > 0); - const T& eps = harris_eps_for_bound(this->m_upper_bounds[j]); if (this->below_bound(this->m_x[j], this->m_upper_bounds[j])) { - limit_theta((this->m_upper_bounds[j] - this->m_x[j] + eps) / m, theta, unlimited); + limit_theta((this->m_upper_bounds[j] - this->m_x[j]) / m, theta, unlimited); if (theta < zero_of_type()) { theta = zero_of_type(); unlimited = false; @@ -603,8 +576,7 @@ class lp_primal_core_solver:public lp_core_solver_base { void limit_theta_on_basis_column_for_feas_case_m_pos_no_check(unsigned j, const T & m, X & theta, bool & unlimited ) { lp_assert(m > 0); - const X& eps = harris_eps_for_bound(this->m_upper_bounds[j]); - limit_theta( (this->m_upper_bounds[j] - this->m_x[j] + eps) / m, theta, unlimited); + limit_theta( (this->m_upper_bounds[j] - this->m_x[j]) / m, theta, unlimited); if (theta < zero_of_type()) { theta = zero_of_type(); } @@ -612,9 +584,9 @@ class lp_primal_core_solver:public lp_core_solver_base { // j is a basic column or the entering, in any case x[j] has to stay feasible. // m is the multiplier. updating t in a way that holds the following - // x[j] + t * m >= this->m_lower_bounds[j]- harris_feasibility_tolerance ( if m < 0 ) + // x[j] + t * m >= this->m_lower_bounds[j]( if m < 0 ) // or - // x[j] + t * m <= this->m_upper_bounds[j] + harris_feasibility_tolerance ( if m > 0) + // x[j] + t * m <= this->m_upper_bounds[j] ( if m > 0) void limit_theta_on_basis_column(unsigned j, T m, X & theta, bool & unlimited) { switch (this->m_column_types[j]) { case column_type::free_column: break; @@ -679,7 +651,6 @@ class lp_primal_core_solver:public lp_core_solver_base { bool column_is_benefitial_for_entering_basis(unsigned j) const; bool column_is_benefitial_for_entering_basis_precise(unsigned j) const; bool can_enter_basis(unsigned j); - bool done(); void init_infeasibility_costs(); void init_infeasibility_cost_for_column(unsigned j); @@ -694,90 +665,8 @@ class lp_primal_core_solver:public lp_core_solver_base { return (a > zero_of_type() && m_sign_of_entering_delta > 0) || (a < zero_of_type() && m_sign_of_entering_delta < 0); } - - bool lower_bounds_are_set() const override { return true; } - void print_bound_info_and_x(unsigned j, std::ostream & out); - - void init_infeasibility_after_update_x_if_inf(unsigned leaving) { - if (this->using_infeas_costs()) { - init_infeasibility_costs_for_changed_basis_only(); - this->m_costs[leaving] = zero_of_type(); - this->remove_column_from_inf_set(leaving); - } - } - void init_inf_set() { - this->clear_inf_set(); - for (unsigned j = 0; j < this->m_n(); j++) { - if (this->m_basis_heading[j] < 0) - continue; - if (!this->column_is_feasible(j)) - this->insert_column_into_inf_set(j); - } - } - - int get_column_out_of_bounds_delta_sign(unsigned j) { - switch (this->m_column_types[j]) { - case column_type::fixed: - case column_type::boxed: - if (this->x_below_low_bound(j)) - return -1; - if (this->x_above_upper_bound(j)) - return 1; - break; - case column_type::lower_bound: - if (this->x_below_low_bound(j)) - return -1; - break; - case column_type::upper_bound: - if (this->x_above_upper_bound(j)) - return 1; - break; - case column_type::free_column: - return 0; - default: - lp_assert(false); - } - return 0; - } - - void init_column_row_non_zeroes() { - this->m_columns_nz.resize(this->m_A.column_count()); - this->m_rows_nz.resize(this->m_A.row_count()); - for (unsigned i = 0; i < this->m_A.column_count(); i++) { - if (this->m_columns_nz[i] == 0) - this->m_columns_nz[i] = this->m_A.m_columns[i].size(); - } - for (unsigned i = 0; i < this->m_A.row_count(); i++) { - if (this->m_rows_nz[i] == 0) - this->m_rows_nz[i] = this->m_A.m_rows[i].size(); - } - } - - - int x_at_bound_sign(unsigned j) { - switch (this->m_column_types[j]) { - case column_type::fixed: - return 0; - case column_type::boxed: - if (this->x_is_at_lower_bound(j)) - return 1; - return -1; - break; - case column_type::lower_bound: - return 1; - break; - case column_type::upper_bound: - return -1; - break; - default: - lp_assert(false); - } - return 0; - - } - unsigned solve_with_tableau(); bool basis_column_is_set_correctly(unsigned j) const { @@ -844,10 +733,6 @@ class lp_primal_core_solver:public lp_core_solver_base { m_beta(A.row_count()), m_epsilon_of_reduced_cost(T(1)/T(10000000)), m_bland_mode_threshold(1000) { - - - m_converted_harris_eps = zero_of_type(); - this->set_status(lp_status::UNKNOWN); } diff --git a/src/math/lp/lp_primal_core_solver_def.h b/src/math/lp/lp_primal_core_solver_def.h index cfd8fc947d6..04c32d41c51 100644 --- a/src/math/lp/lp_primal_core_solver_def.h +++ b/src/math/lp/lp_primal_core_solver_def.h @@ -142,53 +142,6 @@ int lp_primal_core_solver::choose_entering_column(unsigned number_of_benef return choose_entering_column_presize(number_of_benefitial_columns_to_go_over); } - -template bool lp_primal_core_solver::get_harris_theta(X & theta) { - bool unlimited = true; - for (unsigned i : this->m_ed.m_index) { - if (this->m_settings.abs_val_is_smaller_than_pivot_tolerance(this->m_ed[i])) continue; - limit_theta_on_basis_column(this->m_basis[i], - this->m_ed[i] * m_sign_of_entering_delta, theta, unlimited); - if (!unlimited && is_zero(theta)) break; - } - return unlimited; -} - - -template int lp_primal_core_solver:: -find_leaving_on_harris_theta(X const & harris_theta, X & t) { - int leaving = -1; - T pivot_abs_max = zero_of_type(); - // we know already that there is no bound flip on entering - // we also know that harris_theta is limited, so we will find a leaving - zero_harris_eps(); - unsigned steps = this->m_ed.m_index.size(); - unsigned k = this->m_settings.random_next() % steps; - unsigned initial_k = k; - do { - unsigned i = this->m_ed.m_index[k]; - const T & ed = this->m_ed[i]; - if (this->m_settings.abs_val_is_smaller_than_pivot_tolerance(ed)) { - if (++k == steps) - k = 0; - continue; - } - X ratio; - unsigned j = this->m_basis[i]; - bool unlimited = true; - limit_theta_on_basis_column(j, - ed * m_sign_of_entering_delta, ratio, unlimited); - if ((!unlimited) && ratio <= harris_theta) { - if (leaving == -1 || abs(ed) > pivot_abs_max) { - t = ratio; - leaving = j; - pivot_abs_max = abs(ed); - } - } - if (++k == steps) k = 0; - } while (k != initial_k); - return leaving; -} - - template bool lp_primal_core_solver::try_jump_to_another_bound_on_entering(unsigned entering, const X & theta, X & t, @@ -246,68 +199,6 @@ try_jump_to_another_bound_on_entering_unlimited(unsigned entering, X & t ) { return true; } -template int lp_primal_core_solver::find_leaving_and_t_precise(unsigned entering, X & t) { - bool unlimited = true; - unsigned steps = this->m_ed.m_index.size(); - unsigned k = this->m_settings.random_next() % steps; - unsigned initial_k = k; - unsigned row_min_nz = this->m_n() + 1; - m_leaving_candidates.clear(); - do { - unsigned i = this->m_ed.m_index[k]; - const T & ed = this->m_ed[i]; - lp_assert(!numeric_traits::is_zero(ed)); - unsigned j = this->m_basis[i]; - limit_theta_on_basis_column(j, - ed * m_sign_of_entering_delta, t, unlimited); - if (!unlimited) { - m_leaving_candidates.push_back(j); - row_min_nz = this->m_rows_nz[i]; - } - if (++k == steps) k = 0; - } while (unlimited && k != initial_k); - if (unlimited) { - if (try_jump_to_another_bound_on_entering_unlimited(entering, t)) - return entering; - return -1; - } - - X ratio; - while (k != initial_k) { - unsigned i = this->m_ed.m_index[k]; - const T & ed = this->m_ed[i]; - lp_assert(!numeric_traits::is_zero(ed)); - unsigned j = this->m_basis[i]; - unlimited = true; - limit_theta_on_basis_column(j, -ed * m_sign_of_entering_delta, ratio, unlimited); - if (unlimited) { - if (++k == steps) k = 0; - continue; - } - unsigned i_nz = this->m_rows_nz[i]; - if (ratio < t) { - t = ratio; - m_leaving_candidates.clear(); - m_leaving_candidates.push_back(j); - row_min_nz = this->m_rows_nz[i]; - } else if (ratio == t && i_nz < row_min_nz) { - m_leaving_candidates.clear(); - m_leaving_candidates.push_back(j); - row_min_nz = this->m_rows_nz[i]; - } else if (ratio == t && i_nz == row_min_nz) { - m_leaving_candidates.push_back(j); - } - if (++k == steps) k = 0; - } - - ratio = t; - unlimited = false; - if (try_jump_to_another_bound_on_entering(entering, t, ratio, unlimited)) { - t = ratio; - return entering; - } - k = this->m_settings.random_next() % m_leaving_candidates.size(); - return m_leaving_candidates[k]; -} @@ -499,17 +390,6 @@ template void lp_primal_core_solver::one_iteratio -template bool lp_primal_core_solver::done() { - if (this->get_status() == lp_status::OPTIMAL) return true; - if (this->get_status() == lp_status::INFEASIBLE) { - return true; - } - if (this->m_iters_with_no_cost_growing >= this->m_settings.max_number_of_iterations_with_no_improvements) { - this->set_status(lp_status::CANCELLED); - return true; - } - return false; -} template void lp_primal_core_solver::init_infeasibility_costs_for_changed_basis_only() { diff --git a/src/math/lp/lp_primal_core_solver_tableau_def.h b/src/math/lp/lp_primal_core_solver_tableau_def.h index 256853be722..d4355fa800e 100644 --- a/src/math/lp/lp_primal_core_solver_tableau_def.h +++ b/src/math/lp/lp_primal_core_solver_tableau_def.h @@ -43,6 +43,7 @@ template void lp_primal_core_solver::advance_on_e } advance_on_entering_and_leaving_tableau(entering, leaving, t); } + template int lp_primal_core_solver::choose_entering_column_tableau() { //this moment m_y = cB * B(-1) unsigned number_of_benefitial_columns_to_go_over = get_number_of_non_basic_column_to_try_for_enter(); @@ -85,9 +86,6 @@ template void lp_primal_core_solver::advance_on_e } - - - template unsigned lp_primal_core_solver::solve_with_tableau() { init_run_tableau(); diff --git a/src/test/lp/lp.cpp b/src/test/lp/lp.cpp index c993dbd7ddf..24aa8767a37 100644 --- a/src/test/lp/lp.cpp +++ b/src/test/lp/lp.cpp @@ -56,7 +56,6 @@ #include "math/lp/int_cube.h" #include "math/lp/emonics.h" #include "math/lp/static_matrix.h" -#include "math/lp/dense_matrix.h" bool my_white_space(const char & a) { return a == ' ' || a == '\t'; @@ -431,205 +430,8 @@ void change_basis(unsigned entering, unsigned leaving, vector& basis, nbasis[place_in_non_basis] = leaving; } - - -#ifdef Z3DEBUG -void test_small_lu(lp_settings & settings) { - -} - -#endif - - -void fill_long_row(static_matrix &m, int i) { - int n = m.column_count(); - for (int j = 0; j < n; j ++) { - m (i, (j + i) % n) = j * j; - } -} - - - -void fill_long_row_exp(static_matrix &m, int i) { - int n = m.column_count(); - - for (int j = 0; j < n; j ++) { - m(i, j) = my_random() % 20; - } -} - - - - - - int perm_id = 0; - - - - - - - -void init_b(vector & b, static_matrix & m, vector & x) { - for (unsigned i = 0; i < m.row_count(); i++) { - b.push_back(m.dot_product_with_row(i, x)); - } -} - - -void test_lp_0() { - std::cout << " test_lp_0 " << std::endl; - static_matrix m_(3, 7); - m_(0, 0) = 3; m_(0, 1) = 2; m_(0, 2) = 1; m_(0, 3) = 2; m_(0, 4) = 1; - m_(1, 0) = 1; m_(1, 1) = 1; m_(1, 2) = 1; m_(1, 3) = 1; m_(1, 5) = 1; - m_(2, 0) = 4; m_(2, 1) = 3; m_(2, 2) = 3; m_(2, 3) = 4; m_(2, 6) = 1; - vector x_star(7); - x_star[0] = 225; x_star[1] = 117; x_star[2] = 420; - x_star[3] = x_star[4] = x_star[5] = x_star[6] = 0; - vector b; - init_b(b, m_, x_star); - vector basis(3); - basis[0] = 0; basis[1] = 1; basis[2] = 2; - vector costs(7); - costs[0] = 19; - costs[1] = 13; - costs[2] = 12; - costs[3] = 17; - costs[4] = 0; - costs[5] = 0; - costs[6] = 0; - - vector column_types(7, column_type::lower_bound); - vector upper_bound_values; - lp_settings settings; - simple_column_namer cn; - vector nbasis; - vector heading; - - lp_primal_core_solver lpsolver(m_, b, x_star, basis, nbasis, heading, costs, column_types, upper_bound_values, settings, cn); - - lpsolver.solve(); -} - -void test_lp_1() { - std::cout << " test_lp_1 " << std::endl; - static_matrix m(4, 7); - m(0, 0) = 1; m(0, 1) = 3; m(0, 2) = 1; m(0, 3) = 1; - m(1, 0) = -1; m(1, 2) = 3; m(1, 4) = 1; - m(2, 0) = 2; m(2, 1) = -1; m(2, 2) = 2; m(2, 5) = 1; - m(3, 0) = 2; m(3, 1) = 3; m(3, 2) = -1; m(3, 6) = 1; -#ifdef Z3DEBUG - //print_matrix(m, std::cout); -#endif - vector x_star(7); - x_star[0] = 0; x_star[1] = 0; x_star[2] = 0; - x_star[3] = 3; x_star[4] = 2; x_star[5] = 4; x_star[6] = 2; - - vector basis(4); - basis[0] = 3; basis[1] = 4; basis[2] = 5; basis[3] = 6; - - vector b; - b.push_back(3); - b.push_back(2); - b.push_back(4); - b.push_back(2); - - vector costs(7); - costs[0] = 5; - costs[1] = 5; - costs[2] = 3; - costs[3] = 0; - costs[4] = 0; - costs[5] = 0; - costs[6] = 0; - - - - vector column_types(7, column_type::lower_bound); - vector upper_bound_values; - - std::cout << "calling lp\n"; - lp_settings settings; - simple_column_namer cn; - - vector nbasis; - vector heading; - - lp_primal_core_solver lpsolver(m, b, - x_star, - basis, - nbasis, heading, - costs, - column_types, upper_bound_values, settings, cn); - - lpsolver.solve(); -} - - -void test_lp_primal_core_solver() { - test_lp_0(); - test_lp_1(); -} - - - - -#ifdef Z3DEBUG - -void fill_uniformly(dense_matrix & m, unsigned dim) { - int v = 0; - for (unsigned i = 0; i < dim; i++) { - for (unsigned j = 0; j < dim; j++) { - m.set_elem(i, j, v++); - } - } -} - - - - -void test_dense_matrix() { - dense_matrix d(3, 2); - d.set_elem(0, 0, 1); - d.set_elem(1, 1, 2); - d.set_elem(2, 0, 3); - // print_matrix(d); - - dense_matrix unit(2, 2); - d.set_elem(0, 0, 1); - d.set_elem(1, 1, 1); - - dense_matrix c = d * unit; - - // print_matrix(d); - - dense_matrix perm(3, 3); - perm.set_elem(0, 1, 1); - perm.set_elem(1, 0, 1); - perm.set_elem(2, 2, 1); - auto c1 = perm * d; - // print_matrix(c1); - - - dense_matrix p2(2, 2); - p2.set_elem(0, 1, 1); - p2.set_elem(1, 0, 1); - auto c2 = d * p2; -} - -#endif - - - - - -void lp_solver_test() { - // lp_revised_solver lp_revised; - // lp_revised.get_minimal_solution(); -} - bool get_int_from_args_parser(const char * option, argument_parser & args_parser, unsigned & n) { std::string s = args_parser.get_option_value(option); if (!s.empty()) { @@ -650,9 +452,6 @@ bool get_double_from_args_parser(const char * option, argument_parser & args_par -bool values_are_one_percent_close(double a, double b); - - void get_time_limit_and_max_iters_from_parser(argument_parser & args_parser, unsigned & time_limit); // forward definition @@ -939,10 +738,6 @@ void setup_args_parser(argument_parser & parser) { parser.add_option_with_help_string("-gomory", "gomory"); parser.add_option_with_help_string("-intd", "test integer_domain"); parser.add_option_with_help_string("-xyz_sample", "run a small interactive scenario"); - parser.add_option_with_after_string_with_help("--density", "the percentage of non-zeroes in the matrix below which it is not dense"); - parser.add_option_with_after_string_with_help("--harris_toler", "harris tolerance"); - parser.add_option_with_after_string_with_help("--checklu", "the file name for lu checking"); - parser.add_option_with_after_string_with_help("--partial_pivot", "the partial pivot constant, a number somewhere between 10 and 100"); parser.add_option_with_after_string_with_help("--percent_for_enter", "which percent of columns check for entering column"); parser.add_option_with_help_string("--totalinf", "minimizes the total infeasibility instead of diminishing infeasibility of the rows"); parser.add_option_with_after_string_with_help("--rep_frq", "the report frequency, in how many iterations print the cost and other info "); @@ -1202,71 +997,6 @@ void get_matrix_dimensions(std::ifstream & f, unsigned & m, unsigned & n) { n = atoi(r[1].c_str()); } -void read_row_cols(unsigned i, static_matrix& A, std::ifstream & f) { - do { - std::string line; - getline(f, line); - if (line== "row_end") - break; - auto r = split_and_trim(line); - lp_assert(r.size() == 4); - unsigned j = atoi(r[1].c_str()); - double v = atof(r[3].c_str()); - A.set(i, j, v); - } while (true); -} - -bool read_row(static_matrix & A, std::ifstream & f) { - std::string line; - getline(f, line); - if (static_cast(line.find("row")) == -1) - return false; - auto r = split_and_trim(line); - if (r[0] != "row") - std::cout << "wrong row line" << line << std::endl; - unsigned i = atoi(r[1].c_str()); - read_row_cols(i, A, f); - return true; -} - -void read_rows(static_matrix& A, std::ifstream & f) { - while (read_row(A, f)) {} -} - -void read_basis(vector & basis, std::ifstream & f) { - std::cout << "reading basis" << std::endl; - std::string line; - getline(f, line); - lp_assert(line == "basis_start"); - do { - getline(f, line); - if (line == "basis_end") - break; - unsigned j = atoi(line.c_str()); - basis.push_back(j); - } while (true); -} - -void read_indexed_vector(indexed_vector & v, std::ifstream & f) { - std::string line; - getline(f, line); - lp_assert(line == "vector_start"); - do { - getline(f, line); - if (line == "vector_end") break; - auto r = split_and_trim(line); - unsigned i = atoi(r[0].c_str()); - double val = atof(r[1].c_str()); - v.set_value(val, i); - std::cout << "setting value " << i << " = " << val << std::endl; - } while (true); -} - -void check_lu_from_file(std::string lufile_name) { - lp_assert(false); -} - - void print_st(lp_status status) { std::cout << lp_status_to_string(status) << std::endl; @@ -2298,15 +2028,6 @@ void test_lp_local(int argn, char**argv) { test_bound_propagation(); return finalize(0); } - - - std::string lufile = args_parser.get_option_value("--checklu"); - if (!lufile.empty()) { - check_lu_from_file(lufile); - return finalize(0); - } - - if (args_parser.option_is_used("-tbq")) {