use sparse matrix

dev/main.cpp

@@ -35,6 +35,7 @@
 #include "idol/optimizers/robust-optimization/column-and-constraint-generation/Optimizers_ColumnAndConstraintGeneration.h"
 #include "idol/optimizers/robust-optimization/column-and-constraint-generation/stabilizers/NoStabilization.h"
 #include "idol/optimizers/mixed-integer-optimization/wrappers/Mosek/Mosek.h"
+#include "idol/modeling/new-matrix/SparseMatrix.h"
 
 #include <iostream>
 #include <OsiCpxSolverInterface.hpp>
@@ -120,93 +121,16 @@ count_variables_and_constraints(const Model& t_model,
 
 int main(int t_argc, char** t_argv) {
 
-    Env env;
-
-    Model model(env);
-
-    const auto x = model.add_vars(Dim<1>(10), -1., 1., Continuous, "x");
-
-    model.add_ctr(x[0] == 0, "c1");
-    auto c = model.add_ctr(idol_Sum(i, Range(10), i * x[i]) <= 5., "c2");
-
-    model.use(Gurobi());
-
-    model.write("model_before.lp");
-
-    model.set_ctr_row(c, Row(idol_Sum(i, Range(10), (i % 2) * x[i]), 4));
+    SparseMatrix matrix(4);
 
-    model.write("model_after.lp");
+    matrix.add_major({ 0, 1 }, { 2, 3 });
+    matrix.add_major({ 0, 3 }, { 8, 7 });
+    matrix.add_major({ 1, 2 }, { .5, .25 });
+    matrix.add_major({ 2, 3 }, { 1, 1 });
 
-    /*
-
-    if (t_argc != 4) {
-        throw Exception("Expected arguments: <path_to_file> <stabilization=0|1> <time_limit>");
-    }
-
-    const std::string instance_file = t_argv[1];
-    const bool use_stabilization = std::stoi(t_argv[2]);
-    const double time_limit = std::stod(t_argv[3]);
-
-    Env env;
-
-    auto [model,
-          var_annotation,
-          ctr_annotation,
-          lower_level_objective] = Bilevel::read_from_file<Gurobi>(env, instance_file);
-
-    const auto [
-            n_lower_level_vars,
-            n_upper_level_vars,
-            n_lower_level_ctrs,
-            n_upper_level_ctrs
-            ] = count_variables_and_constraints(model, var_annotation, ctr_annotation);
-
-    std::unique_ptr<Bilevel::CCGStabilizer> stabilization;
-    if (use_stabilization) {
-        stabilization.reset(Bilevel::CCGStabilizers::TrustRegion().clone());
-    } else {
-        stabilization.reset(Bilevel::CCGStabilizers::NoStabilization().clone());
-    }
+    std::cout << matrix << std::endl;
 
-    model.use(
-                Bilevel::ColumnAndConstraintGeneration(var_annotation,
-                                                      ctr_annotation,
-                                                      lower_level_objective)
-                    .with_master_optimizer(Gurobi())
-                    .with_lower_level_optimizer(Gurobi())
-                    .with_stabilization(*stabilization)
-                    .with_time_limit(time_limit)
-                    .with_logs(true)
-            );
-
-    model.update();
-
-    model.optimize();
-
-    const auto& optimizer = model.optimizer().as<Optimizers::Bilevel::ColumnAndConstraintGeneration>();
-    const auto& ro_optimizer = optimizer.two_stage_robust_model().optimizer().as<Optimizers::Robust::ColumnAndConstraintGeneration>();
-
-    std::cout
-            << "result,"
-            << instance_file << ','
-            << model.vars().size() << ','
-            << model.ctrs().size() << ','
-            << n_lower_level_vars << ','
-            << n_upper_level_vars << ','
-            << n_lower_level_ctrs << ','
-            << n_upper_level_ctrs << ','
-            << ro_optimizer.uncertainty_set().vars().size() << ','
-            << ro_optimizer.uncertainty_set().ctrs().size() << ','
-            << use_stabilization << ','
-            << model.optimizer().time().count() << ','
-            << model.get_status() << ','
-            << model.get_reason() << ','
-            << model.get_best_bound() << ','
-            << model.get_best_obj() << ','
-            << optimizer.n_iterations() << ','
-            << std::endl;
-
-     */
+    std::cout << matrix.get_coefficient(2, 2) << std::endl;
 
     return 0;
 }

lib/CMakeLists.txt

@@ -248,6 +248,7 @@ add_library(idol STATIC
         include/idol/modeling/robust-optimization/StageDescription.h
         include/idol/optimizers/robust-optimization/convexification/Convexification.h
         src/optimizers/robust-optimization/convexification/Convexification.cpp
+        include/idol/modeling/new-matrix/SparseMatrix.h
 )
 
 find_package(OpenMP REQUIRED)

lib/include/idol/modeling/new-matrix/SparseMatrix.h

@@ -0,0 +1,195 @@
+//
+// Created by henri on 06.09.24.
+//
+
+#ifndef IDOL_SPARSEMATRIX_H
+#define IDOL_SPARSEMATRIX_H
+
+#include <vector>
+#include <ostream>
+
+namespace idol {
+    class SparseMatrix;
+    class SparseConstant;
+
+    std::ostream& operator<<(std::ostream& t_out, const SparseMatrix& t_matrix);
+}
+
+class idol::SparseConstant {
+    double m_value;
+public:
+    [[nodiscard]] double numerical() const { return m_value; }
+};
+
+class idol::SparseMatrix {
+    unsigned int m_minor_dimension = 0;
+    std::vector<double> m_coefficients;
+    std::vector<unsigned int> m_indices;
+    std::vector<unsigned int> m_start_indices;
+public:
+    SparseMatrix() = default;
+
+    SparseMatrix(unsigned int t_minor_dimension) : m_minor_dimension(t_minor_dimension) {}
+
+    [[nodiscard]] unsigned int minor_dimension() const { return m_minor_dimension; }
+
+    [[nodiscard]] unsigned int major_dimension() const { return m_start_indices.size(); }
+
+    void add_major(const std::vector<unsigned int>& t_indices, const std::vector<double>& t_values) {
+        m_start_indices.emplace_back(m_indices.size());
+        m_indices.insert(m_indices.end(), t_indices.begin(), t_indices.end());
+        m_coefficients.insert(m_coefficients.end(), t_values.begin(), t_values.end());
+    }
+
+    const std::vector<double> coefficients() const { return m_coefficients; }
+
+    const std::vector<unsigned int> indices() const { return m_indices; }
+
+    const std::vector<unsigned int> start_indices() const { return m_start_indices; }
+
+    class MajorView {
+        const SparseMatrix *m_matrix;
+        unsigned int m_major = 0;
+    public:
+        MajorView(const SparseMatrix *t_matrix, unsigned int t_major) : m_matrix(t_matrix), m_major(t_major) {}
+
+        class iterator {
+            const SparseMatrix* m_matrix;
+            unsigned int m_index = 0;
+        public:
+            iterator(const SparseMatrix* t_matrix, unsigned int t_index) : m_matrix(t_matrix), m_index(t_index) {}
+
+            iterator(const iterator&) = default;
+            iterator(iterator&&) noexcept = default;
+
+            iterator& operator=(const iterator&) = default;
+            iterator& operator=(iterator&&) noexcept = default;
+
+            bool operator!=(const iterator& t_other) const {
+                return m_index != t_other.m_index;
+            }
+
+            bool operator==(const iterator& t_other) const {
+                return m_index == t_other.m_index;
+            }
+
+            std::pair<unsigned int, double> operator*() const {
+                return {m_matrix->m_indices[m_index], m_matrix->m_coefficients[m_index] };
+            }
+
+            iterator& operator++() {
+                ++m_index;
+                return *this;
+            }
+
+            iterator& operator--() {
+                --m_index;
+                return *this;
+            }
+
+            iterator& operator+=(unsigned int t_offset) {
+                m_index += t_offset;
+                return *this;
+            }
+
+            iterator& operator-=(unsigned int t_offset) {
+                m_index -= t_offset;
+                return *this;
+            }
+
+            iterator operator+(unsigned int t_offset) const {
+                return {
+                    m_matrix,
+                    m_index + t_offset
+                };
+            }
+
+            iterator operator-(unsigned int t_offset) const {
+                return {
+                    m_matrix,
+                    m_index - t_offset
+                };
+            }
+
+            iterator operator-(const iterator& t_other) const {
+                return {
+                    m_matrix,
+                    m_index - t_other.m_index
+                };
+            }
+
+        };
+
+        iterator begin() const {
+            return {
+                m_matrix,
+                m_matrix->m_start_indices[m_major]
+            };
+        }
+
+        iterator end() const {
+
+            if (m_major + 1 == m_matrix->major_dimension()) {
+                return {
+                    m_matrix,
+                    (unsigned int) m_matrix->m_coefficients.size()
+                };
+            }
+
+            return {
+                m_matrix,
+                m_matrix->m_start_indices[m_major + 1]
+            };
+        }
+
+    };
+
+    auto get_major(unsigned int t_major) const {
+        return MajorView(this, t_major);
+    }
+
+    double get_coefficient(unsigned int t_major, unsigned int t_minor) const;
+};
+
+template<>
+struct std::iterator_traits<idol::SparseMatrix::MajorView::iterator> {
+    using iterator_category = std::forward_iterator_tag;
+    using value_type = std::pair<unsigned int, double>;
+    using difference_type = unsigned int;
+    using pointer = value_type*;
+    using reference = value_type&;
+};
+
+double idol::SparseMatrix::get_coefficient(unsigned int t_major, unsigned int t_minor) const {
+
+    const auto& major = get_major(t_major);
+    const auto it = std::lower_bound(
+            major.begin(),
+            major.end(),
+            std::pair<unsigned int, double>(t_minor, .0)
+    );
+
+    if (it == major.end() || (*it).first != t_minor) {
+        return 0.;
+    }
+
+    return (*it).second;
+}
+
+std::ostream& idol::operator<<(std::ostream& t_out, const SparseMatrix& t_matrix) {
+
+    const auto& indices = t_matrix.indices();
+    const auto& start_indices = t_matrix.start_indices();
+    const auto& coefficients = t_matrix.coefficients();
+
+    for (unsigned int major = 0, major_dim = t_matrix.major_dimension() ; major < major_dim ; ++major ) {
+        for (const auto& [index, coefficient] : t_matrix.get_major(major)) {
+            t_out << coefficient << " " << "x_" << index << " + ";
+        }
+        t_out << '\n';
+    }
+
+    return t_out;
+}
+
+#endif //IDOL_SPARSEMATRIX_H