diff options
Diffstat (limited to 'src/measurements.cpp')
-rw-r--r-- | src/measurements.cpp | 70 |
1 files changed, 35 insertions, 35 deletions
diff --git a/src/measurements.cpp b/src/measurements.cpp index ffff06d..17fbdee 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, double Lx, double Ly, double beta) { +void update_distribution_file(std::string id, const std::vector<uint64_t>& data, unsigned 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); @@ -11,7 +11,7 @@ void update_distribution_file(std::string id, const std::vector<uint64_t>& data, if (file.is_open()) { file >> N_old; - for (unsigned int i = 0; i < data.size(); i++) { + for (unsigned i = 0; i < data.size(); i++) { uint64_t num; file >> num; data_old[i] = num; @@ -23,7 +23,7 @@ void update_distribution_file(std::string id, const std::vector<uint64_t>& data, std::ofstream file_out(filename); file_out <<std::fixed<< N_old + N << "\n"; - for (unsigned int i = 0; i < data.size(); i++) { + for (unsigned i = 0; i < data.size(); i++) { file_out <<std::fixed<< data_old[i] + data[i] << " "; } @@ -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, double Lx, double Ly, double beta, unsigned int Mx, unsigned int My) { +void update_field_file(std::string id, const std::vector<T>& data, unsigned N, double Lx, double Ly, double beta, unsigned Mx, unsigned 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); @@ -44,7 +44,7 @@ void update_field_file(std::string id, const std::vector<T>& data, unsigned int if (file.is_open()) { file >> N_old; for (unsigned j = 0; j < data.size(); j++) { - for (unsigned int i = 0; i < data[j].size(); i++) { + for (unsigned i = 0; i < data[j].size(); i++) { file >> data_old[j][i]; } } @@ -55,7 +55,7 @@ void update_field_file(std::string id, const std::vector<T>& data, unsigned int file_out <<std::fixed<< N_old + N << "\n"; for (unsigned j = 0; j < data.size(); j++) { - for (unsigned int i = 0; i < data[j].size(); i++) { + for (unsigned i = 0; i < data[j].size(); i++) { file_out << data_old[j][i] + data[j][i] << " "; } file_out << "\n"; @@ -65,8 +65,8 @@ 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 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++) { +std::vector<fftw_complex> data_transform(unsigned Mx, unsigned My, const std::vector<T>& data, fftw_plan forward_plan, double *fftw_forward_in, fftw_complex *fftw_forward_out) { + for (unsigned i = 0; i < Mx * My; i++) { fftw_forward_in[i] = (double)data[i]; } @@ -74,7 +74,7 @@ std::vector<fftw_complex> data_transform(unsigned int Mx, unsigned int My, const std::vector<fftw_complex> output(Mx * (My / 2 + 1)); - for (unsigned int i = 0; i < Mx * (My / 2 + 1); i++) { + for (unsigned i = 0; i < Mx * (My / 2 + 1); i++) { output[i][0] = fftw_forward_out[i][0]; output[i][1] = fftw_forward_out[i][1]; } @@ -83,24 +83,24 @@ std::vector<fftw_complex> data_transform(unsigned int Mx, unsigned int My, const } template <class T> -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++) { +void correlation(unsigned Mx, unsigned 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 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 < (Mx / 2 + 1) * (My / 2 + 1); i++) { + for (unsigned 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 Mx, unsigned int My, std::vector<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 Mx, unsigned My, std::vector<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 < My * (Mx / 2 + 1); i++) { + for (unsigned 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; } @@ -108,23 +108,23 @@ void autocorrelation(unsigned int Mx, unsigned int My, std::vector<std::vector<T fftw_execute(reverse_plan); for (unsigned j = 0; j < out_data.size(); j++) { - for (unsigned int i = 0; i < (Mx / 2 + 1) * (My / 2 + 1); i++) { + for (unsigned i = 0; i < (Mx / 2 + 1) * (My / 2 + 1); i++) { out_data[j][i] += (T)pow(fftw_reverse_out[Mx * (i / (Mx / 2 + 1)) + i % (Mx / 2 + 1)] / (Mx * My), j + 1); } } } -unsigned int edge_r_to_ind(graph::coordinate r, double Lx, double Ly, unsigned int Mx, unsigned int My) { +unsigned edge_r_to_ind(graph::coordinate r, double Lx, double Ly, unsigned Mx, unsigned My) { return floor((Mx * r.x) / Lx) + Mx * floor((My * r.y) / Ly); } -ma::ma(double Lx, double Ly, unsigned int Mx, unsigned int My, double beta, unsigned Ncum) : - 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), - sC(2 * (unsigned int)ceil(Lx * Ly / 2), 0), - sA(2 * (unsigned int)ceil(Lx * Ly / 2), 0), - sd(3 * (unsigned int)ceil(Lx * Ly / 2), 0), +ma::ma(double Lx, double Ly, unsigned Mx, unsigned My, double beta, unsigned Ncum) : + Lx(Lx), Ly(Ly), Mx(Mx), My(My), beta(beta), G(2 * (unsigned)ceil(Lx * Ly / 2)), + sc(2 * (unsigned)ceil(Lx * Ly / 2), 0), + sa(2 * (unsigned)ceil(Lx * Ly / 2), 0), + sC(2 * (unsigned)ceil(Lx * Ly / 2), 0), + sA(2 * (unsigned)ceil(Lx * Ly / 2), 0), + sd(3 * (unsigned)ceil(Lx * Ly / 2), 0), Ccc(Ncum), Css(Ncum), Cmm(Ncum), @@ -195,7 +195,7 @@ void ma::pre_fracture(const network&) { boost::remove_edge_if(trivial, G); } -void ma::bond_broken(const network& net, const current_info& cur, unsigned int i) { +void ma::bond_broken(const network& net, const current_info& cur, unsigned i) { long double c = logl(cur.conductivity / fabs(cur.currents[i])) + net.thresholds[i]; if (c > lv && avalanches.back().size() > 0) { sa[avalanches.back().size() - 1]++; @@ -221,10 +221,10 @@ void ma::bond_broken(const network& net, const current_info& cur, unsigned int i } void ma::post_fracture(network &n) { - std::vector<unsigned int> component(boost::num_vertices(G)); - unsigned int num = boost::connected_components(G, &component[0]); + std::vector<unsigned> component(boost::num_vertices(G)); + unsigned num = boost::connected_components(G, &component[0]); - std::list<unsigned int> crack = find_minimal_crack(G, n); + std::list<unsigned> crack = find_minimal_crack(G, n); // crack surface correlations std::fill_n(fftw_forward_in, Mx * My, 0.0); @@ -236,17 +236,17 @@ void ma::post_fracture(network &n) { autocorrelation(Mx, My, Css, forward_plan, fftw_forward_in, fftw_forward_out, reverse_plan, fftw_reverse_in, fftw_reverse_out); - unsigned int crack_component = component[n.G.dual_edges[crack.front()].v[0]]; + unsigned crack_component = component[n.G.dual_edges[crack.front()].v[0]]; - std::vector<std::list<unsigned int>> components(num); + std::vector<std::list<unsigned>> components(num); - for (unsigned int i = 0; i < n.G.dual_vertices.size(); i++) { + for (unsigned i = 0; i < n.G.dual_vertices.size(); i++) { components[component[i]].push_back(i); } // non-spanning clusters std::fill_n(fftw_forward_in, Mx * My, 0.0); - for (unsigned int i = 0; i < num; i++) { + for (unsigned i = 0; i < num; i++) { 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; @@ -275,7 +275,7 @@ void ma::post_fracture(network &n) { autocorrelation(Mx, My, Cmm, forward_plan, fftw_forward_in, fftw_forward_out, reverse_plan, fftw_reverse_in, fftw_reverse_out); - std::function<bool(unsigned int)> inCrack = [&](unsigned int i) -> bool { + std::function<bool(unsigned)> inCrack = [&](unsigned i) -> bool { return component[n.G.dual_edges[i].v[0]] == crack_component; }; @@ -303,11 +303,11 @@ void ma::post_fracture(network &n) { autocorrelation(Mx, My, Cll, forward_plan, fftw_forward_in, fftw_forward_out, reverse_plan, fftw_reverse_in, fftw_reverse_out); // damage size distribution - unsigned int total_broken = 0; + unsigned total_broken = 0; std::fill_n(fftw_forward_in, Mx * My, 0.0); - for (unsigned int i = 0; i < n.G.edges.size(); i++) { + for (unsigned i = 0; i < n.G.edges.size(); i++) { if (n.fuses[i]) { total_broken++; fftw_forward_in[edge_r_to_ind(n.G.edges[i].r, Lx, Ly, Mx, My)] += 1.0; @@ -322,7 +322,7 @@ void ma::post_fracture(network &n) { unsigned ii = avalanches.back().front(); long double c = logl(cur.conductivity / fabs(cur.currents[ii])) + n.thresholds[ii]; - for (unsigned int i = 0; i < n.G.edges.size(); i++) { + for (unsigned i = 0; i < n.G.edges.size(); i++) { double x; if (!n.fuses[i]) { stress_file << std::scientific << fabs(cur.currents[i]) << " "; |