diff options
Diffstat (limited to 'src/measurements.cpp')
-rw-r--r-- | src/measurements.cpp | 41 |
1 files changed, 21 insertions, 20 deletions
diff --git a/src/measurements.cpp b/src/measurements.cpp index 43483be..27f239d 100644 --- a/src/measurements.cpp +++ b/src/measurements.cpp @@ -2,7 +2,7 @@ #include "measurements.hpp" #include <iostream> -void update_distribution_file(std::string id, const std::vector<uint64_t>& data, unsigned int N, unsigned int Lx, unsigned int Ly, double beta) { +void update_distribution_file(std::string id, const std::vector<uint64_t>& data, unsigned int N, double Lx, double Ly, double beta) { std::string filename = "fracture_" + std::to_string(Lx) + "_" + std::to_string(Ly) + "_" + std::to_string(beta) + "_" + id + ".dat"; std::ifstream file(filename); @@ -31,7 +31,7 @@ void update_distribution_file(std::string id, const std::vector<uint64_t>& data, } template <class T> -void update_field_file(std::string id, const std::vector<T>& data, unsigned int N, unsigned int Lx, unsigned int Ly, double beta, unsigned int Mx, unsigned int My) { +void update_field_file(std::string id, const std::vector<T>& data, unsigned int N, double Lx, double Ly, double beta, unsigned int Mx, unsigned int My) { std::string filename = "fracture_" + std::to_string(Lx) + "_" + std::to_string(Ly) + "_" + std::to_string(beta) + "_" + id + "_" + std::to_string(Mx) + "_" + std::to_string(My) + ".dat"; std::ifstream file(filename); @@ -57,16 +57,16 @@ void update_field_file(std::string id, const std::vector<T>& data, unsigned int } template <class T> -std::vector<fftw_complex> data_transform(unsigned int Lx, unsigned int Ly, const std::vector<T>& data, fftw_plan forward_plan, double *fftw_forward_in, fftw_complex *fftw_forward_out) { - for (unsigned int i = 0; i < Lx * Ly; i++) { +std::vector<fftw_complex> data_transform(unsigned int Mx, unsigned int My, const std::vector<T>& data, fftw_plan forward_plan, double *fftw_forward_in, fftw_complex *fftw_forward_out) { + for (unsigned int i = 0; i < Mx * My; i++) { fftw_forward_in[i] = (double)data[i]; } fftw_execute(forward_plan); - std::vector<fftw_complex> output(Lx * (Ly / 2 + 1)); + std::vector<fftw_complex> output(Mx * (My / 2 + 1)); - for (unsigned int i = 0; i < Lx * (Ly / 2 + 1); i++) { + for (unsigned int i = 0; i < Mx * (My / 2 + 1); i++) { output[i][0] = fftw_forward_out[i][0]; output[i][1] = fftw_forward_out[i][1]; } @@ -75,43 +75,44 @@ std::vector<fftw_complex> data_transform(unsigned int Lx, unsigned int Ly, const } template <class T> -void correlation(unsigned int Lx, unsigned int Ly, std::vector<T>& data, const std::vector<fftw_complex>& tx1, const std::vector<fftw_complex>& tx2, fftw_plan reverse_plan, fftw_complex *fftw_reverse_in, double *fftw_reverse_out) { - for (unsigned int i = 0; i < Lx * (Ly / 2 + 1); i++) { +void correlation(unsigned int Mx, unsigned int My, std::vector<T>& data, const std::vector<fftw_complex>& tx1, const std::vector<fftw_complex>& tx2, fftw_plan reverse_plan, fftw_complex *fftw_reverse_in, double *fftw_reverse_out) { + for (unsigned int i = 0; i < Mx * (My / 2 + 1); i++) { fftw_reverse_in[i][0] = tx1[i][0] * tx2[i][0] + tx1[i][1] * tx2[i][1]; fftw_reverse_in[i][1] = tx1[i][0] * tx2[i][1] - tx1[i][1] * tx2[i][0]; } fftw_execute(reverse_plan); - for (unsigned int i = 0; i < (Lx / 2 + 1) * (Ly / 2 + 1); i++) { - data[i] += (T)(fftw_reverse_out[Lx * (i / (Lx / 2 + 1)) + i % (Lx / 2 + 1)] / (Lx * Ly)); + for (unsigned int i = 0; i < (Mx / 2 + 1) * (My / 2 + 1); i++) { + data[i] += (T)(fftw_reverse_out[Mx * (i / (Mx / 2 + 1)) + i % (Mx / 2 + 1)] / (Mx * My)); } } template <class T> -void autocorrelation(unsigned int Lx, unsigned int Ly, std::vector<T>& out_data, fftw_plan forward_plan, double *fftw_forward_in, fftw_complex *fftw_forward_out, fftw_plan reverse_plan, fftw_complex *fftw_reverse_in, double *fftw_reverse_out) { +void autocorrelation(unsigned int Mx, unsigned int My, std::vector<T>& out_data, fftw_plan forward_plan, double *fftw_forward_in, fftw_complex *fftw_forward_out, fftw_plan reverse_plan, fftw_complex *fftw_reverse_in, double *fftw_reverse_out) { fftw_execute(forward_plan); - for (unsigned int i = 0; i < Ly * (Lx / 2 + 1); i++) { + for (unsigned int i = 0; i < My * (Mx / 2 + 1); i++) { fftw_reverse_in[i][0] = pow(fftw_forward_out[i][0], 2) + pow(fftw_forward_out[i][1], 2); fftw_reverse_in[i][1] = 0.0; } fftw_execute(reverse_plan); - for (unsigned int i = 0; i < (Lx / 2 + 1) * (Ly / 2 + 1); i++) { - out_data[i] += (T)(fftw_reverse_out[Lx * (i / (Lx / 2 + 1)) + i % (Lx / 2 + 1)] / (Lx * Ly)); + for (unsigned int i = 0; i < (Mx / 2 + 1) * (My / 2 + 1); i++) { + out_data[i] += (T)(fftw_reverse_out[Mx * (i / (Mx / 2 + 1)) + i % (Mx / 2 + 1)] / (Mx * My)); } } -unsigned int edge_r_to_ind(graph::coordinate r, unsigned int Lx, unsigned int Ly, unsigned int Mx, unsigned int My) { +unsigned int edge_r_to_ind(graph::coordinate r, double Lx, double Ly, unsigned int Mx, unsigned int My) { return floor((Mx * r.x) / Lx) + Mx * floor((My * r.y) / Ly); } -ma::ma(unsigned int Lx, unsigned int Ly, unsigned int Mx, unsigned int My, double beta) : Lx(Lx), Ly(Ly), Mx(Mx), My(My), beta(beta), G(Lx * Ly), - sc(Lx * Ly, 0), - sa(Lx * Ly, 0), - sd(Lx * Ly, 0), +ma::ma(double Lx, double Ly, unsigned int Mx, unsigned int My, double beta) : + Lx(Lx), Ly(Ly), Mx(Mx), My(My), beta(beta), G(2 * (unsigned int)ceil(Lx * Ly / 2)), + sc(2 * (unsigned int)ceil(Lx * Ly / 2), 0), + sa(2 * (unsigned int)ceil(Lx * Ly / 2), 0), + sd(2 * (unsigned int)ceil(Lx * Ly / 2), 0), sb(log2(Mx < My ? Mx : My) + 1, 0), Ccc((Mx / 2 + 1) * (My / 2 + 1), 0), Css((Mx / 2 + 1) * (My / 2 + 1), 0), @@ -217,7 +218,7 @@ void ma::post_fracture(network &n) { // non-spanning clusters std::fill_n(fftw_forward_in, Mx * My, 0.0); for (unsigned int i = 0; i < num; i++) { - if (i != crack_component) { + if (i != crack_component && components[i].size() > 0) { for (auto it = components[i].begin(); it != components[i].end(); it++) { fftw_forward_in[edge_r_to_ind(n.G.dual_vertices[*it].r, Lx, Ly, Mx, My)] += 1.0; } |