diff options
| -rw-r--r-- | .gitmodules | 6 | ||||
| -rw-r--r-- | metastable.cpp | 172 | ||||
| m--------- | randutils | 0 | ||||
| m--------- | wolff | 0 |
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 |