halfway through getting end-finding working

1 files changed, 36 insertions, 11 deletions

diff --git a/src/measurements.cpp b/src/measurements.cpp
index 2e45168..16f1225 100644
--- a/src/measurements.cpp
+++ b/src/measurements.cpp
@@ -185,6 +185,7 @@ ma::ma(unsigned n, double a, double beta, double weight)
   Na = 0;
   NA = 0;
   Nq = 0;
+  NS = 0;
 }
 
 ma::ma(unsigned Lx, unsigned Ly, double beta, double weight)
@@ -324,7 +325,7 @@ void ma::post_fracture(network& n) {
     std::vector<bool> vertex_in(n.G.vertices.size());
 
     for (unsigned i = 0; i < n.G.edges.size(); i++) {
-      if (!n.backbone[i]) {
+      if (!last_backbone[i]) {
         vertex_in[n.G.edges[i].v[0]] = true;
         vertex_in[n.G.edges[i].v[1]] = true;
       }
@@ -376,11 +377,13 @@ void ma::post_fracture(network& n) {
     std::map<unsigned, std::list<unsigned>> fragments;
 
     for (unsigned e : crack.second) {
-      unsigned root = last_clusters.findroot(n.G.dual_edges[e].v[0]);
-      if (fragments.contains(root)) {
-        fragments[root].push_back(e);
-      } else {
-        fragments[root] = {e};
+      if (last_backbone[e]) {
+        unsigned root = last_clusters.findroot(n.G.dual_edges[e].v[0]);
+        if (fragments.contains(root)) {
+          fragments[root].push_back(e);
+        } else {
+          fragments[root] = {e};
+        }
       }
     }
 
@@ -389,6 +392,23 @@ void ma::post_fracture(network& n) {
 
       sp[biggest_fragment.size()]++;
 
+      std::map<unsigned, std::list<unsigned>> end_finder;
+
+      for (unsigned e : biggest_fragment) {
+        unsigned v0 = n.G.dual_edges[e].v[0];
+        unsigned v1 = n.G.dual_edges[e].v[1];
+        if (end_finder.contains(v0)) {
+          end_finder[v0].push_back(v0);
+        } else {
+          end_finder[v0] = {v0};
+        }
+        if (end_finder.contains(v1)) {
+          end_finder[v1].push_back(v1);
+        } else {
+          end_finder[v1] = {v1};
+        }
+      }
+
       unsigned left_end, right_end;
 
       bool comp;
@@ -460,7 +480,12 @@ void ma::post_fracture(network& n) {
         right_end = biggest_fragment.back();
       }
 
-      std::cout << n.G.dual_edges[left_end].r.x << " " << n.G.dual_edges[left_end].r.y  << "    " << n.G.dual_edges[right_end].r.x << "   " << n.G.dual_edges[right_end].r.y << "   "<< crack.first[0] % 2 << "     " << Δr.x << " " << Δr.y << "    " << biggest_fragment.size() << "\n";
+      std::cout << crack.first[0] << " " << crack.first[1] << "\n";
+      std::cout << left_end << " " << right_end << "\n";
+      for (unsigned e : biggest_fragment) {
+        std::cout << e << " ";
+      }
+        std::cout << "\n";
 
       for (unsigned i = 0; i < n.G.dual_edges.size(); i++) {
         if (!last_backbone[i]) {
@@ -471,25 +496,25 @@ void ma::post_fracture(network& n) {
             Δxl = Δx_tmpl < 0 ? Δx_tmpl + n.G.L.x : Δx_tmpl;
 
             Δy_tmpl = fabs(n.G.dual_edges[left_end].r.y - n.G.dual_edges[i].r.y);
-            Δyl = Δy_tmpl > n.G.L.y / 2 ? Δy_tmpl - n.G.L.y / 2 : Δy_tmpl;
+            Δyl = Δy_tmpl > n.G.L.y / 2 ? n.G.L.y - Δy_tmpl : Δy_tmpl;
 
             Δx_tmpr = n.G.dual_edges[i].r.x - n.G.dual_edges[left_end].r.x;
             Δxr = Δx_tmpr < 0 ? Δx_tmpr + n.G.L.x : Δx_tmpr;
 
             Δy_tmpr = fabs(n.G.dual_edges[left_end].r.y - n.G.dual_edges[i].r.y);
-            Δyr = Δy_tmpr > n.G.L.y / 2 ? Δy_tmpr - n.G.L.y / 2 : Δy_tmpr;
+            Δyr = Δy_tmpr > n.G.L.y / 2 ? n.G.L.y - Δy_tmpr : Δy_tmpr;
           } else {
             Δx_tmpl = n.G.dual_edges[left_end].r.y - n.G.dual_edges[i].r.y;
             Δxl = Δx_tmpl < 0 ? Δx_tmpl + n.G.L.y : Δx_tmpl;
 
             Δy_tmpl = fabs(n.G.dual_edges[left_end].r.x - n.G.dual_edges[i].r.x);
-            Δyl = Δy_tmpl > n.G.L.x / 2 ? Δy_tmpl - n.G.L.x / 2 : Δy_tmpl;
+            Δyl = Δy_tmpl > n.G.L.x / 2 ? n.G.L.x - Δy_tmpl : Δy_tmpl;
 
             Δx_tmpr = n.G.dual_edges[i].r.y - n.G.dual_edges[left_end].r.y;
             Δxr = Δx_tmpr < 0 ? Δx_tmpr + n.G.L.y : Δx_tmpr;
 
             Δy_tmpr = fabs(n.G.dual_edges[left_end].r.x - n.G.dual_edges[i].r.x);
-            Δyr = Δy_tmpr > n.G.L.x / 2 ? Δy_tmpr - n.G.L.x / 2 : Δy_tmpr;
+            Δyr = Δy_tmpr > n.G.L.x / 2 ? n.G.L.x - Δy_tmpr : Δy_tmpr;
           }
 
           if (Δxl < Δxr) {