1 files changed, 9 insertions, 14 deletions

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());
   }