Created
September 6, 2018 14:55
-
-
Save sbischoff-ai/2903f8049199504de95f8a53ff53d965 to your computer and use it in GitHub Desktop.
Simple Implementation of 2D Ising model. Initialize system with edge length and temperature, sweep until equilibrated. Plot to taste.
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
#include "Ising.h" | |
#include <cmath> | |
#ifndef _CSTDLIB | |
#define _CSTDLIB | |
#include <cstdlib> | |
#endif | |
Ising::Ising(unsigned int L_, double temp_){ | |
L = L_; | |
Lsq = L*L; | |
temp = temp_; | |
p2 = exp(-4.0/temp); | |
p4 = exp(-8.0/temp); | |
grid = new int[L*L]; | |
unsigned int i; | |
for (i = 0; i < Lsq; i++) grid[i] = 1; | |
} | |
Ising::~Ising(){ | |
delete[] grid; | |
}; | |
// Mersenne Twister for better random numbers | |
double random(){ | |
#define N 624 | |
#define M 397 | |
#define HI 0x80000000 | |
#define LO 0x7fffffff | |
#define MAX 0xFFFFFFFF | |
static const uint32_t seed = 5489UL; | |
static const uint32_t A[2] = { 0, 0x9908b0df }; | |
static uint32_t y[N]; | |
static int i = N + 1; | |
if (i > N) { | |
/* Initialize with pseudo random number */ | |
y[0] = seed; | |
for (i = 1; i < N; ++i) { | |
y[i] = (1812433253UL * (y[i - 1] ^ (y[i - 1] >> 30)) + i); | |
} | |
} | |
if (i == N) { | |
/* calculate new state vector */ | |
uint32_t h; | |
for (i = 0; i < N - M; ++i) { | |
h = (y[i] & HI) | (y[i + 1] & LO); | |
y[i] = y[i + M] ^ (h >> 1) ^ A[h & 1]; | |
} | |
for (; i < N - 1; ++i) { | |
h = (y[i] & HI) | (y[i + 1] & LO); | |
y[i] = y[i + (M - N)] ^ (h >> 1) ^ A[h & 1]; | |
} | |
h = (y[N - 1] & HI) | (y[0] & LO); | |
y[N - 1] = y[M - 1] ^ (h >> 1) ^ A[h & 1]; | |
i = 0; | |
} | |
uint32_t e = y[i++]; | |
/* tempering */ | |
e ^= (e >> 11); | |
e ^= (e << 7) & 0x9d2c5680; | |
e ^= (e << 15) & 0xefc60000; | |
e ^= (e >> 18); | |
return (double)e / MAX; | |
#undef N | |
#undef M | |
#undef HI | |
#undef LO | |
} | |
void Ising::sweep(){ | |
unsigned int i, xy, x, y; | |
int deltaE; | |
for(i = 0; i < Lsq; i++){ | |
xy = Lsq*random(); | |
y = xy%L; | |
x = xy/L; | |
deltaE = grid[xy]*(grid[L*((x-1)%L)+y]+grid[L*((x+1)%L)+y]+grid[L*x+(y-1)%L]+grid[L*x+(y+1)%L]); // calculates half of deltaE | |
switch(deltaE){ | |
case 2: | |
if(random() < p2) grid[xy] *= -1; | |
break; | |
case 4: | |
if(random() < p4) grid[xy] *= -1; | |
break; | |
default: grid[xy] *= -1; | |
} | |
} | |
} |
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
#ifndef ISING_H | |
#define ISING_H | |
class Ising{ | |
public: | |
Ising(unsigned int, double); | |
~Ising(); | |
void sweep(); | |
int* grid; | |
private: | |
double temp, p2, p4; | |
unsigned int L; | |
unsigned int Lsq; | |
}; | |
#endif //ISING_H |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment