summaryrefslogtreecommitdiff
diff options
context:
space:
mode:
-rw-r--r--.gitmodules6
-rw-r--r--metastable.cpp172
m---------randutils0
m---------wolff0
4 files changed, 178 insertions, 0 deletions
diff --git a/.gitmodules b/.gitmodules
new file mode 100644
index 0000000..a509f25
--- /dev/null
+++ b/.gitmodules
@@ -0,0 +1,6 @@
+[submodule "randutils"]
+ path = randutils
+ url = https://gist.github.com/imneme/540829265469e673d045
+[submodule "wolff"]
+ path = wolff
+ url = git:research/wolff/code/c++
diff --git a/metastable.cpp b/metastable.cpp
new file mode 100644
index 0000000..cb47ffc
--- /dev/null
+++ b/metastable.cpp
@@ -0,0 +1,172 @@
+
+#define WOLFF_USE_FINITE_STATES
+
+#include <iostream>
+#include <cstring>
+#include <fstream>
+
+
+#include "randutils/randutils.hpp"
+#include "wolff/lib/wolff_models/ising.hpp"
+
+using namespace wolff;
+
+void update_distribution_file(std::string id, const std::vector<uint64_t>& data, std::string model_string) {
+ std::string filename = model_string + id + ".dat";
+ std::ifstream file(filename);
+
+ std::vector<uint64_t> data_old(data.size(), 0);
+
+ if (file.is_open()) {
+ for (unsigned i = 0; i < data.size(); i++) {
+ uint64_t num;
+ file >> num;
+ data_old[i] = num;
+ }
+
+ file.close();
+ }
+
+ std::ofstream file_out(filename);
+
+ for (unsigned i = 0; i < data.size(); i++) {
+ file_out <<std::fixed<< data_old[i] + data[i] << " ";
+ }
+
+ file_out.close();
+}
+
+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();
+}
+
+class meas : public measurement<ising_t, ising_t, graph<>> {
+private:
+ std::vector<uint64_t> mag_dist;
+ std::vector<std::vector<uint64_t>> energy_mag_dist;
+
+ signed E;
+ signed M;
+
+ std::string model_string;
+
+public:
+ meas(const system<ising_t, ising_t, graph<>>& S, double H) : mag_dist(S.nv + 1, 0), energy_mag_dist(S.nv + 1) {
+ M = S.nv * S.s[0];
+ E = S.ne;
+ for (std::vector<uint64_t>& d : energy_mag_dist) {
+ d.resize(S.ne + 1, 0);
+ }
+ model_string = "metastable_" + std::to_string(S.G.D) + "_" + std::to_string(S.G.L) + "_" + std::to_string(S.T) + "_" + std::to_string(H) + "_";
+ }
+
+ ~meas() {
+ update_distribution_file("m", mag_dist, model_string);
+ update_distribution_file("e", energy_mag_dist, model_string);
+ }
+
+ void plain_bond_visited(const system<ising_t, ising_t, graph<>>&, const typename graph<>::halfedge&, const ising_t&,
+ double dE) override {
+ if (dE > 0) {
+ E -= 2;
+ } else {
+ E += 2;
+ }
+ }
+
+ void ghost_bond_visited(const system<ising_t, ising_t, graph<>>&, const typename graph<>::vertex&,
+ const ising_t& s_old, const ising_t& s_new, double dE) override {
+ M += s_new - s_old;
+ }
+
+ void post_cluster(unsigned, unsigned, const system<ising_t, ising_t, graph<>>& S) override {
+ mag_dist[(S.nv + M) / 2]++;
+ energy_mag_dist[(S.nv + M) / 2][(E + S.ne) / 2]++;
+ }
+};
+
+int main(int argc, char* argv[]) {
+ unsigned N = (unsigned)1e4;
+ unsigned D = 2;
+ unsigned L = 128;
+ double T = 2.26918531421;
+ double H = 0.0;
+
+ int opt;
+
+ while ((opt = getopt(argc, argv, "N:D:L:T:H:")) != -1) {
+ switch (opt) {
+ case 'N': // number of steps
+ N = (unsigned)atof(optarg);
+ break;
+ case 'D': // dimension
+ D = atoi(optarg);
+ break;
+ case 'L': // linear size
+ L = atoi(optarg);
+ break;
+ case 'T': // temperature
+ T = atof(optarg);
+ break;
+ case 'H': // external field
+ H = atof(optarg);
+ break;
+ default:
+ exit(EXIT_FAILURE);
+ }
+ }
+
+ // define the spin-spin coupling
+ std::function<double(const ising_t&, const ising_t&)> Z =
+ [](const ising_t& s1, const ising_t& s2) -> double { return (double)(s1 * s2); };
+
+ // define the spin-field coupling
+ std::function<double(const ising_t&)> B = [=](const ising_t& s) -> double { return H * s; };
+
+ // initialize the lattice
+ graph<> G(D, L);
+
+ // initialize the system
+ system<ising_t, ising_t, graph<>> S(G, T, Z, B);
+
+ randutils::auto_seed_128 seeds;
+ std::mt19937 rng(seeds);
+
+ meas A(S, H);
+
+ // run wolff N times
+ S.run_wolff(N, gen_ising<graph<>>, A, rng);
+
+ return 0;
+}
+
diff --git a/randutils b/randutils
new file mode 160000
+Subproject 8486a610a954a8248c12485fb4cfc390a5f5f85
diff --git a/wolff b/wolff
new file mode 160000
+Subproject 4991418a416b75680737ec8f1e47b355b219728