summaryrefslogtreecommitdiff
path: root/examples/ising_animation.cpp
blob: ab105859dbd6cf61f963e97b1f91c4bb67dc3126 (plain)
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;
}