started migration away from boost algorithms is the backbone trimming code

2 files changed, 69 insertions, 32 deletions

diff --git a/lib/include/clusters.hpp b/lib/include/clusters.hpp
new file mode 100644
index 0000000..4cc9073
--- /dev/null
+++ b/lib/include/clusters.hpp
@@ -0,0 +1,52 @@
+
+#include <vector>
+#include "graph.hpp"
+
+class ClusterTree {
+private:
+  std::vector<signed> p;
+  std::vector<unsigned> o;
+
+public:
+  std::vector<unsigned> c;
+
+  ClusterTree(unsigned n) : p(n, -1), o(n, 1), c(n, 0) { c[0] = n; }
+
+  unsigned findroot(unsigned i) {
+    if (p[i] < 0)
+      return i;
+    else
+      return p[i] = this->findroot(p[i]);
+  }
+
+  void join(unsigned i, unsigned j) {
+    c[o[i] - 1]--;
+    c[o[j] - 1]--;
+
+    if (o[i] < o[j]) {
+      p[i] = j;
+      o[j] += o[i];
+      c[o[j] - 1]++;
+    } else {
+      p[j] = i;
+      o[i] += o[j];
+      c[o[i] - 1]++;
+    }
+  }
+
+  void add_bond(const graph::edge& b) {
+    unsigned i = this->findroot(b.v[0]);
+    unsigned j = this->findroot(b.v[1]);
+
+    if (i != j) {
+      this->join(i, j);
+    }
+  }
+
+  bool same_component(unsigned v0, unsigned v1) {
+    unsigned i = this->findroot(v0);
+    unsigned j = this->findroot(v1);
+
+    return i == j;
+  }
+};
diff --git a/lib/include/network.hpp b/lib/include/network.hpp
index 8211cfc..bf9015c 100644
--- a/lib/include/network.hpp
+++ b/lib/include/network.hpp
@@ -14,53 +14,38 @@
 #include "hooks.hpp"
 #include "current_info.hpp"
 #include "array_hash.hpp"
+#include "clusters.hpp"
 
 #define CURRENT_CUTOFF 1e-10
 
-#include <boost/graph/adjacency_list.hpp>
-#include <boost/graph/connected_components.hpp>
-#include <boost/graph/depth_first_search.hpp>
-#include <boost/range/combine.hpp>
-#include <boost/foreach.hpp>
-#include <boost/graph/graph_utility.hpp>
-#include <boost/graph/incremental_components.hpp>
-#include <boost/pending/disjoint_sets.hpp>
-
-struct EdgeProperties {
-  unsigned index;
-};
+class network {
+  private:
+    problem px;
+    problem py;
+    std::unordered_map<std::array<unsigned, 2>, bool> seen_pairs;
 
-typedef boost::adjacency_list<boost::listS, boost::vecS, boost::undirectedS, boost::no_property, EdgeProperties> Graph;
-typedef boost::graph_traits<Graph>::vertex_descriptor Vertex;
-typedef boost::graph_traits<Graph>::vertices_size_type VertexIndex;
+    void update_backbone(const std::vector<double>& c);
+    void break_edge(unsigned, bool unbreak = false);
+    void get_cycle_edges_helper(std::set<unsigned>&, std::set<unsigned>&, unsigned, unsigned) const;
+    std::set<unsigned> get_cycle_edges(unsigned) const;
+    std::array<unsigned, 2> find_cycle(const std::set<unsigned>&, unsigned, unsigned) const;
 
-class network {
   public:
     const graph& G;
 
-    Graph bG;
-    Graph dG;
-    std::vector<VertexIndex> rank;
-    std::vector<Vertex> parent;
-    boost::disjoint_sets<VertexIndex*, Vertex*> ds;
-    std::unordered_map<std::array<unsigned, 2>, bool> seen_pairs;
+    ClusterTree C;
 
     std::vector<bool> fuses;
     std::vector<bool> backbone;
     std::vector<long double> thresholds;
 
-    problem px;
-    problem py;
-
     network(const graph&, cholmod_common*);
     network(const network&);
 
     void set_thresholds(double, std::mt19937&);
-    void fracture(hooks&, bool one_axis = false);
-    void update_backbone(const std::vector<double>& c);
+    void fracture(hooks&, bool one_axis = true);
 
-    void break_edge(unsigned, bool unbreak = false);
-    virtual current_info get_current_info() {current_info empty; return empty;};
+    virtual current_info get_current_info() const {current_info empty; return empty;};
 };
 
 class fuse_network : public network {
@@ -71,7 +56,7 @@ class fuse_network : public network {
     fuse_network(const graph&, cholmod_common*, double weight = 1.0);
     fuse_network(const fuse_network&);
 
-    current_info get_current_info() override;
+    current_info get_current_info() const override;
 };
 
 class elastic_network : public network {
@@ -82,7 +67,7 @@ class elastic_network : public network {
     elastic_network(const graph&, cholmod_common*, double weight = 0.5);
     elastic_network(const elastic_network&);
 
-    current_info get_current_info() override;
+    current_info get_current_info() const override;
 };
 
 class percolation_network : public network {
@@ -90,6 +75,6 @@ class percolation_network : public network {
     percolation_network(const graph&, cholmod_common*);
     percolation_network(const percolation_network&);
 
-    current_info get_current_info() override;
+    current_info get_current_info() const override;
 };