Skip to content

Commit

Permalink
Remove non feasible costs (#6653)
Browse files Browse the repository at this point in the history
* before rm lu

Signed-off-by: Lev Nachmanson <levnach@hotmail.com>

* rm get_column_in_lu_mode

* rm lu

Signed-off-by: Lev Nachmanson <levnach@hotmail.com>

* rm lu

Signed-off-by: Lev Nachmanson <levnach@hotmail.com>

* rm_lp

Signed-off-by: Lev Nachmanson <levnach@hotmail.com>

* rm_lu

* rm lu

* rm lu

Signed-off-by: Lev Nachmanson <levnach@hotmail.com>

* rm lu

Signed-off-by: Lev Nachmanson <levnach@hotmail.com>

* rm lu

* rm lu

Signed-off-by: Lev Nachmanson <levnach@hotmail.com>

* rm lu

Signed-off-by: Lev Nachmanson <levnach@hotmail.com>

* rm lu

* cleanup

* rm breakpoints

* rm dealing with doubles

* Revert "rm dealing with doubles"

This reverts commit 547254a.

* rm lu

Signed-off-by: Lev Nachmanson <levnach@hotmail.com>

* rm lu

* rm lu

* rm scaler

* rm square_sparse_matrix

* more cleanup

* rm dead code

* rp precise

* remove many methods dealing with double

Signed-off-by: Lev Nachmanson <levnach@hotmail.com>

* rm lu related fields from lp_core_solver_base.h

* remove dead code

Signed-off-by: Lev Nachmanson <levnach@hotmail.com>

* more dead code removal

* remove more dead code

* more dead code

* rm dead code

* more dead code

* fix lp_tst

* more dead code

* replace lp_assert(false) with UNREACHABLE

* rm non feas costs

Signed-off-by: Lev Nachmanson <levnach@hotmail.com>

* fix the build

Signed-off-by: Lev Nachmanson <levnach@hotmail.com>

---------

Signed-off-by: Lev Nachmanson <levnach@hotmail.com>
  • Loading branch information
levnach committed Mar 28, 2023
1 parent fe348b8 commit 130400d
Show file tree
Hide file tree
Showing 9 changed files with 21 additions and 299 deletions.
1 change: 0 additions & 1 deletion src/math/lp/lar_core_solver_def.h
Original file line number Diff line number Diff line change
Expand Up @@ -40,7 +40,6 @@ void lar_core_solver::prefix_r() {
if (m_r_solver.m_settings.simplex_strategy() != simplex_strategy_enum::tableau_rows) {
m_r_solver.m_costs.resize(m_r_solver.m_n());
m_r_solver.m_d.resize(m_r_solver.m_n());
m_r_solver.set_using_infeas_costs(true);
}
}

Expand Down
23 changes: 4 additions & 19 deletions src/math/lp/lar_solver.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -191,9 +191,11 @@ namespace lp {
stats().m_max_rows = A_r().row_count();
if (strategy_is_undecided())
decide_on_strategy_and_adjust_initial_state();

auto strategy_was = settings().simplex_strategy();
settings().set_simplex_strategy(simplex_strategy_enum::tableau_rows);
m_mpq_lar_core_solver.m_r_solver.m_look_for_feasible_solution_only = true;
auto ret = solve();
settings().set_simplex_strategy(strategy_was);
return ret;
}

Expand Down Expand Up @@ -355,7 +357,6 @@ namespace lp {
m_basic_columns_with_changed_cost.resize(m_mpq_lar_core_solver.m_r_x.size());
move_non_basic_columns_to_bounds(false);
auto& rslv = m_mpq_lar_core_solver.m_r_solver;
rslv.set_using_infeas_costs(false);
lp_assert(costs_are_zeros_for_r_solver());
lp_assert(reduced_costs_are_zeroes_for_r_solver());
rslv.m_costs.resize(A_r().column_count(), zero_of_type<mpq>());
Expand Down Expand Up @@ -490,11 +491,11 @@ namespace lp {
lp_status lar_solver::maximize_term(unsigned j_or_term,
impq& term_max) {
TRACE("lar_solver", print_values(tout););

lar_term term = get_term_to_maximize(j_or_term);
if (term.is_empty()) {
return lp_status::UNBOUNDED;
}

impq prev_value;
auto backup = m_mpq_lar_core_solver.m_r_x;
if (m_mpq_lar_core_solver.m_r_solver.calc_current_x_is_feasible_include_non_basis()) {
Expand Down Expand Up @@ -710,14 +711,6 @@ namespace lp {
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);

if (tableau_with_costs()) {
if (m_mpq_lar_core_solver.m_r_solver.using_infeas_costs()) {
for (unsigned j : m_basic_columns_with_changed_cost)
m_mpq_lar_core_solver.m_r_solver.update_inf_cost_for_column_tableau(j);
lp_assert(m_mpq_lar_core_solver.m_r_solver.reduced_costs_are_correct_tableau());
}
}
}


Expand Down Expand Up @@ -1344,14 +1337,6 @@ namespace lp {
for (unsigned j : became_feas)
m_mpq_lar_core_solver.m_r_solver.remove_column_from_inf_set(j);


if (use_tableau_costs()) {
for (unsigned j : became_feas)
m_mpq_lar_core_solver.m_r_solver.update_inf_cost_for_column_tableau(j);
for (unsigned j : basic_columns_with_changed_cost)
m_mpq_lar_core_solver.m_r_solver.update_inf_cost_for_column_tableau(j);
lp_assert(m_mpq_lar_core_solver.m_r_solver.reduced_costs_are_correct_tableau());
}
}

bool lar_solver::model_is_int_feasible() const {
Expand Down
2 changes: 0 additions & 2 deletions src/math/lp/lp_core_solver_base.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -65,8 +65,6 @@ template bool lp::lp_core_solver_base<lp::mpq, lp::mpq>::pivot_column_tableau(un
template void lp::lp_core_solver_base<lp::mpq, lp::numeric_pair<lp::mpq> >::transpose_rows_tableau(unsigned int, unsigned int);
template bool lp::lp_core_solver_base<lp::mpq, lp::numeric_pair<lp::mpq> >::inf_set_is_correct() const;
template bool lp::lp_core_solver_base<lp::mpq, lp::mpq>::inf_set_is_correct() const;
template bool lp::lp_core_solver_base<lp::mpq, lp::numeric_pair<lp::mpq> >::infeasibility_costs_are_correct() const;
template bool lp::lp_core_solver_base<lp::mpq, lp::mpq >::infeasibility_costs_are_correct() const;
template bool lp::lp_core_solver_base<lp::mpq, lp::numeric_pair<lp::mpq> >::remove_from_basis(unsigned int);


13 changes: 2 additions & 11 deletions src/math/lp/lp_core_solver_base.h
Original file line number Diff line number Diff line change
Expand Up @@ -63,15 +63,12 @@ class lp_core_solver_base {
bool current_x_is_infeasible() const { return m_inf_set.size() != 0; }
private:
u_set m_inf_set;
bool m_using_infeas_costs;
public:
const u_set& inf_set() const { return m_inf_set; }
u_set& inf_set() { return m_inf_set; }
void inf_set_increase_size_by_one() { m_inf_set.increase_size_by_one(); }
bool inf_set_contains(unsigned j) const { return m_inf_set.contains(j); }
unsigned inf_set_size() const { return m_inf_set.size(); }
bool using_infeas_costs() const { return m_using_infeas_costs; }
void set_using_infeas_costs(bool val) { m_using_infeas_costs = val; }
unsigned inf_set_size() const { return m_inf_set.size(); }
indexed_vector<T> m_pivot_row; // this is the real pivot row of the simplex tableu
static_matrix<T, X> & m_A; // the matrix A
// vector<X> const & m_b; // the right side
Expand Down Expand Up @@ -198,11 +195,7 @@ class lp_core_solver_base {
if (m_settings.simplex_strategy() == simplex_strategy_enum::tableau_rows)
return true;
CASSERT("check_static_matrix", m_A.is_correct());
if (m_using_infeas_costs) {
if (infeasibility_costs_are_correct() == false) {
return false;
}
}


unsigned n = m_A.column_count();
for (unsigned j = 0; j < n; j++) {
Expand Down Expand Up @@ -236,8 +229,6 @@ class lp_core_solver_base {
return below_bound(m_x[p], m_lower_bounds[p]);
}

bool infeasibility_costs_are_correct() const;
bool infeasibility_cost_is_correct_for_column(unsigned j) const;

bool x_above_lower_bound(unsigned p) const {
return above_bound(m_x[p], m_lower_bounds[p]);
Expand Down
56 changes: 0 additions & 56 deletions src/math/lp/lp_core_solver_base_def.h
Original file line number Diff line number Diff line change
Expand Up @@ -43,7 +43,6 @@ lp_core_solver_base(static_matrix<T, X> & A,
m_iters_with_no_cost_growing(0),
m_status(lp_status::FEASIBLE),
m_inf_set(A.column_count()),
m_using_infeas_costs(false),
m_pivot_row(A.column_count()),
m_A(A),
m_basis(basis),
Expand Down Expand Up @@ -94,8 +93,6 @@ pivot_to_reduced_costs_tableau(unsigned i, unsigned j) {





// template <typename T, typename X> void lp_core_solver_base<T, X>::
// update_index_of_ed() {
// m_index_of_ed.clear();
Expand Down Expand Up @@ -434,57 +431,4 @@ template <typename T, typename X> bool lp_core_solver_base<T, X>::remove_from_ba
}


template <typename T, typename X> bool
lp_core_solver_base<T, X>::infeasibility_costs_are_correct() const {
if (! this->m_using_infeas_costs)
return true;
lp_assert(costs_on_nbasis_are_zeros());
for (unsigned j :this->m_basis) {
if (!infeasibility_cost_is_correct_for_column(j)) {
TRACE("lar_solver", tout << "incorrect cost for column " << j << std::endl;);
return false;
}
if (!is_zero(m_d[j])) {
TRACE("lar_solver", tout << "non zero inf cost for basis j = " << j << std::endl;);
return false;
}
}
return true;
}

template <typename T, typename X> bool
lp_core_solver_base<T, X>::infeasibility_cost_is_correct_for_column(unsigned j) const {
T r = -one_of_type<T>();

switch (this->m_column_types[j]) {
case column_type::fixed:
case column_type::boxed:
if (this->x_above_upper_bound(j)) {
return (this->m_costs[j] == r);
}
if (this->x_below_low_bound(j)) {
return (this->m_costs[j] == -r);
}
return is_zero(this->m_costs[j]);

case column_type::lower_bound:
if (this->x_below_low_bound(j)) {
return this->m_costs[j] == -r;
}
return is_zero(this->m_costs[j]);

case column_type::upper_bound:
if (this->x_above_upper_bound(j)) {
return this->m_costs[j] == r;
}
return is_zero(this->m_costs[j]);
case column_type::free_column:
return is_zero(this->m_costs[j]);
default:
UNREACHABLE();
return true;
}
}


}
1 change: 0 additions & 1 deletion src/math/lp/lp_primal_core_solver.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -33,6 +33,5 @@ template unsigned lp_primal_core_solver<mpq, mpq>::solve();
template unsigned lp_primal_core_solver<mpq, numeric_pair<mpq> >::solve();
template bool lp::lp_primal_core_solver<lp::mpq, lp::mpq>::update_basis_and_x_tableau(int, int, lp::mpq const&);
template bool lp::lp_primal_core_solver<lp::mpq, lp::numeric_pair<lp::mpq> >::update_basis_and_x_tableau(int, int, lp::numeric_pair<lp::mpq> const&);
template void lp::lp_primal_core_solver<rational, lp::numeric_pair<rational> >::update_inf_cost_for_column_tableau(unsigned);

}
36 changes: 6 additions & 30 deletions src/math/lp/lp_primal_core_solver.h
Original file line number Diff line number Diff line change
Expand Up @@ -334,28 +334,12 @@ class lp_primal_core_solver : public lp_core_solver_base<T, X> {
void update_x_tableau_rows(unsigned entering, unsigned leaving,
const X &delta) {
this->add_delta_to_x(entering, delta);
if (!this->using_infeas_costs()) {
for (const auto &c : this->m_A.m_columns[entering]) {
if (leaving != this->m_basis[c.var()]) {
this->add_delta_to_x_and_track_feasibility(
this->m_basis[c.var()], -delta * this->m_A.get_val(c));
}
}
} else { // using_infeas_costs() == true
lp_assert(this->column_is_feasible(entering));
lp_assert(this->m_costs[entering] == zero_of_type<T>());
// m_d[entering] can change because of the cost change for basic columns.
for (const auto &c : this->m_A.m_columns[entering]) {
unsigned j = this->m_basis[c.var()];
if (j != leaving)
this->add_delta_to_x(j, -delta * this->m_A.get_val(c));
update_inf_cost_for_column_tableau(j);
if (is_zero(this->m_costs[j]))
this->remove_column_from_inf_set(j);
else
this->insert_column_into_inf_set(j);
}
}
for (const auto &c : this->m_A.m_columns[entering]) {
if (leaving != this->m_basis[c.var()]) {
this->add_delta_to_x_and_track_feasibility(
this->m_basis[c.var()], -delta * this->m_A.get_val(c));
}
}
}

void update_basis_and_x_tableau_rows(int entering, int leaving, X const &tt) {
Expand Down Expand Up @@ -437,7 +421,6 @@ class lp_primal_core_solver : public lp_core_solver_base<T, X> {
}

void decide_on_status_when_cannot_find_entering() {
lp_assert(!need_to_switch_costs());
this->set_status(this->current_x_is_feasible() ? lp_status::OPTIMAL
: lp_status::INFEASIBLE);
}
Expand Down Expand Up @@ -628,11 +611,6 @@ class lp_primal_core_solver : public lp_core_solver_base<T, X> {

bool column_is_benefitial_for_entering_basis(unsigned j) const;
void init_infeasibility_costs();

void init_infeasibility_cost_for_column(unsigned j);
T get_infeasibility_cost_for_column(unsigned j) const;
void init_infeasibility_costs_for_changed_basis_only();

void print_column(unsigned j, std::ostream &out);

void print_bound_info_and_x(unsigned j, std::ostream &out);
Expand All @@ -652,8 +630,6 @@ class lp_primal_core_solver : public lp_core_solver_base<T, X> {

void init_run_tableau();
void update_x_tableau(unsigned entering, const X &delta);
void update_inf_cost_for_column_tableau(unsigned j);

// the delta is between the old and the new cost (old - new)
void update_reduced_cost_for_basic_column_cost_change(const T &delta,
unsigned j) {
Expand Down
112 changes: 0 additions & 112 deletions src/math/lp/lp_primal_core_solver_def.h
Original file line number Diff line number Diff line change
Expand Up @@ -258,122 +258,10 @@ template <typename T, typename X> void lp_primal_core_solver<T, X>::find_feas
solve();
}

template <typename T, typename X>
void lp_primal_core_solver<T, X>::init_infeasibility_costs_for_changed_basis_only() {
for (unsigned i : this->m_ed.m_index)
init_infeasibility_cost_for_column(this->m_basis[i]);
this->set_using_infeas_costs(true);
}


template <typename T, typename X>
void lp_primal_core_solver<T, X>::init_infeasibility_costs() {
lp_assert(this->m_x.size() >= this->m_n());
lp_assert(this->m_column_types.size() >= this->m_n());
for (unsigned j = this->m_n(); j--;)
init_infeasibility_cost_for_column(j);
this->set_using_infeas_costs(true);
}

template <typename T, typename X> T
lp_primal_core_solver<T, X>::get_infeasibility_cost_for_column(unsigned j) const {
if (this->m_basis_heading[j] < 0) {
return zero_of_type<T>();
}
T ret;
// j is a basis column
switch (this->m_column_types[j]) {
case column_type::fixed:
case column_type::boxed:
if (this->x_above_upper_bound(j)) {
ret = 1;
} else if (this->x_below_low_bound(j)) {
ret = -1;
} else {
ret = numeric_traits<T>::zero();
}
break;
case column_type::lower_bound:
if (this->x_below_low_bound(j)) {
ret = -1;
} else {
ret = numeric_traits<T>::zero();
}
break;
case column_type::upper_bound:
if (this->x_above_upper_bound(j)) {
ret = 1;
} else {
ret = numeric_traits<T>::zero();
}
break;
case column_type::free_column:
ret = numeric_traits<T>::zero();
break;
default:
UNREACHABLE();
ret = numeric_traits<T>::zero(); // does not matter
break;
}

ret = - ret;

return ret;
}


// changed m_inf_set too!
template <typename T, typename X> void
lp_primal_core_solver<T, X>::init_infeasibility_cost_for_column(unsigned j) {

if (this->m_basis_heading[j] < 0) {
this->m_costs[j] = numeric_traits<T>::zero();
this->remove_column_from_inf_set(j);
return;
}
// j is a basis column
switch (this->m_column_types[j]) {
case column_type::fixed:
case column_type::boxed:
if (this->x_above_upper_bound(j)) {
this->m_costs[j] = 1;
} else if (this->x_below_low_bound(j)) {
this->m_costs[j] = -1;
} else {
this->m_costs[j] = numeric_traits<T>::zero();
}
break;
case column_type::lower_bound:
if (this->x_below_low_bound(j)) {
this->m_costs[j] = -1;
} else {
this->m_costs[j] = numeric_traits<T>::zero();
}
break;
case column_type::upper_bound:
if (this->x_above_upper_bound(j)) {
this->m_costs[j] = 1;
} else {
this->m_costs[j] = numeric_traits<T>::zero();
}
break;
case column_type::free_column:
this->m_costs[j] = numeric_traits<T>::zero();
break;
default:
UNREACHABLE();
break;
}

if (numeric_traits<T>::is_zero(this->m_costs[j])) {
this->remove_column_from_inf_set(j);
} else {
this->insert_column_into_inf_set(j);
}
this->m_costs[j] = - this->m_costs[j];

}


template <typename T, typename X> void lp_primal_core_solver<T, X>::print_column(unsigned j, std::ostream & out) {
out << this->column_name(j) << " " << j << " " << column_type_to_string(this->m_column_type[j]) << " x = " << this->m_x[j] << " " << "c = " << this->m_costs[j];;
Expand Down
Loading

0 comments on commit 130400d

Please sign in to comment.