Skip to content

Commit

Permalink
max term
Browse files Browse the repository at this point in the history
Signed-off-by: Lev Nachmanson <levnach@hotmail.com>
  • Loading branch information
levnach committed Apr 24, 2020
1 parent 8921ed5 commit 8c5993d
Show file tree
Hide file tree
Showing 6 changed files with 42 additions and 37 deletions.
2 changes: 1 addition & 1 deletion src/math/lp/lar_core_solver_def.h
Original file line number Diff line number Diff line change
Expand Up @@ -199,7 +199,7 @@ void lar_core_solver::calculate_pivot_row(unsigned i) {
m_r_solver.m_breakpoint_indices_queue.resize(m_r_solver.m_n());
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.m_using_infeas_costs = true;
m_r_solver.set_using_infeas_costs(true);
}
}

Expand Down
17 changes: 10 additions & 7 deletions src/math/lp/lar_solver.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -448,7 +448,7 @@ void lar_solver::prepare_costs_for_r_solver(const lar_term & term) {
if (move_non_basic_columns_to_bounds())
find_feasible_solution();
auto & rslv = m_mpq_lar_core_solver.m_r_solver;
rslv.m_using_infeas_costs = false;
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 @@ -574,16 +574,19 @@ lar_term lar_solver::get_term_to_maximize(unsigned j_or_term) const {
lp_status lar_solver::maximize_term(unsigned j_or_term,
impq &term_max) {
TRACE("lar_solver", print_values(tout););
bool was_feasible = m_mpq_lar_core_solver.m_r_solver.calc_current_x_is_feasible_include_non_basis();
impq prev_value;
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 (was_feasible) {
if (m_mpq_lar_core_solver.m_r_solver.calc_current_x_is_feasible_include_non_basis()) {
prev_value = term.apply(m_mpq_lar_core_solver.m_r_x);
} else {
m_mpq_lar_core_solver.m_r_solver.m_look_for_feasible_solution_only = false;
if (solve() != lp_status::OPTIMAL)
return lp_status::UNBOUNDED;
}

m_mpq_lar_core_solver.m_r_solver.m_look_for_feasible_solution_only = false;
Expand Down Expand Up @@ -618,7 +621,7 @@ lp_status lar_solver::maximize_term(unsigned j_or_term,
if (change) {
term_max = term.apply(m_mpq_lar_core_solver.m_r_x);
}
if (was_feasible && term_max < prev_value) {
if (term_max < prev_value) {
term_max = prev_value;
m_mpq_lar_core_solver.m_r_x = backup;
}
Expand Down Expand Up @@ -822,7 +825,7 @@ void lar_solver::update_x_and_inf_costs_for_columns_with_changed_bounds_tableau(
update_x_and_inf_costs_for_column_with_changed_bounds(j);

if (tableau_with_costs()) {
if (m_mpq_lar_core_solver.m_r_solver.m_using_infeas_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
6 changes: 4 additions & 2 deletions src/math/lp/lp_core_solver_base.h
Original file line number Diff line number Diff line change
Expand Up @@ -52,9 +52,11 @@ class lp_core_solver_base {
}
bool current_x_is_infeasible() const { return m_inf_set.size() != 0; }
u_set m_inf_set;
private:
bool m_using_infeas_costs;


public:
bool using_infeas_costs() const { return m_using_infeas_costs; }
void set_using_infeas_costs(bool val) { m_using_infeas_costs = val; }
vector<unsigned> m_columns_nz; // m_columns_nz[i] keeps an approximate value of non zeroes the i-th column
vector<unsigned> m_rows_nz; // m_rows_nz[i] keeps an approximate value of non zeroes in the i-th row
indexed_vector<T> m_pivot_row_of_B_1; // the pivot row of the reverse of B
Expand Down
8 changes: 4 additions & 4 deletions src/math/lp/lp_primal_core_solver.h
Original file line number Diff line number Diff line change
Expand Up @@ -403,7 +403,7 @@ class lp_primal_core_solver:public lp_core_solver_base<T, X> {
if (this->m_settings.simplex_strategy() == simplex_strategy_enum::tableau_rows)
return false;
// lp_assert(calc_current_x_is_feasible() == current_x_is_feasible());
return this->current_x_is_feasible() == this->m_using_infeas_costs;
return this->current_x_is_feasible() == this->using_infeas_costs();
}


Expand Down Expand Up @@ -445,13 +445,13 @@ class lp_primal_core_solver:public lp_core_solver_base<T, X> {
// this version assumes that the leaving already has the right value, and does not update it
void update_x_tableau_rows(unsigned entering, unsigned leaving, const X& delta) {
this->add_delta_to_x(entering, delta);
if (!this->m_using_infeas_costs) {
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 { // m_using_infeas_costs == true
} 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.
Expand Down Expand Up @@ -818,7 +818,7 @@ class lp_primal_core_solver:public lp_core_solver_base<T, X> {
void print_bound_info_and_x(unsigned j, std::ostream & out);

void init_infeasibility_after_update_x_if_inf(unsigned leaving) {
if (this->m_using_infeas_costs) {
if (this->using_infeas_costs()) {
init_infeasibility_costs_for_changed_basis_only();
this->m_costs[leaving] = zero_of_type<T>();
this->m_inf_set.erase(leaving);
Expand Down
28 changes: 14 additions & 14 deletions src/math/lp/lp_primal_core_solver_def.h
Original file line number Diff line number Diff line change
Expand Up @@ -116,7 +116,7 @@ template <typename T, typename X>
bool lp_primal_core_solver<T, X>::column_is_benefitial_for_entering_basis(unsigned j) const {
if (numeric_traits<T>::precise())
return column_is_benefitial_for_entering_basis_precise(j);
if (this->m_using_infeas_costs && this->m_settings.use_breakpoints_in_feasibility_search)
if (this->using_infeas_costs() && this->m_settings.use_breakpoints_in_feasibility_search)
return column_is_benefitial_for_entering_on_breakpoints(j);
const T& dj = this->m_d[j];
switch (this->m_column_types[j]) {
Expand Down Expand Up @@ -150,7 +150,7 @@ bool lp_primal_core_solver<T, X>::column_is_benefitial_for_entering_basis(unsign
template <typename T, typename X>
bool lp_primal_core_solver<T, X>::column_is_benefitial_for_entering_basis_precise(unsigned j) const {
lp_assert (numeric_traits<T>::precise());
if (this->m_using_infeas_costs && this->m_settings.use_breakpoints_in_feasibility_search)
if (this->using_infeas_costs() && this->m_settings.use_breakpoints_in_feasibility_search)
return column_is_benefitial_for_entering_on_breakpoints(j);
const T& dj = this->m_d[j];
TRACE("lar_solver", tout << "dj=" << dj << "\n";);
Expand Down Expand Up @@ -223,7 +223,7 @@ int lp_primal_core_solver<T, X>::choose_entering_column_presize(unsigned number_
return -1;
unsigned entering = *entering_iter;
m_sign_of_entering_delta = this->m_d[entering] > 0 ? 1 : -1;
if (this->m_using_infeas_costs && this->m_settings.use_breakpoints_in_feasibility_search)
if (this->using_infeas_costs() && this->m_settings.use_breakpoints_in_feasibility_search)
m_sign_of_entering_delta = -m_sign_of_entering_delta;
m_non_basis_list.erase(entering_iter);
m_non_basis_list.push_back(entering);
Expand Down Expand Up @@ -263,7 +263,7 @@ int lp_primal_core_solver<T, X>::choose_entering_column(unsigned number_of_benef
if (entering_iter != m_non_basis_list.end()) {
unsigned entering = *entering_iter;
m_sign_of_entering_delta = this->m_d[entering] > 0? 1 : -1;
if (this->m_using_infeas_costs && this->m_settings.use_breakpoints_in_feasibility_search)
if (this->using_infeas_costs() && this->m_settings.use_breakpoints_in_feasibility_search)
m_sign_of_entering_delta = - m_sign_of_entering_delta;
m_non_basis_list.erase(entering_iter);
m_non_basis_list.push_back(entering);
Expand Down Expand Up @@ -648,7 +648,7 @@ template <typename T, typename X> void lp_primal_core_solver<T, X>::init_run(
init_inf_set();
if (this->current_x_is_feasible() && this->m_look_for_feasible_solution_only)
return;
this->m_using_infeas_costs = false;
this->set_using_infeas_costs(false);
if (this->m_settings.backup_costs)
backup_and_normalize_costs();
m_epsilon_of_reduced_cost = numeric_traits<X>::precise()? zero_of_type<T>(): T(1)/T(10000000);
Expand Down Expand Up @@ -685,7 +685,7 @@ void lp_primal_core_solver<T, X>::advance_on_entering_equal_leaving(int entering
return;
}
}
if (this->m_using_infeas_costs) {
if (this->using_infeas_costs()) {
lp_assert(is_zero(this->m_costs[entering]));
init_infeasibility_costs_for_changed_basis_only();
}
Expand All @@ -700,7 +700,7 @@ void lp_primal_core_solver<T, X>::advance_on_entering_equal_leaving(int entering

template <typename T, typename X>void lp_primal_core_solver<T, X>::advance_on_entering_and_leaving(int entering, int leaving, X & t) {
lp_assert(entering >= 0 && m_non_basis_list.back() == static_cast<unsigned>(entering));
lp_assert(this->m_using_infeas_costs || t >= zero_of_type<X>());
lp_assert(this->using_infeas_costs() || t >= zero_of_type<X>());
lp_assert(leaving >= 0 && entering >= 0);
lp_assert(entering != leaving || !is_zero(t)); // otherwise nothing changes
if (entering == leaving) {
Expand Down Expand Up @@ -879,13 +879,13 @@ template <typename T, typename X> unsigned lp_primal_core_solver<T, X>::solve()
return 0;
}
do {
if (this->print_statistics_with_iterations_and_nonzeroes_and_cost_and_check_that_the_time_is_over((this->m_using_infeas_costs? "inf" : "feas"), * this->m_settings.get_message_ostream())) {
if (this->print_statistics_with_iterations_and_nonzeroes_and_cost_and_check_that_the_time_is_over((this->using_infeas_costs()? "inf" : "feas"), * this->m_settings.get_message_ostream())) {
return this->total_iterations();
}
one_iteration();

TRACE("lar_solver", tout << "one iteration: " << this->get_status() << "\n";);
lp_assert(!this->m_using_infeas_costs || this->costs_on_nbasis_are_zeros());
lp_assert(!this->using_infeas_costs() || this->costs_on_nbasis_are_zeros());
switch (this->get_status()) {
case lp_status::OPTIMAL: // double check that we are at optimum
case lp_status::INFEASIBLE:
Expand Down Expand Up @@ -1109,7 +1109,7 @@ 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->m_using_infeas_costs = true;
this->set_using_infeas_costs(true);
}


Expand All @@ -1119,7 +1119,7 @@ void lp_primal_core_solver<T, X>::init_infeasibility_costs() {
lp_assert(this->m_column_types.size() >= this->m_n());
for (unsigned j = this->m_n(); j--;)
init_infeasibility_cost_for_column(j);
this->m_using_infeas_costs = true;
this->set_using_infeas_costs(true);
}

template <typename T, typename X> T
Expand Down Expand Up @@ -1297,13 +1297,13 @@ template <typename T, typename X> void lp_primal_core_solver<T, X>::print_breakp
template <typename T, typename X>
void lp_primal_core_solver<T, X>::init_reduced_costs() {
lp_assert(!this->use_tableau());
if (this->current_x_is_infeasible() && !this->m_using_infeas_costs) {
if (this->current_x_is_infeasible() && !this->using_infeas_costs()) {
init_infeasibility_costs();
} else if (this->current_x_is_feasible() && this->m_using_infeas_costs) {
} else if (this->current_x_is_feasible() && this->using_infeas_costs()) {
if (this->m_look_for_feasible_solution_only)
return;
this->m_costs = m_costs_backup;
this->m_using_infeas_costs = false;
this->set_using_infeas_costs(false);
}

this->init_reduced_costs_for_one_iteration();
Expand Down
18 changes: 9 additions & 9 deletions src/math/lp/lp_primal_core_solver_tableau_def.h
Original file line number Diff line number Diff line change
Expand Up @@ -86,7 +86,7 @@ template <typename T, typename X> int lp_primal_core_solver<T, X>::choose_enteri
return -1;
unsigned entering = *entering_iter;
m_sign_of_entering_delta = this->m_d[entering] > 0 ? 1 : -1;
if (this->m_using_infeas_costs && this->m_settings.use_breakpoints_in_feasibility_search)
if (this->using_infeas_costs() && this->m_settings.use_breakpoints_in_feasibility_search)
m_sign_of_entering_delta = -m_sign_of_entering_delta;
m_non_basis_list.erase(entering_iter);
m_non_basis_list.push_back(entering);
Expand All @@ -110,7 +110,7 @@ unsigned lp_primal_core_solver<T, X>::solve_with_tableau() {
return 0;
}
do {
if (this->print_statistics_with_iterations_and_nonzeroes_and_cost_and_check_that_the_time_is_over((this->m_using_infeas_costs? "inf t" : "feas t"), * this->m_settings.get_message_ostream())) {
if (this->print_statistics_with_iterations_and_nonzeroes_and_cost_and_check_that_the_time_is_over((this->using_infeas_costs()? "inf t" : "feas t"), * this->m_settings.get_message_ostream())) {
return this->total_iterations();
}
if (this->m_settings.use_tableau_rows()) {
Expand Down Expand Up @@ -216,7 +216,7 @@ template <typename T, typename X>void lp_primal_core_solver<T, X>::advance_on_en
lp_assert((this->m_settings.simplex_strategy() ==
simplex_strategy_enum::tableau_rows) ||
m_non_basis_list.back() == static_cast<unsigned>(entering));
lp_assert(this->m_using_infeas_costs || !is_neg(t));
lp_assert(this->using_infeas_costs() || !is_neg(t));
lp_assert(entering != leaving || !is_zero(t)); // otherwise nothing changes
if (entering == leaving) {
advance_on_entering_equal_leaving_tableau(entering, t);
Expand Down Expand Up @@ -359,12 +359,12 @@ update_basis_and_x_tableau(int entering, int leaving, X const & tt) {
template <typename T, typename X> void lp_primal_core_solver<T, X>::
update_x_tableau(unsigned entering, const X& delta) {
this->add_delta_to_x(entering, delta);
if (!this->m_using_infeas_costs) {
if (!this->using_infeas_costs()) {
for (const auto & c : this->m_A.m_columns[entering]) {
unsigned i = c.var();
this->add_delta_to_x_and_track_feasibility(this->m_basis[i], - delta * this->m_A.get_val(c));
}
} else { // m_using_infeas_costs == true
} 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.
Expand All @@ -386,7 +386,7 @@ template <typename T, typename X> void lp_primal_core_solver<T, X>::
update_inf_cost_for_column_tableau(unsigned j) {
lp_assert(this->m_settings.simplex_strategy() != simplex_strategy_enum::tableau_rows);

lp_assert(this->m_using_infeas_costs);
lp_assert(this->using_infeas_costs());

T new_cost = get_infeasibility_cost_for_column(j);
T delta = this->m_costs[j] - new_cost;
Expand All @@ -397,13 +397,13 @@ update_inf_cost_for_column_tableau(unsigned j) {
}

template <typename T, typename X> void lp_primal_core_solver<T, X>::init_reduced_costs_tableau() {
if (this->current_x_is_infeasible() && !this->m_using_infeas_costs) {
if (this->current_x_is_infeasible() && !this->using_infeas_costs()) {
init_infeasibility_costs();
} else if (this->current_x_is_feasible() && this->m_using_infeas_costs) {
} else if (this->current_x_is_feasible() && this->using_infeas_costs()) {
if (this->m_look_for_feasible_solution_only)
return;
this->m_costs = m_costs_backup;
this->m_using_infeas_costs = false;
this->set_using_infeas_costs(false);
}
unsigned size = this->m_basis_heading.size();
for (unsigned j = 0; j < size; j++) {
Expand Down

0 comments on commit 8c5993d

Please sign in to comment.