From 643baba78eb15a685d959aae718ee3eeade2f806 Mon Sep 17 00:00:00 2001 From: Jaron Kent-Dobias Date: Fri, 19 Oct 2018 01:26:14 -0400 Subject: big library overhual, fixed all examples, added documentation with sphinx --- lib/include/wolff.hpp | 21 +++++++----- lib/include/wolff/cluster.hpp | 57 +++++++++++++++++---------------- lib/include/wolff/finite_states.hpp | 6 +++- lib/include/wolff/graph.hpp | 21 ++++++------ lib/include/wolff/measurement.hpp | 18 +++++++---- lib/include/wolff/models/ising.hpp | 11 ++++++- lib/include/wolff/models/orthogonal.hpp | 15 +++++---- lib/include/wolff/models/potts.hpp | 2 +- lib/include/wolff/system.hpp | 30 ++++++++++------- lib/include/wolff/types.h | 2 +- 10 files changed, 109 insertions(+), 74 deletions(-) (limited to 'lib/include') diff --git a/lib/include/wolff.hpp b/lib/include/wolff.hpp index b78cd4b..7e909c8 100644 --- a/lib/include/wolff.hpp +++ b/lib/include/wolff.hpp @@ -3,25 +3,30 @@ #define WOLFF_H #include "wolff/cluster.hpp" +#include "wolff/measurement.hpp" + +namespace wolff{ template -void wolff(N_t N, wolff_system& S, - std::function r_gen, - wolff_measurement& A, std::mt19937& rng) { +void system::run_wolff(N_t N, + std::function &, v_t)> r_gen, + measurement& A, std::mt19937& rng) { - std::uniform_int_distribution dist(0, S.nv - 1); + std::uniform_int_distribution dist(0, nv - 1); for (N_t n = 0; n < N; n++) { v_t i0 = dist(rng); - R_t r = r_gen(rng, S.s[i0]); + R_t r = r_gen(rng, *this, i0); - A.pre_cluster(n, N, S, i0, r); + A.pre_cluster(n, N, *this, i0, r); - wolff_cluster_flip(S, i0, r, rng, A); + this->flip_cluster(i0, r, rng, A); - A.post_cluster(n, N, S); + A.post_cluster(n, N, *this); } } +} + #endif diff --git a/lib/include/wolff/cluster.hpp b/lib/include/wolff/cluster.hpp index 055cdf3..b66f367 100644 --- a/lib/include/wolff/cluster.hpp +++ b/lib/include/wolff/cluster.hpp @@ -7,20 +7,19 @@ #include #include -#include "types.h" #include "system.hpp" -#include "graph.hpp" -#include "measurement.hpp" + +namespace wolff { template -void wolff_cluster_flip(wolff_system& S, v_t i0, const R_t& r, - std::mt19937& rng, wolff_measurement& A) { +void system::flip_cluster(v_t i0, const R_t& r, + std::mt19937& rng, measurement& A) { std::uniform_real_distribution dist(0.0, 1.0); std::queue queue; queue.push(i0); - std::vector marks(S.G.nv, false); + std::vector marks(G.nv, false); while (!queue.empty()) { v_t i = queue.front(); @@ -33,21 +32,21 @@ void wolff_cluster_flip(wolff_system& S, v_t i0, const R_t& r, #ifndef WOLFF_NO_FIELD R_t s0_new; - bool we_are_ghost = (i == S.nv); + bool we_are_ghost = (i == nv); if (we_are_ghost) { - s0_new = r.act(S.s0); + s0_new = r.act(s0); } else #endif { - si_new = r.act(S.s[i]); + si_new = r.act(s[i]); } - for (const v_t &j : S.G.adj[i]) { + for (const v_t &j : G.adjacency[i]) { double dE, p; #ifndef WOLFF_NO_FIELD - bool neighbor_is_ghost = (j == S.nv); + bool neighbor_is_ghost = (j == nv); if (we_are_ghost || neighbor_is_ghost) { X_t s0s_old, s0s_new; @@ -55,18 +54,18 @@ void wolff_cluster_flip(wolff_system& S, v_t i0, const R_t& r, if (neighbor_is_ghost) { non_ghost = i; - s0s_old = S.s0.act_inverse(S.s[i]); - s0s_new = S.s0.act_inverse(si_new); + s0s_old = s0.act_inverse(s[i]); + s0s_new = s0.act_inverse(si_new); } else { non_ghost = j; - s0s_old = S.s0.act_inverse(S.s[j]); - s0s_new = s0_new.act_inverse(S.s[j]); + s0s_old = s0.act_inverse(s[j]); + s0s_new = s0_new.act_inverse(s[j]); } #ifdef WOLFF_SITE_DEPENDENCE - dE = S.B(non_ghost, s0s_old) - S.B(non_ghost, s0s_new); + dE = B(non_ghost, s0s_old) - B(non_ghost, s0s_new); #else - dE = S.B(s0s_old) - S.B(s0s_new); + dE = B(s0s_old) - B(s0s_new); #endif #ifdef WOLFF_FINITE_STATES @@ -75,28 +74,28 @@ void wolff_cluster_flip(wolff_system& S, v_t i0, const R_t& r, #endif // run measurement hooks for encountering a ghost bond - A.ghost_bond_visited(S, non_ghost, s0s_old, s0s_new, dE); + A.ghost_bond_visited(*this, non_ghost, s0s_old, s0s_new, dE); } else // this is a perfectly normal bond! #endif { #ifdef WOLFF_BOND_DEPENDENCE - dE = S.Z(i, S.s[i], j, S.s[j]) - S.Z(i, si_new, j, S.s[j]); + dE = Z(i, s[i], j, s[j]) - Z(i, si_new, j, s[j]); #else - dE = S.Z(S.s[i], S.s[j]) - S.Z(si_new, S.s[j]); + dE = Z(s[i], s[j]) - Z(si_new, s[j]); #endif #ifdef WOLFF_FINITE_STATES - p = finite_states_Zp[finite_states_enum(S.s[i])] + p = finite_states_Zp[finite_states_enum(s[i])] [finite_states_enum(si_new)] - [finite_states_enum(S.s[j])]; + [finite_states_enum(s[j])]; #endif // run measurement hooks for encountering a plain bond - A.plain_bond_visited(S, i, si_new, j, dE); + A.plain_bond_visited(*this, i, si_new, j, dE); } #ifndef FINITE_STATES - p = 1.0 - exp(-dE / S.T); + p = 1.0 - exp(-dE / T); #endif if (dist(rng) < p) { @@ -106,17 +105,19 @@ void wolff_cluster_flip(wolff_system& S, v_t i0, const R_t& r, #ifndef WOLFF_NO_FIELD if (we_are_ghost) { - A.ghost_site_transformed(S, s0_new); - S.s0 = s0_new; + A.ghost_site_transformed(*this, s0_new); + s0 = s0_new; } else #endif { - A.plain_site_transformed(S, i, si_new); - S.s[i] = si_new; + A.plain_site_transformed(*this, i, si_new); + s[i] = si_new; } } } } +} + #endif diff --git a/lib/include/wolff/finite_states.hpp b/lib/include/wolff/finite_states.hpp index 4cca69e..7e54eff 100644 --- a/lib/include/wolff/finite_states.hpp +++ b/lib/include/wolff/finite_states.hpp @@ -8,13 +8,15 @@ #include "system.hpp" +namespace wolff { + std::array, WOLFF_FINITE_STATES_N>, WOLFF_FINITE_STATES_N> finite_states_Zp; #ifndef WOLFF_NO_FIELD std::array, WOLFF_FINITE_STATES_N> finite_states_Bp; #endif template -void finite_states_init(const wolff_system& S) { +void finite_states_init(const system& S) { #ifndef WOLFF_NO_FIELD for (q_t i = 0; i < WOLFF_FINITE_STATES_N; i++) { for (q_t j = 0; j < WOLFF_FINITE_STATES_N; j++) { @@ -31,5 +33,7 @@ void finite_states_init(const wolff_system& S) { } } +} + #endif diff --git a/lib/include/wolff/graph.hpp b/lib/include/wolff/graph.hpp index fcd73fa..65b3941 100644 --- a/lib/include/wolff/graph.hpp +++ b/lib/include/wolff/graph.hpp @@ -5,23 +5,26 @@ #include #include -#include "types.h" +namespace wolff { -typedef enum lattice_t { - SQUARE_LATTICE, - DIAGONAL_LATTICE -} lattice_t; +#include "types.h" -class graph_t { +class graph { public: + D_t D; + L_t L; v_t ne; v_t nv; - std::vector> adj; + std::vector> adjacency; std::vector> coordinate; - graph_t(D_t D, L_t L, lattice_t lat = SQUARE_LATTICE); - void add_ext(); + graph(); + graph(D_t D, L_t L); + + void add_ghost(); }; +} + #endif diff --git a/lib/include/wolff/measurement.hpp b/lib/include/wolff/measurement.hpp index a6a5f79..604eaf5 100644 --- a/lib/include/wolff/measurement.hpp +++ b/lib/include/wolff/measurement.hpp @@ -4,21 +4,25 @@ #include "system.hpp" +namespace wolff { + template -class wolff_measurement { +class measurement { public: - virtual void pre_cluster(N_t, N_t, const wolff_system&, v_t, const R_t&) = 0; + virtual void pre_cluster(N_t, N_t, const system&, v_t, const R_t&) = 0; - virtual void plain_bond_visited(const wolff_system&, v_t, const X_t&, v_t, double) = 0; - virtual void plain_site_transformed(const wolff_system&, v_t, const X_t&) = 0; + virtual void plain_bond_visited(const system&, v_t, const X_t&, v_t, double) = 0; + virtual void plain_site_transformed(const system&, v_t, const X_t&) = 0; #ifndef WOLFF_NO_FIELD - virtual void ghost_bond_visited(const wolff_system&, v_t, const X_t&, const X_t&, double) = 0; - virtual void ghost_site_transformed(const wolff_system&, const R_t&) = 0; + virtual void ghost_bond_visited(const system&, v_t, const X_t&, const X_t&, double) = 0; + virtual void ghost_site_transformed(const system&, const R_t&) = 0; #endif - virtual void post_cluster(N_t, N_t, const wolff_system&) = 0; + virtual void post_cluster(N_t, N_t, const system&) = 0; }; +} + #endif diff --git a/lib/include/wolff/models/ising.hpp b/lib/include/wolff/models/ising.hpp index c6a3a47..824c063 100644 --- a/lib/include/wolff/models/ising.hpp +++ b/lib/include/wolff/models/ising.hpp @@ -3,6 +3,9 @@ #define WOLFF_MODELS_ISING #include "../types.h" +#include "../system.hpp" + +namespace wolff { class ising_t { public: @@ -73,9 +76,15 @@ inline ising_t::F_t operator*(double a, const ising_t& s) { return s * a; } +ising_t gen_ising(std::mt19937&, const system&, v_t) { + return ising_t(true); +}; + #define WOLFF_FINITE_STATES_N 2 const ising_t finite_states_possible[2] = {ising_t(0), ising_t(1)}; -q_t finite_states_enum(ising_t state) { return (q_t)state.x; } +q_t finite_states_enum(const ising_t& state) { return (q_t)state.x; } + +} #endif diff --git a/lib/include/wolff/models/orthogonal.hpp b/lib/include/wolff/models/orthogonal.hpp index 9dd5ddd..58203b7 100644 --- a/lib/include/wolff/models/orthogonal.hpp +++ b/lib/include/wolff/models/orthogonal.hpp @@ -4,7 +4,8 @@ #include #include -#include +#include "../types.h" +#include "../system.hpp" #include "vector.hpp" template @@ -113,7 +114,7 @@ class orthogonal_t : public std::array, q> { }; template -orthogonal_t generate_rotation_uniform (std::mt19937& r, const vector_t & v) { +orthogonal_t generate_rotation_uniform (std::mt19937& r, const system, vector_t>&, v_t) { std::normal_distribution dist(0.0,1.0); orthogonal_t ptr; ptr.is_reflection = true; @@ -135,7 +136,7 @@ orthogonal_t generate_rotation_uniform (std::mt19937& r, const vecto } template -orthogonal_t generate_rotation_perturbation (std::mt19937& r, const vector_t & v0, double epsilon, unsigned int n) { +orthogonal_t generate_rotation_perturbation (std::mt19937& r, const system, vector_t>& S, v_t i0, double epsilon, unsigned int n) { std::normal_distribution dist(0.0,1.0); orthogonal_t m; m.is_reflection = true; @@ -149,14 +150,14 @@ orthogonal_t generate_rotation_perturbation (std::mt19937& r, const double cosr = cos(2 * M_PI * rotation / (double)n / 2.0); double sinr = sin(2 * M_PI * rotation / (double)n / 2.0); - v[0] = v0[0] * cosr - v0[1] * sinr; - v[1] = v0[1] * cosr + v0[0] * sinr; + v[0] = S.s[i0][0] * cosr - S.s[i0][1] * sinr; + v[1] = S.s[i0][1] * cosr + S.s[i0][0] * sinr; for (q_t i = 2; i < q; i++) { - v[i] = v0[i]; + v[i] = S.s[i0][i]; } } else { - v = v0; + v = S.s[i0]; } double m_dot_v = 0; diff --git a/lib/include/wolff/models/potts.hpp b/lib/include/wolff/models/potts.hpp index 903f25f..10b440e 100644 --- a/lib/include/wolff/models/potts.hpp +++ b/lib/include/wolff/models/potts.hpp @@ -46,5 +46,5 @@ class potts_t { const potts_t states[256] = {{0}, {1}, {2}, {3}, {4}, {5}, {6}, {7}, {8}, {9}, {11}, {12}, {13}, {14}, {15}, {16}, {17}, {18}, {19}, {20}, {21}, {22}, {23}, {24}, {25}, {26}, {27}, {28}, {29}, {30}, {31}, {32}, {33}, {34}, {35}, {36}, {37}, {38}, {39}, {40}, {41}, {42}, {43}, {44}, {45}, {46}, {47}, {48}, {49}, {50}, {51}, {52}, {53}, {54}, {55}, {56}, {57}, {58}, {59}, {60}, {61}, {62}, {63}, {64}, {65}, {66}, {67}, {68}, {69}, {70}, {71}, {72}, {73}, {74}, {75}, {76}, {77}, {78}, {79}, {80}, {81}, {82}, {83}, {84}, {85}, {86}, {87}, {88}, {89}, {90}, {91}, {92}, {93}, {94}, {95}, {96}, {97}, {98}, {99}, {100}, {101}, {102}, {103}, {104}, {105}, {106}, {107}, {108}, {109}, {110}, {111}, {112}, {113}, {114}, {115}, {116}, {117}, {118}, {119}, {120}, {121}, {122}, {123}, {124}, {125}, {126}, {127}, {128}, {129}, {130}, {131}, {132}, {133}, {134}, {135}, {136}, {137}, {138}, {139}, {140}, {141}, {142}, {143}, {144}, {145}, {146}, {147}, {148}, {149}, {150}, {151}, {152}, {153}, {154}, {155}, {156}, {157}, {158}, {159}, {160}, {161}, {162}, {163}, {164}, {165}, {166}, {167}, {168}, {169}, {170}, {171}, {172}, {173}, {174}, {175}, {176}, {177}, {178}, {179}, {180}, {181}, {182}, {183}, {184}, {185}, {186}, {187}, {188}, {189}, {190}, {191}, {192}, {193}, {194}, {195}, {196}, {197}, {198}, {199}, {200}, {201}, {202}, {203}, {204}, {205}, {206}, {207}, {208}, {209}, {210}, {211}, {212}, {213}, {214}, {215}, {216}, {217}, {218}, {219}, {220}, {221}, {222}, {223}, {224}, {225}, {226}, {227}, {228}, {229}, {230}, {231}, {232}, {233}, {234}, {235}, {236}, {237}, {238}, {239}, {240}, {241}, {242}, {243}, {244}, {245}, {246}, {247}, {248}, {249}, {250}, {251}, {252}, {253}, {254}, {255}}; template -q_t finite_states_enum(potts_t state) { return (q_t)state.x; } +q_t finite_states_enum(const potts_t& state) { return (q_t)state.x; } diff --git a/lib/include/wolff/system.hpp b/lib/include/wolff/system.hpp index ac82f44..32ad38e 100644 --- a/lib/include/wolff/system.hpp +++ b/lib/include/wolff/system.hpp @@ -4,23 +4,25 @@ #include #include +#include -#include "types.h" #include "graph.hpp" +namespace wolff { + +#include "types.h" + +template +class measurement; + template -class wolff_system { +class system { public: - D_t D; // dimension - L_t L; // linear size v_t nv; // number of vertices v_t ne; // number of edges - graph_t G; // the graph defining the lattice with ghost + graph G; // the graph defining the lattice with ghost double T; // the temperature std::vector s; // the state of the ordinary spins -#ifndef WOLFF_NO_FIELD - R_t s0; // the current state of the ghost site -#endif #ifdef WOLFF_BOND_DEPENDENCE std::function Z; // coupling between sites @@ -29,6 +31,7 @@ class wolff_system { #endif #ifndef WOLFF_NO_FIELD + R_t s0; // the current state of the ghost site #ifdef WOLFF_SITE_DEPENDENCE std::function B; // coupling with the external field #else @@ -36,7 +39,7 @@ class wolff_system { #endif #endif - wolff_system(D_t D, L_t L, double T, + system(graph g, double T, #ifdef WOLFF_BOND_DEPENDENCE std::function Z #else @@ -49,7 +52,7 @@ class wolff_system { , std::function B #endif #endif - , lattice_t lat = SQUARE_LATTICE) : D(D), L(L), G(D, L, lat), T(T), Z(Z) + ) : G(g), T(T), Z(Z) #ifndef WOLFF_NO_FIELD , s0(), B(B) #endif @@ -58,13 +61,18 @@ class wolff_system { ne = G.ne; s.resize(nv); #ifndef WOLFF_NO_FIELD - G.add_ext(); + G.add_ghost(); #endif #ifdef WOLFF_FINITE_STATES finite_states_init(*this); #endif } + + void flip_cluster(v_t, const R_t&, std::mt19937&, measurement&); + void run_wolff(N_t, std::function &, v_t)> r_gen, measurement& A, std::mt19937& rng); }; +} + #endif diff --git a/lib/include/wolff/types.h b/lib/include/wolff/types.h index 68bd09e..1bc8c4d 100644 --- a/lib/include/wolff/types.h +++ b/lib/include/wolff/types.h @@ -8,7 +8,7 @@ typedef uint_fast16_t L_t; // linear size typedef uint_fast64_t N_t; // cycle iterator #define MAX_v UINT_FAST32_MAX -#define MAX_Q UINT_FAST8_MAX +#define MAX_q UINT_FAST8_MAX #define MAX_D UINT_FAST8_MAX #define MAX_L UINT_FAST16_MAX #define MAX_N UINT_FAST64_MAX -- cgit v1.2.3-70-g09d2