summaryrefslogtreecommitdiff
path: root/space_wolff.hpp
diff options
context:
space:
mode:
Diffstat (limited to 'space_wolff.hpp')
-rw-r--r--space_wolff.hpp113
1 files changed, 41 insertions, 72 deletions
diff --git a/space_wolff.hpp b/space_wolff.hpp
index a78116d..02c544d 100644
--- a/space_wolff.hpp
+++ b/space_wolff.hpp
@@ -30,16 +30,13 @@ const std::array<std::array<unsigned, 16>, 16> smiley = {
{{0, 0, 1, 1, 1, 0, 0, 0, 0, 0, 0, 1, 1, 1, 0, 0}},
{{0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0}}}};
-template <class T, int D>
-using Vector = Eigen::Matrix<T, D, 1>;
+template <class T, int D> using Vector = Eigen::Matrix<T, D, 1>;
-template <class T, int D>
-using Matrix = Eigen::Matrix<T, D, D>;
+template <class T, int D> using Matrix = Eigen::Matrix<T, D, D>;
/** brief diff - periodic subtraction of vectors
*/
-template <class T, int D>
-Vector<T, D> diff(T L, Vector<T, D> v1, Vector<T, D> v2) {
+template <class T, int D> Vector<T, D> diff(T L, Vector<T, D> v1, Vector<T, D> v2) {
Vector<T, D> v;
for (unsigned i = 0; i < D; i++) {
@@ -52,19 +49,17 @@ Vector<T, D> diff(T L, Vector<T, D> v1, Vector<T, D> v2) {
return v;
}
-template <class T, int D, class S>
-class Spin {
- public:
+template <class T, int D, class S> class Spin {
+public:
Vector<T, D> x;
S s;
};
-template <class T, int D>
-class Euclidean {
- private:
+template <class T, int D> class Euclidean {
+private:
T L;
- public:
+public:
Vector<T, D> t;
Matrix<T, D> r;
@@ -85,8 +80,7 @@ class Euclidean {
r = r0;
}
- template <class S>
- Spin<T, D, S> act(const Spin<T, D, S>& s) const {
+ template <class S> Spin<T, D, S> act(const Spin<T, D, S>& s) const {
Spin<T, D, S> s_new;
s_new.x = t + r * s.x;
@@ -122,14 +116,13 @@ class Euclidean {
}
};
-template <class T, int D, class S>
-class Octree {
- private:
+template <class T, int D, class S> class Octree {
+private:
T L;
unsigned N;
std::vector<std::set<Spin<T, D, S>*>> data;
- public:
+public:
Octree(T L, unsigned depth) : L(L), N(pow(2, depth)), data(pow(N, D)){};
unsigned ind(Vector<T, D> x) const {
@@ -148,16 +141,14 @@ class Octree {
void remove(Spin<T, D, S>* s) { data[ind(s->x)].erase(s); };
- std::set<Spin<T, D, S>*> neighbors(const Vector<T, D>& x,
- unsigned depth) const {
+ std::set<Spin<T, D, S>*> neighbors(const Vector<T, D>& x, unsigned depth) const {
std::set<Spin<T, D, S>*> ns;
std::set<unsigned> seen;
nearest_neighbors_of(ind(x), depth, ns, seen);
return ns;
};
- void nearest_neighbors_of(unsigned ind, unsigned depth,
- std::set<Spin<T, D, S>*>& ns,
+ void nearest_neighbors_of(unsigned ind, unsigned depth, std::set<Spin<T, D, S>*>& ns,
std::set<unsigned>& seen) const {
ns.insert(data[ind].begin(), data[ind].end());
seen.insert(ind);
@@ -165,9 +156,8 @@ class Octree {
if (depth > 0) {
for (unsigned i = 0; i < D; i++) {
for (signed j : {-1, 1}) {
- unsigned nind =
- pow(N, i + 1) * (ind / ((unsigned)pow(N, i + 1))) +
- fmod(pow(N, i + 1) + ind + j * pow(N, i), pow(N, i + 1));
+ unsigned nind = pow(N, i + 1) * (ind / ((unsigned)pow(N, i + 1))) +
+ fmod(pow(N, i + 1) + ind + j * pow(N, i), pow(N, i + 1));
if (!seen.contains(nind)) {
seen.insert(nind);
@@ -180,13 +170,13 @@ class Octree {
};
class quantity {
- private:
+private:
double total;
double total2;
std::vector<double> C;
unsigned wait;
- public:
+public:
unsigned n;
std::list<double> hist;
@@ -246,8 +236,7 @@ class quantity {
}
double σ() const {
- return 2.0 / (n - wait) * this->τ()[0] *
- (C[0] / (n - wait) - pow(this->avg(), 2));
+ return 2.0 / (n - wait) * this->τ()[0] * (C[0] / (n - wait) - pow(this->avg(), 2));
}
double serr() const { return sqrt(this->σ()); }
@@ -255,32 +244,25 @@ class quantity {
unsigned num_added() const { return n - wait; }
};
-template <class U, int D, class S>
-class Model;
+template <class U, int D, class S> class Model;
-template <class U, int D, class S>
-class measurement {
- public:
- virtual void pre_cluster(const Model<U, D, S>&, unsigned,
- const Euclidean<U, D>&){};
- virtual void plain_bond_visited(const Model<U, D, S>&, const Spin<U, D, S>*,
- const Spin<U, D, S>*, const Spin<U, D, S>&,
- double){};
- virtual void plain_site_transformed(const Model<U, D, S>&,
- const Spin<U, D, S>*,
+template <class U, int D, class S> class measurement {
+public:
+ virtual void pre_cluster(const Model<U, D, S>&, unsigned, const Euclidean<U, D>&){};
+ virtual void plain_bond_visited(const Model<U, D, S>&, const Spin<U, D, S>*, const Spin<U, D, S>*,
+ const Spin<U, D, S>&, double){};
+ virtual void plain_site_transformed(const Model<U, D, S>&, const Spin<U, D, S>*,
const Spin<U, D, S>&){};
- virtual void ghost_bond_visited(const Model<U, D, S>&, const Spin<U, D, S>&,
- const Spin<U, D, S>&, double){};
- virtual void ghost_site_transformed(const Model<U, D, S>&,
- const Euclidean<U, D>&){};
+ virtual void ghost_bond_visited(const Model<U, D, S>&, const Spin<U, D, S>&, const Spin<U, D, S>&,
+ double){};
+ virtual void ghost_site_transformed(const Model<U, D, S>&, const Euclidean<U, D>&){};
virtual void post_cluster(const Model<U, D, S>&){};
};
-template <class U, int D, class S>
-class Model {
- public:
+template <class U, int D, class S> class Model {
+public:
U L;
Euclidean<U, D> s0;
std::vector<Spin<U, D, S>> s;
@@ -294,8 +276,7 @@ class Model {
double E;
measurement<U, D, S>& A;
- void one_sequences(std::list<std::array<double, D>>& sequences,
- unsigned level) {
+ void one_sequences(std::list<std::array<double, D>>& sequences, unsigned level) {
if (level > 0) {
unsigned new_level = level - 1;
unsigned old_length = sequences.size();
@@ -308,18 +289,10 @@ class Model {
}
}
- Model(U L,
- std::function<double(const Spin<U, D, S>&, const Spin<U, D, S>&)> Z,
- std::function<double(const Spin<U, D, S>&)> B, unsigned dDepth,
- unsigned nDepth, measurement<U, D, S>& A)
- : L(L),
- s0(L),
- dict(L, dDepth),
- Z(Z),
- B(B),
- dDepth(dDepth),
- nDepth(nDepth),
- A(A) {
+ Model(U L, std::function<double(const Spin<U, D, S>&, const Spin<U, D, S>&)> Z,
+ std::function<double(const Spin<U, D, S>&)> B, unsigned dDepth, unsigned nDepth,
+ measurement<U, D, S>& A)
+ : L(L), s0(L), dict(L, dDepth), Z(Z), B(B), dDepth(dDepth), nDepth(nDepth), A(A) {
std::array<double, D> ini_sequence;
ini_sequence.fill(1);
std::list<std::array<double, D>> sequences;
@@ -327,7 +300,7 @@ class Model {
one_sequences(sequences, D);
- sequences.pop_front(); // don't want the identity matrix!
+ sequences.pop_front(); // don't want the identity matrix!
for (std::array<double, D> sequence : sequences) {
Matrix<U, D> m;
@@ -421,13 +394,10 @@ class Model {
} else {
Spin<U, D, S> si_new = r.act(*si);
std::set<Spin<U, D, S>*> all_neighbors;
- std::set<Spin<U, D, S>*> current_neighbors =
- dict.neighbors(si->x, nDepth);
- std::set<Spin<U, D, S>*> new_neighbors =
- dict.neighbors(si_new.x, nDepth);
+ std::set<Spin<U, D, S>*> current_neighbors = dict.neighbors(si->x, nDepth);
+ std::set<Spin<U, D, S>*> new_neighbors = dict.neighbors(si_new.x, nDepth);
- all_neighbors.insert(current_neighbors.begin(),
- current_neighbors.end());
+ all_neighbors.insert(current_neighbors.begin(), current_neighbors.end());
all_neighbors.insert(new_neighbors.begin(), new_neighbors.end());
all_neighbors.insert(NULL);
for (Spin<U, D, S>* sj : all_neighbors) {
@@ -461,8 +431,7 @@ class Model {
std::uniform_real_distribution<double> t_dist(0, L);
std::uniform_int_distribution<unsigned> r_dist(0, D - 1);
std::uniform_int_distribution<unsigned> ind_dist(0, s.size() - 1);
- std::uniform_int_distribution<unsigned> coin(
- 0, mats.size() + steps.size() - 1);
+ std::uniform_int_distribution<unsigned> coin(0, mats.size() + steps.size() - 1);
for (unsigned i = 0; i < N; i++) {
Vector<U, D> t;