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
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
|
#include <getopt.h>
#include <iostream>
#include <chrono>
#include <GL/glut.h>
#define WOLFF_USE_FINITE_STATES
#include <wolff/models/ising.hpp>
#include <wolff.hpp>
using namespace wolff;
class draw_ising : public measurement<ising_t, ising_t> {
private:
unsigned int frame_skip;
v_t C;
public:
draw_ising(const system<ising_t, ising_t>& S, unsigned int window_size, unsigned int frame_skip, int argc, char *argv[]) : frame_skip(frame_skip){
glutInit(&argc, argv);
glutInitDisplayMode(GLUT_SINGLE | GLUT_RGB);
glutInitWindowSize(window_size, window_size);
glutCreateWindow("wolff");
glClearColor(0.0,0.0,0.0,0.0);
glMatrixMode(GL_PROJECTION);
glLoadIdentity();
gluOrtho2D(0.0, S.G.L, 0.0, S.G.L);
}
void pre_cluster(N_t, N_t, const system<ising_t, ising_t>& S, v_t, const ising_t&) {
glClear(GL_COLOR_BUFFER_BIT);
for (v_t i = 0; i < pow(S.G.L, 2); i++) {
if (S.s[i].x == S.s0.x) {
glColor3f(0.0, 0.0, 0.0);
} else {
glColor3f(1.0, 1.0, 1.0);
}
glRecti(i / S.G.L, i % S.G.L, (i / S.G.L) + 1, (i % S.G.L) + 1);
}
glFlush();
C = 0;
}
void plain_bond_visited(const system<ising_t, ising_t>&, v_t, const ising_t&, v_t, double dE) {}
void ghost_bond_visited(const system<ising_t, ising_t>&, v_t, const ising_t& s_old, const ising_t& s_new, double dE) {}
void plain_site_transformed(const system<ising_t, ising_t>& S, v_t i, const ising_t&) {
glColor3f(1.0, 0.0, 0.0);
glRecti(i / S.G.L, i % S.G.L, (i / S.G.L) + 1, (i % S.G.L) + 1);
C++;
if (C % frame_skip == 0) {
glFlush();
}
}
void ghost_site_transformed(const system<ising_t, ising_t>&, const ising_t&) {}
void post_cluster(N_t, N_t, const system<ising_t, ising_t>&) {}
};
int main(int argc, char *argv[]) {
// set defaults
N_t N = (N_t)1e4;
D_t D = 2;
L_t L = 128;
double T = 2.26918531421;
double H = 0.0;
unsigned int window_size = 512;
unsigned int frame_skip = 1;
int opt;
// take command line arguments
while ((opt = getopt(argc, argv, "N:D:L:T:H:w:f:")) != -1) {
switch (opt) {
case 'N': // number of steps
N = (N_t)atof(optarg);
break;
case 'D': // dimension
D = atoi(optarg);
break;
case 'L': // linear size
L = atoi(optarg);
break;
case 'T': // temperature
T = atof(optarg);
break;
case 'H': // external field
H = atof(optarg);
break;
case 'w':
window_size = atoi(optarg);
break;
case 'f':
frame_skip = atoi(optarg);
break;
default:
exit(EXIT_FAILURE);
}
}
// define the spin-spin coupling
std::function <double(const ising_t&, const ising_t&)> Z = [] (const ising_t& s1, const ising_t& s2) -> double {
return (double)(s1 * s2);
};
// define the spin-field coupling
std::function <double(const ising_t&)> B = [=] (const ising_t& s) -> double {
return H * s;
};
// initialize the lattice
graph G(D, L);
// initialize the system
system<ising_t, ising_t> S(G, T, Z, B);
// define function that generates self-inverse rotations
std::function <ising_t(std::mt19937&, const system<ising_t, ising_t>&, v_t)> gen_R = [] (std::mt19937&, const system<ising_t, ising_t>&, v_t) -> ising_t {
return ising_t(true);
};
// initailze the measurement object
draw_ising A(S, window_size, frame_skip, argc, argv);
// initialize the random number generator
auto seed = std::chrono::high_resolution_clock::now().time_since_epoch().count();
std::mt19937 rng{seed};
// run wolff N times
S.run_wolff(N, gen_R, A, rng);
// exit
return 0;
}
|