mirror of
https://github.com/Z3Prover/z3
synced 2025-04-23 17:15:31 +00:00
rm lp_solver
This commit is contained in:
Lev Nachmanson
parent
5e4bca3d26
commit
ff1dc0424c
8 changed files with 2 additions and 929 deletions
|
|
@ -43,15 +43,4 @@ public:
|
|||
m_value = val;
|
||||
}
|
||||
};
|
||||
#ifdef Z3DEBUG
|
||||
template <typename X>
|
||||
bool check_vector_for_small_values(indexed_vector<X> & w, lp_settings & settings) {
|
||||
for (unsigned i : w.m_index) {
|
||||
const X & v = w[i];
|
||||
if ((!is_zero(v)) && settings.abs_val_is_smaller_than_drop_tolerance(v))
|
||||
return false;
|
||||
}
|
||||
return true;
|
||||
}
|
||||
#endif
|
||||
}
|
||||
|
|
|
|||
|
|
@ -651,35 +651,7 @@ public:
|
|||
}
|
||||
}
|
||||
|
||||
void scale_problem_for_doubles(
|
||||
static_matrix<double, double>& A,
|
||||
vector<double> & lower_bounds,
|
||||
vector<double> & upper_bounds) {
|
||||
vector<double> column_scale_vector;
|
||||
vector<double> right_side_vector(A.column_count());
|
||||
settings().reps_in_scaler = 5;
|
||||
scaler<double, double > scaler(right_side_vector,
|
||||
A,
|
||||
settings().scaling_minimum,
|
||||
settings().scaling_maximum,
|
||||
column_scale_vector,
|
||||
settings());
|
||||
if (! scaler.scale()) {
|
||||
// the scale did not succeed, unscaling
|
||||
A.clear();
|
||||
create_double_matrix(A);
|
||||
} else {
|
||||
for (unsigned j = 0; j < A.column_count(); j++) {
|
||||
if (m_r_solver.column_has_upper_bound(j)) {
|
||||
upper_bounds[j] /= column_scale_vector[j];
|
||||
}
|
||||
if (m_r_solver.column_has_lower_bound(j)) {
|
||||
lower_bounds[j] /= column_scale_vector[j];
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
}
|
||||
|
||||
// returns the trace of basis changes
|
||||
vector<unsigned> find_solution_signature_with_doubles(lar_solution_signature & signature) {
|
||||
if (m_d_solver.m_factorization == nullptr || m_d_solver.m_factorization->get_status() != LU_status::OK) {
|
||||
|
|
|
|||
|
|
@ -30,7 +30,6 @@ Revision History:
|
|||
#include <cstdlib>
|
||||
#include <algorithm>
|
||||
#include "math/lp/lu.h"
|
||||
#include "math/lp/lp_solver.h"
|
||||
#include "math/lp/static_matrix.h"
|
||||
#include "math/lp/core_solver_pretty_printer.h"
|
||||
#include "math/lp/lp_core_solver_base.h"
|
||||
|
|
|
|||
|
|
@ -1,55 +0,0 @@
|
|||
/*++
|
||||
Copyright (c) 2017 Microsoft Corporation
|
||||
|
||||
Module Name:
|
||||
|
||||
<name>
|
||||
|
||||
Abstract:
|
||||
|
||||
<abstract>
|
||||
|
||||
Author:
|
||||
|
||||
Lev Nachmanson (levnach)
|
||||
|
||||
Revision History:
|
||||
|
||||
|
||||
--*/
|
||||
#include <string>
|
||||
#include "math/lp/lp_solver_def.h"
|
||||
template void lp::lp_solver<double, double>::add_constraint(lp::lp_relation, double, unsigned int);
|
||||
template void lp::lp_solver<double, double>::cleanup();
|
||||
template void lp::lp_solver<double, double>::count_slacks_and_artificials();
|
||||
template void lp::lp_solver<double, double>::fill_m_b();
|
||||
template void lp::lp_solver<double, double>::fill_matrix_A_and_init_right_side();
|
||||
template void lp::lp_solver<double, double>::flip_costs();
|
||||
template double lp::lp_solver<double, double>::get_column_cost_value(unsigned int, lp::column_info<double>*) const;
|
||||
template int lp::lp_solver<double, double>::get_column_index_by_name(std::string) const;
|
||||
template double lp::lp_solver<double, double>::get_column_value_with_core_solver(unsigned int, lp::lp_core_solver_base<double, double>*) const;
|
||||
template lp::column_info<double>* lp::lp_solver<double, double>::get_or_create_column_info(unsigned int);
|
||||
template void lp::lp_solver<double, double>::give_symbolic_name_to_column(std::string, unsigned int);
|
||||
template void lp::lp_solver<double, double>::print_statistics_on_A(std::ostream & out);
|
||||
template bool lp::lp_solver<double, double>::problem_is_empty();
|
||||
template void lp::lp_solver<double, double>::scale();
|
||||
template void lp::lp_solver<double, double>::set_scaled_cost(unsigned int);
|
||||
template lp::lp_solver<double, double>::~lp_solver();
|
||||
template void lp::lp_solver<lp::mpq, lp::mpq>::add_constraint(lp::lp_relation, lp::mpq, unsigned int);
|
||||
template void lp::lp_solver<lp::mpq, lp::mpq>::cleanup();
|
||||
template void lp::lp_solver<lp::mpq, lp::mpq>::count_slacks_and_artificials();
|
||||
template void lp::lp_solver<lp::mpq, lp::mpq>::fill_m_b();
|
||||
template void lp::lp_solver<lp::mpq, lp::mpq>::fill_matrix_A_and_init_right_side();
|
||||
template void lp::lp_solver<lp::mpq, lp::mpq>::flip_costs();
|
||||
template lp::mpq lp::lp_solver<lp::mpq, lp::mpq>::get_column_cost_value(unsigned int, lp::column_info<lp::mpq>*) const;
|
||||
template int lp::lp_solver<lp::mpq, lp::mpq>::get_column_index_by_name(std::string) const;
|
||||
template lp::mpq lp::lp_solver<lp::mpq, lp::mpq>::get_column_value_by_name(std::string) const;
|
||||
template lp::mpq lp::lp_solver<lp::mpq, lp::mpq>::get_column_value_with_core_solver(unsigned int, lp::lp_core_solver_base<lp::mpq, lp::mpq>*) const;
|
||||
template lp::column_info<lp::mpq>* lp::lp_solver<lp::mpq, lp::mpq>::get_or_create_column_info(unsigned int);
|
||||
template void lp::lp_solver<lp::mpq, lp::mpq>::give_symbolic_name_to_column(std::string, unsigned int);
|
||||
template void lp::lp_solver<lp::mpq, lp::mpq>::print_statistics_on_A(std::ostream & out);
|
||||
template bool lp::lp_solver<lp::mpq, lp::mpq>::problem_is_empty();
|
||||
template void lp::lp_solver<lp::mpq, lp::mpq>::scale();
|
||||
template void lp::lp_solver<lp::mpq, lp::mpq>::set_scaled_cost(unsigned int);
|
||||
template lp::lp_solver<lp::mpq, lp::mpq>::~lp_solver();
|
||||
template double lp::lp_solver<double, double>::get_column_value_by_name(std::string) const;
|
||||
|
|
@ -1,260 +0,0 @@
|
|||
/*++
|
||||
Copyright (c) 2017 Microsoft Corporation
|
||||
|
||||
Module Name:
|
||||
|
||||
<name>
|
||||
|
||||
Abstract:
|
||||
|
||||
<abstract>
|
||||
|
||||
Author:
|
||||
|
||||
Lev Nachmanson (levnach)
|
||||
|
||||
Revision History:
|
||||
|
||||
|
||||
--*/
|
||||
|
||||
#pragma once
|
||||
#include <string>
|
||||
#include <unordered_map>
|
||||
#include <algorithm>
|
||||
#include "util/vector.h"
|
||||
#include "math/lp/lp_settings.h"
|
||||
#include "math/lp/column_info.h"
|
||||
#include "math/lp/static_matrix.h"
|
||||
#include "math/lp/lp_core_solver_base.h"
|
||||
#include "math/lp/scaler.h"
|
||||
#include "math/lp/bound_analyzer_on_row.h"
|
||||
namespace lp {
|
||||
enum lp_relation {
|
||||
Less_or_equal,
|
||||
Equal,
|
||||
Greater_or_equal
|
||||
};
|
||||
|
||||
template <typename T, typename X>
|
||||
struct lp_constraint {
|
||||
X m_rs; // right side of the constraint
|
||||
lp_relation m_relation;
|
||||
lp_constraint() {} // empty constructor
|
||||
lp_constraint(T rs, lp_relation relation): m_rs(rs), m_relation(relation) {}
|
||||
};
|
||||
|
||||
|
||||
template <typename T, typename X>
|
||||
class lp_solver : public column_namer {
|
||||
column_info<T> * get_or_create_column_info(unsigned column);
|
||||
|
||||
protected:
|
||||
T get_column_cost_value(unsigned j, column_info<T> * ci) const;
|
||||
public:
|
||||
unsigned m_total_iterations;
|
||||
static_matrix<T, X>* m_A; // this is the matrix of constraints
|
||||
vector<T> m_b; // the right side vector
|
||||
unsigned m_first_stage_iterations;
|
||||
unsigned m_second_stage_iterations;
|
||||
std::unordered_map<unsigned, lp_constraint<T, X>> m_constraints;
|
||||
std::unordered_map<var_index, column_info<T>*> m_map_from_var_index_to_column_info;
|
||||
std::unordered_map<unsigned, std::unordered_map<unsigned, T> > m_A_values;
|
||||
std::unordered_map<std::string, unsigned> m_names_to_columns; // don't have to use it
|
||||
std::unordered_map<unsigned, unsigned> m_external_rows_to_core_solver_rows;
|
||||
std::unordered_map<unsigned, unsigned> m_core_solver_rows_to_external_rows;
|
||||
std::unordered_map<unsigned, unsigned> m_core_solver_columns_to_external_columns;
|
||||
vector<T> m_column_scale;
|
||||
std::unordered_map<unsigned, std::string> m_name_map;
|
||||
unsigned m_artificials;
|
||||
unsigned m_slacks;
|
||||
vector<column_type> m_column_types;
|
||||
vector<T> m_costs;
|
||||
vector<T> m_x;
|
||||
vector<T> m_upper_bounds;
|
||||
vector<unsigned> m_basis;
|
||||
vector<unsigned> m_nbasis;
|
||||
vector<int> m_heading;
|
||||
|
||||
|
||||
lp_status m_status;
|
||||
|
||||
lp_settings m_settings;
|
||||
lp_solver():
|
||||
m_A(nullptr), // this is the matrix of constraints
|
||||
m_first_stage_iterations (0),
|
||||
m_second_stage_iterations (0),
|
||||
m_artificials (0),
|
||||
m_slacks (0),
|
||||
m_status(lp_status::UNKNOWN)
|
||||
{}
|
||||
|
||||
unsigned row_count() const { return this->m_A->row_count(); }
|
||||
|
||||
void add_constraint(lp_relation relation, T right_side, unsigned row_index);
|
||||
|
||||
void set_cost_for_column(unsigned column, T column_cost) {
|
||||
get_or_create_column_info(column)->set_cost(column_cost);
|
||||
}
|
||||
std::string get_variable_name(unsigned j) const override;
|
||||
|
||||
void set_row_column_coefficient(unsigned row, unsigned column, T const & val) {
|
||||
m_A_values[row][column] = val;
|
||||
}
|
||||
// returns the current cost
|
||||
virtual T get_current_cost() const = 0;
|
||||
// do not have to call it
|
||||
void give_symbolic_name_to_column(std::string name, unsigned column);
|
||||
|
||||
virtual T get_column_value(unsigned column) const = 0;
|
||||
|
||||
T get_column_value_by_name(std::string name) const;
|
||||
|
||||
// returns -1 if not found
|
||||
virtual int get_column_index_by_name(std::string name) const;
|
||||
|
||||
void set_lower_bound(unsigned i, T bound) {
|
||||
column_info<T> *ci = get_or_create_column_info(i);
|
||||
ci->set_lower_bound(bound);
|
||||
}
|
||||
|
||||
void set_upper_bound(unsigned i, T bound) {
|
||||
column_info<T> *ci = get_or_create_column_info(i);
|
||||
ci->set_upper_bound(bound);
|
||||
}
|
||||
|
||||
void unset_lower_bound(unsigned i) {
|
||||
get_or_create_column_info(i)->unset_lower_bound();
|
||||
}
|
||||
|
||||
void unset_upper_bound(unsigned i) {
|
||||
get_or_create_column_info(i)->unset_upper_bound();
|
||||
}
|
||||
|
||||
void set_fixed_value(unsigned i, T val) {
|
||||
column_info<T> *ci = get_or_create_column_info(i);
|
||||
ci->set_fixed_value(val);
|
||||
}
|
||||
|
||||
void unset_fixed_value(unsigned i) {
|
||||
get_or_create_column_info(i)->unset_fixed();
|
||||
}
|
||||
|
||||
lp_status get_status() const {
|
||||
return m_status;
|
||||
}
|
||||
|
||||
void set_status(lp_status st) {
|
||||
m_status = st;
|
||||
}
|
||||
|
||||
|
||||
~lp_solver() override;
|
||||
|
||||
void flip_costs();
|
||||
|
||||
virtual void find_maximal_solution() = 0;
|
||||
void set_time_limit(unsigned time_limit_in_seconds) {
|
||||
m_settings.time_limit = time_limit_in_seconds;
|
||||
}
|
||||
|
||||
|
||||
protected:
|
||||
bool problem_is_empty();
|
||||
|
||||
void scale();
|
||||
|
||||
|
||||
void print_rows_scale_stats(std::ostream & out);
|
||||
|
||||
void print_columns_scale_stats(std::ostream & out);
|
||||
|
||||
void print_row_scale_stats(unsigned i, std::ostream & out);
|
||||
|
||||
void print_column_scale_stats(unsigned j, std::ostream & out);
|
||||
|
||||
void print_scale_stats(std::ostream & out);
|
||||
|
||||
void get_max_abs_in_row(std::unordered_map<unsigned, T> & row_map);
|
||||
|
||||
void pin_vars_down_on_row(std::unordered_map<unsigned, T> & row) {
|
||||
pin_vars_on_row_with_sign(row, - numeric_traits<T>::one());
|
||||
}
|
||||
|
||||
void pin_vars_up_on_row(std::unordered_map<unsigned, T> & row) {
|
||||
pin_vars_on_row_with_sign(row, numeric_traits<T>::one());
|
||||
}
|
||||
|
||||
void pin_vars_on_row_with_sign(std::unordered_map<unsigned, T> & row, T sign );
|
||||
|
||||
bool get_minimal_row_value(std::unordered_map<unsigned, T> & row, T & lower_bound);
|
||||
|
||||
bool get_maximal_row_value(std::unordered_map<unsigned, T> & row, T & lower_bound);
|
||||
|
||||
bool row_is_zero(std::unordered_map<unsigned, T> & row);
|
||||
|
||||
bool row_e_is_obsolete(std::unordered_map<unsigned, T> & row, unsigned row_index);
|
||||
|
||||
bool row_ge_is_obsolete(std::unordered_map<unsigned, T> & row, unsigned row_index);
|
||||
|
||||
bool row_le_is_obsolete(std::unordered_map<unsigned, T> & row, unsigned row_index);
|
||||
|
||||
// analyse possible max and min values that are derived from var boundaries
|
||||
// Let us say that the we have a "ge" constraint, and the min value is equal to the rs.
|
||||
// Then we know what values of the variables are. For each positive coeff of the row it has to be
|
||||
// the low boundary of the var and for a negative - the upper.
|
||||
|
||||
// this routing also pins the variables to the boundaries
|
||||
bool row_is_obsolete(std::unordered_map<unsigned, T> & row, unsigned row_index );
|
||||
|
||||
void remove_fixed_or_zero_columns();
|
||||
|
||||
void remove_fixed_or_zero_columns_from_row(unsigned i, std::unordered_map<unsigned, T> & row);
|
||||
|
||||
unsigned try_to_remove_some_rows();
|
||||
|
||||
void cleanup();
|
||||
|
||||
void map_external_rows_to_core_solver_rows();
|
||||
|
||||
void map_external_columns_to_core_solver_columns();
|
||||
|
||||
unsigned number_of_core_structurals() {
|
||||
return static_cast<unsigned>(m_core_solver_columns_to_external_columns.size());
|
||||
}
|
||||
|
||||
void restore_column_scales_to_one() {
|
||||
for (unsigned i = 0; i < m_column_scale.size(); i++) m_column_scale[i] = numeric_traits<T>::one();
|
||||
}
|
||||
|
||||
void unscale();
|
||||
|
||||
void fill_A_from_A_values();
|
||||
|
||||
void fill_matrix_A_and_init_right_side();
|
||||
|
||||
void count_slacks_and_artificials();
|
||||
|
||||
void count_slacks_and_artificials_for_row(unsigned i);
|
||||
|
||||
T lower_bound_shift_for_row(unsigned i);
|
||||
|
||||
void fill_m_b();
|
||||
|
||||
T get_column_value_with_core_solver(unsigned column, lp_core_solver_base<T, X> * core_solver) const;
|
||||
void set_scaled_cost(unsigned j);
|
||||
void print_statistics_on_A(std::ostream & out) {
|
||||
out << "extended A[" << this->m_A->row_count() << "," << this->m_A->column_count() << "]" << std::endl;
|
||||
}
|
||||
|
||||
public:
|
||||
lp_settings & settings() { return m_settings;}
|
||||
void print_model(std::ostream & s) const {
|
||||
s << "objective = " << get_current_cost() << std::endl;
|
||||
s << "column values\n";
|
||||
for (auto & it : m_names_to_columns) {
|
||||
s << it.first << " = " << get_column_value(it.second) << std::endl;
|
||||
}
|
||||
}
|
||||
};
|
||||
}
|
||||
|
|
@ -1,571 +0,0 @@
|
|||
/*++
|
||||
Copyright (c) 2017 Microsoft Corporation
|
||||
|
||||
Module Name:
|
||||
|
||||
<name>
|
||||
|
||||
Abstract:
|
||||
|
||||
<abstract>
|
||||
|
||||
Author:
|
||||
|
||||
Lev Nachmanson (levnach)
|
||||
|
||||
Revision History:
|
||||
|
||||
|
||||
--*/
|
||||
#pragma once
|
||||
|
||||
#include <string>
|
||||
#include <algorithm>
|
||||
#include "util/vector.h"
|
||||
#include "math/lp/lp_solver.h"
|
||||
namespace lp {
|
||||
template <typename T, typename X> column_info<T> * lp_solver<T, X>::get_or_create_column_info(unsigned column) {
|
||||
auto it = m_map_from_var_index_to_column_info.find(column);
|
||||
return (it == m_map_from_var_index_to_column_info.end())? (m_map_from_var_index_to_column_info[column] = new column_info<T>()) : it->second;
|
||||
}
|
||||
|
||||
template <typename T, typename X>
|
||||
std::string lp_solver<T, X>::get_variable_name(unsigned j) const { // j here is the core solver index
|
||||
if (!m_settings.print_external_var_name())
|
||||
return std::string("j")+T_to_string(j);
|
||||
auto it = this->m_core_solver_columns_to_external_columns.find(j);
|
||||
if (it == this->m_core_solver_columns_to_external_columns.end())
|
||||
return std::string("x")+T_to_string(j);
|
||||
unsigned external_j = it->second;
|
||||
auto t = this->m_map_from_var_index_to_column_info.find(external_j);
|
||||
if (t == this->m_map_from_var_index_to_column_info.end()) {
|
||||
return std::string("x") +T_to_string(external_j);
|
||||
}
|
||||
return t->second->get_name();
|
||||
}
|
||||
|
||||
template <typename T, typename X> T lp_solver<T, X>::get_column_cost_value(unsigned j, column_info<T> * ci) const {
|
||||
if (ci->is_fixed()) {
|
||||
return ci->get_cost() * ci->get_fixed_value();
|
||||
}
|
||||
return ci->get_cost() * get_column_value(j);
|
||||
}
|
||||
template <typename T, typename X> void lp_solver<T, X>::add_constraint(lp_relation relation, T right_side, unsigned row_index) {
|
||||
lp_assert(m_constraints.find(row_index) == m_constraints.end());
|
||||
lp_constraint<T, X> cs(right_side, relation);
|
||||
m_constraints[row_index] = cs;
|
||||
}
|
||||
|
||||
template <typename T, typename X> void lp_solver<T, X>::give_symbolic_name_to_column(std::string name, unsigned column) {
|
||||
auto it = m_map_from_var_index_to_column_info.find(column);
|
||||
column_info<T> *ci;
|
||||
if (it == m_map_from_var_index_to_column_info.end()){
|
||||
m_map_from_var_index_to_column_info[column] = ci = new column_info<T>;
|
||||
} else {
|
||||
ci = it->second;
|
||||
}
|
||||
ci->set_name(name);
|
||||
m_names_to_columns[name] = column;
|
||||
}
|
||||
|
||||
|
||||
template <typename T, typename X> T lp_solver<T, X>::get_column_value_by_name(std::string name) const {
|
||||
auto it = m_names_to_columns.find(name);
|
||||
if (it == m_names_to_columns.end()) {
|
||||
std::stringstream s;
|
||||
s << "get_column_value_by_name " << name;
|
||||
throw_exception(s.str());
|
||||
}
|
||||
return get_column_value(it -> second);
|
||||
}
|
||||
|
||||
// returns -1 if not found
|
||||
template <typename T, typename X> int lp_solver<T, X>::get_column_index_by_name(std::string name) const {
|
||||
auto t = m_names_to_columns.find(name);
|
||||
if (t == m_names_to_columns.end()) {
|
||||
return -1;
|
||||
}
|
||||
return t->second;
|
||||
}
|
||||
|
||||
|
||||
template <typename T, typename X> lp_solver<T, X>::~lp_solver(){
|
||||
delete m_A;
|
||||
for (auto t : m_map_from_var_index_to_column_info) {
|
||||
delete t.second;
|
||||
}
|
||||
}
|
||||
|
||||
template <typename T, typename X> void lp_solver<T, X>::flip_costs() {
|
||||
for (auto t : m_map_from_var_index_to_column_info) {
|
||||
column_info<T> *ci = t.second;
|
||||
ci->set_cost(-ci->get_cost());
|
||||
}
|
||||
}
|
||||
|
||||
template <typename T, typename X> bool lp_solver<T, X>::problem_is_empty() {
|
||||
for (auto & c : m_A_values)
|
||||
if (!c.second.empty())
|
||||
return false;
|
||||
return true;
|
||||
}
|
||||
|
||||
template <typename T, typename X> void lp_solver<T, X>::scale() {
|
||||
if (numeric_traits<T>::precise() || m_settings.use_scaling == false) {
|
||||
m_column_scale.clear();
|
||||
m_column_scale.resize(m_A->column_count(), one_of_type<T>());
|
||||
return;
|
||||
}
|
||||
|
||||
T smin = T(m_settings.scaling_minimum);
|
||||
T smax = T(m_settings.scaling_maximum);
|
||||
|
||||
scaler<T, X> scaler(m_b, *m_A, smin, smax, m_column_scale, this->m_settings);
|
||||
if (!scaler.scale()) {
|
||||
unscale();
|
||||
}
|
||||
}
|
||||
|
||||
|
||||
template <typename T, typename X> void lp_solver<T, X>::print_rows_scale_stats(std::ostream & out) {
|
||||
out << "rows max" << std::endl;
|
||||
for (unsigned i = 0; i < m_A->row_count(); i++) {
|
||||
print_row_scale_stats(i, out);
|
||||
}
|
||||
out << std::endl;
|
||||
}
|
||||
|
||||
template <typename T, typename X> void lp_solver<T, X>::print_columns_scale_stats(std::ostream & out) {
|
||||
out << "columns max" << std::endl;
|
||||
for (unsigned i = 0; i < m_A->column_count(); i++) {
|
||||
print_column_scale_stats(i, out);
|
||||
}
|
||||
out << std::endl;
|
||||
}
|
||||
|
||||
template <typename T, typename X> void lp_solver<T, X>::print_row_scale_stats(unsigned i, std::ostream & out) {
|
||||
out << "(" << std::min(m_A->get_min_abs_in_row(i), abs(m_b[i])) << " ";
|
||||
out << std::max(m_A->get_max_abs_in_row(i), abs(m_b[i])) << ")";
|
||||
}
|
||||
|
||||
template <typename T, typename X> void lp_solver<T, X>::print_column_scale_stats(unsigned j, std::ostream & out) {
|
||||
out << "(" << m_A->get_min_abs_in_row(j) << " ";
|
||||
out << m_A->get_max_abs_in_column(j) << ")";
|
||||
}
|
||||
|
||||
template <typename T, typename X> void lp_solver<T, X>::print_scale_stats(std::ostream & out) {
|
||||
print_rows_scale_stats(out);
|
||||
print_columns_scale_stats(out);
|
||||
}
|
||||
|
||||
template <typename T, typename X> void lp_solver<T, X>::get_max_abs_in_row(std::unordered_map<unsigned, T> & row_map) {
|
||||
T ret = numeric_traits<T>::zero();
|
||||
for (auto jp : row_map) {
|
||||
T ac = numeric_traits<T>::abs(jp->second);
|
||||
if (ac > ret) {
|
||||
ret = ac;
|
||||
}
|
||||
}
|
||||
return ret;
|
||||
}
|
||||
|
||||
template <typename T, typename X> void lp_solver<T, X>::pin_vars_on_row_with_sign(std::unordered_map<unsigned, T> & row, T sign ) {
|
||||
for (auto t : row) {
|
||||
unsigned j = t.first;
|
||||
column_info<T> * ci = m_map_from_var_index_to_column_info[j];
|
||||
T a = t.second;
|
||||
if (a * sign > numeric_traits<T>::zero()) {
|
||||
lp_assert(ci->upper_bound_is_set());
|
||||
ci->set_fixed_value(ci->get_upper_bound());
|
||||
} else {
|
||||
lp_assert(ci->lower_bound_is_set());
|
||||
ci->set_fixed_value(ci->get_lower_bound());
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
template <typename T, typename X> bool lp_solver<T, X>::get_minimal_row_value(std::unordered_map<unsigned, T> & row, T & lower_bound) {
|
||||
lower_bound = numeric_traits<T>::zero();
|
||||
for (auto & t : row) {
|
||||
T a = t.second;
|
||||
column_info<T> * ci = m_map_from_var_index_to_column_info[t.first];
|
||||
if (a > numeric_traits<T>::zero()) {
|
||||
if (ci->lower_bound_is_set()) {
|
||||
lower_bound += ci->get_lower_bound() * a;
|
||||
} else {
|
||||
return false;
|
||||
}
|
||||
} else {
|
||||
if (ci->upper_bound_is_set()) {
|
||||
lower_bound += ci->get_upper_bound() * a;
|
||||
} else {
|
||||
return false;
|
||||
}
|
||||
}
|
||||
}
|
||||
return true;
|
||||
}
|
||||
|
||||
template <typename T, typename X> bool lp_solver<T, X>::get_maximal_row_value(std::unordered_map<unsigned, T> & row, T & lower_bound) {
|
||||
lower_bound = numeric_traits<T>::zero();
|
||||
for (auto & t : row) {
|
||||
T a = t.second;
|
||||
column_info<T> * ci = m_map_from_var_index_to_column_info[t.first];
|
||||
if (a < numeric_traits<T>::zero()) {
|
||||
if (ci->lower_bound_is_set()) {
|
||||
lower_bound += ci->get_lower_bound() * a;
|
||||
} else {
|
||||
return false;
|
||||
}
|
||||
} else {
|
||||
if (ci->upper_bound_is_set()) {
|
||||
lower_bound += ci->get_upper_bound() * a;
|
||||
} else {
|
||||
return false;
|
||||
}
|
||||
}
|
||||
}
|
||||
return true;
|
||||
}
|
||||
|
||||
template <typename T, typename X> bool lp_solver<T, X>::row_is_zero(std::unordered_map<unsigned, T> & row) {
|
||||
for (auto & t : row) {
|
||||
if (!is_zero(t.second))
|
||||
return false;
|
||||
}
|
||||
return true;
|
||||
}
|
||||
|
||||
template <typename T, typename X> bool lp_solver<T, X>::row_e_is_obsolete(std::unordered_map<unsigned, T> & row, unsigned row_index) {
|
||||
T rs = m_constraints[row_index].m_rs;
|
||||
if (row_is_zero(row)) {
|
||||
if (!is_zero(rs))
|
||||
m_status = lp_status::INFEASIBLE;
|
||||
return true;
|
||||
}
|
||||
|
||||
T lower_bound;
|
||||
bool lb = get_minimal_row_value(row, lower_bound);
|
||||
if (lb) {
|
||||
T diff = lower_bound - rs;
|
||||
if (!val_is_smaller_than_eps(diff, m_settings.refactor_tolerance)){
|
||||
// lower_bound > rs + m_settings.refactor_epsilon
|
||||
m_status = lp_status::INFEASIBLE;
|
||||
return true;
|
||||
}
|
||||
if (val_is_smaller_than_eps(-diff, m_settings.refactor_tolerance)){
|
||||
pin_vars_down_on_row(row);
|
||||
return true;
|
||||
}
|
||||
}
|
||||
|
||||
T upper_bound;
|
||||
bool ub = get_maximal_row_value(row, upper_bound);
|
||||
if (ub) {
|
||||
T diff = rs - upper_bound;
|
||||
if (!val_is_smaller_than_eps(diff, m_settings.refactor_tolerance)) {
|
||||
// upper_bound < rs - m_settings.refactor_tolerance
|
||||
m_status = lp_status::INFEASIBLE;
|
||||
return true;
|
||||
}
|
||||
if (val_is_smaller_than_eps(-diff, m_settings.refactor_tolerance)){
|
||||
pin_vars_up_on_row(row);
|
||||
return true;
|
||||
}
|
||||
}
|
||||
|
||||
return false;
|
||||
}
|
||||
|
||||
template <typename T, typename X> bool lp_solver<T, X>::row_ge_is_obsolete(std::unordered_map<unsigned, T> & row, unsigned row_index) {
|
||||
T rs = m_constraints[row_index].m_rs;
|
||||
if (row_is_zero(row)) {
|
||||
if (rs > zero_of_type<X>())
|
||||
m_status = lp_status::INFEASIBLE;
|
||||
return true;
|
||||
}
|
||||
|
||||
T upper_bound;
|
||||
if (get_maximal_row_value(row, upper_bound)) {
|
||||
T diff = rs - upper_bound;
|
||||
if (!val_is_smaller_than_eps(diff, m_settings.refactor_tolerance)) {
|
||||
// upper_bound < rs - m_settings.refactor_tolerance
|
||||
m_status = lp_status::INFEASIBLE;
|
||||
return true;
|
||||
}
|
||||
if (val_is_smaller_than_eps(-diff, m_settings.refactor_tolerance)){
|
||||
pin_vars_up_on_row(row);
|
||||
return true;
|
||||
}
|
||||
}
|
||||
|
||||
return false;
|
||||
}
|
||||
|
||||
template <typename T, typename X> bool lp_solver<T, X>::row_le_is_obsolete(std::unordered_map<unsigned, T> & row, unsigned row_index) {
|
||||
T lower_bound;
|
||||
T rs = m_constraints[row_index].m_rs;
|
||||
if (row_is_zero(row)) {
|
||||
if (rs < zero_of_type<X>())
|
||||
m_status = lp_status::INFEASIBLE;
|
||||
return true;
|
||||
}
|
||||
|
||||
if (get_minimal_row_value(row, lower_bound)) {
|
||||
T diff = lower_bound - rs;
|
||||
if (!val_is_smaller_than_eps(diff, m_settings.refactor_tolerance)){
|
||||
// lower_bound > rs + m_settings.refactor_tolerance
|
||||
m_status = lp_status::INFEASIBLE;
|
||||
return true;
|
||||
}
|
||||
if (val_is_smaller_than_eps(-diff, m_settings.refactor_tolerance)){
|
||||
pin_vars_down_on_row(row);
|
||||
return true;
|
||||
}
|
||||
}
|
||||
|
||||
return false;
|
||||
}
|
||||
|
||||
// analyse possible max and min values that are derived from var boundaries
|
||||
// Let us say that the we have a "ge" constraint, and the min value is equal to the rs.
|
||||
// Then we know what values of the variables are. For each positive coeff of the row it has to be
|
||||
// the low boundary of the var and for a negative - the upper.
|
||||
|
||||
// this routing also pins the variables to the boundaries
|
||||
template <typename T, typename X> bool lp_solver<T, X>::row_is_obsolete(std::unordered_map<unsigned, T> & row, unsigned row_index ) {
|
||||
auto & constraint = m_constraints[row_index];
|
||||
switch (constraint.m_relation) {
|
||||
case lp_relation::Equal:
|
||||
return row_e_is_obsolete(row, row_index);
|
||||
|
||||
case lp_relation::Greater_or_equal:
|
||||
return row_ge_is_obsolete(row, row_index);
|
||||
|
||||
case lp_relation::Less_or_equal:
|
||||
return row_le_is_obsolete(row, row_index);
|
||||
}
|
||||
lp_unreachable();
|
||||
return false; // it is unreachable
|
||||
}
|
||||
|
||||
template <typename T, typename X> void lp_solver<T, X>::remove_fixed_or_zero_columns() {
|
||||
for (auto & i_row : m_A_values) {
|
||||
remove_fixed_or_zero_columns_from_row(i_row.first, i_row.second);
|
||||
}
|
||||
}
|
||||
|
||||
template <typename T, typename X> void lp_solver<T, X>::remove_fixed_or_zero_columns_from_row(unsigned i, std::unordered_map<unsigned, T> & row) {
|
||||
auto & constraint = m_constraints[i];
|
||||
vector<unsigned> removed;
|
||||
for (auto & col : row) {
|
||||
unsigned j = col.first;
|
||||
lp_assert(m_map_from_var_index_to_column_info.find(j) != m_map_from_var_index_to_column_info.end());
|
||||
column_info<T> * ci = m_map_from_var_index_to_column_info[j];
|
||||
if (ci->is_fixed()) {
|
||||
removed.push_back(j);
|
||||
T aj = col.second;
|
||||
constraint.m_rs -= aj * ci->get_fixed_value();
|
||||
} else {
|
||||
if (numeric_traits<T>::is_zero(col.second)){
|
||||
removed.push_back(j);
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
for (auto j : removed) {
|
||||
row.erase(j);
|
||||
}
|
||||
}
|
||||
|
||||
template <typename T, typename X> unsigned lp_solver<T, X>::try_to_remove_some_rows() {
|
||||
vector<unsigned> rows_to_delete;
|
||||
for (auto & t : m_A_values) {
|
||||
if (row_is_obsolete(t.second, t.first)) {
|
||||
rows_to_delete.push_back(t.first);
|
||||
}
|
||||
|
||||
if (m_status == lp_status::INFEASIBLE) {
|
||||
return 0;
|
||||
}
|
||||
}
|
||||
if (!rows_to_delete.empty()) {
|
||||
for (unsigned k : rows_to_delete) {
|
||||
m_A_values.erase(k);
|
||||
}
|
||||
}
|
||||
remove_fixed_or_zero_columns();
|
||||
return static_cast<unsigned>(rows_to_delete.size());
|
||||
}
|
||||
|
||||
template <typename T, typename X> void lp_solver<T, X>::cleanup() {
|
||||
int n = 0; // number of deleted rows
|
||||
int d;
|
||||
while ((d = try_to_remove_some_rows()) > 0)
|
||||
n += d;
|
||||
|
||||
if (n == 1) {
|
||||
LP_OUT(m_settings, "deleted one row" << std::endl);
|
||||
} else if (n) {
|
||||
LP_OUT(m_settings, "deleted " << n << " rows" << std::endl);
|
||||
}
|
||||
}
|
||||
|
||||
template <typename T, typename X> void lp_solver<T, X>::map_external_rows_to_core_solver_rows() {
|
||||
unsigned size = 0;
|
||||
for (auto & row : m_A_values) {
|
||||
m_external_rows_to_core_solver_rows[row.first] = size;
|
||||
m_core_solver_rows_to_external_rows[size] = row.first;
|
||||
size++;
|
||||
}
|
||||
}
|
||||
|
||||
template <typename T, typename X> void lp_solver<T, X>::map_external_columns_to_core_solver_columns() {
|
||||
unsigned size = 0;
|
||||
for (auto & row : m_A_values) {
|
||||
for (auto & col : row.second) {
|
||||
if (col.second == numeric_traits<T>::zero() || m_map_from_var_index_to_column_info[col.first]->is_fixed()) {
|
||||
throw_exception("found fixed column");
|
||||
}
|
||||
unsigned j = col.first;
|
||||
auto column_info_it = m_map_from_var_index_to_column_info.find(j);
|
||||
lp_assert(column_info_it != m_map_from_var_index_to_column_info.end());
|
||||
|
||||
auto j_column = column_info_it->second->get_column_index();
|
||||
if (!is_valid(j_column)) { // j is a newcomer
|
||||
m_map_from_var_index_to_column_info[j]->set_column_index(size);
|
||||
m_core_solver_columns_to_external_columns[size++] = j;
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
template <typename T, typename X> void lp_solver<T, X>::unscale() {
|
||||
delete m_A;
|
||||
m_A = nullptr;
|
||||
fill_A_from_A_values();
|
||||
restore_column_scales_to_one();
|
||||
fill_m_b();
|
||||
}
|
||||
|
||||
template <typename T, typename X> void lp_solver<T, X>::fill_A_from_A_values() {
|
||||
m_A = new static_matrix<T, X>(static_cast<unsigned>(m_A_values.size()), number_of_core_structurals());
|
||||
for (auto & t : m_A_values) {
|
||||
auto row_it = m_external_rows_to_core_solver_rows.find(t.first);
|
||||
lp_assert(row_it != m_external_rows_to_core_solver_rows.end());
|
||||
unsigned row = row_it->second;
|
||||
for (auto k : t.second) {
|
||||
auto column_info_it = m_map_from_var_index_to_column_info.find(k.first);
|
||||
lp_assert(column_info_it != m_map_from_var_index_to_column_info.end());
|
||||
column_info<T> *ci = column_info_it->second;
|
||||
unsigned col = ci->get_column_index();
|
||||
lp_assert(is_valid(col));
|
||||
bool col_is_flipped = m_map_from_var_index_to_column_info[k.first]->is_flipped();
|
||||
if (!col_is_flipped) {
|
||||
(*m_A)(row, col) = k.second;
|
||||
} else {
|
||||
(*m_A)(row, col) = - k.second;
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
template <typename T, typename X> void lp_solver<T, X>::fill_matrix_A_and_init_right_side() {
|
||||
map_external_rows_to_core_solver_rows();
|
||||
map_external_columns_to_core_solver_columns();
|
||||
lp_assert(m_A == nullptr);
|
||||
fill_A_from_A_values();
|
||||
m_b.resize(m_A->row_count());
|
||||
}
|
||||
|
||||
template <typename T, typename X> void lp_solver<T, X>::count_slacks_and_artificials() {
|
||||
for (int i = row_count() - 1; i >= 0; i--) {
|
||||
count_slacks_and_artificials_for_row(i);
|
||||
}
|
||||
}
|
||||
|
||||
template <typename T, typename X> void lp_solver<T, X>::count_slacks_and_artificials_for_row(unsigned i) {
|
||||
lp_assert(this->m_constraints.find(this->m_core_solver_rows_to_external_rows[i]) != this->m_constraints.end());
|
||||
auto & constraint = this->m_constraints[this->m_core_solver_rows_to_external_rows[i]];
|
||||
switch (constraint.m_relation) {
|
||||
case Equal:
|
||||
m_artificials++;
|
||||
break;
|
||||
case Greater_or_equal:
|
||||
m_slacks++;
|
||||
if (this->m_b[i] > 0) {
|
||||
m_artificials++;
|
||||
}
|
||||
break;
|
||||
case Less_or_equal:
|
||||
m_slacks++;
|
||||
if (this->m_b[i] < 0) {
|
||||
m_artificials++;
|
||||
}
|
||||
break;
|
||||
}
|
||||
}
|
||||
|
||||
template <typename T, typename X> T lp_solver<T, X>::lower_bound_shift_for_row(unsigned i) {
|
||||
T ret = numeric_traits<T>::zero();
|
||||
|
||||
auto row = this->m_A_values.find(i);
|
||||
if (row == this->m_A_values.end()) {
|
||||
throw_exception("cannot find row");
|
||||
}
|
||||
for (auto col : row->second) {
|
||||
ret += col.second * this->m_map_from_var_index_to_column_info[col.first]->get_shift();
|
||||
}
|
||||
return ret;
|
||||
}
|
||||
|
||||
template <typename T, typename X> void lp_solver<T, X>::fill_m_b() {
|
||||
for (int i = this->row_count() - 1; i >= 0; i--) {
|
||||
lp_assert(this->m_constraints.find(this->m_core_solver_rows_to_external_rows[i]) != this->m_constraints.end());
|
||||
unsigned external_i = this->m_core_solver_rows_to_external_rows[i];
|
||||
auto & constraint = this->m_constraints[external_i];
|
||||
this->m_b[i] = constraint.m_rs - lower_bound_shift_for_row(external_i);
|
||||
}
|
||||
}
|
||||
|
||||
template <typename T, typename X> T lp_solver<T, X>::get_column_value_with_core_solver(unsigned column, lp_core_solver_base<T, X> * core_solver) const {
|
||||
auto cit = this->m_map_from_var_index_to_column_info.find(column);
|
||||
if (cit == this->m_map_from_var_index_to_column_info.end()) {
|
||||
return numeric_traits<T>::zero();
|
||||
}
|
||||
|
||||
column_info<T> * ci = cit->second;
|
||||
|
||||
if (ci->is_fixed()) {
|
||||
return ci->get_fixed_value();
|
||||
}
|
||||
|
||||
unsigned cj = ci->get_column_index();
|
||||
if (cj != static_cast<unsigned>(-1)) {
|
||||
T v = core_solver->get_var_value(cj) * this->m_column_scale[cj];
|
||||
if (ci->is_free()) {
|
||||
return v;
|
||||
}
|
||||
if (!ci->is_flipped()) {
|
||||
return v + ci->get_lower_bound();
|
||||
}
|
||||
|
||||
// the flipped case when there is only upper bound
|
||||
return -v + ci->get_upper_bound(); //
|
||||
}
|
||||
|
||||
return numeric_traits<T>::zero(); // returns zero for out of boundary columns
|
||||
}
|
||||
|
||||
template <typename T, typename X> void lp_solver<T, X>::set_scaled_cost(unsigned j) {
|
||||
// grab original costs but modify it with the column scales
|
||||
lp_assert(j < this->m_column_scale.size());
|
||||
column_info<T> * ci = this->m_map_from_var_index_to_column_info[this->m_core_solver_columns_to_external_columns[j]];
|
||||
T cost = ci->get_cost();
|
||||
if (ci->is_flipped()){
|
||||
cost *= T(-1);
|
||||
}
|
||||
lp_assert(ci->is_fixed() == false);
|
||||
this->m_costs[j] = cost * this->m_column_scale[j];
|
||||
}
|
||||
}
|
||||
Loading…
Add table
Add a link
Reference in a new issue