Skip to content

Commit

Permalink
allow for quadratic terms in Constant
Browse files Browse the repository at this point in the history
  • Loading branch information
@hlefebvr
hlefebvr committed Jul 31, 2023
1 parent 70c34ce commit 7b412d6
Show file tree
Hide file tree
Showing 4 changed files with 108 additions and 37 deletions.
49 changes: 35 additions & 14 deletions lib/include/modeling/expressions/Constant.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -8,19 +8,28 @@
#include "../../containers/Map.h"
#include "../parameters/Param.h"
#include "../numericals.h"
#include "../../containers/IteratorForward.h"

namespace idol {

namespace Solution {
class Primal;

class Dual;
}

class Algorithm;
class Constant;
class QuadParam;
}

struct idol::QuadParam {
const Param& key1;
const Param& key2;
const Constant& constant;
QuadParam(const std::pair<Param, Param>& t_pair, const Constant& t_constant)
: key1(t_pair.first), key2(t_pair.second), constant(t_constant) {}
};

/**
* Constant term modeling object.
*
Expand All @@ -38,10 +47,12 @@ namespace idol {
* 2 * xi where xi is a Param will yield a Constant.
*/
class idol::Constant {
Map<Param, double> m_products;
Map<Param, double> m_linear_terms;
Map<std::pair<Param, Param>, double, idol::impl::symmetric_pair_hash, idol::impl::symmetric_pair_equal_to> m_quadratic_terms;
double m_constant = 0.;

void insert_or_add(const Param& t_param, double t_value);
void insert_or_add(const Param& t_param_1, const Param& t_param_2, double t_value);
public:
/**
* Create a new constant term equal to zero.
Expand Down Expand Up @@ -76,6 +87,16 @@ class idol::Constant {
*/
void set(const Param& t_param, double t_value);

/**
* Sets the factor for a given pair of parameters in the constant term.
*
* If t_factor is zero, the entry for t_param is removed (i.e., not stored).
* @param t_param_1 The parameter 1.
* @param t_param_2 The parameter 3.
* @param t_value The factor associated to the pair of parameters.
*/
void set(const Param& t_param_1, const Param& t_param_2, double t_value);

/**
* Returns the factor associated to the parameter given as argument.
*
Expand All @@ -84,6 +105,8 @@ class idol::Constant {
*/
double get(const Param& t_param) const;

double get(const Param& t_param_1, const Param& t_param_2) const;

/**
* Sets the numerical term equal to the constant given as parameter.
* @param t_constant The numerical term.
Expand All @@ -103,22 +126,20 @@ class idol::Constant {
/**
* Returns the number of Param-double pairs stored inside the Constant.
*/
unsigned int size() const { return m_products.size(); }
unsigned int size() const { return m_linear_terms.size(); }

/**
* Returns true if the Constant does not contain any Param-double pair.
*/
bool is_numerical() const;

using iterator = Map<Param, double>::iterator;
using const_iterator = Map<Param, double>::const_iterator;
auto linear() { return IteratorForward(m_linear_terms); }

auto linear() const { return ConstIteratorForward(m_linear_terms); }

auto quadratic() { return IteratorForward(m_quadratic_terms); }

iterator begin() { return m_products.begin(); }
iterator end() { return m_products.end(); }
const_iterator begin() const { return m_products.begin(); }
const_iterator end() const { return m_products.end(); }
const_iterator cbegin() const { return m_products.begin(); }
const_iterator cend() const { return m_products.end(); }
auto quadratic() const { return ConstIteratorForward(m_quadratic_terms); }

/**
* Multiplies the Constant by t_factor (i.e., every Param-double pairs and the numerical term are multiplied).
Expand All @@ -138,7 +159,7 @@ class idol::Constant {
* Resulting zero entries are removed.
* @param t_term The parameter to add.
*/
Constant& operator+=(Param t_term);
Constant& operator+=(const Param& t_term);

/**
* Adds up another constant term (i.e., both the numerical terms and every Param-double pairs are added up).
Expand Down Expand Up @@ -172,8 +193,8 @@ namespace idol {

const double constant = t_coefficient.numerical();

auto it = t_coefficient.begin();
const auto end = t_coefficient.end();
auto it = t_coefficient.linear().begin();
const auto end = t_coefficient.linear().end();

if (it == end) {
return t_os << constant;
Expand Down
84 changes: 67 additions & 17 deletions lib/src/modeling/expressions/Constant.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -8,9 +8,9 @@

idol::Constant idol::Constant::Zero;

idol::Constant::Constant(const Param &t_param, double t_value) : m_products({ { t_param, t_value } }) {
idol::Constant::Constant(const Param &t_param, double t_value) : m_linear_terms({{t_param, t_value } }) {
if (equals(t_value, 0., ToleranceForSparsity)) {
m_products.clear();
m_linear_terms.clear();
}
}

Expand All @@ -21,31 +21,40 @@ idol::Constant::Constant(double t_constant) : m_constant(t_constant) {
void idol::Constant::set(const Param &t_param, double t_value) {

if (equals(t_value, 0., ToleranceForSparsity)) {
m_products.erase(t_param);
m_linear_terms.erase(t_param);
return;
}

auto [it, success] = m_products.emplace(t_param, t_value);
auto [it, success] = m_linear_terms.emplace(t_param, t_value);
if (!success) {
it->second = t_value;
}
}

double idol::Constant::get(const Param &t_param) const {
auto it = m_products.find(t_param);
return it == m_products.end() ? 0. : it->second;
auto it = m_linear_terms.find(t_param);
return it == m_linear_terms.end() ? 0. : it->second;
}

double idol::Constant::get(const idol::Param &t_param_1, const idol::Param &t_param_2) const {
auto it = m_quadratic_terms.find(std::make_pair(t_param_1, t_param_2));
return it == m_quadratic_terms.end() ? 0. : it->second;
}

idol::Constant &idol::Constant::operator*=(double t_factor) {

if (equals(t_factor, 0., ToleranceForSparsity)) {
m_constant = 0;
m_products.clear();
m_linear_terms.clear();
m_quadratic_terms.clear();
return *this;
}

m_constant *= t_factor;
for (auto& [param, value] : m_products) {
for (auto& [param, value] : m_linear_terms) {
value *= t_factor;
}
for (auto& [pair, value] : m_quadratic_terms) {
value *= t_factor;
}

Expand All @@ -57,16 +66,19 @@ idol::Constant &idol::Constant::operator+=(double t_term) {
return *this;
}

idol::Constant &idol::Constant::operator+=(Param t_term) {
idol::Constant &idol::Constant::operator+=(const Param& t_term) {
insert_or_add(t_term, 1.);
return *this;
}

idol::Constant &idol::Constant::operator+=(const Constant &t_term) {
m_constant += t_term.m_constant;
for (auto [param, value] : t_term) {
for (auto [param, value] : t_term.linear()) {
insert_or_add(param, value);
}
for (auto [pair, value] : t_term.quadratic()) {
insert_or_add(pair.first, pair.second, value);
}
return *this;
}

Expand All @@ -82,9 +94,12 @@ idol::Constant &idol::Constant::operator-=(Param t_term) {

idol::Constant &idol::Constant::operator-=(const Constant &t_term) {
m_constant -= t_term.m_constant;
for (auto [param, value] : t_term) {
for (auto [param, value] : t_term.linear()) {
insert_or_add(param, -value);
}
for (auto [pair, value] : t_term.quadratic()) {
insert_or_add(pair.first, pair.second, -value);
}
return *this;
}

Expand All @@ -93,35 +108,70 @@ void idol::Constant::insert_or_add(const Param &t_param, double t_value) {
return;
}

auto [it, success] = m_products.emplace(t_param, t_value);
auto [it, success] = m_linear_terms.emplace(t_param, t_value);
if (!success) {
it->second += t_value;
if (equals(it->second, 0., ToleranceForSparsity)) {
m_products.erase(it);
m_linear_terms.erase(it);
}
}
}

void idol::Constant::insert_or_add(const idol::Param &t_param_1, const idol::Param &t_param_2, double t_value) {
if (equals(t_value, 0., ToleranceForSparsity)) {
return;
}

auto [it, success] = m_quadratic_terms.emplace(std::make_pair(t_param_1, t_param_2), t_value);
if (!success) {
it->second += t_value;
if (equals(it->second, 0., ToleranceForSparsity)) {
m_quadratic_terms.erase(it);
}
}
}

bool idol::Constant::is_zero() const {
return m_products.empty() && equals(m_constant, 0., ToleranceForSparsity);
return m_linear_terms.empty() && m_quadratic_terms.empty() && equals(m_constant, 0., ToleranceForSparsity);
}

bool idol::Constant::is_numerical() const {
return m_products.empty();
return m_linear_terms.empty() && m_quadratic_terms.empty();
}

double idol::Constant::fix(const Solution::Primal &t_primals) const {
double result = m_constant;
for (const auto& [param, coeff] : m_products) {
for (const auto& [param, coeff] : m_linear_terms) {
result += coeff * t_primals.get(param.as<Var>());
}
for (const auto& [pair, coeff] : m_quadratic_terms) {
const auto [param1, param2] = pair;
result += coeff * t_primals.get(param1.as<Var>()) * t_primals.get(param2.as<Var>());
}
return result;
}

double idol::Constant::fix(const Solution::Dual &t_duals) const {
double result = m_constant;
for (const auto& [param, coeff] : m_products) {
for (const auto& [param, coeff] : m_linear_terms) {
result += coeff * t_duals.get(param.as<Ctr>());
}
for (const auto& [pair, coeff] : m_quadratic_terms) {
const auto [param1, param2] = pair;
result += coeff * t_duals.get(param1.as<Ctr>()) * t_duals.get(param2.as<Ctr>());
}
return result;
}

void idol::Constant::set(const idol::Param &t_param_1, const idol::Param &t_param_2, double t_value) {

if (equals(t_value, 0., ToleranceForSparsity)) {
m_quadratic_terms.erase(std::make_pair( t_param_1, t_param_2 ));
return;
}

auto [it, success] = m_quadratic_terms.emplace(std::make_pair(t_param_1, t_param_2), t_value);
if (!success) {
it->second = t_value;
}
}
4 changes: 2 additions & 2 deletions lib/src/optimizers/branch-and-bound/callbacks/CutSeparation.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -34,13 +34,13 @@ void idol::impl::CutSeparation::operator()(CallbackEvent t_event) {

::idol::Expr objective = row.rhs().numerical();

for (const auto& [param, coeff] : row.rhs()) {
for (const auto& [param, coeff] : row.rhs().linear()) {
objective += coeff * param.as<Var>();
}

for (const auto& [var, constant] : row.linear()) {
::idol::Expr term = -constant.numerical();
for (const auto& [param, coeff] : constant) {
for (const auto& [param, coeff] : constant.linear()) {
term += -coeff * param.as<Var>();
}
objective += term * current_solution.get(var);
Expand Down
8 changes: 4 additions & 4 deletions lib/src/optimizers/column-generation/Optimizers_ColumnGeneration.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -695,13 +695,13 @@ void idol::Optimizers::ColumnGeneration::Subproblem::update_objective(bool t_far

for (const auto &[ctr, constant] : m_generation_pattern.linear()) {
objective += constant.numerical() * -t_duals.get(ctr);
for (const auto &[param, coeff]: constant) {
for (const auto &[param, coeff] : constant.linear()) {
objective += -t_duals.get(ctr) * coeff * param.as<Var>();
}
}

if (!t_farkas_pricing) {
for (const auto &[param, coeff] : m_generation_pattern.obj()) {
for (const auto &[param, coeff] : m_generation_pattern.obj().linear()) {
objective += coeff * param.as<Var>();
}
}
Expand Down Expand Up @@ -757,12 +757,12 @@ double idol::Optimizers::ColumnGeneration::Subproblem::compute_reduced_cost(cons

for (const auto &[ctr, constant] : m_generation_pattern.linear()) {
result += constant.numerical() * -t_duals.get(ctr);
for (const auto &[param, coeff]: constant) {
for (const auto &[param, coeff] : constant.linear()) {
result += -t_duals.get(ctr) * coeff * primals.get(param.as<Var>());
}
}

for (const auto &[param, coeff] : m_generation_pattern.obj()) {
for (const auto &[param, coeff] : m_generation_pattern.obj().linear()) {
result += coeff * primals.get(param.as<Var>());
}

Expand Down

0 comments on commit 7b412d6

Please sign in to comment.