1
2
3
4
5
6
7
8
9
10
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
|
#pragma once
#include <array>
#include <list>
#include <vector>
#include "matrix.hpp"
#include "vector.hpp"
#include "euclidean.hpp"
#include "transformation.hpp"
template <int D> 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();
for (std::array<double, D>& sequence : sequences) {
std::array<double, D> new_sequence = sequence;
new_sequence[new_level] = -1;
sequences.push_front(new_sequence);
}
one_sequences<D>(sequences, new_level);
}
}
template <class U, int D> std::vector<Vector<U, D>> torus_vecs(U L) {
std::vector<Vector<U, D>> vecs;
std::array<double, D> ini_sequence;
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) {
Vector<U, D> v;
for (unsigned i = 0; i < D; i++) {
if (sequence[i] == 1) {
v(i) = 0;
} else {
v(i) = L / 2;
}
}
vecs.push_back(v);
}
return vecs;
}
template <class U, int D> std::vector<Matrix<U, D>> torus_mats() {
std::vector<Matrix<U, D>> mats;
std::array<double, D> ini_sequence;
ini_sequence.fill(1);
std::list<std::array<double, D>> sequences;
sequences.push_back(ini_sequence);
one_sequences<D>(sequences, D);
sequences.pop_back(); // don't want the identity matrix!
for (std::array<double, D> sequence : sequences) {
Matrix<U, D> m;
for (unsigned i = 0; i < D; i++) {
for (unsigned j = 0; j < D; j++) {
if (i == j) {
m(i, j) = sequence[i];
} else {
m(i, j) = 0;
}
}
}
mats.push_back(m);
}
for (unsigned i = 0; i < D; i++) {
for (unsigned j = 0; j < D; j++) {
if (i != j) {
Matrix<U, D> m;
for (unsigned k = 0; k < D; k++) {
for (unsigned l = 0; l < D; l++) {
if ((k == i && l == j) || (k == j && l == i)) {
if (i < j) {
m(k, l) = 1;
} else {
m(k, l) = -1;
}
} else {
m(k, l) = 0;
}
}
}
mats.push_back(m);
}
}
}
return mats;
}
template <int D, class S>
Gen<double, D, TorusGroup<double, D>, S> uniformGenTorus(double L) {
std::vector<Vector<double, D>> torusVectors = torus_vecs<double, D>(L);
std::vector<Matrix<double, D>> torusMatrices = torus_mats<double, D>();
return [L, torusVectors,
torusMatrices](Model<double, D, TorusGroup<double, D>, S>& M,
Rng& r) -> Transformation<double, D, TorusGroup<double, D>, S>* {
Matrix<double, D> m = r.pick(torusMatrices);
Vector<double, D> t;
for (unsigned i = 0; i < D; i++) {
t(i) = r.uniform<double>(0, L);
}
TorusGroup<double, D> g = TorusGroup<double, D>({(double)L, t - m * t, m});
Spin<double, D, S>* ss = r.pick(M.s);
return new SpinFlip<double, D, TorusGroup<double, D>, S>(M, g, ss);
};
}
|
