computational physicscorrelation-functionscritical-phenomenaising-model
I'm simulating the Ising Model in 2D up to 5D and I want to calculate the spin-spin correlation, correlation length, and critical exponent of the system. What is a good way to go about doing this?
Generally, I know I need:
$$<s_os_r> \cong \frac{1}{N}\sum_i^N s_os_r $$
Where $N$ is the number of Monte Carlo trials run.
What is a good way to actually calculate $s_os_r$? And how does this help me with calculating the correlation length $\xi$ and the critical exponent $\nu$
That is not a definition of correlation length. (It is a definition of the critical exponent.)
The correlation length is defined in terms of the 2-point correlation function of spin observables. Pick points $x$ and $y$ on the lattice, and consider the expectation value $\langle s(x) s(y) \rangle$ of the product of the spin observable at $x$ and the spin observable at $y$. This quantity tells you how strongly correlated the spin at $x$ and the spin at $y$ are, as a function of the temperature, coupling constant, and the distance between $x$ and $y$. If $T > T_c$, then the correlation function dies off exponentially fast in $|x-y|$.
$\langle s(x) s(y) \rangle \sim e^{-\frac{|x-y|}{\xi(T)}}$
The correlation length is, by definition, the constant (in $x$ and $y$, but not in $T$) which tells you how fast the correlation function vanishes.
I would argue that this maybe due to the way you calculate your autocorrelation. An autocorrelation like that straight line is the result of a large square signal.
The Ising model has a phase transition at the critical temperature. Above it, it's disordered; below it, it becomes ordered, which means that the magnetization stops flipping back and forth. This was shown analytically by L. Osanger in his article Crystal Statistics. I. A Two-Dimensional Model with an Order-Disorder Transition. Now, I assume that if you are using The Metropolis algorithm you are using a finite lattice. This just causes the transition to become less sharp (even if you use periodic boundary conditions), but it's still there, as can be seen in this graph which also uses the Metropolis algorithm, in a grid of 100 spins:
So you can see that it's not unexpected that below the critical temperature, all the spins align and you just get a constant magnetization. Now, a constant signal should really give you a constant autocorrelation, but if your integration is done over a finite domain, which I assume it is, you would get an slopping autocorrelation like that. This picture should help seeing why:
The value of the green area will decrease linearly with T.
Best Answer
I will explain how I measured the spin-spin correlation function for the 2d Ising model. Generalization to more than 2 dimensions should be straightforward as long as you have hypercubic lattices.
Just to get the notation straight: Let's use the name $\sigma_{(i,j)}$ for the spin at position $\vec r_{(i,j)} = \begin{pmatrix} i \\ j \end{pmatrix}$. Let's assume our Ising model has $L \times L$ spins, so we have $1 \leq i,j \leq L$.
The first thing to note is that I assume that the spin-spin correlation function - lets call it $\chi_{(i,j)(k,l)}$ - only depends on the absolute distance between the two spins $\sigma_{(i,j)}$ and $\sigma_{(k,l)}$.
$$\chi_{(i,j)(k,l)} = \chi( | \vec r_{(i,j)} - \vec r_{(k,l)} | ) = \chi(r)$$
I'm not sure if this is true in general, but there is a proof that it is true for the 2d Ising model at the critical temperature. (A rotationally symmetric spin-spin correlation function is proven, and that is just another way of saying that it depends on the distance only.) Don't take my word for it, but for now I wouldn't worry about this assumption ...
Knowing this, it is easy to see that at every step in your Markov chain you don't have to look at all possible $L^2$ pairs of spins. Instead you can only look at a subset of pairs, where you have all distances $r$ you are interested in within this subset.
For my implementation I chose to use all spins $\sigma_{(i,i)}$ along the diagonal of the lattice and to only measure the correlation along the two directions of the lattice vectors. So I only consider the correlation between $\sigma_{(i,i)}$ and all $\sigma_{(i,j)}$ (horizontally) and $\sigma_{(j,i)}$ (vertically). Remember that the spin-spin correlation function is rotationally symmetric, so you don't gain anything by considering pairs of spins along any other direction. Measuring along the axes makes your life much easier as you will only have integer distances $r$, meaning that you can use a plain $L-1$ element array with the distance as an index to accumulate the terms of your sum. This is what you should calculate at every Monte Carlo step:
$$ \text{sum}(r) = \sum_{i=1}^L \sum_{j=1}^L \left( \sigma_{(i,i)} \sigma_{(i,j)} + \sigma_{(i,i)} \sigma_{(j,i)} \right) \; \text{ with } \; r = |i-j|$$
$$ \text{samples}(r) = \sum_{i=1}^L \sum_{j=1}^L 2 \; \text{ with } \; r = |i-j|$$
$$ \chi_n(r) = \frac{\text{sum}(r)}{\text{samples}(r)}$$
This is exactly what this function IsingModel2d::ss_corr() in my code does. (Keeping track of the samples is only needed because you make a different number of measurements for different values of $r$. Consider as an example the extreme cases: All but two spins $\sigma_{(i,i)}$ on the diagonal have 4 nearest neighbors at distance $r=1$, but only $\sigma_{(1,1)}$ and $\sigma_{(L,L)}$ have two partners each at distance $r = L-1$. So you have to compensate for the fact that you measure the correlation between close spins more often than between far away spins.)
Note that you will get a different $\chi_n(r)$ at every Monte Carlo step, depending on your current configuration. In the end you should average them (just like you do with any other observable) to get your final result for the spin-spin correlation function:
$$\chi(r) = \frac{1}{N-1} \sum_{n=1}^N \chi_n(r)$$
Once you have this, you can fit $\chi(r) = r^{-\beta}$ where you should hopefully find $\beta = 0.25$. Check out this other question on details regarding the fit.
Hope this helps!