diff options
Diffstat (limited to 'lib/include')
| -rw-r--r-- | lib/include/wolff.hpp | 15 | ||||
| -rw-r--r-- | lib/include/wolff/cluster.hpp | 20 | ||||
| -rw-r--r-- | lib/include/wolff/graph.hpp | 9 | ||||
| -rw-r--r-- | lib/include/wolff/meas.h | 19 | ||||
| -rw-r--r-- | lib/include/wolff/state.hpp | 81 | 
5 files changed, 46 insertions, 98 deletions
diff --git a/lib/include/wolff.hpp b/lib/include/wolff.hpp index c10a211..b730c8d 100644 --- a/lib/include/wolff.hpp +++ b/lib/include/wolff.hpp @@ -3,7 +3,7 @@  #include "wolff/state.hpp"  template <class R_t, class X_t> -void wolff(count_t N, state_t <R_t, X_t>& s, std::function <R_t(std::mt19937&, X_t)> gen_R, std::function <void(const state_t <R_t, X_t>&)> measurements, std::mt19937& r, bool silent) { +void wolff(count_t N, state_t <R_t, X_t>& s, std::function <R_t(std::mt19937&, X_t)> gen_R, wolff_measurement<R_t, X_t>& m, std::mt19937& r) {  #ifdef FINITE_STATES  #ifdef NOFIELD @@ -15,21 +15,16 @@ void wolff(count_t N, state_t <R_t, X_t>& s, std::function <R_t(std::mt19937&, X    std::uniform_int_distribution<v_t> dist(0, s.nv); -  if (!silent) printf("\n");    for (count_t steps = 0; steps < N; steps++) { -    if (!silent) printf("\033[F\033[JWOLFF: step %" PRIu64 " / %" PRIu64 ": E = %.2f, S = %" PRIv "\n", steps, N, s.E, s.last_cluster_size); -      v_t v0 = dist(r);      R_t step = gen_R(r, s.spins[v0]); -    flip_cluster <R_t, X_t> (s, v0, step, r); -    measurements(s); -  } +    m.pre_cluster(s, steps, N, v0, step); + +    flip_cluster<R_t, X_t>(s, v0, step, r, m); -  if (!silent) { -    printf("\033[F\033[J"); +    m.post_cluster(s, steps, N);    } -  printf("WOLFF: step %" PRIu64 " / %" PRIu64 ": E = %.2f, S = %" PRIv "\n", N, N, s.E, s.last_cluster_size);  } diff --git a/lib/include/wolff/cluster.hpp b/lib/include/wolff/cluster.hpp index e9dff7b..805e2c3 100644 --- a/lib/include/wolff/cluster.hpp +++ b/lib/include/wolff/cluster.hpp @@ -9,11 +9,11 @@  #include "types.h"  #include "state.hpp"  #include "graph.hpp" +#include "meas.h"  template <class R_t, class X_t> -void flip_cluster(state_t<R_t, X_t>& s, v_t v0, const R_t& r, std::mt19937& rand) { -  std::uniform_real_distribution<double> dist(0.0,1.0); -  v_t nv = 0; +void flip_cluster(state_t<R_t, X_t>& s, v_t v0, const R_t& r, std::mt19937& rand, wolff_measurement<R_t, X_t>& m) { +  std::uniform_real_distribution<double> dist(0.0, 1.0);    std::stack<v_t> stack;    stack.push(v0); @@ -75,8 +75,8 @@ void flip_cluster(state_t<R_t, X_t>& s, v_t v0, const R_t& r, std::mt19937& rand            prob = H_probs[state_to_ind(rs_old)][state_to_ind(rs_new)];  #endif -          s.update_magnetization(rs_old, rs_new); -          s.update_fourierZero(non_ghost, rs_old, rs_new); +          // run measurement hooks for encountering a ghost bond +          m.ghost_bond_added(non_ghost, rs_old, rs_new, dE);          } else // this is a perfectly normal bond!  #endif          { @@ -90,9 +90,10 @@ void flip_cluster(state_t<R_t, X_t>& s, v_t v0, const R_t& r, std::mt19937& rand  #ifdef FINITE_STATES            prob = J_probs[state_to_ind(s.spins[v])][state_to_ind(si_new)][state_to_ind(s.spins[vn])];  #endif -        } -        s.update_energy(dE); +          // run measurement hooks for encountering a plain bond +          m.plain_bond_added(v, s.spins[v], si_new, vn, s.spins[vn], dE); +        }  #ifndef FINITE_STATES          prob = 1.0 - exp(-dE / s.T); @@ -105,16 +106,15 @@ void flip_cluster(state_t<R_t, X_t>& s, v_t v0, const R_t& r, std::mt19937& rand  #ifndef NOFIELD        if (v_is_ghost) { +        m.ghost_site_transformed(s.R, R_new);          s.R = R_new;        } else  #endif        { +        m.plain_site_transformed(v, s.spins[v], si_new);          s.spins[v] = si_new; -        nv++;        }      }    } - -  s.last_cluster_size = nv;  } diff --git a/lib/include/wolff/graph.hpp b/lib/include/wolff/graph.hpp index 165544a..0aeb6af 100644 --- a/lib/include/wolff/graph.hpp +++ b/lib/include/wolff/graph.hpp @@ -1,13 +1,16 @@  #pragma once -#include <inttypes.h>  #include <cmath> -#include <stdlib.h>  #include <vector>  #include "types.h" +typedef enum lattice_t { +  SQUARE_LATTICE, +  DIAGONAL_LATTICE +} lattice_t; +  class graph_t {    public:      v_t ne; @@ -15,7 +18,7 @@ class graph_t {      std::vector<std::vector<v_t>> v_adj;      std::vector<std::vector<double>> coordinate; -    graph_t(D_t D, L_t L); +    graph_t(D_t D, L_t L, lattice_t lat = SQUARE_LATTICE);      void add_ext();  }; diff --git a/lib/include/wolff/meas.h b/lib/include/wolff/meas.h new file mode 100644 index 0000000..4895eac --- /dev/null +++ b/lib/include/wolff/meas.h @@ -0,0 +1,19 @@ + +#pragma once + +#include "types.h" + +template <class R_t, class X_t> +class wolff_measurement { +  public: +    virtual void pre_cluster(const state_t<R_t, X_t>&, count_t, count_t, v_t, const R_t&) = 0; + +    virtual void plain_bond_added(v_t, const X_t&, const X_t&, v_t, const X_t&, double) = 0; +    virtual void ghost_bond_added(v_t, const X_t&, const X_t&, double) = 0; + +    virtual void plain_site_transformed(v_t, const X_t&, const X_t&) = 0; +    virtual void ghost_site_transformed(const R_t&, const R_t&) = 0; + +    virtual void post_cluster(const state_t<R_t, X_t>&, count_t, count_t) = 0; +}; + diff --git a/lib/include/wolff/state.hpp b/lib/include/wolff/state.hpp index e7c0ac3..4909881 100644 --- a/lib/include/wolff/state.hpp +++ b/lib/include/wolff/state.hpp @@ -9,26 +9,17 @@  template <class R_t, class X_t>  class state_t { -  private: -    // updating fourier terms F requires many cos and sin calls, faster to do it beforehand. -    std::vector<std::vector<double>> precomputed_cos; -    std::vector<std::vector<double>> precomputed_sin;    public: -    D_t D; -    L_t L; -    v_t nv; // the number of vertices in the lattice -    v_t ne; // the number of edges in the lattice -    graph_t g; // the graph defining the lattice without ghost +    D_t D; // the dimension of the system +    L_t L; // the linear size of the lattice +    v_t nv; // the number of vertices in the original lattice +    v_t ne; // the number of edges in the original lattice +    graph_t g; // the graph defining the lattice with ghost      double T; // the temperature      std::vector<X_t> spins; // the state of the ordinary spins  #ifndef NOFIELD      R_t R; // the current state of the ghost site  #endif -    double E; // the system's total energy -    typename X_t::M_t M; // the "sum" of the spins, like the total magnetization -    v_t last_cluster_size; // the size of the last cluster -    std::vector<typename X_t::F_t> ReF; -    std::vector<typename X_t::F_t> ImF;  #ifdef BOND_DEPENDENCE      std::function <double(v_t, const X_t&, v_t, const X_t&)> J; // coupling between sites @@ -57,7 +48,7 @@ class state_t {          , std::function <double(const X_t&)> H  #endif  #endif -           ) : D(D), L(L), g(D, L), T(T), +        , lattice_t lat = SQUARE_LATTICE) : D(D), L(L), g(D, L, lat), T(T),  #ifndef NOFIELD                R(),  #endif @@ -69,69 +60,9 @@ class state_t {        nv = g.nv;        ne = g.ne;        spins.resize(nv); -#ifdef BOND_DEPENDENCE -      E = 0; -      for (v_t v = 0; v < nv; v++) { -        for (const v_t &vn : g.v_adj[v]) { -          if (v < vn) { -            E -= J(v, spins[v], vn, spins[vn]); -          } -        } -      } -#else -      E = - (double)ne * J(spins[0], spins[0]); -#endif -  #ifndef NOFIELD        g.add_ext(); -#ifdef SITE_DEPENDENCE -      for (v_t i = 0; i < nv; i++) { -        E -= H(i, spins[i]); -      } -#else -      E -= (double)nv * H(spins[0]); -#endif -#endif - -      M = spins[0] * nv; -      last_cluster_size = 0; -      ReF.resize(D); -      ImF.resize(D); -      for (D_t i = 0; i < D; i++) { -        ReF[i] = spins[0] * 0.0; -        ImF[i] = spins[0] * 0.0; -      } -      precomputed_cos.resize(nv); -      precomputed_sin.resize(nv); -      for (v_t i = 0; i < nv; i++) { -        precomputed_cos[i].resize(D); -        precomputed_sin[i].resize(D); -        for (D_t j = 0; j < D; j++) { -          precomputed_cos[i][j] = cos(2 * M_PI * g.coordinate[i][j] / (double)L); -          precomputed_sin[i][j] = sin(2 * M_PI * g.coordinate[i][j] / (double)L); -        } -      } -    } - -    void update_magnetization(const X_t& s_old, const X_t& s_new) { -      M += s_new - s_old; -    } - -    void update_energy(const double& dE) { -      E += dE; -    } - -    void update_fourierZero(v_t v, const X_t& s_old, const X_t& s_new) { -#ifdef DIMENSION -      for (D_t i = 0; i < DIMENSION; i++) { -#else -      for (D_t i = 0; i < D; i++) {  #endif -        ReF[i] += (s_new - s_old) * precomputed_cos[v][i]; -        ImF[i] += (s_new - s_old) * precomputed_sin[v][i]; -      }      }  }; - -  | 
