summaryrefslogtreecommitdiff
diff options
context:
space:
mode:
-rw-r--r--.gitignore2
-rw-r--r--.gitmodules3
-rw-r--r--percolation.cpp126
m---------randutils0
4 files changed, 131 insertions, 0 deletions
diff --git a/.gitignore b/.gitignore
new file mode 100644
index 0000000..9f92fb3
--- /dev/null
+++ b/.gitignore
@@ -0,0 +1,2 @@
+*.dat
+*.processed
diff --git a/.gitmodules b/.gitmodules
new file mode 100644
index 0000000..94345fd
--- /dev/null
+++ b/.gitmodules
@@ -0,0 +1,3 @@
+[submodule "randutils"]
+ path = randutils
+ url = https://gist.github.com/imneme/540829265469e673d045
diff --git a/percolation.cpp b/percolation.cpp
new file mode 100644
index 0000000..f0b0610
--- /dev/null
+++ b/percolation.cpp
@@ -0,0 +1,126 @@
+
+#include <fstream>
+#include <random>
+#include <cinttypes>
+#include <vector>
+#include <string>
+#include <algorithm>
+
+#include "randutils/randutils.hpp"
+
+#include <boost/graph/adjacency_list.hpp>
+#include <boost/graph/incremental_components.hpp>
+#include <boost/pending/disjoint_sets.hpp>
+#include <boost/foreach.hpp>
+
+typedef boost::adjacency_list<boost::listS, boost::vecS, boost::undirectedS> Graph;
+typedef boost::graph_traits<Graph>::vertex_descriptor Vertex;
+typedef boost::graph_traits<Graph>::vertices_size_type VertexIndex;
+typedef VertexIndex* Rank;
+typedef Vertex* Parent;
+
+void update_distribution_file(std::string id, const std::vector<std::vector<uint64_t>>& data, std::string model_string) {
+ std::string filename = model_string + id + ".dat";
+ std::ifstream file(filename);
+
+ std::vector<std::vector<uint64_t>> data_old(data.size());
+
+ for (unsigned i = 0; i < data.size(); i++) {
+ data_old[i].resize(data[0].size());
+ }
+
+ if (file.is_open()) {
+ for (unsigned i = 0; i < data.size(); i++) {
+ for (unsigned j = 0; j < data[0].size(); j++) {
+ uint64_t num;
+ file >> num;
+ data_old[i][j] = num;
+ }
+ }
+
+ file.close();
+ }
+
+ std::ofstream file_out(filename);
+
+ for (unsigned i = 0; i < data.size(); i++) {
+ for (unsigned j = 0; j < data[0].size(); j++) {
+ file_out << std::fixed << data_old[i][j] + data[i][j] << " ";
+ }
+ file_out << "\n";
+ }
+
+ file_out.close();
+}
+
+int main(int argc, char* argv[]) {
+ unsigned L = 16;
+ unsigned N = 100;
+
+ int opt;
+
+ while ((opt = getopt(argc, argv, "L:N:")) != -1) {
+ switch (opt) {
+ case 'L':
+ L = (unsigned)atof(optarg);
+ break;
+ case 'N':
+ N = (unsigned)atof(optarg);
+ break;
+ default:
+ exit(1);
+ }
+ }
+
+ randutils::auto_seed_128 seeds;
+ std::mt19937 rng{seeds};
+
+ unsigned n_vertices = pow(L, 2);;
+ unsigned n_bonds = 2 * n_vertices;
+
+ std::vector<std::array<unsigned, 2>> bonds(n_bonds);
+
+ for (unsigned v_i = 0; v_i < n_vertices; v_i++) {
+ bonds[2 * v_i + 0] = {v_i, L * (v_i / L) + (v_i + 1) % L};
+ bonds[2 * v_i + 1] = {v_i, L * (((v_i / L) + 1) % L) + v_i % L};
+ }
+
+ std::vector<std::vector<uint64_t>> sn(n_bonds);
+
+ for (std::vector<uint64_t>& sni : sn) {
+ sni.resize(n_vertices);
+ }
+
+ for (unsigned trial = 0; trial < N; trial++) {
+ Graph G(n_vertices);
+ std::vector<VertexIndex> rank(n_vertices);
+ std::vector<Vertex> parent(n_vertices);
+
+ boost::disjoint_sets<VertexIndex*, Vertex*> ds(&rank[0], &parent[0]);
+
+ initialize_incremental_components(G, ds);
+ incremental_components(G, ds);
+
+ std::shuffle(bonds.begin(), bonds.end(), rng);
+
+ for (unsigned p = 0; p < n_bonds; p++) {
+ boost::add_edge(bonds[p][0], bonds[p][1], G);
+ ds.union_set(bonds[p][0], bonds[p][1]);
+
+ boost::component_index<VertexIndex> components(parent.begin(), parent.end());
+ BOOST_FOREACH(VertexIndex current_index, components) {
+ unsigned comp_size = 0;
+ BOOST_FOREACH(VertexIndex child_index, components[current_index]) {
+ comp_size++;
+ }
+ sn[p][comp_size - 1]++;
+ }
+ }
+ }
+
+ update_distribution_file("sn", sn, "percolation_" + std::to_string(L) + "_");
+
+ return 0;
+}
+
+
diff --git a/randutils b/randutils
new file mode 160000
+Subproject 8486a610a954a8248c12485fb4cfc390a5f5f85