diff --git a/dev/main.cpp b/dev/main.cpp index 2d70b0990..488e363d8 100644 --- a/dev/main.cpp +++ b/dev/main.cpp @@ -3,35 +3,39 @@ #include "idol/mixed-integer/modeling/variables/Var.h" #include "idol/mixed-integer/modeling/models/Model.h" #include "idol/mixed-integer/modeling/objects/Env.h" +#include "idol/mixed-integer/modeling/expressions/operations/operators.h" using namespace idol; int main(int t_argc, const char** t_argv) { Env env; - Model model(env); - const auto x = model.add_vars(Dim<1>(10), 1., 0., Continuous, "x"); + Model model(env, Model::Storage::RowOriented); + const auto x = model.add_vars(Dim<1>(10), 0., 1., Continuous, "x"); - SparseVector v1; - SparseVector v2; + Column col; + col.set_obj(1.); + model.add_var(0, 10, Continuous, col); - for (unsigned int i = 0 ; i < 10 ; ++i) { - v1.push_back(x[i], i); - //v1.push_back(x[i], 10 - i); + for (unsigned int i = 0; i < 9; ++i) { + model.add_ctr(x[i] + x[i+1] >= 0.5); } - for (unsigned int i = 0 ; i < 10 ; ++i) { - v2.push_back(x[i], i); - v2.push_back(x[i], 10 - i); - } - std::cout << v1.is_reduced() << std::endl; - std::cout << v2.is_reduced() << std::endl; - std::cout << v1.is_sorted_by_index() << std::endl; - std::cout << v2.is_sorted_by_index() << std::endl; - auto sum = v1; - sum += v2; + //std::cout << model << std::endl; + + const auto column = model.get_var_column(x[3]); + + model.dump(); + + model.set_storage(Model::Storage::ColumnOriented, true); + + const auto& row = model.get_ctr_row(model.get_ctr_by_index(2)); + + model.dump(); + + model.set_storage(Model::Storage::Both); - std::cout << sum << std::endl; + model.dump(); return 0; } diff --git a/lib/include/idol/general/numericals.h b/lib/include/idol/general/numericals.h index b21258513..670aba729 100644 --- a/lib/include/idol/general/numericals.h +++ b/lib/include/idol/general/numericals.h @@ -129,6 +129,10 @@ namespace idol { return std::abs(t_ub - t_lb); } + static bool is_zero(double t_value, double t_tolerance) { + return std::abs(t_value) <= t_tolerance; + } + static bool is_integer(double t_value, double t_tolerance) { return std::abs(t_value - std::round(t_value)) <= t_tolerance; } @@ -138,27 +142,6 @@ namespace idol { return std::round(t_value * multiplier) / multiplier; } - static double multiply_with_precision(double t_a, double t_b, unsigned int t_n_digits) { - const auto n_a = std::min(std::floor( std::log10(t_a) ) + 1, t_n_digits); - const auto n_b = std::min(std::floor( std::log10(t_b) ) + 1, t_n_digits); - const long scaled_a = long(long(t_a * std::pow(10., t_n_digits - n_a)) * std::pow(10, n_a)); - const long scaled_b = long(long(t_b * std::pow(10., t_n_digits - n_b)) * std::pow(10, n_b)); - const long exact_product = scaled_a * scaled_b; - const double result = exact_product / std::pow(10, 2 * t_n_digits); - return result; - } - - static double multiply_with_precision_by_power_of_10(double t_x, unsigned int t_exponent, unsigned int t_n_digits) { - - if (equals(t_x, 0., Tolerance::Sparsity)) { - return 0.; - } - - const auto n_x = std::min(std::floor( std::log10(std::abs(t_x)) ) + 1, t_n_digits); - const double result = long(t_x * std::pow(10., t_n_digits - n_x)) * std::pow(10, n_x + t_exponent - t_n_digits); - return result; - } - } #endif //OPTIMIZE_NUMERICALS_H diff --git a/lib/include/idol/general/utils/SparseVector.h b/lib/include/idol/general/utils/SparseVector.h index ebf5f0b01..a6a1167b5 100644 --- a/lib/include/idol/general/utils/SparseVector.h +++ b/lib/include/idol/general/utils/SparseVector.h @@ -10,6 +10,7 @@ #include #include "sort.h" #include "idol/general/utils/exceptions/Exception.h" +#include "idol/general/numericals.h" namespace idol { template @@ -19,8 +20,6 @@ namespace idol { template, idol::identity, idol::get_id>> class idol::SparseVector { public: - //static_assert(std::is_default_constructible_v); - enum class SortingCriteria { Index, Value, @@ -60,6 +59,8 @@ class idol::SparseVector { [[nodiscard]] const ValueT& value_at(unsigned int t_index) const { return m_values[t_index]; } + [[nodiscard]] bool has_index(const IndexT& t_index) const; + [[nodiscard]] ValueT get(const IndexT& t_index1) const; /** @@ -135,9 +136,32 @@ class idol::SparseVector { [[nodiscard]] const_iterator cend() const { return const_iterator(size(), *this); } + void sparsify(); + SparseVector& merge_without_conflict(const SparseVector& t_vec); }; +template +void idol::SparseVector::sparsify() { + + unsigned int i = 0, j = 0; + const unsigned int n = m_indices.size(); + while (i < n) { + if (!is_zero(m_values[i], Tolerance::Sparsity)) { + m_indices[j] = std::move(m_indices[i]); + m_values[j] = std::move(m_values[i]); + ++j; + } + ++i; + } + + for (unsigned int k = j; k < n; ++k) { + m_indices.pop_back(); + m_values.pop_back(); + } + +} + template idol::SparseVector & idol::SparseVector::merge_without_conflict(const SparseVector &t_vec) { @@ -179,6 +203,18 @@ idol::SparseVector::binary_operation_on_sorted_ return *this; } + if (empty()) { + reserve(t_vec2.size()); + for (unsigned int i = 0, n = t_vec2.size() ; i < n ; ++i) { + push_back(t_vec2.index_at(i), t_operation(ValueT{}, t_vec2.value_at(i))); + } + return *this; + } + + if (t_vec2.empty()) { + return *this; + } + SparseVector result; if (IndexExtractorT()(m_indices.back()) < IndexExtractorT()(t_vec2.m_indices.front())) { @@ -609,6 +645,41 @@ void idol::SparseVector::reduce() { } +template +bool idol::SparseVector::has_index(const IndexT &t_index) const { + + if (!m_is_reduced) { + for (unsigned int i = 0, n = m_indices.size() ; i < n ; ++i) { + if (IndexExtractorT()(m_indices[i]) == IndexExtractorT()(t_index)) { + return true; + } + } + return false; + } + + if (m_sorting_criteria != SortingCriteria::Index) { + for (unsigned int i = 0, n = m_indices.size() ; i < n ; ++i) { + if (IndexExtractorT()(m_indices[i]) == IndexExtractorT()(t_index)) { + return true; + } + } + return false; + } + + const auto it = std::lower_bound( + m_indices.begin(), + m_indices.end(), + t_index, + [](const IndexT& t_index1, const IndexT& t_index2) { return IndexExtractorT()(t_index1) < IndexExtractorT()(t_index2); } + ); + + if (it == m_indices.end() || IndexExtractorT()(*it) != IndexExtractorT()(t_index)) { + return false; + } + + return true; +} + template ValueT idol::SparseVector::get(const IndexT &t_index) const { diff --git a/lib/include/idol/mixed-integer/modeling/constraints/CtrVersion.h b/lib/include/idol/mixed-integer/modeling/constraints/CtrVersion.h index f63195ac7..9e724bf54 100644 --- a/lib/include/idol/mixed-integer/modeling/constraints/CtrVersion.h +++ b/lib/include/idol/mixed-integer/modeling/constraints/CtrVersion.h @@ -16,6 +16,10 @@ class idol::CtrVersion : public Version, public TempCtr { public: CtrVersion(unsigned int t_index, TempCtr&& t_temp_ctr) : Version(t_index), TempCtr(std::move(t_temp_ctr)) {} CtrVersion(unsigned int t_index, const TempCtr& t_temp_ctr) : Version(t_index), TempCtr(t_temp_ctr) {} + + using TempCtr::has_row; + using TempCtr::reset_row; + using TempCtr::set_row; }; #endif //IDOL_CTRVERSION_H diff --git a/lib/include/idol/mixed-integer/modeling/constraints/TempCtr.h b/lib/include/idol/mixed-integer/modeling/constraints/TempCtr.h index 018d46cc8..faae7cf00 100644 --- a/lib/include/idol/mixed-integer/modeling/constraints/TempCtr.h +++ b/lib/include/idol/mixed-integer/modeling/constraints/TempCtr.h @@ -37,8 +37,12 @@ namespace idol { * ``` */ class idol::TempCtr { - Row m_row; + std::unique_ptr m_row; CtrType m_type = LessOrEqual; +protected: + void set_row(Row&& t_row) { m_row = std::make_unique(std::move(t_row)); } + bool has_row() const { return (bool) m_row; } + void reset_row() { m_row.reset(); } public: /** * Default constructor. @@ -54,13 +58,13 @@ class idol::TempCtr { * @param t_row The desired row. * @param t_type The desired constraint type. */ - TempCtr(Row&& t_row, CtrType t_type) : m_row(std::move(t_row)), m_type(t_type) {} + TempCtr(Row&& t_row, CtrType t_type) : m_row(std::make_unique(std::move(t_row))), m_type(t_type) {} /** * Copy constructor. * @param t_src The object to copy. */ - TempCtr(const TempCtr& t_src) = default; + TempCtr(const TempCtr& t_src) : m_row(t_src.m_row ? std::make_unique(*t_src.m_row) : std::unique_ptr()), m_type(t_src.m_type) {} /** * Move constructor. @@ -84,13 +88,13 @@ class idol::TempCtr { * Returns the row of the temporary constraint (see Row). * @return The row of the temporary constraint. */ - [[nodiscard]] const Row& row() const { return m_row; } + [[nodiscard]] const Row& row() const { return *m_row; } /** * Returns the row of the temporary constraint (see Row). * @return The row of the temporary constraint. */ - Row& row() { return m_row; } + Row& row() { return *m_row; } /** * Returns the temporary constraint type. diff --git a/lib/include/idol/mixed-integer/modeling/expressions/Expr.h b/lib/include/idol/mixed-integer/modeling/expressions/Expr.h index 376ab43d6..e8685cf56 100644 --- a/lib/include/idol/mixed-integer/modeling/expressions/Expr.h +++ b/lib/include/idol/mixed-integer/modeling/expressions/Expr.h @@ -125,14 +125,16 @@ idol::impl::Expr::Expr(const QuadExpr &t_expr) : m_quadratic(t template idol::impl::Expr::Expr(const LinExpr &t_lin_expr, const QuadExpr &t_quad_expr, double t_constant) : m_linear(t_lin_expr), - m_quadratic(t_quad_expr) { + m_quadratic(t_quad_expr), + m_constant(t_constant) { } template idol::impl::Expr::Expr(const Expr &t_src) : m_linear(t_src.m_linear), - m_quadratic(t_src.m_quadratic) { + m_quadratic(t_src.m_quadratic), + m_constant(t_src.m_constant) { } diff --git a/lib/include/idol/mixed-integer/modeling/expressions/LinExpr.h b/lib/include/idol/mixed-integer/modeling/expressions/LinExpr.h index 6d7f3373e..c9f9cffee 100644 --- a/lib/include/idol/mixed-integer/modeling/expressions/LinExpr.h +++ b/lib/include/idol/mixed-integer/modeling/expressions/LinExpr.h @@ -40,4 +40,23 @@ template idol::LinExpr::LinExpr(double t_factor, const Key &t_key) : SparseVector({ t_key }, { t_factor }, SparseVector::SortingCriteria::Index, true) { } +namespace idol { + template + std::ostream& operator<<(std::ostream& t_os, const LinExpr& t_expr) { + + if (t_expr.empty()) { + return t_os << "0"; + } + + unsigned int i = 0; + t_os << t_expr.value_at(i) << " " << t_expr.index_at(i); + ++i; + for (unsigned int n = t_expr.size() ; i < n ; ++i) { + t_os << " + " << t_expr.value_at(i) << " " << t_expr.index_at(i); + } + + return t_os; + } +} + #endif //OPTIMIZE_EXPR_H diff --git a/lib/include/idol/mixed-integer/modeling/expressions/QuadExpr.h b/lib/include/idol/mixed-integer/modeling/expressions/QuadExpr.h index 5a8162cf9..82a3275f5 100644 --- a/lib/include/idol/mixed-integer/modeling/expressions/QuadExpr.h +++ b/lib/include/idol/mixed-integer/modeling/expressions/QuadExpr.h @@ -63,4 +63,23 @@ double idol::QuadExpr::get(const Key1 &t_a, const Key2 &t_b) const { return SparseVector, double>::get({ t_a, t_b }); } +namespace idol { + template + std::ostream& operator<<(std::ostream& t_os, const QuadExpr& t_expr) { + + if (t_expr.empty()) { + return t_os << "0"; + } + + unsigned int i = 0; + t_os << t_expr.value_at(i) << " " << t_expr.index_at(i) << " + "; + ++i; + for (unsigned int n = t_expr.size() ; i < n ; ++i) { + t_os << t_expr.value_at(i) << " " << t_expr.index_at(i) << " + "; + } + + return t_os; + } +} + #endif //IDOL_QUADEXPR_H diff --git a/lib/include/idol/mixed-integer/modeling/expressions/operations/operators_Var.h b/lib/include/idol/mixed-integer/modeling/expressions/operations/operators_Var.h index a9aa36cf0..bdd2edbdc 100644 --- a/lib/include/idol/mixed-integer/modeling/expressions/operations/operators_Var.h +++ b/lib/include/idol/mixed-integer/modeling/expressions/operations/operators_Var.h @@ -61,6 +61,10 @@ namespace idol { idol::LinExpr operator+(const idol::Var &t_a, const idol::Var &t_b); + idol::Expr operator+(double t_term, idol::LinExpr &&t_lin_expr); + + idol::Expr operator+(idol::LinExpr &&t_lin_expr, double t_term); + idol::LinExpr operator+(idol::LinExpr &&t_lin_expr, const idol::Var &t_var); idol::LinExpr operator+(const idol::Var &t_var, idol::LinExpr &&t_lin_expr); @@ -77,6 +81,10 @@ namespace idol { idol::LinExpr operator+(const idol::LinExpr &t_a, const idol::LinExpr &t_b); + idol::Expr operator+(double t_term, idol::QuadExpr &&t_quad_expr); + + idol::Expr operator+(idol::QuadExpr &&t_quad_expr, double t_term); + idol::QuadExpr operator+(idol::QuadExpr &&t_a, const idol::QuadExpr &t_b); idol::QuadExpr operator+(const idol::QuadExpr &t_a, idol::QuadExpr &&t_b); diff --git a/lib/include/idol/mixed-integer/modeling/models/Model.h b/lib/include/idol/mixed-integer/modeling/models/Model.h index 53e3ad33d..fce7c4279 100644 --- a/lib/include/idol/mixed-integer/modeling/models/Model.h +++ b/lib/include/idol/mixed-integer/modeling/models/Model.h @@ -40,16 +40,21 @@ namespace idol { * This class is used to represent a mathematical optimization model. */ class idol::Model { +public: + enum Storage { ColumnOriented, RowOriented, Both }; +private: Env& m_env; const unsigned int m_id; bool m_has_been_moved = false; ObjectiveSense m_sense = Minimize; - std::unique_ptr> m_objective; - std::unique_ptr> m_rhs; + Expr m_objective; + LinExpr m_rhs; std::vector m_variables; std::vector m_constraints; - bool m_is_row_oriented = true; + + Storage m_storage = RowOriented; + mutable bool m_has_minor_representation = false; std::unique_ptr m_optimizer; std::unique_ptr m_optimizer_factory; @@ -57,6 +62,15 @@ class idol::Model { void throw_if_no_optimizer() const { if (!m_optimizer) { throw Exception("No optimizer was found."); } } Model(const Model& t_src); + + template void throw_if_unknown_object(const LinExpr& t_expr); + template void throw_if_unknown_object(const QuadExpr& t_expr); + void add_column_to_rows(const Var& t_var); + void add_row_to_columns(const Ctr& t_ctr); + void build_row(const Ctr& t_ctr); + void build_column(const Var& t_var); + void build_rows(); + void build_columns(); public: /** * Creates a new model for a mathematical optimization problem. @@ -68,9 +82,7 @@ class idol::Model { * ``` * @param t_env the optimization environment which will store the model */ - explicit Model(Env& t_env); - - Model(Env& t_env, ObjectiveSense t_sense); + explicit Model(Env& t_env, Storage t_storage = RowOriented); Model(Model&&) noexcept; @@ -1185,6 +1197,13 @@ class idol::Model { */ void set_solution_index(unsigned int t_index); + void dump(std::ostream& t_os = std::cout) const; + + Storage storage() const { return m_storage; } + + void set_storage(Storage t_storage, bool t_reset_minor_representation = false); + + void reset_minor_representation(); }; template @@ -1336,6 +1355,11 @@ namespace idol { using namespace idol; + if (t_model.storage() == Model::Storage::ColumnOriented) { + std::cerr << "Warning: print a column-oriented model leads to the generation of a minor representation. " + "This operation can be slow and memory-consuming." << std::endl; + } + LimitedWidthStream stream(t_os, 120); if (t_model.get_obj_sense() == Minimize) { @@ -1356,11 +1380,11 @@ namespace idol { stream << '\t' << ctr << ": "; if (linear.empty()) { - stream << "ERROR CANNOT WRITE DOWN QUADRATIC TERMS"; // quadratic; + stream << quadratic; } else { stream << linear; if (!quadratic.empty()) { - stream << " + " << "ERROR CANNOT WRITE DOWN QUADRATIC TERMS"; //quadratic; + stream << " + " << quadratic; } } diff --git a/lib/include/idol/mixed-integer/modeling/objects/Env.h b/lib/include/idol/mixed-integer/modeling/objects/Env.h index bf6df4377..a26d3cf26 100644 --- a/lib/include/idol/mixed-integer/modeling/objects/Env.h +++ b/lib/include/idol/mixed-integer/modeling/objects/Env.h @@ -96,11 +96,11 @@ class idol::impl::Env { } ObjectId create_var(std::string t_name, TempVar&& t_temp_var) { - return create(std::move(t_name), "Var", m_variables, std::move(t_temp_var)); + return create(std::move(t_name), "x", m_variables, std::move(t_temp_var)); } ObjectId create_ctr(std::string t_name, TempCtr&& t_temp_ctr) { - return create(std::move(t_name), "Ctr", m_constraints, std::move(t_temp_ctr)); + return create(std::move(t_name), "c", m_constraints, std::move(t_temp_ctr)); } public: diff --git a/lib/include/idol/mixed-integer/modeling/variables/TempVar.h b/lib/include/idol/mixed-integer/modeling/variables/TempVar.h index 29fd774e5..74b690947 100644 --- a/lib/include/idol/mixed-integer/modeling/variables/TempVar.h +++ b/lib/include/idol/mixed-integer/modeling/variables/TempVar.h @@ -30,7 +30,11 @@ class idol::TempVar { double m_lb = 0.; double m_ub = Inf; VarType m_type = Continuous; - Column m_column; + std::unique_ptr m_column; +protected: + void set_column(Column&& t_column) { m_column = std::make_unique(std::move(t_column)); } + bool has_column() const { return (bool) m_column; } + void reset_column() { m_column.reset(); } public: /** * Default constructor. @@ -48,7 +52,7 @@ class idol::TempVar { * @param t_type The desired variable type. * @param t_column The desired column. */ - TempVar(double t_lb, double t_ub, VarType t_type, Column&& t_column) : m_lb(t_lb), m_ub(t_ub), m_type(t_type), m_column(std::move(t_column)) {} + TempVar(double t_lb, double t_ub, VarType t_type, Column&& t_column) : m_lb(t_lb), m_ub(t_ub), m_type(t_type), m_column(std::make_unique(std::move(t_column))) {} /** * Copy constructor. @@ -60,7 +64,7 @@ class idol::TempVar { * Move constructor. * @param t_src The object to move. */ - TempVar(const TempVar& t_src) = default; + TempVar(const TempVar& t_src) : m_lb(t_src.m_lb), m_ub(t_src.m_ub), m_type(t_src.m_type), m_column(t_src.m_column ? std::make_unique(*t_src.m_column) : std::unique_ptr()) {} /** * Copy-assignment operator. @@ -78,13 +82,13 @@ class idol::TempVar { * Returns the column of the temporary variable (see `Column`). * @return The column of the temporary variable. */ - [[nodiscard]] const Column& column() const { return m_column; } + [[nodiscard]] const Column& column() const { return *m_column; } /** * Returns the column of the temporary variable (see `Column`). * @return The column of the temporary variable. */ - Column& column() { return m_column; } + Column& column() { return *m_column; } /** * Returns the lower bound of the temporary variable. diff --git a/lib/include/idol/mixed-integer/modeling/variables/VarVersion.h b/lib/include/idol/mixed-integer/modeling/variables/VarVersion.h index 7bd079362..d34c08ed9 100644 --- a/lib/include/idol/mixed-integer/modeling/variables/VarVersion.h +++ b/lib/include/idol/mixed-integer/modeling/variables/VarVersion.h @@ -18,6 +18,10 @@ class idol::VarVersion : public Version, public TempVar { VarVersion(unsigned int t_index, double t_lb, double t_ub, VarType t_type, Column&& t_column) : Version(t_index), TempVar(t_lb, t_ub, t_type, std::move(t_column)) {} VarVersion(unsigned int t_index, TempVar&& t_temp_var) : Version(t_index), TempVar(std::move(t_temp_var)) {} VarVersion(unsigned int t_index, const TempVar& t_temp_var) : Version(t_index), TempVar(t_temp_var) {} + + using TempVar::has_column; + using TempVar::reset_column; + using TempVar::set_column; }; #endif //IDOL_VARVERSION_H diff --git a/lib/src/mixed-integer/modeling/constraints/TempCtr.cpp b/lib/src/mixed-integer/modeling/constraints/TempCtr.cpp index ffdb311ae..765ddfff5 100644 --- a/lib/src/mixed-integer/modeling/constraints/TempCtr.cpp +++ b/lib/src/mixed-integer/modeling/constraints/TempCtr.cpp @@ -32,9 +32,9 @@ TempCtr operator==(Expr&& t_lhs, const Expr& t_rhs) { return std::move TempCtr operator==(const Expr& t_lhs, const Expr& t_rhs) { return Expr(t_lhs) == Expr(t_rhs); } bool TempCtr::is_violated(const PrimalPoint &t_solution) const { - const double rhs = m_row.rhs(); + const double rhs = m_row->rhs(); double lhs = 0.; - for (const auto& [var, coeff] : m_row.linear()) { + for (const auto& [var, coeff] : m_row->linear()) { lhs += coeff * t_solution.get(var); } switch (m_type) { diff --git a/lib/src/mixed-integer/modeling/expressions/Constant.cpp b/lib/src/mixed-integer/modeling/expressions/Constant.cpp index 85b57f596..40a7eb90b 100644 --- a/lib/src/mixed-integer/modeling/expressions/Constant.cpp +++ b/lib/src/mixed-integer/modeling/expressions/Constant.cpp @@ -325,49 +325,6 @@ void idol::Constant::round() { } -idol::Constant &idol::Constant::multiply_with_precision(double t_factor, unsigned int t_n_digits) { - - m_constant = ::idol::multiply_with_precision(m_constant, t_factor, t_n_digits); - - if (m_linear_terms) { - for (auto &[param, value]: *m_linear_terms) { - value = ::idol::multiply_with_precision(value, t_factor, t_n_digits); - } - } - - if (m_quadratic_terms) { - - for (auto &[params, value]: *m_quadratic_terms) { - value = ::idol::multiply_with_precision(value, t_factor, t_n_digits); - } - - } - - return *this; -} - -idol::Constant & -idol::Constant::multiply_with_precision_by_power_of_10(unsigned int t_exponent, unsigned int t_n_digits) { - - m_constant = ::idol::multiply_with_precision_by_power_of_10(m_constant, t_exponent, t_n_digits); - - if (m_linear_terms) { - for (auto &[param, value]: *m_linear_terms) { - value = ::idol::multiply_with_precision_by_power_of_10(value, t_exponent, t_n_digits); - } - } - - if (m_quadratic_terms) { - - for (auto& [params, value] : *m_quadratic_terms) { - value = ::idol::multiply_with_precision_by_power_of_10(value, t_exponent, t_n_digits); - } - - } - - return *this; -} - idol::Constant &idol::Constant::operator=(const idol::Constant &t_rhs) { if (this == &t_rhs) { diff --git a/lib/src/mixed-integer/modeling/expressions/operations/operators_Var.cpp b/lib/src/mixed-integer/modeling/expressions/operations/operators_Var.cpp index 3c18fbbc0..b223390dd 100644 --- a/lib/src/mixed-integer/modeling/expressions/operations/operators_Var.cpp +++ b/lib/src/mixed-integer/modeling/expressions/operations/operators_Var.cpp @@ -103,6 +103,18 @@ LinExpr idol::operator+(const Var& t_a, const Var& t_b) { return result; } +idol::Expr idol::operator+(double t_term, idol::LinExpr &&t_lin_expr) { + Expr result(std::move(t_lin_expr)); + result.constant() += t_term; + return result; +} + +idol::Expr idol::operator+(idol::LinExpr &&t_lin_expr, double t_term) { + Expr result(std::move(t_lin_expr)); + result.constant() += t_term; + return result; +} + LinExpr idol::operator+(LinExpr&& t_lin_expr, const Var& t_var) { LinExpr result(std::move(t_lin_expr)); result += LinExpr(t_var); @@ -139,6 +151,18 @@ LinExpr idol::operator+(const LinExpr& t_a, const LinExpr& t_b) { return LinExpr(t_a) + t_b; } +idol::Expr idol::operator+(double t_term, idol::QuadExpr &&t_quad_expr) { + Expr result(std::move(t_quad_expr)); + result.constant() += t_term; + return result; +} + +idol::Expr idol::operator+(idol::QuadExpr &&t_quad_expr, double t_term) { + Expr result(std::move(t_quad_expr)); + result.constant() += t_term; + return result; +} + QuadExpr idol::operator+(QuadExpr&& t_a, const QuadExpr& t_b) { QuadExpr result(std::move(t_a)); result += t_b; diff --git a/lib/src/mixed-integer/modeling/models/Model.cpp b/lib/src/mixed-integer/modeling/models/Model.cpp index d834f3018..a5a8655d0 100644 --- a/lib/src/mixed-integer/modeling/models/Model.cpp +++ b/lib/src/mixed-integer/modeling/models/Model.cpp @@ -3,13 +3,9 @@ // #include "idol/mixed-integer/modeling/models/Model.h" #include "idol/mixed-integer/modeling/objects/Env.h" +#include "idol/mixed-integer/modeling/expressions/operations/operators.h" -idol::Model::Model(Env &t_env) : m_env(t_env), m_id(t_env.create_model_id()) {} - - -idol::Model::Model(idol::Env &t_env, idol::ObjectiveSense t_sense) : Model(t_env) { - set_obj_sense(t_sense); -} +idol::Model::Model(Env &t_env, Storage t_storage) : m_env(t_env), m_id(t_env.create_model_id()), m_storage(t_storage) {} idol::Model::Model(idol::Model && t_src) noexcept : m_env(t_src.m_env), @@ -47,6 +43,13 @@ idol::Model::~Model() { void idol::Model::add(const Var &t_var, TempVar t_temp_var) { + + // TODO sort and reduce by index in the model + + throw_if_unknown_object(t_temp_var.column().obj_quadratic()); + throw_if_unknown_object(t_temp_var.column().linear()); + throw_if_unknown_object(t_temp_var.column().quadratic()); + m_env.create_version(*this, t_var, m_variables.size(), @@ -55,14 +58,48 @@ void idol::Model::add(const Var &t_var, TempVar t_temp_var) { t_temp_var.type(), std::move(t_temp_var.column()) ); + m_variables.emplace_back(t_var); if (has_optimizer()) { optimizer().add(t_var); } - // add_column_to_rows(t_var); - throw Exception("Model. Work in progress..."); + const auto& column = get_var_column(t_var); + m_objective += column.obj() * t_var + t_var * column.obj_quadratic(); + + if (m_storage == ColumnOriented && !m_has_minor_representation) { + return; + } + + add_column_to_rows(t_var); + + if (m_storage == RowOriented) { + m_env.version(*this, t_var).reset_column(); + } + +} + +void idol::Model::add_column_to_rows(const idol::Var &t_var) { + + const auto& column = get_var_column(t_var); + + for (const auto& [ctr, constant] : column.linear()) { + auto& version = m_env.version(*this, ctr); + if (version.has_row()) { + version.row().linear().push_back(t_var, constant); + } + } + + for (const auto& [pair, constant] : column.quadratic()) { + + auto& version = m_env.version(*this, pair.first); + if (version.has_row()) { + version.row().quadratic().push_back({ pair.second, t_var }, constant); + } + + } + } void idol::Model::add(const Var &t_var) { @@ -106,6 +143,10 @@ void idol::Model::remove(const Ctr &t_ctr) { } void idol::Model::add(const Ctr &t_ctr, TempCtr t_temp_ctr) { + + throw_if_unknown_object(t_temp_ctr.row().linear()); + throw_if_unknown_object(t_temp_ctr.row().quadratic()); + m_env.create_version(*this, t_ctr, m_constraints.size(), std::move(t_temp_ctr)); m_constraints.emplace_back(t_ctr); @@ -113,8 +154,46 @@ void idol::Model::add(const Ctr &t_ctr, TempCtr t_temp_ctr) { optimizer().add(t_ctr); } - //add_row_to_columns(t_ctr); - throw Exception("Model. Work in progress..."); + const auto& row = get_ctr_row(t_ctr); + m_rhs += row.rhs() * t_ctr; + + if (m_storage == RowOriented && !m_has_minor_representation) { + return; + } + + add_row_to_columns(t_ctr); + + if (m_storage == ColumnOriented) { + m_env.version(*this, t_ctr).reset_row(); + } + +} + +void idol::Model::add_row_to_columns(const idol::Ctr &t_ctr) { + + const auto& row = get_ctr_row(t_ctr); + + for (const auto& [var, constant] : row.linear()) { + auto& version = m_env.version(*this, var); + if (version.has_column()) { + version.column().linear().push_back(t_ctr, constant); + } + } + + for (const auto& [pair, constant] : row.quadratic()) { + + auto& version1 = m_env.version(*this, pair.first); + if (version1.has_column()) { + version1.column().quadratic().push_back({ t_ctr, pair.second }, constant); + } + + auto& version2 = m_env.version(*this, pair.second); + if (version2.has_column()) { + version2.column().quadratic().push_back({ t_ctr, pair.first }, constant); + } + + } + } void idol::Model::add(const Ctr &t_ctr) { @@ -130,17 +209,11 @@ idol::ObjectiveSense idol::Model::get_obj_sense() const { } const idol::Expr &idol::Model::get_obj_expr() const { - if (!m_objective) { - throw Exception("Objective expression is not available."); - } - return *m_objective; + return m_objective; } const idol::LinExpr &idol::Model::get_rhs_expr() const { - if (!m_rhs) { - throw Exception("Right hand side expression is not available."); - } - return *m_rhs; + return m_rhs; } double idol::Model::get_mat_coeff(const Ctr &t_ctr, const Var &t_var) const { @@ -148,7 +221,14 @@ double idol::Model::get_mat_coeff(const Ctr &t_ctr, const Var &t_var) const { } const idol::Row &idol::Model::get_ctr_row(const Ctr &t_ctr) const { - return m_env.version(*this, t_ctr).row(); + auto& version = m_env.version(*this, t_ctr); + + if (!version.has_row()) { + const_cast(this)->build_row(t_ctr); + m_has_minor_representation = true; + } + + return version.row(); } idol::CtrType idol::Model::get_ctr_type(const Ctr &t_ctr) const { @@ -168,7 +248,14 @@ double idol::Model::get_var_ub(const Var &t_var) const { } const idol::Column &idol::Model::get_var_column(const Var &t_var) const { - return m_env.version(*this, t_var).column(); + auto& version = m_env.version(*this, t_var); + + if (!version.has_column()) { + const_cast(this)->build_column(t_var); + m_has_minor_representation = true; + } + + return version.column(); } unsigned int idol::Model::get_var_index(const Var &t_var) const { @@ -193,6 +280,9 @@ idol::Model* idol::Model::clone() const { idol::Model::Model(const Model& t_src) : Model(t_src.m_env) { + reserve_vars(t_src.vars().size()); + reserve_ctrs(t_src.ctrs().size()); + for (const auto& var : t_src.vars()) { add(var, TempVar( t_src.get_var_lb(var), @@ -212,6 +302,9 @@ idol::Model::Model(const Model& t_src) : Model(t_src.m_env) { set_obj_sense(t_src.get_obj_sense()); set_obj_expr(t_src.get_obj_expr()); + // TODO set storage to initial value + throw Exception("Model. Work in progress..."); + if (t_src.m_optimizer_factory) { use(*t_src.m_optimizer_factory); } @@ -323,7 +416,7 @@ void idol::Model::set_rhs_expr(LinExpr &&t_rhs) { void idol::Model::set_obj_const(double t_constant) { - m_objective->constant() = t_constant; + m_objective.constant() = t_constant; if (has_optimizer()) { optimizer().update_obj_constant(); @@ -526,3 +619,214 @@ void idol::Model::reserve_ctrs(unsigned int t_size) { double idol::Model::get_var_reduced_cost(const idol::Var &t_var) const { return optimizer().get_var_reduced_cost(t_var); } + +void idol::Model::build_row(const idol::Ctr &t_ctr) { + + Row row; + for (const auto& var : vars()) { + const auto& column = get_var_column(var); + const double value = column.linear().get(t_ctr); + row.linear().push_back(var, value); + if (!row.quadratic().empty()) { + throw Exception("Building a row with quadratic terms is not implemented."); + } + } + row.set_rhs(m_rhs.get(t_ctr)); + + row.linear().sparsify(); + row.quadratic().sparsify(); + + m_env.version(*this, t_ctr).set_row(std::move(row)); + +} + +void idol::Model::build_column(const idol::Var &t_var) { + + Column column; + for (const auto& ctr : ctrs()) { + const auto& row = get_ctr_row(ctr); + column.linear().push_back(ctr, row.linear().get(t_var)); + if (!column.quadratic().empty()) { + throw Exception("Building a column with quadratic terms is not implemented."); + } + } + column.set_obj(m_objective.linear().get(t_var)); + if (!column.obj_quadratic().empty()) { + throw Exception("Building a column with quadratic terms is not implemented."); + } + + column.linear().sparsify(); + column.quadratic().sparsify(); + + m_env.version(*this, t_var).set_column(std::move(column)); + +} + +template +void idol::Model::throw_if_unknown_object(const idol::QuadExpr &t_expr) { + for (const auto& [pair, constant] : t_expr) { + if (!has(pair.first)) { + throw Exception("Object " + pair.first.name() + " is not part of the model."); + } + if (!has(pair.second)) { + throw Exception("Object " + pair.second.name() + " is not part of the model."); + } + } +} + +template +void idol::Model::throw_if_unknown_object(const idol::LinExpr &t_expr) { + for (const auto& [obj, constant] : t_expr) { + if (!has(obj)) { + throw Exception("Object " + obj.name() + " is not part of the model."); + } + } +} + +void idol::Model::dump(std::ostream &t_os) const { + + t_os << "Objective Function:\n"; + t_os << "-------------------\n"; + + t_os << m_objective << '\n'; + + t_os << "RHS:\n"; + t_os << "----\n"; + + t_os << m_rhs << '\n'; + const int cell_width = 6; // Fixed cell width for consistent alignment + + // Column-Oriented Representation + t_os << "Column Oriented Representation:\n"; + t_os << "-------------------------------\n"; + + // Print variable names with fixed width for header alignment + t_os << std::setw(cell_width) << " "; + for (const auto& var : vars()) { + t_os << std::setw(cell_width) << var; + } + t_os << '\n'; + + // Print each row with fixed-width formatting + for (const auto& ctr : ctrs()) { + t_os << std::setw(cell_width) << ctr; + for (const auto& var : vars()) { + auto& version = m_env.version(*this, var); + if (version.has_column()) { + const auto& linear = version.column().linear(); + if (linear.has_index(ctr)) { + t_os << std::setw(cell_width) << linear.get(ctr); + } else { + t_os << std::setw(cell_width) << '_'; + } + } else { + t_os << std::setw(cell_width) << 'X'; + } + } + t_os << '\n'; + } + + // Row-Oriented Representation + t_os << "\nRow Oriented Representation:\n"; + t_os << "----------------------------\n"; + + t_os << std::setw(cell_width) << " "; + for (const auto& var : vars()) { + t_os << std::setw(cell_width) << var; + } + t_os << '\n'; + + for (const auto& ctr : ctrs()) { + t_os << std::setw(cell_width) << ctr; + for (const auto& var : vars()) { + auto& version = m_env.version(*this, ctr); + if (version.has_row()) { + const auto& linear = version.row().linear(); + if (linear.has_index(var)) { + t_os << std::setw(cell_width) << linear.get(var); + } else { + t_os << std::setw(cell_width) << '_'; + } + } else { + t_os << std::setw(cell_width) << 'X'; + } + } + t_os << '\n'; + } + +} + +void idol::Model::set_storage(idol::Model::Storage t_storage, bool t_reset_minor_representation) { + + if (t_storage == m_storage) { + if (t_reset_minor_representation) { + reset_minor_representation(); + } + return; + } + + if (m_storage == Both) { + m_storage = t_storage; + if (t_reset_minor_representation) { + reset_minor_representation(); + } + return; + } + + if (t_storage == Both) { + if (m_storage == ColumnOriented) { + build_rows(); + } else { + build_columns(); + } + } else if (t_storage == ColumnOriented) { + build_columns(); + } else if (t_storage == RowOriented) { + build_rows(); + } else { [[unlikely]] + throw Exception("Unsupported storage type."); + } + + m_storage = t_storage; + + if (t_reset_minor_representation) { + reset_minor_representation(); + } + +} + +void idol::Model::reset_minor_representation() { + + m_has_minor_representation = false; + + if (m_storage == ColumnOriented) { + for (const auto& ctr : ctrs()) { + m_env.version(*this, ctr).reset_row(); + } + return; + } + + if (m_storage == RowOriented) { + for (const auto& var : vars()) { + m_env.version(*this, var).reset_column(); + } + return; + } + +} + +void idol::Model::build_rows() { + for (const auto& ctr : ctrs()) { + if (!m_env.version(*this, ctr).has_row()) { + build_row(ctr); + } + } +} + +void idol::Model::build_columns() { + for (const auto& var : vars()) { + if (!m_env.version(*this, var).has_column()) { + build_column(var); + } + } +} diff --git a/tests/modeling/Expr.operators-type.cpp b/tests/modeling/Expr.operators-type.cpp index f17aa6ab9..86ffda26e 100644 --- a/tests/modeling/Expr.operators-type.cpp +++ b/tests/modeling/Expr.operators-type.cpp @@ -61,19 +61,19 @@ TEST_CASE("Expr: operators type deduction", "[unit][modeling-old][Expr]") { CHECK(std::is_same_v>); CHECK(std::is_same_v>); CHECK(std::is_same_v>); - CHECK(std::is_same_v>); - CHECK(std::is_same_v>); + //CHECK(std::is_same_v>); + //CHECK(std::is_same_v>); CHECK(std::is_same_v>); CHECK(std::is_same_v>); //CHECK(std::is_same_v>); //CHECK(std::is_same_v>); CHECK(std::is_same_v>); - CHECK(std::is_same_v>); + //CHECK(std::is_same_v>); //CHECK(std::is_same_v>); //CHECK(std::is_same_v>); - CHECK(std::is_same_v>); - CHECK(std::is_same_v>); + //CHECK(std::is_same_v>); + //CHECK(std::is_same_v>); } @@ -123,19 +123,19 @@ TEST_CASE("Expr: operators type deduction", "[unit][modeling-old][Expr]") { CHECK(std::is_same_v>); CHECK(std::is_same_v>); CHECK(std::is_same_v>); - CHECK(std::is_same_v>); - CHECK(std::is_same_v>); + //CHECK(std::is_same_v>); + //CHECK(std::is_same_v>); CHECK(std::is_same_v>); CHECK(std::is_same_v>); //CHECK(std::is_same_v>); //CHECK(std::is_same_v>); CHECK(std::is_same_v>); - CHECK(std::is_same_v>); + //CHECK(std::is_same_v>); //CHECK(std::is_same_v>); //CHECK(std::is_same_v>); - CHECK(std::is_same_v>); - CHECK(std::is_same_v>); + //CHECK(std::is_same_v>); + //CHECK(std::is_same_v>); } diff --git a/tests/modeling/Model.add-by-column.cpp b/tests/modeling/Model.add-by-column.cpp index 71e4ede45..29767487d 100644 --- a/tests/modeling/Model.add-by-column.cpp +++ b/tests/modeling/Model.add-by-column.cpp @@ -11,6 +11,8 @@ SCENARIO("Model: Add a variable by column", "[unit][modeling-old][Model]") { Env env; + auto storage = GENERATE(Model::Storage::RowOriented, Model::Storage::ColumnOriented, Model::Storage::Both); + Model model(env); GIVEN("An initial model with some constraints and no variable") { diff --git a/tests/modeling/QuadExpr.rotated-cone.cpp b/tests/modeling/QuadExpr.rotated-cone.cpp index 86f9f5ee0..ca1941108 100644 --- a/tests/modeling/QuadExpr.rotated-cone.cpp +++ b/tests/modeling/QuadExpr.rotated-cone.cpp @@ -22,8 +22,8 @@ double eval(const LinExpr& t_expr, const PrimalPoint& t_primal) { double eval(const QuadExpr& t_expr, const PrimalPoint& t_primal) { double result = 0; - for (const auto& [var1, var2, constant] : t_expr) { - result += constant * t_primal.get(var1) * t_primal.get(var2); + for (const auto& [vars, constant] : t_expr) { + result += constant * t_primal.get(vars.first) * t_primal.get(vars.second); } return result; }