From 0bb21914e67774a87056eb4a12f1c44b817e4d74 Mon Sep 17 00:00:00 2001 From: Lev Nachmanson Date: Thu, 19 May 2016 14:59:32 -0700 Subject: [PATCH] dev(lp): simplify the design of lar_solver Signed-off-by: Lev Nachmanson --- src/util/lp/canonic_left_side.h | 13 +- src/util/lp/column_info.h | 35 +++- src/util/lp/lar_solver.cpp | 275 ++++++++++++++++-------------- src/util/lp/lar_solver.h | 42 +++-- src/util/lp/lp_primal_simplex.cpp | 6 +- src/util/lp/lp_solver.cpp | 25 ++- src/util/lp/lp_solver.h | 7 +- src/util/lp/lp_utils.h | 4 +- 8 files changed, 229 insertions(+), 178 deletions(-) diff --git a/src/util/lp/canonic_left_side.h b/src/util/lp/canonic_left_side.h index 64b5c3cfe..a469dcd91 100644 --- a/src/util/lp/canonic_left_side.h +++ b/src/util/lp/canonic_left_side.h @@ -12,8 +12,7 @@ #include #include "util/lp/column_info.h" namespace lean { -typedef unsigned var_index; -typedef unsigned constraint_index; + enum lconstraint_kind { LE = -2, LT = -1 , GE = 2, GT = 1, EQ = 0 }; @@ -25,10 +24,10 @@ inline bool compare(const std::pair & a, const std::pair(-1); + var_index m_additional_var_index = static_cast(-1); // this is the index of the additional variable created for this constraint std::vector> m_coeffs; - column_info m_column_info; lar_normalized_constraint * m_low_bound_witness = nullptr; lar_normalized_constraint * m_upper_bound_witness = nullptr; @@ -42,10 +41,6 @@ public: normalize(); } - void set_name(std::string name) { - m_column_info.set_name(name); - } - unsigned size() const { return static_cast(m_coeffs.size()); } void normalize() { diff --git a/src/util/lp/column_info.h b/src/util/lp/column_info.h index b7435284d..4a40330a7 100644 --- a/src/util/lp/column_info.h +++ b/src/util/lp/column_info.h @@ -12,6 +12,10 @@ #include #include "util/lp/lp_settings.h" namespace lean { +inline bool is_valid(unsigned j) { return static_cast(j) >= 0;} +typedef unsigned var_index; +typedef unsigned constraint_index; + template class column_info { std::string m_name; @@ -24,13 +28,40 @@ class column_info { T m_cost = numeric_traits::zero(); T m_fixed_value; bool m_is_fixed = false; - + unsigned m_column_index; public: + void set_column_index(unsigned j) { + m_column_index = j; + } + // the default constructor + column_info() {} + + column_info(unsigned column_index) : m_column_index(column_index) { + } + + column_info(const column_info & ci) { + m_name = ci.m_name; + m_low_bound_is_set = ci.m_low_bound_is_set; + m_low_bound_is_strict = ci.m_low_bound_is_strict; + m_upper_bound_is_set = ci.m_upper_bound_is_set; + m_upper_bound_is_strict = ci.m_upper_bound_is_strict; + m_low_bound = ci.m_low_bound; + m_upper_bound = ci.m_upper_bound; + m_cost = ci.m_cost; + m_fixed_value = ci.m_fixed_value; + m_is_fixed = ci.m_is_fixed; + m_column_index = ci.m_column_index; + } + + unsigned get_column_index() const { + return m_column_index; + } + column_type get_column_type() const { return m_is_fixed? fixed : (m_low_bound_is_set? (m_upper_bound_is_set? boxed : low_bound) : (m_upper_bound_is_set? upper_bound: free_column)); } - column_type get_column_type_no_flipping() { + column_type get_column_type_no_flipping() { if (m_is_fixed) { return column_type::fixed; } diff --git a/src/util/lp/lar_solver.cpp b/src/util/lp/lar_solver.cpp index cec41b6aa..5547349a8 100644 --- a/src/util/lp/lar_solver.cpp +++ b/src/util/lp/lar_solver.cpp @@ -33,9 +33,15 @@ double conversion_helper ::get_upper_bound(const column_info & ci) canonic_left_side * lar_solver::create_or_fetch_existing_left_side(const buffer>& left_side_par) { auto left_side = new canonic_left_side(left_side_par); lean_assert(left_side->size() > 0); - auto it = m_canonic_left_sides.find(left_side); - if (it == m_canonic_left_sides.end()) { - m_canonic_left_sides.insert(left_side); + auto it = m_set_of_canonic_left_sides.find(left_side); + if (it == m_set_of_canonic_left_sides.end()) { + m_set_of_canonic_left_sides.insert(left_side); + lean_assert(m_map_from_var_index_to_column_info_with_cls.find(m_available_var_index) == m_map_from_var_index_to_column_info_with_cls.end()); + unsigned vj = m_available_var_index; + m_map_from_var_index_to_column_info_with_cls[vj] = column_info_with_cls(left_side); + m_map_from_var_index_to_column_info_with_cls[vj].m_column_info.set_name("_s" + T_to_string(vj)); + left_side->m_additional_var_index = vj; + m_available_var_index++; } else { delete left_side; left_side = *it; @@ -43,14 +49,7 @@ canonic_left_side * lar_solver::create_or_fetch_existing_left_side(const buffer< return left_side; } - -bool lar_solver::var_is_fixed(unsigned j) { - auto it = m_map_from_var_index_to_left_side.find(j); - lean_assert(it != m_map_from_var_index_to_left_side.end()); - return it->second->m_column_info.is_fixed(); -} - -mpq lar_solver::find_ratio(canonic_left_side * ls, const lar_constraint & constraint) { +mpq lar_solver::find_ratio_of_original_constraint_to_normalized(canonic_left_side * ls, const lar_constraint & constraint) { lean_assert(ls->m_coeffs.size() > 0); auto first_pair = ls->m_coeffs[0]; lean_assert(first_pair.first == numeric_traits::one()); @@ -60,81 +59,68 @@ mpq lar_solver::find_ratio(canonic_left_side * ls, const lar_constraint & constr return it->second; } -void lar_solver::add_canonic_left_side_for_var(var_index i, std::string var_name) { - buffer> b; - b.push_back(std::make_pair(numeric_traits::one(), i)); - auto can_ls = new canonic_left_side(b); - can_ls->set_name(var_name); - - lean_assert(m_canonic_left_sides.find(can_ls) == m_canonic_left_sides.end()); - m_canonic_left_sides.insert(can_ls); - m_map_from_var_index_to_left_side[i] = can_ls; -} - -void lar_solver::map_left_side_to_column_of_A(canonic_left_side* left_side, unsigned & j) { - lean_assert(m_column_indices_to_canonic_left_sides.find(j) == m_column_indices_to_canonic_left_sides.end()); - left_side->m_column_index = j; // assigning this index does not change the hash of canonic_left_side - if (left_side->size() > 1) { // if size is one we will not create a row for this left side - left_side->m_row_index = static_cast(m_basis.size()); - m_basis.push_back(j); // j will be a basis column, so we put it into the basis as well - } - m_column_indices_to_canonic_left_sides[j++] = left_side; // pointing from the column to the left side +void lar_solver::map_left_side_to_A_of_core_solver(canonic_left_side* left_side, unsigned j) { + var_index additional_var = left_side->m_additional_var_index; + lean_assert(valid_index(additional_var)); + auto it = m_map_from_var_index_to_column_info_with_cls.find(additional_var); + lean_assert(it != m_map_from_var_index_to_column_info_with_cls.end()); + column_info & ci = it->second.m_column_info; + lean_assert(!is_valid(ci.get_column_index())); + lean_assert(left_side->size() > 0); // if size is zero we have an empty row + left_side->m_row_index = m_basis.size(); + m_basis.push_back(j); // j will be a basis column, so we put it into the basis as well + lean_assert(m_map_from_column_indices_to_var_index.find(j) == m_map_from_column_indices_to_var_index.end()); + ci.set_column_index(j); + m_map_from_column_indices_to_var_index[j] = additional_var; } -void lar_solver::map_left_sides_to_columns_of_A() { - unsigned j = 0; - for (auto it : m_canonic_left_sides) - map_left_side_to_column_of_A(it, j); +void lar_solver::map_left_sides_to_A_of_core_solver() { + unsigned j = m_map_from_column_indices_to_var_index.size(); + for (auto it : m_set_of_canonic_left_sides) + map_left_side_to_A_of_core_solver(it, j++); } -bool valid_index(unsigned j) { return static_cast(j) >= 0;} - - // this adds a row to A template -void lar_solver::add_row_to_A(static_matrix & A, unsigned i, canonic_left_side * ls) { +void lar_solver::fill_row_of_A(static_matrix & A, unsigned i, canonic_left_side * ls) { for (auto & t : ls->m_coeffs) { var_index vi = t.second; - canonic_left_side *var_left_side = m_map_from_var_index_to_left_side[vi]; - unsigned column = var_left_side->m_column_index; - lean_assert(valid_index(column)); + unsigned column = get_column_index_from_var_index(vi); + lean_assert(is_valid(column)); A.set(i, column, convert_struct::convert(t.first)); } - A.set(i, ls->m_column_index, - one_of_type()); + unsigned additional_column = get_column_index_from_var_index(ls->m_additional_var_index); + lean_assert(is_valid(additional_column)); + A.set(i, additional_column, - one_of_type()); +} + +void lar_solver::fill_set_of_active_var_indices() { + for (auto & t : m_set_of_canonic_left_sides) + for (auto & a : t->m_coeffs) + m_set_of_active_var_indices.insert(a.second); } template void lar_solver::create_matrix_A(static_matrix & A) { - unsigned n = static_cast(m_column_indices_to_canonic_left_sides.size()); - unsigned m = static_cast(m_basis.size()); + unsigned m = m_set_of_canonic_left_sides.size(); + fill_set_of_active_var_indices(); + unsigned n = m_set_of_active_var_indices.size() + m; A.init_empty_matrix(m, n); unsigned i = 0; - for (auto & t : m_column_indices_to_canonic_left_sides) - if (valid_index(t.second->m_row_index)) - add_row_to_A(A, i++, t.second); -#ifdef LEAN_DEBUG - // print_matrix(m_A); -#endif + for (auto t : m_set_of_canonic_left_sides) { + lean_assert(t->size() > 0); + if (is_valid(t->m_row_index)) + fill_row_of_A(A, i++, t); + } } -// void lar_solver::fill_column_info_names() { -// for (unsigned j = 0; j < m_A.column_count(); j++) { -// column_info t; -// m_column_infos.push_back(t); -// if (j < m_map_from_var_index_to_name_left_side_pair.size()) { -// m_column_infos.back().set_name(m_map_from_var_index_to_name_left_side_pair[j]); -// } else { -// string pref("_s_"); -// m_column_infos.back().set_name(pref + T_to_string(j)); -// } -// m_column_names -// } -// } void lar_solver::set_upper_bound_for_column_info(lar_normalized_constraint * norm_constr) { const mpq & v = norm_constr->m_right_side; canonic_left_side * ls = norm_constr->m_canonic_left_side; - column_info & ci = ls->m_column_info; + var_index additional_var_index = ls->m_additional_var_index; + lean_assert(is_valid(additional_var_index)); + column_info & ci = get_column_info_from_var_index(additional_var_index); lean_assert(norm_constr->m_kind == LE || norm_constr->m_kind == LT || norm_constr->m_kind == EQ); bool strict = norm_constr->m_kind == LT; if (!ci.upper_bound_is_set()) { @@ -165,7 +151,7 @@ bool lar_solver::try_to_set_fixed(column_info & ci) { void lar_solver::set_low_bound_for_column_info(lar_normalized_constraint * norm_constr) { const mpq & v = norm_constr->m_right_side; canonic_left_side * ls = norm_constr->m_canonic_left_side; - column_info & ci = ls->m_column_info; + column_info & ci = get_column_info_from_var_index(ls->m_additional_var_index); lean_assert(norm_constr->m_kind == GE || norm_constr->m_kind == GT || norm_constr->m_kind == EQ); bool strict = norm_constr->m_kind == GT; if (!ci.low_bound_is_set()) { @@ -221,10 +207,10 @@ column_type lar_solver::get_column_type(column_info & ci) { void lar_solver::fill_column_names() { m_column_names.clear(); - for (auto t : m_canonic_left_sides) { - auto & ci = t->m_column_info; - unsigned j = t->m_column_index; - lean_assert(valid_index(j)); + for (auto & t : m_map_from_var_index_to_column_info_with_cls) { + column_info & ci = t.second.m_column_info; + unsigned j = ci.get_column_index(); + lean_assert(is_valid(j)); std::string name = ci.get_name(); if (name.size() == 0) name = std::string("_s") + T_to_string(j); @@ -234,24 +220,25 @@ void lar_solver::fill_column_names() { void lar_solver::fill_column_types() { m_column_types.clear(); - m_column_types.resize(m_canonic_left_sides.size(), free_column); - for (auto t : m_canonic_left_sides) { - auto & ci = t->m_column_info; - unsigned j = t->m_column_index; - lean_assert(valid_index(j)); - m_column_types[j] = get_column_type(ci); + m_column_types.resize(m_map_from_var_index_to_column_info_with_cls.size(), free_column); + for (auto t : m_set_of_canonic_left_sides) { + var_index additional_vj = t->m_additional_var_index; + unsigned j = get_column_index_from_var_index(additional_vj); + lean_assert(is_valid(j)); + m_column_types[j] = get_column_type(get_column_info_from_var_index(additional_vj)); } } template void lar_solver::fill_bounds_for_core_solver(std::vector & lb, std::vector & ub) { - unsigned n = static_cast(m_canonic_left_sides.size()); // this is the number of columns + unsigned n = static_cast(m_map_from_var_index_to_column_info_with_cls.size()); // this is the number of columns lb.resize(n); ub.resize(n); - for (auto t : m_canonic_left_sides) { - auto & ci = t->m_column_info; - unsigned j = t->m_column_index; - lean_assert(valid_index(j)); + for (auto t : m_set_of_canonic_left_sides) { + auto & ci = get_column_info_from_var_index(t->m_additional_var_index); + unsigned j = ci.get_column_index(); + lean_assert(is_valid(j)); + lean_assert(j < n); if (ci.low_bound_is_set()) lb[j] = conversion_helper::get_low_bound(ci); if (ci.upper_bound_is_set()) @@ -290,31 +277,32 @@ template V lar_solver::get_column_val(std::vector & low_bound, s lar_solver::~lar_solver() { std::vector to_delete; - for (auto it : m_canonic_left_sides) + for (auto it : m_set_of_canonic_left_sides) to_delete.push_back(it); for (auto t : to_delete) delete t; } - - - var_index lar_solver::add_var(std::string s) { - auto got = m_var_names_to_var_index.find(s); - if (got != m_var_names_to_var_index.end()) return got->second; - + auto it = m_var_names_to_var_index.find(s); + if (it != m_var_names_to_var_index.end()) + return it->second; var_index i = m_available_var_index++; + lean_assert(m_map_from_var_index_to_column_info_with_cls.find(i) == m_map_from_var_index_to_column_info_with_cls.end()); + auto ci_with_cls = column_info_with_cls(); + ci_with_cls.m_column_info.set_name(s); + m_map_from_var_index_to_column_info_with_cls[i] = ci_with_cls; m_var_names_to_var_index[s] = i; - add_canonic_left_side_for_var(i, s); return i; } constraint_index lar_solver::add_constraint(const buffer>& left_side, lconstraint_kind kind_par, mpq right_side_par) { lean_assert(left_side.size() > 0); constraint_index i = m_available_constr_index++; + lean_assert(m_normalized_constraints.find(i) == m_normalized_constraints.end()); lar_constraint original_constr(left_side, kind_par, right_side_par, i); canonic_left_side * ls = create_or_fetch_existing_left_side(left_side); - mpq ratio = find_ratio(ls, original_constr); + mpq ratio = find_ratio_of_original_constraint_to_normalized(ls, original_constr); auto kind = ratio.is_neg()? flip_kind(kind_par): kind_par; mpq right_side = right_side_par / ratio; lar_normalized_constraint normalized_constraint(ls, ratio, kind, right_side, original_constr); @@ -455,9 +443,21 @@ mpq lar_solver::sum_of_right_sides_of_evidence(const buffersecond; + auto it2 = m_map_from_var_index_to_column_info_with_cls.find(var_j); + lean_assert(it2 != m_map_from_var_index_to_column_info_with_cls.end()); + canonic_left_side * ls = it2->second.m_canonic_left_side; int adj_sign = coeff.is_pos() ? inf_sign : -inf_sign; lar_normalized_constraint * bound_constr = adj_sign < 0? ls->m_upper_bound_witness : ls->m_low_bound_witness; @@ -577,10 +582,10 @@ void lar_solver::get_infeasibility_evidence_for_inf_sign(buffer::one(); - for (auto t : m_canonic_left_sides) { - auto & ci = t->m_column_info; - unsigned j = t->m_column_index; - lean_assert (valid_index(j)); + for (auto t : m_set_of_canonic_left_sides) { + auto & ci = get_column_info_from_var_index(t->m_additional_var_index); + unsigned j = ci.get_column_index(); + lean_assert (is_valid(j)); if (ci.low_bound_is_set()) restrict_delta_on_low_bound_column(delta, j); if (ci.upper_bound_is_set()) @@ -623,18 +628,21 @@ void lar_solver::restrict_delta_on_upper_bound(mpq& delta, unsigned j) { void lar_solver::get_model(std::unordered_map & variable_values){ lean_assert(m_status == OPTIMAL); mpq delta = find_delta_for_strict_bounds(); - for (auto & it : m_map_from_var_index_to_left_side) { - numeric_pair & rp = m_x[it.second->m_column_index]; + for (auto & it : m_map_from_var_index_to_column_info_with_cls) { + column_info & ci = it.second.m_column_info; + unsigned j = ci.get_column_index(); + numeric_pair & rp = m_x[j]; variable_values[it.first] = rp.x + delta * rp.y; } } -std::string lar_solver::get_variable_name(var_index vi) { - if (m_map_from_var_index_to_left_side.size() <= vi) { +std::string lar_solver::get_variable_name(var_index vi) const { + auto it = m_map_from_var_index_to_column_info_with_cls.find(vi); + if (it == m_map_from_var_index_to_column_info_with_cls.end()) { std::string s = "variable " + T_to_string(vi) + " is not found"; return s; } - return m_map_from_var_index_to_left_side[vi]->m_column_info.get_name(); + return it->second.m_column_info.get_name(); } // ********** print region start @@ -650,7 +658,7 @@ void lar_solver::print_constraint(constraint_index ci, std::ostream & out) { void lar_solver::print_canonic_left_side(const canonic_left_side & c, std::ostream & out) { bool first = true; - for (auto it : c.m_coeffs) { + for (auto & it : c.m_coeffs) { auto val = it.first; if (first) { first = false; @@ -664,13 +672,13 @@ void lar_solver::print_canonic_left_side(const canonic_left_side & c, std::ostre } if (val != numeric_traits::one()) out << T_to_string(val); - out << m_map_from_var_index_to_left_side[it.second]->m_column_info.get_name(); + out << get_variable_name(it.second); } } void lar_solver::print_left_side_of_constraint(const lar_base_constraint * c, std::ostream & out) { bool first = true; - for (auto it : c->get_left_side_coefficients()) { + for (auto & it : c->get_left_side_coefficients()) { auto val = it.first; if (numeric_traits::is_zero(val)) continue; if (first) { @@ -686,29 +694,29 @@ void lar_solver::print_left_side_of_constraint(const lar_base_constraint * c, st if (val != numeric_traits::one()) out << val; - out << m_map_from_var_index_to_left_side[it.second]->m_column_info.get_name(); + out << get_variable_name(it.second); } } -void lar_solver::print_info_on_column(unsigned j, std::ostream & out) { - for (auto ls : m_canonic_left_sides) { - if (static_cast(ls->m_column_index) == j) { - auto & ci = ls->m_column_info; - if (ci.low_bound_is_set()) { - out << "l = " << ci.get_low_bound(); - } - if (ci.upper_bound_is_set()) { - out << "u = " << ci.get_upper_bound(); - } - out << std::endl; - m_mpq_core_solver.print_column_info(j, out); - } - } -} +// void lar_solver::print_info_on_column(unsigned j, std::ostream & out) { +// for (auto ls : m_set_of_canonic_left_sides) { +// if (static_cast(ls->get_ci()) == j) { +// auto & ci = ls->m_column_info; +// if (ci.low_bound_is_set()) { +// out << "l = " << ci.get_low_bound(); +// } +// if (ci.upper_bound_is_set()) { +// out << "u = " << ci.get_upper_bound(); +// } +// out << std::endl; +// m_mpq_core_solver.print_column_info(j, out); +// } +// } +// } mpq lar_solver::get_infeasibility_of_solution(std::unordered_map & solution) { mpq ret = numeric_traits::zero(); - for (auto it : m_normalized_constraints) { + for (auto & it : m_normalized_constraints) { ret += get_infeasibility_of_constraint(it.second, solution); } return ret; @@ -735,12 +743,9 @@ mpq lar_solver::get_infeasibility_of_constraint(const lar_normalized_constraint mpq lar_solver::get_canonic_left_side_val(canonic_left_side * ls, std::unordered_map & solution) { mpq ret = numeric_traits::zero(); - for (auto it : ls->m_coeffs) { - var_index j = it.second; - auto vi = m_map_from_var_index_to_left_side.find(j); - lean_assert(vi != m_map_from_var_index_to_left_side.end()); - canonic_left_side * var_ls = vi->second; - std::string s = var_ls->m_column_info.get_name(); + for (auto & it : ls->m_coeffs) { + var_index vi = it.second; + std::string s = get_variable_name(vi); auto t = solution.find(s); lean_assert(t != solution.end()); ret += it.first * (t->second); @@ -750,7 +755,7 @@ mpq lar_solver::get_canonic_left_side_val(canonic_left_side * ls, std::unordered mpq lar_solver::get_left_side_val(const lar_constraint & cns, const std::unordered_map & var_map) { mpq ret = numeric_traits::zero(); - for (auto it : cns.m_left_side) { + for (auto & it : cns.m_left_side) { var_index j = it.first; auto vi = var_map.find(j); lean_assert(vi != var_map.end()); @@ -763,5 +768,21 @@ void lar_solver::print_constraint(const lar_base_constraint * c, std::ostream & print_left_side_of_constraint(c, out); out <<" " << lconstraint_kind_string(c->m_kind) << " " << c->m_right_side; } + +unsigned lar_solver::get_column_index_from_var_index(var_index vi) const{ + auto it = m_map_from_var_index_to_column_info_with_cls.find(vi); + if (it == m_map_from_var_index_to_column_info_with_cls.end()) + return static_cast(-1); + return it->second.m_column_info.get_column_index(); +} + +column_info & lar_solver::get_column_info_from_var_index(var_index vi) { + auto it = m_map_from_var_index_to_column_info_with_cls.find(vi); + if (it == m_map_from_var_index_to_column_info_with_cls.end()) { + lean_assert(false); + throw 0; + } + return it->second.m_column_info; +} } diff --git a/src/util/lp/lar_solver.h b/src/util/lp/lar_solver.h index 279583bd7..051f01f7c 100644 --- a/src/util/lp/lar_solver.h +++ b/src/util/lp/lar_solver.h @@ -33,6 +33,14 @@ struct conversion_helper { } }; +struct column_info_with_cls { // column_info_with canonic_left_side + canonic_left_side * m_canonic_left_side; + column_info m_column_info; + + column_info_with_cls(): m_canonic_left_side(nullptr), m_column_info(static_cast(-1)) {} + column_info_with_cls(canonic_left_side * cls) : m_canonic_left_side(cls), m_column_info(static_cast(-1)) {} +}; + template<> struct conversion_helper { static double get_low_bound(const column_info & ci); @@ -43,11 +51,12 @@ class lar_solver { unsigned m_available_var_index = 0; unsigned m_available_constr_index = 0; lp_status m_status = UNKNOWN; + std::unordered_set m_set_of_active_var_indices; // a variable is active if it is referenced in a left side std::unordered_map m_var_names_to_var_index; - std::unordered_map m_map_from_var_index_to_left_side; - std::unordered_set m_canonic_left_sides; - std::unordered_map m_column_indices_to_canonic_left_sides; + std::unordered_set m_set_of_canonic_left_sides; + std::unordered_map m_map_from_column_indices_to_var_index; std::unordered_map m_normalized_constraints; + std::unordered_map m_map_from_var_index_to_column_info_with_cls; static_matrix> m_A; lar_core_solver> m_mpq_core_solver; lp_settings m_settings; @@ -62,24 +71,20 @@ class lar_solver { canonic_left_side * m_infeasible_canonic_left_side = nullptr; // such can be found at the initialization step canonic_left_side * create_or_fetch_existing_left_side(const buffer>& left_side_par); - - bool var_is_fixed(unsigned j); - - mpq find_ratio(canonic_left_side * ls, const lar_constraint & constraint); + mpq find_ratio_of_original_constraint_to_normalized(canonic_left_side * ls, const lar_constraint & constraint); void add_canonic_left_side_for_var(var_index i, std::string var_name); - void map_left_side_to_column_of_A(canonic_left_side* left_side, unsigned & j); + void map_left_side_to_A_of_core_solver(canonic_left_side* left_side, var_index vj); - - void map_left_sides_to_columns_of_A(); + void map_left_sides_to_A_of_core_solver(); bool valid_index(unsigned j) { return static_cast(j) >= 0;} // this adds a row to A template - void add_row_to_A(static_matrix & A, unsigned i, canonic_left_side * ls); + void fill_row_of_A(static_matrix & A, unsigned i, canonic_left_side * ls); template void create_matrix_A(static_matrix & A); @@ -123,6 +128,12 @@ class lar_solver { std::vector & upper_bound, const lar_solution_signature & signature); template V get_column_val(std::vector & low_bound, std::vector & upper_bound, non_basic_column_value_position pos_type, unsigned j); + void map_var_indices_to_columns_of_A(); + + void register_in_map(std::unordered_map & coeffs, lar_constraint & cn, const mpq & a); + unsigned get_column_index_from_var_index(var_index vi) const; + column_info & get_column_info_from_var_index(var_index vi); + void fill_set_of_active_var_indices(); public: ~lar_solver(); @@ -149,10 +160,6 @@ public: constraint_index add_constraint(const buffer>& left_side, lconstraint_kind kind_par, mpq right_side_par); - bool is_infeasible(const column_info & ci) { - return ci.low_bound_is_set() && ci.upper_bound_is_set() && ci.get_low_bound() > ci.get_upper_bound(); - } - bool all_constraints_hold(); bool constraint_holds(const lar_constraint & constr, std::unordered_map & var_map); @@ -163,7 +170,6 @@ public: bool the_relations_are_of_same_type(const buffer> & evidence, lconstraint_kind & the_kind_of_sum); - void register_in_map(std::unordered_map & coeffs, lar_constraint & cn, const mpq & a); bool the_left_sides_sum_to_zero(const buffer> & evidence); bool the_righ_sides_do_not_sum_to_zero(const buffer> & evidence); @@ -214,7 +220,7 @@ public: void get_model(std::unordered_map & variable_values); - std::string get_variable_name(var_index vi); + std::string get_variable_name(var_index vi) const; void print_constraint(constraint_index ci, std::ostream & out); @@ -224,7 +230,7 @@ public: numeric_pair get_infeasibility_from_core_solver(std::unordered_map & solution); - void print_info_on_column(unsigned j, std::ostream & out); + // void print_info_on_column(unsigned j, std::ostream & out); mpq get_infeasibility_of_solution(std::unordered_map & solution); diff --git a/src/util/lp/lp_primal_simplex.cpp b/src/util/lp/lp_primal_simplex.cpp index cdd117ed5..85d343566 100644 --- a/src/util/lp/lp_primal_simplex.cpp +++ b/src/util/lp/lp_primal_simplex.cpp @@ -145,9 +145,9 @@ template void lp_primal_simplex::set_core_solver_ this->m_upper_bounds.resize(total_vars); for (auto cit : this->m_columns) { column_info * ci = cit.second; - auto p = this->m_external_columns_to_core_solver_columns.find(cit.first); - if (p == this->m_external_columns_to_core_solver_columns.end()) continue; - unsigned j = p->second; + unsigned j = ci->get_column_index(); + if (!is_valid(j)) + continue; // the variable is not mapped to a column switch (this->m_column_types[j] = ci->get_column_type()){ case fixed: this->m_upper_bounds[j] = numeric_traits::zero(); diff --git a/src/util/lp/lp_solver.cpp b/src/util/lp/lp_solver.cpp index 77e19ec75..d4417a728 100644 --- a/src/util/lp/lp_solver.cpp +++ b/src/util/lp/lp_solver.cpp @@ -401,9 +401,9 @@ template void lp_solver::map_external_columns_to_ throw_exception("found fixed column"); } unsigned j = col.first; - auto j_place = m_external_columns_to_core_solver_columns.find(j); - if (j_place == m_external_columns_to_core_solver_columns.end()) { // j is a newcomer - m_external_columns_to_core_solver_columns[j] = size; + auto j_column = m_columns[j]->get_column_index(); + if (!is_valid(j_column)) { // j is a newcomer + m_columns[j]->set_column_index(size); m_core_solver_columns_to_external_columns[size++] = j; } } @@ -412,10 +412,9 @@ template void lp_solver::map_external_columns_to_ template void lp_solver::fill_column_names_for_core_solver() { for (auto it : this->m_columns) { - auto p = this->m_external_columns_to_core_solver_columns.find(it.first); - if (p != this->m_external_columns_to_core_solver_columns.end()) { - this->m_name_map[p->second] = it.second->get_name(); - } + unsigned j = it.second->get_column_index(); + if (is_valid(j)) + this->m_name_map[j] = it.second->get_name(); } } @@ -434,10 +433,9 @@ template void lp_solver::fill_A_from_A_values() { lean_assert(m_external_rows_to_core_solver_rows.find(t.first) != m_external_rows_to_core_solver_rows.end()); unsigned row = m_external_rows_to_core_solver_rows[t.first]; for (auto k : t.second) { - lean_assert(m_external_columns_to_core_solver_columns.find(k.first) != m_external_columns_to_core_solver_columns.end()); - unsigned col = m_external_columns_to_core_solver_columns[k.first]; + unsigned col = m_columns[k.first]->get_column_index(); + lean_assert(is_valid(col)); bool col_is_flipped = m_columns[k.first]->is_flipped(); - if (!col_is_flipped) { (*m_A)(row, col) = k.second; } else { @@ -505,7 +503,7 @@ template void lp_solver::fill_m_b() { } } -template T lp_solver::get_column_value_with_core_solver(unsigned column, lp_core_solver_base * core_solver) const { +template T lp_solver::get_column_value_with_core_solver(unsigned column, lp_core_solver_base * core_solver) const { auto cit = this->m_columns.find(column); if (cit == this->m_columns.end()) { return numeric_traits::zero(); @@ -517,9 +515,8 @@ template T lp_solver::get_column_value_with_co return ci->get_fixed_value(); } - auto t = this->m_external_columns_to_core_solver_columns.find(column); - if (t != this->m_external_columns_to_core_solver_columns.end()){ - unsigned cj = t->second; + unsigned cj = ci->get_column_index(); + if (cj != static_cast(-1)) { T v = core_solver->get_var_value(cj) * this->m_column_scale[cj]; if (ci->is_free()) { return v; diff --git a/src/util/lp/lp_solver.h b/src/util/lp/lp_solver.h index 63a1e8020..5bccc816d 100644 --- a/src/util/lp/lp_solver.h +++ b/src/util/lp/lp_solver.h @@ -43,12 +43,11 @@ public: unsigned m_first_stage_iterations = 0; unsigned m_second_stage_iterations = 0; std::unordered_map> m_constraints; - std::unordered_map*> m_columns; + std::unordered_map*> m_columns; std::unordered_map > m_A_values; std::unordered_map m_names_to_columns; // don't have to use it std::unordered_map m_external_rows_to_core_solver_rows; std::unordered_map m_core_solver_rows_to_external_rows; - std::unordered_map m_external_columns_to_core_solver_columns; std::unordered_map m_core_solver_columns_to_external_columns; std::vector m_column_scale; std::unordered_map m_name_map; @@ -197,7 +196,9 @@ protected: void fill_column_names_for_core_solver(); - unsigned number_of_core_structurals() { return static_cast(m_external_columns_to_core_solver_columns.size()); } + unsigned number_of_core_structurals() { + return static_cast(m_core_solver_columns_to_external_columns.size()); + } void restore_column_scales_to_one() { diff --git a/src/util/lp/lp_utils.h b/src/util/lp/lp_utils.h index a15624dd4..edbd0bb68 100644 --- a/src/util/lp/lp_utils.h +++ b/src/util/lp/lp_utils.h @@ -18,7 +18,7 @@ namespace lean { inline void lean_assert(bool b) {} #else #define lean_assert(_x_) {} - #endif +#endif inline void lean_unreachable() { lean_assert(false); } template inline X zero_of_type() { return numeric_traits::zero(); } template inline X one_of_type() { return numeric_traits::one(); } @@ -75,7 +75,7 @@ template inline bool is_zero(const X & v) { return numeric_traits inline double get_double(const X & v) { return numeric_traits::get_double(v); } template inline T zero_of_type() {return numeric_traits::zero();} inline void throw_exception(std::string str) { throw exception(str); } -template inline T from_string(std::string const & str) { lean_unreachable();} +template inline T from_string(std::string const & ) { lean_unreachable();} template <> double inline from_string(std::string const & str) { return atof(str.c_str());} template <> mpq inline from_string(std::string const & str) { return mpq(atof(str.c_str()));