Commit 1a695420 authored by jonasseidel's avatar jonasseidel
Browse files

linear programming overhaul (includes closing of memory leaks by implementing...

linear programming overhaul (includes closing of memory leaks by implementing proper destr and copy constr)
parent 7697d3a8
......@@ -6,8 +6,8 @@ Coefficient::Coefficient(){}
Coefficient::Coefficient(double initial_val) : _coeff({{{}, initial_val}}) {}
bool Coefficient::is_non_zero(){
for(std::pair<const Monomial, double>& m : this->_coeff){
bool Coefficient::is_non_zero() const {
for(const std::pair<const Monomial, double>& m : this->_coeff){
if(std::abs(m.second) > 1e-100){
return true;
}
......@@ -15,17 +15,14 @@ bool Coefficient::is_non_zero(){
return false;
}
bool Coefficient::is_unity(){
bool Coefficient::is_unity() const {
auto search_result = this->_coeff.find({});
if(search_result != this->_coeff.end()){
if(std::abs(search_result->second - 1) > 1e-100){
return false;
}else{
auto tmp = search_result->second;
search_result->second = 0;
bool non_unity = this->is_non_zero();
search_result->second = tmp;
if(non_unity){
if( (Coefficient{-1}+(*this)).is_non_zero() ){
return false;
}
}
......@@ -34,7 +31,7 @@ bool Coefficient::is_unity(){
}
double Coefficient::value(){
double Coefficient::value() const {
double curr = 0;
for(std::pair<Monomial, double> pair : this->_coeff){
curr += pair.second*pair.first.value();
......
......@@ -29,11 +29,25 @@ public:
friend std::ostream& operator<<(std::ostream& os, Coefficient coeff);
Coefficient();
Coefficient(const Coefficient& other) : _coeff(other._coeff){}
Coefficient(double initial_val);
bool is_non_zero();
bool is_unity();
double value();
Coefficient(const Coefficient& other, const std::unordered_map<const Variable*, Variable*>& variable_lookup){
for(auto [monom, coeff_sec_order] : other._coeff){
this->_coeff.insert({Monomial(monom, variable_lookup), coeff_sec_order});
}
}
Coefficient(Coefficient&& other) : _coeff(std::move(other._coeff)){}
void operator=(const Coefficient& other){
this->_coeff = other._coeff;
}
void operator=(Coefficient&& other){
this->_coeff = std::move(other._coeff);
}
bool is_non_zero() const ;
bool is_unity() const ;
double value() const ;
static const bool always_add_sign = true;
};
......
#include "Constraint.h"
Constraint::Constraint(const Constraint& other){
*this = other;
}
Constraint::Constraint(Constraint&& other){
*this = std::move(other);
}
Constraint::Constraint(std::string description, relation relation_type, std::unordered_map<Variable*, Coefficient> lhs, Coefficient rhs)
: _description(description), _relation(relation_type), _lhs(lhs), _rhs(rhs){}
Constraint::Constraint(relation relation_type, std::unordered_map<Variable*, Coefficient> lhs, Coefficient rhs)
: _relation(relation_type), _lhs(lhs), _rhs(rhs){}
std::string Constraint::description(){
Constraint::Constraint(const Constraint& other, const std::unordered_map<const Variable*, Variable*>& variable_lookup)
: _description(other._description),
_relation(other._relation),
_rhs(other._rhs, variable_lookup)
{
for(const auto& [var, coeff] : other._lhs){
auto var_search = variable_lookup.find(var);
if(var_search == variable_lookup.end()) throw std::invalid_argument("variable_lookup didn't include necessary entries");
this->_lhs.insert({var_search->second, Coefficient(coeff, variable_lookup)});
}
}
std::string Constraint::description() const {
return this->_description;
}
relation Constraint::is_equality(){
relation Constraint::is_equality() const {
return this->_relation;
}
......@@ -44,12 +63,25 @@ SCIP_CONS* Constraint::computational_con(SCIP* scip, std::unordered_map<Variable
return cons;
}
std::ostream& operator<<(std::ostream& os, Constraint& constraint){
void Constraint::operator=(const Constraint& other){
this->_description = other._description;
this->_relation = other._relation;
this->_lhs = other._lhs;
this->_rhs = other._rhs;
}
void Constraint::operator=(Constraint&& other){
this->_description = std::move(other._description);
this->_relation = std::move(other._relation);
this->_lhs = std::move(other._lhs);
this->_rhs = std::move(other._rhs);
}
std::ostream& operator<<(std::ostream& os, const Constraint& constraint){
os << constraint.description() << ":\t";
bool empty = true;
bool add_plus = false;
for(auto summand : constraint.lhs()){
for(auto summand : constraint._lhs){
if(summand.second.is_non_zero()){
empty = false;
if(add_plus && !Coefficient::always_add_sign){
......@@ -67,12 +99,12 @@ std::ostream& operator<<(std::ostream& os, Constraint& constraint){
}else{
os << " <= ";
}
if(!constraint.rhs().is_non_zero()){
if(!constraint._rhs.is_non_zero()){
os << "0";
}else if(constraint.rhs().is_unity()){
}else if(constraint._rhs.is_unity()){
os << "1";
}else{
os << constraint.rhs();
os << constraint._rhs;
}
return os;
}
......@@ -21,18 +21,27 @@ class Constraint{
std::unordered_map<Variable*, Coefficient> _lhs;
Coefficient _rhs;
public:
Constraint(){}
Constraint(const Constraint& other);
Constraint(Constraint&& other);
Constraint(std::string description, relation relation_type, std::unordered_map<Variable*, Coefficient> lhs, Coefficient rhs);
Constraint(relation relation_type, std::unordered_map<Variable*, Coefficient> lhs, Coefficient rhs);
std::string description();
relation is_equality();
Constraint(const Constraint& other, const std::unordered_map<const Variable*, Variable*>& variable_lookup);
std::string description() const ;
relation is_equality() const ;
std::unordered_map<Variable*, Coefficient> lhs();
Coefficient lhs_coefficient(Variable* var);
Coefficient rhs();
void operator=(const Constraint& other);
void operator=(Constraint&& other);
SCIP_CONS* computational_con(SCIP* scip, std::unordered_map<Variable*, SCIP_VAR*>& variable_lookup);
friend std::ostream& operator<<(std::ostream& os, const Constraint& constraint);
};
std::ostream& operator<<(std::ostream& os, Constraint& constraint);
std::ostream& operator<<(std::ostream& os, const Constraint& constraint);
#endif
#include "Linear_Program.h"
Linear_Program::Linear_Program(Linear_Program& lp){
Linear_Program::Linear_Program(const Linear_Program& lp){
*this = lp;
}
......@@ -10,18 +10,15 @@ Linear_Program::Linear_Program(Linear_Program&& lp){
Linear_Program::Linear_Program(std::string description, bool maximum) : _description(description), _maximum(maximum) {}
Linear_Program::Linear_Program(std::string description, bool maximum, std::unordered_map<Variable*, Coefficient> direction, Polyeder polyeder)
: _description(description), _maximum(std::move(maximum)), _direction(std::move(direction)), _polyeder(std::move(polyeder)) {}
std::string Linear_Program::description(){
std::string Linear_Program::description() const {
return this->_description;
}
bool Linear_Program::is_maximum(){
bool Linear_Program::is_maximum() const {
return this->_maximum;
}
std::unordered_map<Variable*, Coefficient> Linear_Program::direction(){
std::unordered_map<Variable*, Coefficient>& Linear_Program::direction(){
return this->_direction;
}
......@@ -39,6 +36,10 @@ Polyeder& Linear_Program::polyeder(){
return this->_polyeder;
}
const Polyeder& Linear_Program::polyeder() const {
return this->_polyeder;
}
std::pair<SCIP*, std::unordered_map<Variable*, SCIP_VAR*> > Linear_Program::computational_model(){
SCIP* scip;
SCIP_CALL_ABORT( SCIPcreate(&scip));
......@@ -70,41 +71,33 @@ std::pair<SCIP*, std::unordered_map<Variable*, SCIP_VAR*> > Linear_Program::comp
return {scip, variable_lookup};
}
void Linear_Program::operator=(Linear_Program& lp){
void Linear_Program::operator=(const Linear_Program& lp){
this->_description = lp._description;
this->_maximum = lp._maximum;
this->_direction = lp._direction;
this->_polyeder = lp._polyeder;
auto [copy, variable_lookup] = lp._polyeder.copy_ret_lookup();
this->_polyeder = std::move(copy);
this->_direction = {};
for( const auto& [var, coeff] : lp._direction){
auto variable_search = variable_lookup.find(var);
if( variable_search == variable_lookup.end() ) throw std::range_error("lookup doesn't include all variables");
this->_direction.insert({variable_search->second, Coefficient(coeff, variable_lookup)});
}
}
void Linear_Program::operator=(Linear_Program&& lp){
this->_description = lp._description;
this->_maximum = lp._maximum;
this->_direction = std::move(lp._direction);
lp._direction = std::unordered_map<Variable*, Coefficient>();
this->_polyeder = std::move(lp._polyeder);
lp._polyeder = Polyeder();
lp._direction = std::unordered_map<Variable*, Coefficient>();
}
Linear_Program Linear_Program::relaxation_dual(){
std::vector<Constraint> constraints;
constraints.reserve(this->polyeder().vs_dim());
class bimap_counter{
boost::bimap<Variable*, size_t> _bimap;
size_t _counter;
public:
bimap_counter() : _counter(0) {}
void insert(Variable* var){
_bimap.insert({var, _counter});
_counter++;
}
boost::bimap<Variable*, size_t>& unwrap(){
return _bimap;
}
} variables;
std::unordered_map<Variable*, Coefficient> direction;
std::stringstream dual_name;
dual_name << this->description() << "_dual";
Linear_Program dual(dual_name.str(), !this->is_maximum());
Polyeder& p = dual.polyeder();
for(size_t index = 0; index < this->polyeder().number_of_constraints(); ++index){
std::stringstream name;
......@@ -112,9 +105,9 @@ Linear_Program Linear_Program::relaxation_dual(){
if(this->polyeder().constraint(index).is_equality()){
variables.insert({new Variable(name.str(), Continuous)});
p.add_variable(Variable(name.str(), Continuous));
}else{
variables.insert(new Variable(name.str(), Continuous, {true, 0}));
p.add_variable(Variable(name.str(), Continuous, {true, 0}));
}
}
......@@ -123,7 +116,7 @@ Linear_Program Linear_Program::relaxation_dual(){
// gather lhs_coefficient
std::unordered_map<Variable*, Coefficient> lhs;
for(size_t index = 0; index < this->polyeder().number_of_constraints(); ++index){
lhs.insert({variables.unwrap().right.find(index)->second, {-this->polyeder().constraint(index).lhs_coefficient(curr_variable.left)}});
lhs.insert({p.variables().right.find(index)->second, {-this->polyeder().constraint(index).lhs_coefficient(curr_variable.left)}});
}
......@@ -137,46 +130,40 @@ Linear_Program Linear_Program::relaxation_dual(){
if(curr_variable.left->upper_bound().first == true){
std::stringstream ub_name;
ub_name << "upper_bound_" << curr_variable.left->description();
Variable* upper_bound_var = new Variable(ub_name.str(), Continuous, {true, 0});
variables.insert(upper_bound_var);
Variable* upper_bound_var = p.add_variable(Variable(ub_name.str(), Continuous, {true, 0})).first;
lhs.insert({upper_bound_var, {-1}});
direction.insert({upper_bound_var, curr_variable.left->upper_bound().second});
dual.add_direction_coefficient({upper_bound_var, curr_variable.left->upper_bound().second});
}
}
{
if(curr_variable.left->lower_bound().first == true){
if(curr_variable.left->lower_bound().second == 0){
constraints.push_back(Constraint(name.str(), Inequality, lhs, -this->direction_coefficient(curr_variable.left)));
p.add_constraint(Constraint(name.str(), Inequality, lhs, -this->direction_coefficient(curr_variable.left)));
}else{
std::stringstream lb_name;
lb_name << "lower_bound_" << curr_variable.left->description();
Variable* lower_bound_var = new Variable(lb_name.str(), Continuous, {true, 0});
variables.insert(lower_bound_var);
Variable* lower_bound_var = p.add_variable(Variable(lb_name.str(), Continuous, {true, 0})).first;
lhs.insert({lower_bound_var, {1}});
constraints.push_back(Constraint(name.str(), Equality, lhs, -this->direction_coefficient(curr_variable.left)));
direction.insert({lower_bound_var, {-1}});
p.add_constraint(Constraint(name.str(), Equality, lhs, -this->direction_coefficient(curr_variable.left)));
dual.add_direction_coefficient({lower_bound_var, {-1}});
}
}else{
constraints.push_back(Constraint(name.str(), Equality, lhs, -this->direction_coefficient(curr_variable.left)));
p.add_constraint(Constraint(name.str(), Equality, lhs, -this->direction_coefficient(curr_variable.left)));
}
}
}
Polyeder polyeder (constraints, variables.unwrap());
for(size_t index = 0; index < this->polyeder().number_of_constraints(); ++index){
direction.insert({variables.unwrap().right.find(index)->second, this->polyeder().constraint(index).rhs()});
dual.add_direction_coefficient({p.variables().right.find(index)->second, this->polyeder().constraint(index).rhs()});
}
std::stringstream dual_name;
dual_name << this->description() << "_dual";
return Linear_Program(dual_name.str(), !this->is_maximum(), direction, polyeder);
return dual;
}
std::ostream& operator<<(std::ostream& os, Linear_Program& linear_program){
std::ostream& operator<<(std::ostream& os, const Linear_Program& linear_program){
if(linear_program.is_maximum()){
os << "Maximize\n";
}else{
......@@ -184,9 +171,9 @@ std::ostream& operator<<(std::ostream& os, Linear_Program& linear_program){
}
os << "\tobj: ";
bool empty = true;
if(linear_program.direction().size() > 0){
if(linear_program._direction.size() > 0){
bool add_plus = false;
for(std::pair<Variable*, Coefficient> direction_summand : linear_program.direction()){
for(std::pair<Variable*, Coefficient> direction_summand : linear_program._direction){
if(direction_summand.second.is_non_zero()){
empty = false;
if(add_plus && !Coefficient::always_add_sign){
......@@ -202,7 +189,7 @@ std::ostream& operator<<(std::ostream& os, Linear_Program& linear_program){
}
os << "\n";
os << linear_program.polyeder() << "\n";
os << linear_program._polyeder << "\n";
os << "End";
return os;
}
......@@ -21,31 +21,33 @@
class Linear_Program{
bool _maximum;
std::string _description;
std::unordered_map<Variable*, Coefficient> _direction;
Polyeder _polyeder;
std::unordered_map<Variable*, Coefficient> _direction;
public:
Linear_Program(Linear_Program& lp);
Linear_Program(const Linear_Program& lp);
Linear_Program(Linear_Program&& lp);
Linear_Program(std::string description, bool maximum);
Linear_Program(std::string description, bool maximum, std::unordered_map<Variable*, Coefficient> direction, Polyeder polyeder);
std::string description();
std::string description() const ;
bool is_maximum();
std::unordered_map<Variable*, Coefficient> direction();
bool is_maximum() const ;
std::unordered_map<Variable*, Coefficient>& direction();
void add_direction_coefficient(std::pair<Variable*, Coefficient> summand);
Coefficient direction_coefficient(Variable* index);
Polyeder& polyeder();
const Polyeder& polyeder() const ;
std::pair<SCIP*, std::unordered_map<Variable*, SCIP_VAR*> > computational_model();
void operator=(Linear_Program& lp);
void operator=(const Linear_Program& lp);
void operator=(Linear_Program&& lp);
Linear_Program relaxation_dual();
Linear_Program Danzig_Wolfe_IPM_complete_optimal_inequality_relaxation();
Linear_Program Benders_Reformulierung_integral_relaxation();
friend std::ostream& operator<<(std::ostream& os, const Linear_Program& linear_program);
};
std::ostream& operator<<(std::ostream& os, Linear_Program& linear_program);
std::ostream& operator<<(std::ostream& os, const Linear_Program& linear_program);
#endif
#include "Monomial.h"
Monomial::Monomial(const Monomial& other) : _vars(other._vars) {}
Monomial::Monomial(Monomial&& other) : _vars(std::move(other._vars)) {}
Monomial::Monomial(std::set<Variable*> vars) : _vars(vars) {}
Monomial::Monomial(const Monomial& other, const std::unordered_map<const Variable*, Variable*>& variable_lookup){
for(Variable* var : other._vars){
auto variable_search = variable_lookup.find(var);
if(variable_search == variable_lookup.end()) throw std::invalid_argument("variable_lookup didn't include all relevant variables");
this->_vars.insert(variable_search->second);
}
}
double Monomial::value(){
double curr = 1;
for(Variable* var : this->_vars){
......@@ -11,7 +21,6 @@ double Monomial::value(){
}
return curr;
}
bool Monomial::less::operator()(const Monomial a, const Monomial b) const {
std::stringstream key_a;
for(auto var_ptr : a._vars){
......
......@@ -2,6 +2,8 @@
#define MONOMIAL_H
#include <set>
#include <unordered_set>
#include <unordered_map>
#include <sstream>
#include "Variable.h"
......@@ -10,7 +12,9 @@ class Monomial{
std::set<Variable*> _vars;
public:
Monomial(const Monomial& other);
Monomial(Monomial&& other);
Monomial(std::set<Variable*> vars = {});
Monomial(const Monomial& other, const std::unordered_map<const Variable*, Variable*>& variable_lookup);
double value();
......
......@@ -6,15 +6,30 @@ Polyeder::Polyeder(const Polyeder& p){
*this = p;
}
Polyeder::Polyeder(std::vector<Constraint> constraints, boost::bimap<Variable*, size_t> variables) : _constraints(constraints), _variables(variables) {
// TODO: check that variables used in constraints are registered in variables?
Polyeder::Polyeder(Polyeder&& p){
*this = std::move(p);
}
std::pair<Polyeder, std::unordered_map<const Variable*, Variable*>> Polyeder::copy_ret_lookup() const {
Polyeder copy;
std::unordered_map<const Variable*, Variable*> variable_lookup;
for(auto& [var, id] : this->_variables){
auto [new_var, new_id] = copy.add_variable(*var);
variable_lookup.insert({var, new_var});
}
for(const Constraint& cons : this->_constraints){
copy.add_constraint(cons, variable_lookup);
}
return {std::move(copy), std::move(variable_lookup)};
}
boost::bimap<Variable*, size_t>& Polyeder::variables(){
return this->_variables;
}
size_t Polyeder::vs_dim(){
size_t Polyeder::vs_dim() const {
return this->_variables.size();
}
......@@ -24,20 +39,26 @@ Variable& Polyeder::variable(size_t index){
return *search_result->second;
}
size_t Polyeder::variable(Variable& var){
const Variable& Polyeder::variable(size_t index) const {
auto search_result = this->_variables.right.find(index);
if(search_result == this->_variables.right.end()) throw std::range_error("no such variable!");
return *search_result->second;
}
size_t Polyeder::variable(Variable& var) const {
auto search_result = this->_variables.left.find(&var);
if(search_result == this->_variables.left.end()) throw std::range_error("no such variable!");
return search_result->second;
}
std::pair<Variable*, size_t> Polyeder::add_variable(Variable var){
std::pair<Variable*, size_t> Polyeder::add_variable(const Variable& var){
Variable* new_var = new Variable(var);
size_t position = this->_variables.size();
this->_variables.insert({new_var, position});
return {new_var, position};
}
std::string Polyeder::variable_identifier(size_t index){
std::string Polyeder::variable_identifier(size_t index) const {
std::stringstream name;
if(Variable::VERBOSE_IDENT){
name << this->variable(index).description();
......@@ -47,11 +68,15 @@ std::string Polyeder::variable_identifier(size_t index){
return name.str();
}
std::vector<Constraint> Polyeder::constraints(){
std::vector<Constraint>& Polyeder::constraints(){
return this->_constraints;
}
size_t Polyeder::number_of_constraints(){
const std::vector<Constraint>& Polyeder::constraints() const {
return this->_constraints;
}
size_t Polyeder::number_of_constraints() const {
return this->_constraints.size();
}
......@@ -59,28 +84,40 @@ Constraint& Polyeder::constraint(size_t index){
return this->_constraints[index];
}
void Polyeder::add_constraint(Constraint constraint){
void Polyeder::add_constraint(const Constraint& constraint){
// TODO: check that variables used in constraint are registered in polyeder?
this->_constraints.push_back(constraint);
}
void Polyeder::add_constraint(const Constraint& constraint, const std::unordered_map<const Variable*, Variable*>& variable_lookup){
// TODO: check that variables used in constraint are registered in polyeder?
this->_constraints.push_back(Constraint(constraint, variable_lookup));
}
Polyeder::~Polyeder(){
// for(auto var : this->_variables){
// delete var.left;
// }
for(auto& [var, is] : this->_variables){
delete var;
}
}
void Polyeder::operator=(const Polyeder& p){
this->_constraints = p._constraints;
this->_variables = p._variables;
*this = std::move(this->copy_ret_lookup().first);
}
void Polyeder::operator=(Polyeder&& p){
this->_constraints = std::move(p._constraints);
this->_variables = std::move(p._variables);
if(this != &p){
for(auto& [vars, index] : this->_variables){
delete vars;
}
}
this->_variables = p._variables;