Keep only unique positions in module.

Other minor changes:
1. Remove bi-grades and unnecessary IntPoint class.
2. Remove ModulePresentation ctor from file.

4 files changed, 59 insertions, 74 deletions

diff --git a/matching/include/common_util.h b/matching/include/common_util.h
index b2ee22d..2d8dcb0 100644
--- a/matching/include/common_util.h
+++ b/matching/include/common_util.h
@@ -9,6 +9,7 @@
 #include <sstream>
 #include <string>
 #include <map>
+#include <functional>
 
 #include "common_defs.h"
 #include "phat/helpers/misc.h"
@@ -26,36 +27,6 @@ namespace md {
 
     using Column = std::vector<Index>;
 
-    struct IntPoint {
-        Index x;
-        Index y;
-
-        IntPoint(Index x, Index y) :x(x), y(y) { }
-
-        IntPoint() :x(0), y(0) { }
-
-        inline bool operator==(const IntPoint& v) const
-        {
-            return this->x == v.x && this->y == v.y;
-        }
-
-        // compare both coordinates, as needed in persistence
-        // do not overload operator<, requirements are not satisfied
-        inline bool is_less(const IntPoint& other, bool strict = false) const
-        {
-            if (x <= other.x and y <= other.y) {
-                if (strict) {
-                    return x != other.x or y != other.y;
-                }
-                else {
-                    return true;
-                }
-            }
-            return false;
-        }
-    };
-
-
     struct Point {
         Real x;
         Real y;
@@ -99,6 +70,7 @@ namespace md {
         }
     };
 
+
     using PointVec = std::vector<Point>;
 
     Point operator+(const Point& u, const Point& v);
@@ -201,11 +173,9 @@ namespace md {
 
     Rational midpoint(Rational a, Rational b);
 
+    // return true, if s is empty or starts with # (commented out line)
+    // whitespaces in the beginning of s are ignored
     bool ignore_line(const std::string& s);
-//    bool is_less(Rational a, Rational b);
-//    bool is_leq(Rational a, Rational b);
-//    bool is_greater(Rational a, Rational b);
-//    bool is_geq(Rational a, Rational b);
 
     // split string by delimeter
     template<typename Out>
@@ -222,8 +192,20 @@ namespace md {
         split_by_delim(s, delim, std::back_inserter(elems));
         return elems;
     }
+}
 
+namespace std {
+    template<>
+    struct hash<md::Point>
+    {
+        std::size_t operator()(const md::Point& p) const
+        {
+            auto hx = std::hash<decltype(p.x)>()(p.x);
+            auto hy = std::hash<decltype(p.y)>()(p.y);
+            return hx ^ (hy << 1);
+        }
+    };
+};
 
-}
 
 #endif //MATCHING_DISTANCE_COMMON_UTIL_H
diff --git a/matching/include/persistence_module.h b/matching/include/persistence_module.h
index ea9cfd0..a1fc67e 100644
--- a/matching/include/persistence_module.h
+++ b/matching/include/persistence_module.h
@@ -12,28 +12,43 @@
 
 namespace md {
 
+    /* ModulePresentation contains only information needed for matching
+     * distance computation over Z/2.
+     * Generators are represented as points (critical values),
+     * id i of generator g_i = its index in * vector generators_.
+     *
+     * Relations are represented by struct Relation, which has two members:
+     *   position_ is a point at which relation appears,
+     *   components_ contains indices of generators that sum to zero:
+     *     if components_ = [i, ..., j], then the relation is g_i +...+ g_j = 0.
+     *
+     *  ModulePresentation has member positions_ that contains all
+     *  distinct positions of generators and relations;
+     *  this member simplifies computing local linear bound.
+     */
+
 
     class ModulePresentation {
     public:
 
-        struct Relation {
+        enum Format { rivet_firep };
 
+        struct Relation {
             Point position_;
-            IntPoint bigrade_;
+            IndexVec components_;
 
             Relation() {}
-            Relation(const Point& _pos) : position_(_pos) { }
+            Relation(const Point& _pos, const IndexVec& _components);
 
             Real get_x() const { return position_.x; }
             Real get_y() const { return position_.y; }
-
-            IndexVec components_;
-
         };
 
-        ModulePresentation() { }
+        using RelVec = std::vector<Relation>;
+
+        ModulePresentation() {}
 
-        ModulePresentation(const std::string& fname); // read from file
+        ModulePresentation(const PointVec& _generators, const RelVec& _relations);
 
         Diagram weighted_slice_diagram(const DualPoint& line) const;
 
@@ -41,7 +56,7 @@ namespace md {
         void translate(Real a);
 
         // return minimal value of x- and y-coordinates
-        Real minimal_coordinate() const { return std::min(min_x(), min_x()); }
+        Real minimal_coordinate() const { return std::min(min_x(), min_y()); }
 
         // return box that contains all positions of all simplices
         Box bounding_box() const;
@@ -60,30 +75,22 @@ namespace md {
 
     private:
 
+        PointVec generators_;
+        std::vector<Relation> relations_;
+        PointVec positions_;
+
+
         Real max_x_ { std::numeric_limits<Real>::max() };
         Real max_y_ { std::numeric_limits<Real>::max() };
         Real min_x_ { -std::numeric_limits<Real>::max() };
         Real min_y_ { -std::numeric_limits<Real>::max() };
         Box bounding_box_;
 
-        PointVec generators_;
-        std::vector<Relation> relations_;
-
-        PointVec positions_;
-
         void init_boundaries();
         void project_generators(const DualPoint& slice, IndexVec& sorted_indices, RealVec& projections) const;
         void project_relations(const DualPoint& slice, IndexVec& sorted_indices, RealVec& projections) const;
     };
-
-//    template<typename R = double>
-//    R matching_distance(const PersistenceModule&, const PersistenceModule&, R)
-//    {
-//        return 0.0;
-//    };
-
 } // namespace md
 
 
-
 #endif //MATCHING_DISTANCE_PERSISTENCE_MODULE_H
diff --git a/matching/src/bifiltration.cpp b/matching/src/bifiltration.cpp
index 4131e6b..44b12cf 100644
--- a/matching/src/bifiltration.cpp
+++ b/matching/src/bifiltration.cpp
@@ -82,8 +82,6 @@ namespace md {
                 simplices_.emplace_back(index++, s, BifiltrationFormat::rivet);
             }
         }
-
-        
     }
 
     void Bifiltration::phat_like_format_reader(std::ifstream& ifstr)
diff --git a/matching/src/persistence_module.cpp b/matching/src/persistence_module.cpp
index f759a2b..efb20ef 100644
--- a/matching/src/persistence_module.cpp
+++ b/matching/src/persistence_module.cpp
@@ -1,5 +1,6 @@
 #include <numeric>
 #include <algorithm>
+#include <unordered_set>
 
 #include <phat/boundary_matrix.h>
 #include <phat/compute_persistence_pairs.h>
@@ -28,14 +29,13 @@ namespace md {
 
     // helper function to initialize const member positions_ in ModulePresentation
     PointVec
-    concat_gen_and_rel_positions(const PointVec& generators, const std::vector<ModulePresentation::Relation>& relations)
+    concat_gen_and_rel_positions(const PointVec& generators, const ModulePresentation::RelVec& relations)
     {
-        PointVec result(generators);
-        result.reserve(result.size() + relations.size());
+        std::unordered_set<Point> ps(generators.begin(), generators.end());
         for(const auto& rel : relations) {
-            result.push_back(rel.position_);
+            ps.insert(rel.position_);
         }
-        return result;
+        return PointVec(ps.begin(), ps.end());
     }
 
 
@@ -56,6 +56,14 @@ namespace md {
         bounding_box_ = Box(Point(min_x_, min_y_), Point(max_x_, max_y_));
     }
 
+
+    ModulePresentation::ModulePresentation(const PointVec& _generators, const RelVec& _relations) :
+        generators_(_generators),
+        relations_(_relations)
+    {
+        init_boundaries();
+    }
+
     void ModulePresentation::translate(md::Real a)
     {
         for(auto& g : generators_) {
@@ -71,15 +79,6 @@ namespace md {
     }
 
 
-    ModulePresentation::ModulePresentation(const std::string& fname) // read from file
-            :
-            generators_(),
-            relations_(),
-            positions_(concat_gen_and_rel_positions(generators_, relations_))
-    {
-        init_boundaries();
-    }
-
     /**
      *
      * @param slice line on which generators are projected
@@ -170,7 +169,6 @@ namespace md {
         return positions_;
     }
 
-
     Box ModulePresentation::bounding_box() const
     {
         return bounding_box_;