fixed problem with asymmetric correlation functions

1 files changed, 10 insertions, 16 deletions

diff --git a/src/measurements.cpp b/src/measurements.cpp
index 8613c27..24ffe77 100644
--- a/src/measurements.cpp
+++ b/src/measurements.cpp
@@ -96,7 +96,7 @@ tx1[i][1] * tx2[i][1]; fftw_reverse_in[i][1] = tx1[i][0] * tx2[i][1] - tx1[i][1]
 */
 
 void autocorrelation2(double Lx, double Ly, unsigned Mx, unsigned My, std::vector<uint64_t>& data,
-                      const std::list<graph::coordinate>& pos, std::array<unsigned, 2> count, bool debug = false) {
+                      const std::list<graph::coordinate>& pos, std::array<unsigned, 2> count) {
   for (std::list<graph::coordinate>::const_iterator it1 = pos.begin(); it1 != pos.end(); it1++) {
     for (std::list<graph::coordinate>::const_iterator it2 = it1; it2 != pos.end(); it2++) {
       double Δx_tmp = fabs(it1->x - it2->x);
@@ -105,16 +105,10 @@ void autocorrelation2(double Lx, double Ly, unsigned Mx, unsigned My, std::vecto
       double Δy_tmp = fabs(it1->y - it2->y);
       double Δy = Δy_tmp < Ly / 2 ? Δy_tmp : Ly - Δy_tmp;
 
-      if (debug) {
-        if (Δx > 0 || Δy > 0) {
-          std::cout << Δx << " " << Δy << "\n";
-        }
-      }
-
       if (count[0] % 2 == 0) {
-        data[(unsigned)((1 + 2 * (Mx / 2)) * (Δx / Lx)) + (Mx / 2 + 1) * (unsigned)((1 + 2 * (My / 2)) * (Δy / Ly))]++;
+        data[(unsigned)floor((1 + 2 * (Mx / 2)) * (Δx / Lx)) + (Mx / 2 + 1) * (unsigned)floor((1 + 2 * (My / 2)) * (Δy / Ly))]++;
       } else {
-        data[(unsigned)((1 + 2 * (Mx / 2)) * (Δy / Ly)) + (Mx / 2 + 1) * (unsigned)((1 + 2 * (My / 2)) * (Δx / Lx))]++;
+        data[(unsigned)floor((1 + 2 * (Mx / 2)) * (Δy / Ly)) + (Mx / 2 + 1) * (unsigned)floor((1 + 2 * (My / 2)) * (Δx / Lx))]++;
       }
     }
   }
@@ -165,10 +159,10 @@ unsigned edge_r_to_ind(graph::coordinate r, double Lx, double Ly, unsigned Mx, u
 
 ma::ma(unsigned n, double a, double beta, double weight)
     : sc(2 * n), sn(2 * n), ss(2 * n), sm(2 * n), sl(2 * n), sb(n + 1), sd(3 * n), sa(3 * n),
-      sA(3 * n), si(10000), sI(10000), cc(pow(1+(unsigned)sqrt(n), 2)), cn(pow(1+(unsigned)sqrt(n), 2)),
-      cs(pow(1+(unsigned)sqrt(n), 2)), cm(pow(1+(unsigned)sqrt(n), 2)), cl(pow(1+(unsigned)sqrt(n), 2)),
-      cb(pow(1+(unsigned)sqrt(n), 2)), ca(pow(1+(unsigned)sqrt(n), 2)), cA(pow(1+(unsigned)sqrt(n), 2)),
-      cp(pow(1+(unsigned)sqrt(n), 2)), cq(pow(1+(unsigned)sqrt(n), 2)), last_clusters(2 * n) {
+      sA(3 * n), si(10000), sI(10000), cc(pow(1+2*(unsigned)ceil(sqrt(n)), 2)), cn(pow(1+2*(unsigned)ceil(sqrt(n)), 2)),
+      cs(pow(1+2*(unsigned)ceil(sqrt(n)), 2)), cm(pow(1+2*(unsigned)ceil(sqrt(n)), 2)), cl(pow(1+2*(unsigned)ceil(sqrt(n)), 2)),
+      cb(pow(1+2*(unsigned)ceil(sqrt(n)), 2)), ca(pow(1+2*(unsigned)ceil(sqrt(n)), 2)), cA(pow(1+2*(unsigned)ceil(sqrt(n)), 2)),
+      cp(pow(1+2*(unsigned)ceil(sqrt(n)), 2)), cq(pow(1+2*(unsigned)ceil(sqrt(n)), 2)), last_clusters(2 * n) {
   if (beta != 0.0) {
     model_string = "fracture_" + std::to_string(n) + "_" + std::to_string(a) + "_" +
                    std::to_string(beta) + "_" + std::to_string(weight) + "_";
@@ -176,8 +170,8 @@ ma::ma(unsigned n, double a, double beta, double weight)
     model_string = "fracture_" + std::to_string(n) + "_" + std::to_string(a) + "_INF_" +
                    std::to_string(weight) + "_";
   }
-  Mx = 2 * sqrt(n) * a;
-  My = 2 * sqrt(n) / a;
+  Mx = 4 * ceil(sqrt(n)) * a;
+  My = 4 * ceil(sqrt(n)) / a;
 }
 
 ma::ma(unsigned Lx, unsigned Ly, double beta, double weight)
@@ -289,7 +283,7 @@ void ma::post_fracture(network& n) {
           sm[components[i].size() - 1]++;
         }
         autocorrelation2(n.G.L.x, n.G.L.y, Mx, My, cm, components[i], crack.first);
-        autocorrelation2(n.G.L.x, n.G.L.y, Mx, My, cp, components[i], {0, 1}, true);
+        autocorrelation2(n.G.L.x, n.G.L.y, Mx, My, cp, components[i], {0, 1});
       } else {
         ss[components[i].size() - 1]++;
         autocorrelation2(n.G.L.x, n.G.L.y, Mx, My, cs, components[i], crack.first);