2026-10-01 · chapter Yang–Mills and the mass gapPart 5 Lattice gauge theorythe Monte Carlo method

The Monte Carlo method — make configurations that follow the weight, count the correlations, and give the error

The $\beta$ of this chapter is the coupling of the Wilson action (5-01), not an inverse temperature.

A lattice expectation value is an integral over all links (5-01 §06; $dU_\ell$ is the Haar measure).

A four-dimensional lattice of side 6 has 5,184 links, and $SU(2)$ is 3-dimensional per link, so it is an integral in 15,552 dimensions Computed. So we make configurations $U^{(n)}$ that follow the weight $e^{-S_W}/Z$ and take the sample mean.

The error decreases as $1/\sqrt N$ whatever the dimension (1-08) Known. Drawing uniformly does not work. At $\beta=2.3$ on $6^4$, the mean of the 7,776 plaquettes is 0.6026 in equilibrium, while uniform configurations scatter around 0 with standard deviation 0.00567, so the equilibrium value is 106 standard deviations away Computed.

We make a random rule that changes the configuration one step at a time (a Markov chain, with transition probability $P(x,y)$), so that after a long time it follows the weight $\pi$. A sufficient condition for the distribution not to change in one step (stationarity) is detailed balance.

Summing the right-hand side over $x$ gives the left. A chain satisfying the right-hand side is called reversible. If moreover one can get from anywhere to anywhere (irreducible), the time average of a finite-state chain converges to the expectation value (the ergodic theorem). If it is also free of periodicity (aperiodic), the distribution itself approaches $\pi$ (the convergence theorem) Known. For continuous state spaces too, under conditions such as the transitions having a positive density, there are a law of large numbers and a central limit theorem of the same form (Tierney 1994) Known.

The Metropolis method (Metropolis et al. 1953; the first application was the equation of state of two-dimensional hard spheres) proposes a candidate $y$ with a symmetric proposal $q(x,y)=q(y,x)$, accepts it with the following probability, and stays put if it is rejected.

If $\pi(y)\le\pi(x)$, then $\pi(x)q\cdot\pi(y)/\pi(x)=\pi(y)q\cdot1$, so (3) holds, and $Z$ is not needed. Example: two states, $\pi=(0.2,0.8)$, with the proposal "go to the other one", gives $P=\begin{pmatrix}0&1\\0.25&0.75\end{pmatrix}$ with eigenvalues $1,-0.25$, and the distribution from state 0 is $(0.1992,0.8008)$ after 5 steps Computed. The deviation decreases as powers of the eigenvalue of largest absolute value other than 1 (here $-0.25$), and this is the source of the autocorrelation time of §05.

In four dimensions 6 plaquettes contain the link $U_\ell$, and the part of the action that depends on $U_\ell$ has the following form.

A staple is the product of the three edges of a plaquette other than $U_\ell$ (5-01 §03). A unit that goes once through all links is called a sweep.

The Metropolis proposal is near the current value. If the step size $\delta$ is large, rejections increase; if it is small, the chain moves only a little. We look at a single $U(1)$ link (weight $e^{\kappa\cos\theta}$).

ComputedThe numbers in this figure were computed on this machine (the readings are from Table 1; the chain in the figure is an example run on the spot).

Widening $\delta$ from 0.1 to 2.5 halves the acceptance rate and makes $\tau_{\rm int}$ (§05) about 120 times smaller. A high acceptance rate does not indicate efficiency. At $\delta=0.1$ the mean is off from the exact value by 2.9 standard deviations; for a slow chain, the error estimate is itself uncertain.

The heat-bath method redraws a single link directly from the conditional distribution $\propto e^{\frac{\beta}{2}\mathrm{tr}(U_\ell A_\ell)}$ with the rest held fixed. A sum of $SU(2)$ matrices is a real multiple of an $SU(2)$ matrix (the quaternions of 1-04), so we can write $A_\ell=kV$ ($V\in SU(2)$); $X=U_\ell V$ also moves with the Haar measure, and the weight is $e^{\beta kx_0}$ ($x_0=\tfrac12\mathrm{tr}X$). $SU(2)$ is the three-sphere $x_0^2+|\vec x|^2=1$, and the Haar measure is the uniform measure on it. The slice at fixed $x_0$ is a two-sphere of radius $\sqrt{1-x_0^2}$, with area proportional to $1-x_0^2$. When $x_0$ moves by $dx_0$, the width (arc length) on the three-sphere is $dx_0/\sqrt{1-x_0^2}$. The product is $\sqrt{1-x_0^2}\,dx_0$, so

we draw $x_0$ from this, choose the remaining direction uniformly, and set $U_\ell=XV^{-1}$ Known. This is the method said to have been used in the $SU(2)$ calculation of Creutz (1980). Kennedy–Pendleton (1985) is an improvement that draws (6) quickly, and Cabibbo–Marinari (1982) is a method that updates $SU(N)$ one $SU(2)$ subgroup at a time.

Overrelaxation is a move that goes far without changing the action; for $SU(2)$ it is $U_\ell\mapsto V^{-1}U_\ell^{-1}V^{-1}$. With $X'=X^{-1}$ the action does not change, two applications return to the start, and it preserves the Haar measure, so it satisfies detailed balance. On its own it is not ergodic, so it is mixed with the heat bath (Adler 1981, Creutz 1987, Brown–Woch 1987). Note that a sweep in a fixed order is in general not reversible (the chain seen backwards in time is the one updated in reverse order), and only stationarity remains Known.

A chain starts either from a cold start, with all links set to the identity matrix, or from a hot start drawn uniformly, and the first stretch is discarded (thermalisation). For the plaquette at $6^4$, $\beta=2.3$ (equilibrium 0.6026, the heat bath of Table 2), the heat bath enters the equilibrium band from both sides in 20 sweeps, while the hot side of Metropolis lags, at 0.482 after 20 sweeps and 0.593 after 150 Computed. Agreement of the two starting points is only a necessary condition (5-06 §04; the mathematical counterpart is the mixing time of 8-03 §06).

Even in equilibrium the configurations are not independent. We define the autocorrelation function, the correlation of the time series at a separation of $t$ steps divided by the variance, and the integrated autocorrelation time

($\delta O=O-\langle O\rangle$, $N\gg\tau_{\rm int}$). The effectively independent samples number $N/(2\tau_{\rm int})$ (formulas (23) and (25) of Wolff 2004) Known. For the heat bath of Table 2 that is $9{,}500/(2\cdot5.08)\approx935$. The time $\tau_{\exp}$ of the slowest decay is set by $\lambda_\star$, the eigenvalue of the transition of largest absolute value other than 1, as $\tau_{\exp}=-1/\log\lambda_\star$ (in the example of §02, $\lambda_\star=0.25$).

On the mathematical side, slowness is measured by the Poincaré constant $C_P$. $C_P$ is the smallest constant such that $\mathrm{Var}\,f\le C_P\,\mathbb E|\nabla f|^2$ for every function $f$ (1-10 §07; for a discrete chain the right-hand side is replaced by half the mean square change in one step, $\tfrac12\sum_{x,y}\pi(x)P(x,y)(f(y)-f(x))^2$). Langevin dynamics is a continuous-time dynamics that adds noise to motion down the gradient of the action, and there $\tau_{\exp}=C_P$ holds with equality. For a discrete-time chain, the following three are different quantities Known (Levin–Peres–Wilmer §12.2, 8-03 §06).

• $C_P=1/\gamma$: $\gamma=1-\lambda_2$ is the gap between 1 and the second largest eigenvalue $\lambda_2$.
• The relaxation time $t_{\rm rel}=1/\gamma_\star$: $\gamma_\star=1-\lambda_\star$ is the absolute gap.
• $\tau_{\exp}=-1/\log\lambda_\star$.

When $\gamma_\star$ is small, the ratio of the last two is $1+O(\gamma_\star)$. The Polyakov loop in the last column of Table 2 is the quantity obtained by taking, at each spatial point, the trace of the product of links once around the time direction, and averaging over space (here we look at its absolute value). It is an example of a slowly moving quantity spread over the whole lattice.

The sample sizes are 9,500, 9,500 and 4,500. At the same sample size the error of Metropolis is 2.7 times that of the heat bath, consistent with $\sqrt{42/5.1}\approx2.9$. Comparing $\tau$ multiplied by the cost, overrelaxation loses for the plaquette, 8.3 against 5.1, and gains for the Polyakov loop, 3.5 against 6.0. What gains are slow global quantities Computed. For $SU(3)$ too, overrelaxation is reported to have shrunk the $\tau$ of the magnitude of the Polyakov loop from 28 sweeps to 9 (Brown–Woch 1987), and the tendency agrees Physics.

$\tau_{\rm int}$ is estimated from the data. Binning cuts the time series into blocks of length $B$ and treats the block means as independent. One reads off where the error levels off for $B\gg\tau_{\rm int}$ (Flyvbjerg–Petersen 1989). The jackknife gives the error of a non-linear quantity $f$ from the means $\bar b_{(i)}$ with one block left out at a time.

It is used for Creutz ratios and masses (5-02 §06, 5-08). For the heat-bath plaquette of Table 2, (bin error / naive error)² grows as 1.00, 2.85, 6.23, 10.19 for $B=1,5,20,100$, reaching the Γ-method $2\tau_{\rm int}=10.2$. For $B=200$–$1000$ (47 to 9 blocks) it fluctuates downwards, 7.4 to 7.7. With few blocks the error of the error matters (for the square of the ratio the relative error is about $\sqrt{2/(n-1)}$) Computed.

The Γ method estimates the sum (7) directly, cut off at a window $W$.

The right-hand side is formula (42) of Wolff (2004). It is derived by error propagation from the approximate formula for the error of the windowed sum given by Madras–Sokal (1988) (Wolff's formula (40)) Known. If $W$ is short there is bias; if long, noise matters. The $\tau_{\rm int}$ of this chapter use Wolff's (2004) automatic window ($S=1.5$). Wolff asks one to “verify by eye” that $\tau_{\rm int}(W)$ is flat. In the lower panel of Figure 2 (§07) one can check by eye that it is flat around the window $W$ of the vertical line.

As the continuum limit is approached the correlation length $\xi$ (in lattice units) grows (5-05 §01), and the $\tau$ of local updates grows too. Physics writes $\tau\propto\xi^{z}$ and calls $z$ the dynamical exponent (the review of Hohenberg–Halperin 1977, 8-03 §07) Physics. The lattice free field (1-08 §07) is where this can be checked by a formula.

In a checkerboard heat bath, let $M_e$ be the sum over even points and $M_o$ over odd. An odd point neighbours $2d$ even points, so $M_e'=rM_o+$ (noise), with $r=2d/(2d+m^2)$. For the zero momentum $M=M_e+M_o$

From $\cosh(1/\xi)=1+m^2/2$, $\xi\approx1/m$, so $\tau\approx d\,\xi^2$, and for this dynamics $z=2$ Computed. The intuition is a random walk: a local update takes $\xi^2$ steps to carry information a distance $\xi$. Precisely, the $r=1-m^2/(2d+m^2)$ of (11) approaches 1 as $m\to0$, so the fraction that shrinks in one sweep becomes $1-O(\xi^{-2})$.

ComputedThe numbers in this figure were computed on this machine (the curves compute (11) on the spot. If JavaScript is disabled, Table 3).

Each doubling of $\xi$ multiplies $\tau_{\rm int}(M)$ by about 4. The cost of a calculation of the same precision is estimated to grow roughly as $\xi^{d+z}$ (5-05) Physics. For Langevin dynamics (§05) the relaxation time is $1/m^2=C_P$ (8-03; in the free-field normalisation where the field values $\phi_x$ themselves are the coordinates in which gradients are measured; changing the choice of coordinates changes $C_P$ only by a constant factor). Both are proportional to $m^{-2}$ ($\propto\xi^2$). The $\tau$ of the heat bath is measured in sweeps and the $\tau$ of Langevin in Langevin time, so we do not compare the coefficients.

HMC (hybrid Monte Carlo, Duane et al. 1987) moves all links at once. Momenta $P_\ell$ (traceless Hermitian matrices; real linear combinations of the $T^a$ of 1-03) are given normal distributions, and

the classical mechanics is run with the leapfrog method and accepted with probability $\min(1,e^{-\Delta H})$. The leapfrog method can be run backwards in time and preserves phase-space volume, so even with step errors $e^{-S_W}$ is exactly the stationary distribution Known. Estimating the autocorrelation time of HMC by the Langevin $C_P$ is an assumption (the generators are different).

As the continuum limit is approached, the field splits into sectors of topological charge $Q$ (the number counting how the field is twisted, integer-valued for continuum fields; for its definition on the lattice see 3-05 §08), and small steps find it hard to move between them. Schaefer–Sommer–Virotta (2011) reported, for HMC in pure gauge theory, an “effective dynamical critical exponent of about 5” for $Q^2$, and held that Wilson loops can be decoupled from this slow mode Physics. Lüscher–Schaefer (2011) showed a method that lets topological charge flow in and out through open boundary conditions in the time direction Physics. As a metaphor it is crossing passes in a landscape with many valleys. Precisely, if the measure near the boundaries of sectors is small the Poincaré constant becomes large (8-03 §08) Known, but whether that is the cause of freezing has not been established.

On a finite lattice (1) is a probability measure on a compact space (5-01 §06). What can be said as theorems is the convergence of time averages by the ergodic theorem (for continuous state spaces, under the conditions of §02) and the central limit theorem. For a reversible chain, if the asymptotic variance is finite, a central limit theorem of the form (7) holds (Kipnis–Varadhan 1986) Known. What remains observation is the error bar. If there is a slow mode invisible at finite $N$, the error bar comes out too small ($\delta=0.1$ in Table 1, §08). Extrapolation from finite volume and finite $a$ (5-05 §07, 5-09) needs further assumptions.

The glueball masses and string tension of 5-08 are physical results that rely on that extrapolation Physics, not a proof of the existence and the gap of the problem statement (4-11). On the other hand $C_P$ is a mathematical object, and at strong coupling a Poincaré inequality for lattice Yang–Mills has been shown (Shen–Zhu–Zhu 2023, 5-04) Known. For four-dimensional Yang–Mills there is no theorem fixing the $\xi$-dependence of $C_P$ on the side $\xi\to\infty$.

What this chapter can say

Sources and reproduction

Next chapter to read: 5-08 Measuring glueballs and the string tension (reading masses and the string tension from configurations made by the methods of this chapter).

← 5-06 How U(1) differs from the non-abelian caseContents5-08 Measuring glueballs and the string tension →