[Physics] 2D Ising model, understanding autocorrelations, Monte Carlo

correlation-functionsising-modelstatistical mechanics

I've been struggling quite a bit with implementing an autocorrelation code into my current project. The autocorrelation as it is now, is increasing exponentially from 1 at the start of my MC run, and hitting 2 halfway through the MC simulation regardless of how many sweeps I do through the lattice.

The system

10×10 square lattice with no external magnetic field and ferromagnetic coupling. The reason for 10×10 is for fast execution of the code in order to build it.

Here is what I've done so far:

  1. Letting the Metropolis Monte Carlo work until the system is in equilibrium (checking this by running two different initial states with different random seeds).
  2. Then I start sweeping through the lattice, updating the energies and magnetization for each attempt at flipping one spin. When I've done one sweep over the lattice, the last value of the energy and magnetization gets stored. Then it continues to the next sweep and updating the values further. Thus, the energies and magnetization only gets stored once per sweep.
  3. When I've done, say 2000 sweeps, I calculate the autocorrelation for the system according to Newman & Barkema (Eq 3.21 in http://itf.fys.kuleuven.be/~fpspXIII/material/Barkema_FPSPXIII.pdf). The formula reads: $$\chi(t)=\frac{1}{t_{max}-t}\sum_{t'=0}^{t_{max}-t}m(t')m(t'+t)-\frac{1}{t_{max}-t}\sum_{t'=0}^{t_{max}-t}m(t')\times\frac{1}{t_{max}-t}\sum_{t'=0}^{t_{max}-t}m(t'+t)$$ where t defines number of sweeps of the lattice, i.e the displacement/lag from some value.

My problem
From my simulation at $k_{b}T/J$, where $k_{b}$ is set to 1, $J=1$ is the ferromagnetic coupling, and $T=1$ the autocorrelation function grows. I've tried to normalize it by dividing by the first value resulting in a start at 1, but it acts strange as stated. Thus, I started to calculate by hand trying to see if I had coded something wrong. I worked with a system where each spin in the lattice had spin 1 for every $t\Rightarrow m(t')=m(t'+t)=\langle m \rangle=1$. I then cutoff the autocorrelation at $t$=1000 when running 2000 sweeps. The formula then reduces down to $$\chi(t)=\frac{1}{2000-t}\sum_{t'=0}^{2000-t}1-\frac{1}{2000-t}\sum_{t'=0}^{2000-t}1\times\frac{1}{2000-t}\sum_{t'=0}^{2000-t}1$$. Then for some values: $$\chi(0)=\frac{1}{2000}\sum_{t'=0}^{2000}1-\frac{1}{2000}\sum_{t'=0}^{2000}1\times\frac{1}{2000}\sum_{t'=0}^{2000}1$$
$$=\frac{2001}{2000}-\left(\frac{2001}{2000}\right)^{2}\approx-5\cdot10^{-4}$$
$$\chi(500)=\frac{1501}{1500}-\left(\frac{1501}{1500}\right)^{2}\approx-6.6\cdot10^{-4}$$
$$\chi(1000)=\frac{1001}{1000}-\left(\frac{1001}{1000}\right)^{2}\approx-1\cdot10^{-3}$$
As we see, the autocorrelation value has doubled when check for half of the number of sweeps through the lattice.

Questions:

  1. I'd expect that the correlation function would behave as an exponentially decaying function like $e^{-t/\tau}$ where $\tau$ is the correlation time, but rather, the plots show exponential growth with values as calculated above.
  2. The expression for the autocorrelation as stated in this post is a discretization of $$\chi(t)=\int dt'(m(t')-\langle m\rangle)(m(t'+t)-\langle m \rangle)$$ which implies that if all values are 1, the integral should yield some constant, not an increasing function. What went wrong?

Best Answer

Your sum needs to go from $1$ to $2000$ instead of $0$ to $2000$. Essentially you are taking an average so it doesn't make sense to divide by fewer values than you are summing together (i.e. $2001/2000$). If you make this correction, all your calculations above work out to $0$. This leads to an indeterminate auto-correlation function as you are left with $0/0$ when normalizing by the variance.

Related Question