now measuring clusters before last avalanche

2 files changed, 35 insertions, 0 deletions

diff --git a/src/measurements.cpp b/src/measurements.cpp
index 0b570cd..2af1ba3 100644
--- a/src/measurements.cpp
+++ b/src/measurements.cpp
@@ -143,18 +143,21 @@ ma::ma(unsigned n, double a, unsigned Mx, unsigned My, double beta) :
   sm(3 * n, 0),
   sa(3 * n, 0),
   sl(3 * n, 0),
+  sb(3 * n, 0),
   sD(3 * n, 0),
   ccc((Mx / 2) * (My / 2), 0),
   css((Mx / 2) * (My / 2), 0),
   cmm((Mx / 2) * (My / 2), 0),
   caa((Mx / 2) * (My / 2), 0),
   cll((Mx / 2) * (My / 2), 0),
+  cbb((Mx / 2) * (My / 2), 0),
   cDD((Mx / 2) * (My / 2), 0),
   csD((Mx / 2) * (My / 2), 0)
 {
   N = 0;
   Nc = 0;
   Na = 0;
+  Nb = 0;
 
   if (beta != 0.0) {
     model_string = "fracture_" + std::to_string(n) + "_" + std::to_string(a) + "_" + std::to_string(beta) + "_";
@@ -183,18 +186,21 @@ ma::ma(double Lx, double Ly, unsigned Mx, unsigned My, double beta) :
   sm(3 * (unsigned)ceil(Lx * Ly / 2), 0),
   sa(3 * (unsigned)ceil(Lx * Ly / 2), 0),
   sl(3 * (unsigned)ceil(Lx * Ly / 2), 0),
+  sb(3 * (unsigned)ceil(Lx * Ly / 2), 0),
   sD(3 * (unsigned)ceil(Lx * Ly / 2), 0),
   ccc((Mx / 2) * (My / 2), 0),
   css((Mx / 2) * (My / 2), 0),
   cmm((Mx / 2) * (My / 2), 0),
   caa((Mx / 2) * (My / 2), 0),
   cll((Mx / 2) * (My / 2), 0),
+  cbb((Mx / 2) * (My / 2), 0),
   cDD((Mx / 2) * (My / 2), 0),
   csD((Mx / 2) * (My / 2), 0)
 {
   N = 0;
   Nc = 0;
   Na = 0;
+  Nb = 0;
 
   if (beta != 0.0) {
     model_string = "fracture_" + std::to_string(Lx) + "_" + std::to_string(Ly) + "_" + std::to_string(beta) + "_";
@@ -233,6 +239,7 @@ ma::~ma() {
   update_distribution_file("sm", sm, N, model_string);
   update_distribution_file("sa", sa, Na, model_string);
   update_distribution_file("sl", sl, N, model_string);
+  update_distribution_file("sb", sb, Nb, model_string);
   update_distribution_file("sD", sD, N, model_string);
 
   update_field_file("ccc", ccc, Nc, model_string, Mx, My);
@@ -240,6 +247,7 @@ ma::~ma() {
   update_field_file("cmm", cmm, N, model_string, Mx, My);
   update_field_file("caa", caa, Na, model_string, Mx, My);
   update_field_file("cll", cll, N, model_string, Mx, My);
+  update_field_file("cbb", cbb, Nb, model_string, Mx, My);
   update_field_file("cDD", cDD, N, model_string, Mx, My);
   update_field_file("csD", csD, N, model_string, Mx, My);
 
@@ -262,8 +270,10 @@ void ma::bond_broken(const network& net, const current_info& cur, unsigned i) {
 
     lv = c;
     avalanches.push_back({net.G.edges[i].r});
+    last_avalanche = {i};
   } else {
     avalanches.back().push_back(net.G.edges[i].r);
+    last_avalanche.push_back(i);
   }
 
   boost::add_edge(net.G.dual_edges[i].v[0], net.G.dual_edges[i].v[1], {i}, G);
@@ -328,6 +338,27 @@ void ma::post_fracture(network &n) {
 
   autocorrelation2(n.G.L.x, n.G.L.y, Mx, My, cll, avalanches.back());
 
+  for (unsigned e : last_avalanche) {
+    boost::remove_edge(n.G.dual_edges[e].v[0], n.G.dual_edges[e].v[1], G);
+  }
+
+  num = boost::connected_components(G, &component[0]);
+
+  std::vector<std::list<graph::coordinate>> new_components(num);
+
+  for (unsigned i = 0; i < n.G.dual_vertices.size(); i++) {
+    new_components[component[i]].push_back(n.G.dual_vertices[i].r);
+  }
+
+  for (unsigned i = 0; i < num; i++) {
+    if (new_components[i].size() > 0) {
+      sb[new_components[i].size() - 1]++;
+      Nb++;
+
+      autocorrelation2(n.G.L.x, n.G.L.y, Mx, My, cbb, new_components[i]);
+    }
+  }
+
   N++;
 }
 
diff --git a/src/measurements.hpp b/src/measurements.hpp
index 600c6ac..610c266 100644
--- a/src/measurements.hpp
+++ b/src/measurements.hpp
@@ -39,21 +39,25 @@ class ma : public hooks {
     std::vector<uint64_t> sm; // spanning cluster size distribution
     std::vector<uint64_t> sa; // non-final avalanche size distribution
     std::vector<uint64_t> sl; // final avalanche size distribution
+    std::vector<uint64_t> sb; // final avalanche size distribution
     std::vector<uint64_t> sD; // post-fracture damage distribution
     std::vector<uint64_t> ccc; // cluster-cluster correlations
     std::vector<uint64_t> css; // surface-surface correlations
     std::vector<uint64_t> cmm; // surface-surface correlations
     std::vector<uint64_t> caa; // avalanche-avalanche correlations
     std::vector<uint64_t> cll; // damage-damage distribution
+    std::vector<uint64_t> cbb; // damage-damage distribution
     std::vector<uint64_t> cDD; // damage-damage distribution
     std::vector<uint64_t> csD; // damage-damage distribution
     uint64_t Nc;
     uint64_t Na;
+    uint64_t Nb;
 
   public:
     long double lv;
 
     std::list<std::list<graph::coordinate>> avalanches;
+    std::list<unsigned> last_avalanche;
     std::string model_string;
 
     ma(double Lx, double Ly, unsigned Mx, unsigned My, double beta);