# Averages of absolute values in Monte Carlo simulation of Ising Model

Consider the 2D Ising model in $$0$$ field, with Hamiltonian

$$H=J\sum_{\langle i,j\rangle}\sigma_i\sigma_j$$

The magnetization per spin is defined as

$$M=\frac{1}{N}\sum_i \sigma_i$$

Where $$N$$ is the number of spins. Since $$H$$ is symmetric under the flipping of every spins, it seems to me that for all $$i$$, $$\langle\sigma_i\rangle=\langle -\sigma_i\rangle=0$$. Put differently, for every configuration with magnetization $$M$$ there is a configuration with magnetization $$-M$$ with equal probability, hence $$\langle M \rangle=0$$. Wikipedia claims that this is a "bogus" argument that only works on a finite volume, and that is fine by me, but it sure looks true for finite $$N$$.

If I understand correctly when doing finite size Monte Carlo simulations, this problem is solved by computing $$\langle |M|\rangle$$ instead of $$\langle M\rangle$$. But what if I have another observable $$\Phi$$ and I want to calculate the covariance

$$\langle M-\langle M\rangle,\Phi-\langle\Phi\rangle\rangle$$

$$\langle M-\langle |M|\rangle,\Phi-\langle|\Phi|\rangle\rangle$$

In what case should I swap an average for the average of the corresponding absolute value?

The argument is not "bogus", it simply need to be qualified. Indeed, it is at the heart of the phenomena of (Landau) phase transitions, which is built on spontaneous symmetry breaking. You have a system that has some symmetry, yet we claim that it undergoes phase transition which breaks this symmetry - what gives?! In formal terms, we have an analytic expression for the partition function $$Z = \sum \exp(-\beta E_n)$$ yet we are looking for a non-analytical behavior which characterizes a phase transition.

The answer to these problems is indeed that the analytical nature of $$Z$$ holds only for finite sums. At the thermodynamic limit we can get non analytical behavior. In terms of symmetry, under ergodicity we assume that the system travels among all possible configurations with probability proportional to $$\exp(-\beta E_n)$$, but if we now set the magnetization going on one direction, it is clear that there is a huge energy barrier to flipping all the spins in the other direction. When taking the ergodic assumptions we said that averaging over systems is like averaging over long times, but it is clear that for the system to flip all the spins, the time scale is exponentially small in the energy barrier (for a low enough temperature). Even for not very large systems that may be much larger than the age of the universe! So we potentially solved both problems - however it is not yet clear how to address it analytically. As straight-forward averaging over quantities using $$Z$$ for finite system will still assume ergodicity and the symmetry would be manifest.

So there are two ways. In Monte-Carlo simulations, the system picks one magnetization when lowering the temperature. So each simulation by itself would actually give a finite magnetization that will persist over rather long cycles, at low temperatures. Averaging over different simulations might be problematic as different runs will pick different directions, randomly, so averaging over absolute value should solve that. But you can check that your simulation actually works when you lower the temperature and it randomly choose one magnetization and just sticks with it. It's fun!

A more formal way, which is also an analytical approach to phase transition, is to add a small symmetry-breaking term to the Hamiltonian, and then take a limit of large-system-small-term. So we write $$H = J \sum \sigma_i \sigma_j - h \sum \sigma_i$$ and set $$h$$ to be very small. At the limit $$N\to\infty$$, $$h\to 0$$ such that $$Nh\to \infty$$, the phase transition will be evident. And again, this works even for not very large systems when doing averages using different numerical methods. For $$T>T_c$$ the effect of $$h$$ is small enough and the average magnetization is zero, while for $$T the small term $$h$$ helps the system pick the right direction of magnetization, but the absolute value of the magnetization is not meaningfully affected by it. Again - you can try it at home.

• Hi! Thanks for the answer, this clarifies a bunch of things. Can I conclude then that: (1) in most situation not taking the absolute value shouldn't be a problem, because the simulation will choose a side of the configuration space naturally but (2) taking the absolute value shouldn't change anything, if different simulations choose different sides then it solves the problem, otherwise it does nothing, hence it is better to always take the absolute value. Also (3) this is valid for any observable, not just magnetization Commented May 5, 2020 at 10:39
• Also, I'm using a code that uses Wolff cluster updates, and in each simulation if I compute the average magnetization it's actually 0. In this case, does taking the absolute value of any observable yield the value that it would have in the infinite volume case? Commented May 5, 2020 at 10:44
• I would say "yes" to 1 and 2 and cautious "yes" to 3. Regarding the Wolf cluster algorithm - the problem is that the algorithm circumvents the thermodynamic limit by allowing flipping of macroscopic blocks together. So indeed in this case you will need to take the absolute value of the magnetization over time, as in one time step it will have all spins pointing one direction and in the other all spins pointing in the other direction. The second method I suggested will solve this without needing to take absolute value, as the energy cost of flipping together $N$ spins would be $Nh$.
– user245141
Commented May 5, 2020 at 10:56
• thanks! I implemented the small magnetic field thing using the method described here. If I may ask a quick, trivial question: in the case where the whole cluster update is rejected due to the magnetic field, should I still measure the observables (which then will be equal to the ones in the previous iteration)? Commented May 12, 2020 at 7:45
• I don't know exactly how you implement the averaging, but I would say that you should take that into the average process, as the event of a rejected flip is a physical one, and when you average quantities over time a period in which they do not change should be accounted for. Otherwise you will miscalculate the stability of the system.
– user245141
Commented May 12, 2020 at 8:09

I think there are two different questions in this one.

• First of all, from the point of view of statistical averaging (i.e. from the point of view of the partition function), the states with the opposite magnetization are equally probable and therefore the average magnetization is always zero. However, these opposite magnetization states are separated by a very large energy barrier. Thus, once the system is stuck in one of them, it is very unlikely that it will be found in the other. In experiment or a simulation this is due to the finite observation/simulation time, whereas in the theory of phase transitions the same effect is achieved by taking thermodynamic limit or calculating partial statistical sums.

• The main reason for using $$|M|$$ in Monte Carlo simulations is that the simulation is typically run multiple times, e.g., for different values of temperature. In every particular instance the simulation is likely to get stuck in a minimum with a particular orientation of the magnetization, but this minimum is not always the same, which creates difficulties when, e.g., plotting the magnetization dependence on temperature. Thus is why the absolute value is used, which allows to see the magnetization being smoothly reduced to zero rather than jump between $$|M|$$ and $$-|M|$$. If you are coding MC, it takes little effort to keep track of both $$M$$ and $$|M|$$ (and also \$M^2 = |M|^2 for calculating the susceptibility).

I think there are some things to mention more clearly:

• This may be obvious, but Monte Carlo simulations do not give any time evolution of your system. If you look at the time series representing the states of your Marcov chain during simulation, it does not represent the real time evolution.
• In very special cases the Marcov chain DOES represent the time evolution of your system (or is at least similar to it). This is the case for the SINGLE-FLIP algorithm used in many Monte Carlo simulations of the ising model. The reason is that also the real system kind of acts in the same way: It flips spins one by one.

Now, to answer your question: The Monte Carlo simulation for the ising model gives you an average magnetization of zero. And this is TRUE also in reality.

I know this seems counterintuitive, but let’s think about it: Have you tried to run a Monte Carlo simulation for a system with more than 1000 sites? You will see that the autocorrelation time gets huge. If you go to macroscopic systems the autocorrelation time tends to be infinity. And as we pointed out above: The simulation (in this case) represents the time evolution of the real system. You will never see a macroscopic ferromagnet switching sign, but only because there are so many sites (the autocorrelation time is large). In fact, if you wait infinitely long, the magnet will switch polarity and in the average the magnetization is zero (as the simulation correctly predicts).

A famous example for this is superparamagnetism. There, you have very tiny ferromagnets and therefore short autocorrelation times. Also in the real world you can observe these small ferromagnets switching polarity. In fact, this is a problem for producing magnetic hard discs because it limits the density of magnetic sites on the disk.

Edit: To clarify the „bogus“ argument of taking an absolute value: Yes this is bogus. You can not simply introduce an absolute value to get the result you want to have. For small systems and large systems the simulation gives zero average magnetization (and this is correct in reality). Only for infinitely large systems the magnetization is not zero (a kind of spontaneous symmetry breaking). However, if you want to calculate the magnetization of a macroscopic ferromagnet at a point in time, you have to take time evolution into account. The macroscopic magnet does stay for a very long time in one of the magnetizations and does not switch between then a lot (as mentioned above).