A Conventional MC Study of the $Q$ = 10 Potts Model

An impressive example of the power of the WL method appears in the latter of the two original WL papers [23]. Here, the object of study is the Potts model [27]. This is an $L\times L$ lattice of $N$ sites (i.e., $N = L^2$) where each site can have a spin value, $q_i$, between 1 and $Q$. WL considered the ten-state model; $Q$ = 10. The Hamiltonian of the Potts model is

$$\mathscr{H} = -\sum_{\left<ij\right>}\delta\left(q_i,q_j\right)$$ (353)

where the sum is over all nearest neighbor pairs.

Under periodic boundary conditions, the ground state energy of the Potts model is given by

$$E_0 = -2L^2$$ (354)

and the density of states at $E_0$ (i.e., the degeneracy of the ground state) is easily seen to be

$$\Omega(E_0) = Q.$$ (355)

The $N$ energy levels are

$$-E_0,-E_0+8,-E_0+12,\dots,-4,0,4,\dots,E_0-12,E_0-8,E_0$$ (356)

The ten-state Potts model displays a first-order phase transition with a critical temperature of $T_c = 0.701232$ in the thermodynamic limit ( $N\rightarrow\infty$). The two coexisting phases are characterized by an average energy per particle of -0.965 and -1.671, respectively.

First, let's conduct normal NVT sampling of a ten-state, $L$ = 12 Potts model at temperatures well below, near, and well-above the critical temperature. (The dual purpose code wl-w.c can accomplish this.) At precisely the critical temperature, we should see two peaks of equal height in the energy histogram. Our trial move will consist of randomly selecting one of the 144 spins and then randomly selecting a spin value between 1 and 10. Below, I show the histograms of energy per spin after 10$^6$ cycles (one cycle = 144 attempted flips) for various temperatures, all begun from the same initial (randomly assigned) configuration which is initially in the high-energy regime (E/N $\approx$ -0.25). ($T$ = 0.70991 is the critical temperature, $T_c$, for the $L$ = 12, $Q$ = 10 Potts model as calculated by Wang and Landau.)

Histograms of energy per spin, $H(E/N)$, for the $L$ = 12 ten-state Potts model at temperatures 0.5, 0.7, 0.70991, and 1.0, 5.0, and 10.0, after 10$^6$ cycles using conventional Monte Carlo.

Notice that precisely pinpointing $T_c$ using a series of NVT MC simulations (if we assume that 10$^6$ cycles is enough - I have not yet claimed that it is) would be extremely difficult. For $T$ = 0.7, the peaks are not of equal height, but an increase of a mere 0.0991 brings them to the same height (the signature of a critical temperature). How would we know how to zoom in on $T$ = 0.70991? It is conceivable that we could embed NVT MC inside some nonlinear optimization routine whose objective function is a measure of relative peak heights, which must be minimized by allowing $T$ to vary. This could be automated and could in principle provide a very accurate $T_c$ after a finite number of runs. But we don't know ahead of time how many runs would be required, nor if our optimization scheme is efficient enough to allow us to find $T_c$ to some tolerable level of precision in a reasonable time.

Now, imagine that we could in some way compute the density of states, $\Omega (E)$, from a single simulation. We could then easily arrive at an estimate for $T_c$ by simply evaluating canonical energy histograms at various $T$, each constructed from $\Omega (E)$, until we find one for which the peak heights are the same. Consider:

$$H(E/N) \propto \Omega(E)\exp\left(-E/k_BT\right).$$ (357)

Can we compute $\Omega (E)$ from an NVT MC simulation? In principle, yes. We could rearrange Eq. 358:

$$H(E/N)\exp\left(+\beta E\right) \propto \Omega(E)$$ (358)

This certainly suggests that if we run NVT MC, tabulate a histogram of energy states, and the postmultiply it by the ``anti-Boltzmann'' factor $\exp\left(+\beta E\right)$, we will recover $\Omega (E)$. But the data from the 10$^6$-cycle MC runs shown in the above figure kills any hope of being able to do this in practice. You can see that none of the histograms cover the entire domain of accessible energy levels, so we cannot use any single histogram to produce $\Omega (E)$. You can also see that the very highest energy levels are not accessed even at extremely high temperatures. (Interestingly, a negative temperature resolves this part of the energy spectrum quite well; this indicates that the entropy of the Potts model for $E > -0.25$ decreases with increasing energy.) The histogram taken at the lowest temperature appears to cover a broad part of the domain, but most of it is covered with very poor statistical accuracy, as evidenced by the large fluctuations. The histograms taken near the critical temperature of $T$ = 0.70991 cover a relatively broad domain relatively well, so it is at least conceivable that we can produce part of $\Omega (E)$ from these histograms.

The figure below shows the application of Eq. 359 to determine partial densities of states for each histogram shown in the previous figure.

Partial densities of state computed from $\Omega (E) = H(E/N)exp(+E/k_BT)$ and scaled such that $\Omega (E_0) = 10$ for the $L$ = 12 ten-state Potts model at temperatures 0.5, 0.7, 0.70991, 1.0, 5.0, 10.0, and -1.0 after 10$^6$ cycles using conventional Monte Carlo. $\Omega (E)$ for $T$ = 1.0 is scaled to match $\Omega (E)$ for $T$ = 0.70991, and $\Omega (E)$'s at higher $T$'s (and $T$ = 1.0) are matched to those at neighboring lower $T$.

From this data, we can see the trace of the curve defining the true density of states. It appears to have a maximum at $E/M\approx -0.25$ of a whopping 10$^{145}$ states. We also see that the density of states decreases with energy for $E/N > -0.25$; energy levels higher than this are apparently sampled well in an MC run with negative temperature. Now, in order to piece the true $\Omega (E)$ together from these many runs, we took advantage of the fact that $\Omega(E_0) = Q$, and scaled the densities to have them match in the overlapping regions. Notice that the two partial densities of states for $T$ = 0.7 and 0.70991 agree (not surprising), but that the density of states for $T$ = 0.5 appears much too large. All told, it appears that the total curve is adequately represented by the anti-Boltzmann-treated histograms from just three 10$^6$-cycle NVT MC runs: $T \in \left\{0.7,1.0,-1.0\right\}$.

Now, I have to admit that I have cheated somewhat. I knew ahead of time that the density of states has a maximum, so I knew that negative temperature MC would resolve some part of it. I also knew that the critical temperature is 0.70991, but it was nice to see that number supported by standard MC. Let us now learn how to compute $\Omega$ from a single MC run using the Wang-Landau technique.