diff options
-rw-r--r-- | space_wolff.hpp | 20 | ||||
-rw-r--r-- | spheres.cpp | 23 |
2 files changed, 18 insertions, 25 deletions
diff --git a/space_wolff.hpp b/space_wolff.hpp index 1335f07..b542fbb 100644 --- a/space_wolff.hpp +++ b/space_wolff.hpp @@ -271,7 +271,7 @@ void one_sequences(std::list<std::array<double, D>>& sequences, unsigned level) new_sequence[new_level] = -1; sequences.push_front(new_sequence); } - one_sequences(sequences, new_level); + one_sequences<D>(sequences, new_level); } } @@ -282,6 +282,7 @@ std::vector<Vector<U, D>> torus_vecs(U L) { ini_sequence.fill(1); std::list<std::array<double, D>> sequences; sequences.push_back(ini_sequence); + one_sequences<D>(sequences, D); sequences.pop_front(); // don't want the identity matrix! for (std::array<double, D> sequence : sequences) { @@ -309,7 +310,7 @@ std::vector<Matrix<U, D>> torus_mats() { std::list<std::array<double, D>> sequences; sequences.push_back(ini_sequence); - one_sequences(sequences, D); + one_sequences<D>(sequences, D); sequences.pop_front(); // don't want the identity matrix! @@ -375,7 +376,6 @@ public: void step(double T, unsigned ind, Euclidean<U, D> r) { unsigned cluster_size = 0; - std::uniform_real_distribution<double> dist(0.0, 1.0); std::queue<Spin<U, D, S>*> queue; queue.push(&(s[ind])); @@ -396,7 +396,7 @@ public: Spin<U, D, S> s0s_new = s0_new.inverse().act(ss); double p = 1.0 - exp(-(B(s0s_new) - B(s0s_old)) / T); A.ghost_bond_visited(*this, s0s_old, s0s_new, p); - if (dist(rng) < p) { + if (rng.uniform(0.0, 1.0) < p) { queue.push(&ss); } } @@ -423,7 +423,7 @@ public: p = 1.0 - exp(-(Z(*si, *sj) - Z(si_new, *sj)) / T); A.plain_bond_visited(*this, si, sj, si_new, p); } - if (dist(rng) < p) { + if (rng.uniform(0.0, 1.0) < p) { queue.push(sj); } } @@ -439,15 +439,13 @@ public: } void wolff(double T, unsigned N) { - std::uniform_int_distribution<unsigned> ind_dist(0, s.size() - 1); - for (unsigned i = 0; i < N; i++) { - R g = g(*this, rng); - unsigned ind = ind_dist(rng); + R r = g(*this, rng); + unsigned ind = rng.uniform((unsigned)0, (unsigned)(s.size() - 1)); - A.pre_cluster(*this, ind, g); + A.pre_cluster(*this, ind, r); - this->step(T, ind, g, rng); + this->step(T, ind, r); A.post_cluster(*this); } diff --git a/spheres.cpp b/spheres.cpp index 9e67b87..4ed329e 100644 --- a/spheres.cpp +++ b/spheres.cpp @@ -37,17 +37,14 @@ class animation : public measurement<double, 2, Euclidean<double, D>, double> { std::function<Euclidean<double, D>(const model&, randutils::mt19937_rng&)> eGen(const std::vector<Matrix<double, 2>>& mats, const std::vector<Vector<double, 2>>& vecs) { return [&mats, &vecs] (const model& M, randutils::mt19937_rng& rng) -> Euclidean<double, 2> { - std::uniform_real_distribution<double> t_dist(0, M.L); - std::uniform_int_distribution<unsigned> r_dist(0, D - 1); - std::uniform_int_distribution<unsigned> ind_dist(0, M.s.size() - 1); - std::uniform_int_distribution<unsigned> coin(0, mats.size() + vecs.size() - 1); - Vector<double, D> t; Matrix<double, D> m; - unsigned flip = coin(rng); + unsigned flip = rng.uniform((unsigned)0, (unsigned)(mats.size() + vecs.size() - 1)); if (flip < mats.size()) { + unsigned f_ind = rng.uniform((unsigned)0, (unsigned)M.s.size()); + t = M.s[f_ind].x; for (unsigned j = 0; j < D; j++) { - t(j) = (double)t_dist(rng); + t(j) += rng.variate<double, std::normal_distribution>(0.0, 0.1); } m = mats[flip]; } else { @@ -120,7 +117,7 @@ int main(int argc, char* argv[]) { std::function<double(Spin<double, D, double>)> B = [L, H] (Spin<double, D, double> s) -> double { - return H * sin(2 * M_PI * 3 * s.x(0) / L); + return H * s.x(1); }; std::vector<Matrix<double, D>> mats = torus_mats<double, D>(); @@ -129,16 +126,14 @@ int main(int argc, char* argv[]) { animation A(L, 750, argc, argv); model sphere(L, Z, B, g, std::floor(log2(L)), 2, A); - randutils::auto_seed_128 seeds; - std::mt19937 rng{seeds}; - - std::uniform_real_distribution<double> dist(0.0, L); + randutils::mt19937_rng rng; sphere.s.reserve(n); + unsigned nx = floor(sqrt(n)); for (unsigned i = 0; i < n; i++) { - Vector<double, D> pos = {dist(rng), dist(rng)}; - sphere.s.push_back({pos, dist(rng) / L}); + Vector<double, D> pos = {(i / nx) * L / nx + rng.uniform(0.0, L / (4 * nx)), (i % nx) * L / nx + rng.uniform(0.0, L / (4 * nx))}; + sphere.s.push_back({pos, 0.5}); sphere.dict.insert(&sphere.s.back()); } |