text:
# Static quark potential and string tension for compact U(1) in (2+1) dimensions.

## I Introduction

Classical Monte Carlo simulations of the path integral in Euclidean lattice gauge theory have been very successful, and this is currently the preferred method for ab initio calculations in quantum chromodynamics (QCD) in the low energy regime. Monte Carlo approaches to the Hamiltonian version of QCD propounded by Kogut and Susskind have been less successful, however, and lag at least ten years behind the Euclidean calculations. Our aim in this paper is to see whether useful results can be obtained for the Hamiltonian version by using the standard Euclidean Monte Carlo methods for anisotropic lattices, and extrapolating to the Hamiltonian limit in which the time variable becomes continuous, i.e. the lattice spacing in the time direction goes to zero. The Hamiltonian version of lattice gauge theory is less popular than the Euclidean version, but is still worthy of study. It can provide a valuable check of the universality of the Euclidean results, and it allows the application of many techniques imported from quantum many-body theory and condensed matter physics, such as strong coupling expansions, the t-expansion, the coupled-cluster method and more recently the density matrix renormalization group (DMRG). None of these techniques has proved as useful as Monte Carlo in (3+1) dimensions; but in lower dimensions they are more competitive.

A number of Quantum Monte Carlo methods have been applied to Hamiltonian lattice gauge theory in the past, with somewhat mixed results. A “Projector Monte Carlo” approach using a strong coupling representation for the gauge fields and a Greens Function Monte Carlo (GFMC) approach , which makes use of a weak coupling representation of the gauge fields. However both the approaches run into difficulties in that the former approach runs into difficulties for non-Abelian models, and the later requires the use of a “trial wave function” to guide random walkers in the ensemble towards the preferred regions of configuration space. This introduces a variational element into the procedure, in that the results may exhibit a systematic dependence on the trial wave function . For this reason, we are forced to look yet again for an alternative approach.

As mentioned above our aim in this paper is to use standard Euclidean path integral Monte Carlo (PIMC) techniques for anisotropic lattices, and see whether useful results can be obtained in the Hamiltonian limit. Morningstar and Peardon showed some time ago that the use of anisotropic lattices can be advantageous in any case, particularly for the measurement of glueball masses. We use a number of their techniques in what follows.

As a first trial of this approach, we treat the U(1) gauge model in (2+1)D, which is one of the simplest models with dynamical gauge degrees of freedom, and has also been studied extensively by other means (see Section II). Path integral Monte Carlo (PIMC) methods were applied to this model a long time ago by Hey and collaborators, but the techniques used at that time were not very sophisticated, and the results were rather qualitative. Very little has been done since then using this approach on this particular model.

The rest of this paper is organised as follows. In section II we discuss the U(1) model in (2+1) dimensions in its lattice formulation, and outline some of the work done on it previously. The details of the simulations, including the generation of the gauge-field configurations, the construction of the Wilson loop operators, and the extraction of the potential and string tension estimates are described in section III. In section IV we present our main results for the mean plaquette, static quark potential, string tension. The static quark potential has not previously been exhibited for this model, as far as we are aware. Finally we make an extrapolation to the Hamiltonian limit, and comparisons are made with estimates obtained by other means in that limit. (Sec. 1V). We find that indeed the PIMC method can give better results than other methods, even in the Hamiltonian limit. Our conclusions are summarized in Sec. VI.

## II THE U(1) MODEL

Consider the isotropic Abelian U(1) lattice gauge theory in three dimensions. The theory is defined by the action

$$S=\beta \underset{r,\mu ,\nu }{}\text{ReTr}P_{\mu \nu }$$
(1)
where

$$\text{P}_{\mu \nu }(r)=[1ReTr\{\text{U}_\mu (r)\text{U}_\nu (r+\widehat{\mu })\text{U}_\mu ^{}(r+\widehat{\nu })\text{U}_\nu ^{}(r)\}]$$
(2)
is the plaquette variable given by the product of the link variables taken around an elementary plaquette. The link variable $`\text{U}_\mu (r)`$ is defined by

$$\text{U}_\mu (r)=\text{exp}[ieaA_\mu (r)]=\text{exp}[i\theta _\mu (r)]$$
(3)
where in the compact form of the model, $`\theta _\mu (r)(=eaA_\mu (r))[0,2\pi ]`$ represents the gauge field on the directed link $`rr+\widehat{\mu }`$. The parameter $`\beta `$ is related to the bare gauge coupling by

$$\beta =\frac{1}{g^2}$$
(4)
where $`g^2=ae^2`$, in (2+1) dimensions.

The lattice U(1) model in (2+1) dimensions has been studied by many authors, and possesses some important similarities with QCD (for a more extensive review, see for example ref. ). If one takes the “naive” continuum limit at a fixed energy scale, one regains the simple continuum theory of non-interacting photons; but if one renormalizes or rescales in the standard way so as to maintain the mass gap constant, then one obtains a confining theory of free massive bosons. Polyakov showed that a linear potential appears between static charges due to instantons in the lattice theory; and Göpfert and Mack proved that in the continuum limit the theory converses to a scalar free field theory of massive bosons. They found that in that limit the mass gap behaves as

$$am_D=\sqrt{\frac{8\pi ^2}{g^2}}\text{exp}\left(\frac{\pi ^2}{g^2}v(0)\right)$$
(5)
while the string tension is bounded by

$$a^2\sigma =c\sqrt{\frac{g^2}{2\pi ^2}}\text{exp}(\frac{\pi ^2}{g^2}v(0))$$
(6)
where $`v(0)`$ is the Coulomb potential at zero distance, and has a value in lattice units

$$v(0)=0.2527$$
(7)
for the isotropic case. They argue that eq. (6) represents the true asymptotic behaviour of the string tension, where the constant $`c`$ is equal to $`8`$ in classical approximation. The theory has a non-vanishing string tension for arbitrary large $`\beta `$, similiar to behaviour expected for the string tension in non-abelian lattice gauge theories in 4-dimensions.

In the small-$`\beta `$ limit the string tension behaves as:

$$a^2\sigma =\text{ln}\beta /2+\mathrm{}.$$
(8)
and the large-$`\beta `$ behaviour is given by eq. (6).

Similar behaviour will apply in the anisotropic case. Generalizing discussions by Banks et al and Ben-Menahem, we find (see the appendix A) that the exponential factor takes exactly the same form in the anisotropic case, with a Coulomb potential $`v(0)`$ changing to $`0.3214`$ in the Hamiltonian limit.

For an anisotropic lattice, the gauge action becomes

$$S=\beta _s\underset{r,i<j}{}P_{ij}(r)+\underset{r,i}{}\beta _tP_{it}(r)$$
(9)
where $`P_{ij}`$ and $`P_{it}`$ are the spatial and temporal plaquette variables respectively. In the classical limit

$`\beta _s`$ $`=`$ $`{\displaystyle \frac{a_t}{e^2a_s^2}}={\displaystyle \frac{1}{g^2}}\mathrm{\Delta }\tau `$ (10)
$`\beta _t`$ $`=`$ $`{\displaystyle \frac{1}{e^2a_t}}={\displaystyle \frac{1}{g^2}}{\displaystyle \frac{1}{\mathrm{\Delta }\tau }}`$ (11)
where $`\mathrm{\Delta }\tau =a_t/a_s`$ is the anisotropy parameter, $`a_s`$ is the lattice spacing in the space direction, and $`a_t`$ is the temporal spacing. The above action can be written as

$`S`$ $`=`$ $`\beta \left[\mathrm{\Delta }\tau {\displaystyle \underset{r}{}}{\displaystyle \underset{i<j}{}}\left(1\mathrm{cos}\theta _{ij}(r)\right)+{\displaystyle \frac{1}{\mathrm{\Delta }\tau }}{\displaystyle \underset{r,i}{}}\left(1\mathrm{cos}\theta _{it}(r)\right)\right]`$ (12)
In the limit $`\mathrm{\Delta }\tau 0`$, the time variable becomes continuous, and we obtain the Hamiltonian limit of the model (modulo a Wick rotation back to Minkowski space). The Hamiltonian version of the model has been studied by many methods: some recent studies include series expansions , finite-lattice techniques, the t-expansion , and coupled-cluster techniques, as well as Quantum Monte Carlo methods. Quite accurate estimates have been obtained for the string tension and the mass gaps, which can be used as comparison for our present results. The finite-size scaling properties of the model can be predicted using an effective Lagrangian approach combined with a weak-coupling expansion , and the predictions agree very well with finite-lattice data .

Some analytic results for the model in the weak-coupling approximation appear in the appendix.

## III METHODS

### A Path Integral Monte Carlo algorithm

We perform standard path integral Monte Carlo simulations on a finite lattice of size $`N_s^2\times N_\tau `$, where $`N_s`$ is the number of lattice sites in the space direction and $`N_\tau `$ in the temporal direction, with spacing ratio $`\mathrm{\Delta }\tau =a_t/a_s`$. By varying $`\mathrm{\Delta }\tau `$ it is possible to change $`a_t`$, while keeping the spacing in the spatial direction fixed. The simulations were performed on lattices with $`N_s=16`$ sites in each of the two spatial directions and $`N_t=1664`$ in the temporal direction for a range of couplings $`\beta =13`$.

The ensembles of field configurations were generated by using a Metropolis algorithm. Starting from an arbitrary initial configuration of link angles, we successively update link angles $`\theta _\mu (\stackrel{}{r},\tau )`$ at positions $`(\stackrel{}{r},\tau )`$ which are chosen randomly each time. We propose a change $`\mathrm{\Delta }\theta `$ to this link angle, which is randomly drawn from a uniform distribution on $`[\mathrm{\Delta },\mathrm{\Delta }]`$, where $`\mathrm{\Delta }`$ is adjusted for each set of parameters to give an acceptable “hit rate” around 70-80%. The change is accepted or rejected according to the standard Metropolis procedure.

For high anisotropy ($`\mathrm{\Delta }\tau <<1`$), any change in a time-like plaquette will produce a large change in the action, whereas changes to the space-like plaquettes will cause a much smaller change in the action. This makes the system very “stiff” against variations in the time-like plaquettes, and therefore very slow to equilibrate, with long autocorrelation times. To alleviate this problem, we used a Fourier update procedure. Here proposed changes are made to space-like links which are designed to alter space-like plaquette values much more than the time-like plaquette values. At randomly chosen intervals and random locations, we propose a non-local change $`\mathrm{\Delta }\theta (\stackrel{}{r},\tau )=X\mathrm{sin}k(\tau \tau _0)`$ on a “ladder” of space-like links extending half a wavelength $`(\lambda =2\pi /k)`$ in the time direction, where both $`k`$ and $`X`$ are randomly chosen at each update from uniform distributions in $`[0,\pi /2]`$. We replaced approximately $`30\%`$ of the ordinary Metropolis updates with Fourier updates for anisotropy $`\mathrm{\Delta }\tau >0.444`$ and $`50\%`$ for highly anisotropic cases ($`\mathrm{\Delta }\tau <0.444`$). These moves satisfy the requirements of detailed balance and ergodicity for the algorithm.

A single sweep involves attempting $`N`$ changes to randomly chosen links of the lattice, where $`N(=3N_tN_s^2)`$ is the total number of links on the lattice. The first several thousand sweeps are discarded to allow the system to relax to equilibrium. Fig. 1 shows measurements of the mean plaquette for $`\beta =2.0`$ and $`\mathrm{\Delta }\tau =1.0`$, and it can be seen that equilibrium is reached after about 50,000 thousand sweeps, with the measurements fluctuating about the equilibrum value thereafter. For highly anisotropic cases the system was much slower to equilibrate, despite the Fourier acceleration, and in the worst case the equilibration time was of order 100,000 thousand sweeps.

After discarding the initial sweeps, the configurations were stored every 250 sweeps thereafter for later analysis. Ensembles of about 1,000 configurations were stored to measure the static quark potential and glueball masses at each $`\beta `$ for $`\mathrm{\Delta }\tau 0.4`$, and 1,400 configurations for $`\mathrm{\Delta }\tau 0.333`$. Measurements made on these stored configurations were grouped into 5 blocks, and then the mean and standard deviation of the final quantities were estimated by averaging over the “block averages”, treated as independent measurements. Each block average thus comprised 50,000 - 70,000 sweeps.

### B Interquark Potential

The static quark-antiquark potential, $`V(𝐫)`$ for various spatial separations $`𝐫`$ is extracted from the expectation values of the Wilson loops. The timelike Wilson loops are expected to behave as:

$$W(𝐫,\tau )Z(𝐫)\text{exp}[\tau V(𝐫)]+\text{(excited state contributions)}$$
(13)
We have averaged only over loops $`(𝐱_\mathrm{𝟎},\tau _0)(𝐱_\mathrm{𝟎}+𝐫,\tau _0)(𝐱_\mathrm{𝟎}+𝐫,\tau _0+\tau )(𝐱_\mathrm{𝟎},\tau _0+\tau )(𝐱_\mathrm{𝟎},\tau _0)`$ which follow either two sides of a rectangle between $`𝐱_\mathrm{𝟎}`$ and $`𝐱_\mathrm{𝟎}+𝐫`$ or a single-step ‘staircase’ route, to estimate $`W(𝐫,\tau )`$. To suppress the excited state contributions, a simple APE smearing technique was used on the space-like variables. In this technique an iterative smearing procedure is used to construct Wilson loop (and glueball) operators with a very high degree of overlap with the lowest-lying state. In our single-link smoothing procedure, we replace every space-like link variable by

$$U_iP\left[\alpha U_i+\frac{(1\alpha )}{2}\underset{s}{}U_s\right]$$
(14)
where the sum over “s” refers to the “staples”, or 3-link paths bracketing the given link on either side in the spatial plane, and P denotes a projection onto the group U(1), achieved by renormalizing the magnitude to unity. We used a smearing parameter $`\alpha =0.7`$ and up to 10 iterations of the smearing process.

To further reduce the statistical errors, the timelike Wilson loops were constructed from “ thermally averaged” temporal links. That is, the temporal links $`U_t`$ in each Wilson loop were replaced by their thermal averages

$$\overline{U_t}=𝑑UU\mathrm{exp}(S[U])/𝑑U\mathrm{exp}(S[U])$$
(15)
where the integration is done over the one link only, and depends on the neighbouring links. For the U(1) model, the result can easily be computed in terms of Bessel functions involving the “staples” adjacent to the link in question. This was done for all temporal links except those adjacent to the spatial legs of the loop, which are not “independent”. The procedure has a dramatic effect in reducing the statistical noise, by up to an order of magnitude 10, which corresponds to a factor 100 in the statistics.

Fig. 2 shows an example of the exponential decrease of the $`4\times 1`$Wilson loop $`W(\stackrel{}{r},\tau )`$ with $`\tau `$. Estimates of the potential can now be found from the ratios of successive loop values:

$$V(𝐫)=\frac{1}{\mathrm{\Delta }\tau }\text{ln}\left[\frac{W(𝐫,\tau +\mathrm{\Delta }\tau )}{W(𝐫,\tau )}\right]$$
(16)
This ratio is expected to be independent of $`\tau `$ for $`\tau >>0`$. Fig.3 shows the effective potential as a function of Euclidean time, $`\tau `$, obtained from $`\mathrm{\Delta }\tau =1.0`$ at $`\beta =2.0`$ for $`r=4`$ simulation. We find that ratio reaches its plateau behavior at very early Euclidean times, which reflects a good optimized smearing. We used $`\tau =1`$ to make our estimates of $`V(𝐫)`$.

## IV Simulation results at finite temperature

Simulations were carried out for lattices of $`N_s^2\times N_t`$ sites, with $`N_s=16`$ and $`N_t`$ ranging from 16 to 48 sites, with periodic boundary conditions. Each run involved 250,000 sweeps (350,000 for high anisotropy) of the lattice, with 60,000 sweeps (100,000 high anisotropy) discarded to allow for equilibrum, and configurations recorded every 250 sweeps thereafter. Coupling values from $`\beta =1.0`$ to 3.0 were explored at anisotropies $`\mathrm{\Delta }\tau `$ ranging from 1 to 1/3. We fixed $`\mathrm{\Delta }\tau =16/N_t`$ in the first pass, so that the lattice size $`T=N_ta_t`$ in the time direction remains constant, corresponding to a fixed (low) temperature.

### A Mean Plaquette

The mean plaquette, or equivalently the $`1\times 1`$ Wilson loop, was measured for various $`16\times 16\times N_t`$ (with $`N_t`$ = 16 up to 48) sizes and for a variety of $`\beta `$ values, ranging from strong coupling to weak coupling regime. Fig. 4 shows the behaviour of the mean spatial plaquette $`<P>`$ for different $`\beta `$, at fixed $`\mathrm{\Delta }\tau =1`$ (isotropic case). At small coupling (large $`\beta `$) the plaquette expectation value is expected to approach 1, and conversely it approaches 0 in the high coupling (small $`\beta `$) limit. The plot also shows the $`O(\beta ^4)`$ strong-coupling  and $`O(1/\beta ^5)`$ weak coupling  expansion. The dashed curve represents the weak-coupling series and strong-coupling series is represented by the solid curve. As can be seen that the data follow closely the weak coupling expansion down to $`\beta =2.0`$ and strong coupling expansion upto $`\beta =1.70`$ and the results agree very well. The variation of $`<P>`$ with coupling is extremely smooth, with no sign of phase transition, as we should expect. The cross-over from strong to weak coupling seems, from the Monte Carlo data, to take place quite raidly in the region $`\beta 1.82.0`$. Horsley and Wolff  investigates the effect of finite size lattice in their weak-coupling expansion calculations. They find that such effects enter at order $`1/\beta ^2`$ as a correction of order $`1/l^D`$, where $`l`$ is the lattice size and $`D`$ is the number of dimensions. Thus for U(1) in (2+1) and $`l>10`$, finite size effects should be essentially negligible.
Fig. 5 shows a plot of our estimates of $`<P>`$ as a function of anisotropy $`\mathrm{\Delta }\tau ^2`$ for the case $`\beta =\sqrt{2}`$. We would like to make contact with previous Hamiltonian studies by showing that the mean plaquette value approaches to previously known values in the Hamiltonian limit $`\mathrm{\Delta }\tau 0`$. The extrapolation was performed using a simple cubic fit in powers of $`\mathrm{\Delta }\tau ^2`$. In this limit our results agree very well with the Hamiltonian estimate obtained by Hamer et al. The estimates of the mean plaquette $`<P>`$ in the Hamiltonian limit are graphed as a function of $`\beta `$ in Fig.6 and compared with the series estimates at strong and weak coupling , and Greens function Monte Carlo results. At large $`\beta `$, the data is in good agreement with predictions of weak-coupling perturbation theory. At small $`\beta `$, the data follow closely the strong-coupling and Greens function Monte Carlo estimates, except at very small $`\beta `$. We attribute the discrepnacy to the finite-temperature effects and a better agreement is expected at zero-temperature .

### B Static quark potential and string tension

The effective potential is obtained from the ratio of successive Wilson loop values:

$$V_\tau (𝐫)=\text{ln}\left(\frac{W(𝐫,\tau +\mathrm{\Delta }\tau )}{W(𝐫,\tau )}\right)$$
(17)
Fig. 7 shows a graph of the static quark potential V(R) as a function of radius R at $`\beta =2.0`$ and $`\mathrm{\Delta }\tau =1.0`$. To extract the string tension, the curve is well fitted by a form

$$V(R)=a+b\mathrm{ln}R+\sigma R,$$
(18)
including a logarithmic Coulomb term as expected for classical QED in (2+1) dimensions which dominates the behaviour at small distances, and a linear term as predicted by Polyakov and Göpfert and Mack dominating the behaviour at large distances. The linear behaviour at large distances is very clear, but the data do not extend to very small distances, so there is no real test of the presumed logarithmic behaviour in this regime.

Figs. 8 shows the behaviour of the fitted value of the string tension $`K(=a^2\sigma )`$ as a function of $`\beta `$ for the isotropic case ( $`\mathrm{\Delta }\tau =1`$), plotted in comparison to results predicted for string tension using the large-$`\beta `$ approximation (eq. 6). We also include the strong-coupling values and show the leading order strong-coupling prediction (eq. 7). As can be seen that our estimates of the string tension, in the weak-coupling region are much larger than the Villain prediction, though slopes are more or less well determined. However one does not expect agreement of precise magnitude, since our estimates are obtained using Wilson action, and not the Villian approximation. But one does expect that the slope of the string tension in both cases should ultimately be the same as one approaches the continuum limit. Mack and Göpfert  argued that eq. (6) represents the true asymptotic behaviour of the string tension, apart from the value of the renormalization constant $`c`$. Such a behaviour was obtained by a semi-classical treatment of the theory with the effective action without the correction terms. Thus it is possible that eq. (6), including the normalization constant, becomes exact as $`1/\beta 0`$, if one estimates the actual value of $`c`$. We estimated the actual value of $`c`$ by graphing the $`\text{e}^{2.494\beta }\times 4.443K`$ against $`1/\beta `$ and extrapolating to $`1/\beta 0`$. Using estimated value of the renormalization constant, $`c=44`$, we plotted the asymptotic form of the string tension in weak-coupling limit in comparison to our estimates. The result agrees very well with our estimates.

Fig. 9 shows the behaviour of the string tension ($`K`$) as a function of the anisotropy $`\mathrm{\Delta }\tau ^2`$, for fixed coupling $`\beta =\sqrt{2}`$. An extrapolation to the Hamiltonian limit $`\mathrm{\Delta }\tau 0`$ is performed by a simple cubic fit. Again the extrapolation agrees rather well with earlier Hamiltonian estimates. Note that this quantity depends rather strongly on $`\mathrm{\Delta }\tau `$: there is a factor of three difference between the values at $`\mathrm{\Delta }\tau =0`$ and $`\mathrm{\Delta }\tau =1`$. Extrapolating to $`\mathrm{\Delta }\tau 0`$, estimates of the string tension in the Hamiltonian limit, for various $`\beta `$ values are obtained. The plot of these estimates, in comparison to the large-$`\beta `$ prediction of eq. (6) is shown in Fig. 10. As can be seen, the weak-coupling results, with estimated value for constant coefficient, agree very well with our estimates. The comparison with the theoretical prediction using Villian form is not bad either, especially in the large-$`\beta `$ region. In this region, upto $`\beta 1.8`$, the results agree well within errors and slopes are very well defined.

Fig. 11 shows the plot of Hamiltonian estimates as a function of $`\beta `$, plotted in comparison to the previous estimates obtained from strong-coupling , Greens function Monte Carlo  and finite-size scaling  The agreement with the series expansion and GFMC results is once again excellent. The estimates of finite-size scaling describe the data at $`\beta =2.0`$.

## V Summary and conclusion

In this paper, we used numerical simulation techniques on anisotropic lattices for U(1) model in (2+1) dimensions to obtain a clear picture of the confinement and improve our knowledge of the string tension. Our first aim was to implement Path Integral Monte Carlo method to the U(1) model, which is one of the simplest models and has been studied by various other methods. It seems to work extremely well for the mean plaquette in the strong-coupling and weak-coupling regimes. The simulation results were extrapolated to $`\mathrm{\Delta }\tau 0`$, and Hamiltonian estimates for mean plaquette were obtained. A comparison was made with the previous estimates obtained in this limit and agreement was excellant. The second aim was the get a better picture of the confinement and the numerical study of the asymptotic behaviour of the string tension in the large-$`\beta `$ limit. It was shown that at large separations the U(1) theory in (2+1) dimensions shows a linear behaviour as predicted by Polyakov. However, the presumed logarithmic term could not be tested as our data do not extend to very small separations. The estimates of string tension have confirmed very well the weak-coupling predictions, provided that one takes the estimated value of the constant coefficient $`c`$ rather than its value in the classical approximation using effective action without correction terms. For $`\beta <2.0`$, the data is in good agreement with $`O(\beta ^2)`$ strong-coupling predictions. Once again our Hamiltonian estimates for the string tension fit very well with the series and GFMC results.

This seems to suggest that PIMC does work better than other methods, as stated by Ceperley. The first clear picture of the static quark potential has been obtained, showing clear evidence of the linear confining potential at large distance predicted by Polyakov. Our results from the PIMC do appear to be in agreement with results of the series expansions.

## Acknowledgments

This work was supported by the Australian Research Council. We are grateful for access to the computing facilities of the Australian Centre for Advanced Computing and Communications (ac3) and the Australian Partnership for Advanced Computing (APAC).

## A Weak-coupling prediction using Villian approximation

Here we extend the discussion of the weak-coupling behaviour of the mass gap to the anisotropic case, following the discussions of Banks, Myerson and Kogut and Ben-Menahem. This can be probed by means of the expectation value of the Wilson loop, $`<W(C)>`$, which for pure gauge theory is given by

$$<W(C)>=\frac{Z(J)}{Z(0)}$$
(A1)
where

$$Z(J)=\underset{r,\mu }{}\left\{𝑑\theta _\mu (r)\right\}\text{exp}\left[Si\underset{r,\mu }{}J_\mu (r)\theta _\mu (r)\right]$$
(A2)
and $`J_\mu (r)`$ is a vector field associated with a rectangular contour C of width R and length T. The extra piece of action associated with the external charge is

$$S_J=d^3rJ_\mu (r)\theta _\mu (r)$$

For a point charge, $`S_J`$ is given by the path integral of the gauge field along the world line of the charge. Therefore we choose

$$J_\mu (r)=\{\begin{array}{cc}1\hfill & \text{if }rr+\widehat{\mu }\text{ is on C}\hfill \\ 1\hfill & \text{if }r+\widehat{\mu }\text{ is on C}\hfill \\ 0\hfill & \text{otherwise}\hfill \end{array}$$
(A3)
and this current is conserved: $`\mathrm{\Delta }_\mu J_\mu (r)=0`$. We will assume a spacelike Wilson loop for this discussion.

It is convenient to work in the temporal gauge

$$\theta _t(\stackrel{}{r})=0$$
(A4)
Then one can separate the angular variable $`\theta _i(r)`$ into longitudinal and transverse parts,

$$\theta _i(r)=\mathrm{}_i\rho (r)ϵ_{ij}\frac{1}{\mathrm{}_1^2+\mathrm{}_2^2}\mathrm{}_j\theta (r)$$
(A5)
whence it follows that $`\theta _{ij}(r)=\theta (r)`$. In the temporal gauge eq. (A2) can be written as

$$Z(J)=\underset{r}{}_{\mathrm{}}^{\mathrm{}}\{d\theta (r)\}\text{exp}\left[\beta \mathrm{\Delta }\tau \underset{r}{}\left(1\mathrm{cos}\theta (r)\right)\frac{\beta }{2\mathrm{\Delta }\tau }\underset{r}{}\theta (r)\frac{\overline{\mathrm{}}_t\mathrm{}_t}{\overline{\mathrm{}}_i\mathrm{}_i}\theta (r)+i\underset{r}{}J_iϵ_{ij}\frac{1}{\mathrm{}_k^2}\mathrm{}_j\theta (r)\right]$$
(A6)
with $`\mathrm{}_i^2=\mathrm{}_1^2+\mathrm{}_2^2`$. $`\overline{\mathrm{}}`$ and $`\mathrm{}`$ are defined are the lattice difference operators. To calculate the path integral we replace the potential term by its Villain form

$$\text{e}^{\beta (1\mathrm{cos}\theta (r))}\underset{l=\mathrm{}}{\overset{\mathrm{}}{}}\text{e}^{il\theta }\text{e}^\beta I_l(\beta ).$$
(A7)
where $`I_ł(\beta )`$ is the modified Bessel function and $`\text{exp}(il\theta )`$ is the character of the plaquette variable in the $`l`$th irreducible representation of the compact U(1) group.

For large $`\beta `$ the modified Bessel function can be approximated by

$$e^\beta I_l(\beta )\frac{1}{\sqrt{2\pi \beta }}\text{e}^{l^2/2\beta }$$
(A8)
The Gaussian integration over $`\theta (r)`$ reduces the path integral to

$`Z(J)`$ $``$ $`{\displaystyle \underset{r}{}}{\displaystyle \underset{l(r)=\mathrm{}}{\overset{\mathrm{}}{}}}\text{exp}[{\displaystyle \frac{1}{2\beta \mathrm{\Delta }\tau }}{\displaystyle \underset{r}{}}l(r)(\mathrm{\Delta }\tau ^2{\displaystyle \frac{\overline{\mathrm{}}_i\mathrm{}_i}{\overline{\mathrm{}}_t\mathrm{}_t}}+1)l(r){\displaystyle \frac{\mathrm{\Delta }\tau }{\beta }}{\displaystyle \underset{r}{}}l(r){\displaystyle \frac{1}{\overline{\mathrm{}}_t\mathrm{}_t}}ϵ_{ij}J_i\mathrm{}_j`$ (A10)
$`{\displaystyle \frac{\mathrm{\Delta }\tau }{2\beta }}{\displaystyle \underset{r}{}}J_i(r){\displaystyle \frac{1}{\overline{\mathrm{}}_t\mathrm{}_t}}J_i(r)]`$
The convergence of the sum over $`l`$ can be improved by making use of the Poisson summation formula,

$$\underset{l=\mathrm{}}{\overset{\mathrm{}}{}}h(l)=\underset{m=\mathrm{}}{\overset{\mathrm{}}{}}𝑑\varphi h(\varphi )\text{e}^{2i\pi mr}$$
(A11)
where $`h`$ is an arbitrary function.

Substituting $`l(r)=\mathrm{}_tL(r)`$ and using eq. (A8), Eq. (A10) becomes

$`Z(J)`$ $``$ $`{\displaystyle \underset{r}{}}{\displaystyle _{\mathrm{}}^{\mathrm{}}}d\varphi (r){\displaystyle \underset{m(r)=\mathrm{}}{\overset{\mathrm{}}{}}}[\text{exp}\{{\displaystyle \frac{1}{2\beta \mathrm{\Delta }\tau }}{\displaystyle \underset{r}{}}\varphi (r)(\mathrm{\Delta }\tau ^2\overline{\mathrm{}}_i\mathrm{}_i+\overline{\mathrm{}}_t\mathrm{}_t)\varphi (r)`$ (A14)
$`+{\displaystyle \frac{\mathrm{\Delta }\tau }{\beta }}{\displaystyle \underset{r}{}}\varphi (r){\displaystyle \frac{1}{\mathrm{}_t}}ϵ_{ij}J_i\mathrm{}_j+2i\pi {\displaystyle \underset{r}{}}m(r)\varphi (r)\}]`$
$`\times \text{exp}\left[{\displaystyle \frac{\mathrm{\Delta }\tau }{2\beta }}{\displaystyle \underset{r}{}}J_i(r){\displaystyle \frac{1}{\overline{\mathrm{}}_t\mathrm{}_t}}J_i(r)\right]`$
Finally doing the Gaussian integral over $`\varphi `$ yields;

$`Z(J)`$ $``$ $`{\displaystyle \underset{r}{}}{\displaystyle \underset{m(r)=\mathrm{}}{\overset{\mathrm{}}{}}}\text{exp}[2\pi ^2\beta \mathrm{\Delta }\tau {\displaystyle \underset{r}{}}m(r)\left({\displaystyle \frac{1}{\mathrm{\Delta }\tau ^2\overline{\mathrm{}}_i\mathrm{}_i+\overline{\mathrm{}}_t\mathrm{}_t}}\right)m(r)`$ (A17)
$`+2\pi i\mathrm{\Delta }\tau ^2{\displaystyle \underset{r}{}}m(r)\left({\displaystyle \frac{1}{\mathrm{\Delta }\tau ^2\overline{\mathrm{}}_i\mathrm{}_i+\overline{\mathrm{}}_t\mathrm{}_t}}\right){\displaystyle \frac{1}{\mathrm{}_t}}ϵ_{ij}J_i\mathrm{}_j]`$
$`\times \text{exp}\left[{\displaystyle \frac{\mathrm{\Delta }\tau }{2\beta }}{\displaystyle \underset{r,i}{}}J_i{\displaystyle \frac{1}{\mathrm{\Delta }\tau ^2\overline{\mathrm{}}_i\mathrm{}_i+\overline{\mathrm{}}_t\mathrm{}_t}}J_i(r)\right]`$
The partition function, $`Z(J=0)`$, is given by

$$Z(J=0)=\underset{r}{}\underset{m(r)=\mathrm{}}{\overset{\mathrm{}}{}}\text{exp}\left[2\pi ^2\beta \underset{r,r^{}}{}m(r)G(rr^{})m(r^{})\right]$$
(A18)
where

$$G(rr^{})=\frac{\mathrm{\Delta }\tau \delta _{rr^{}}}{\mathrm{\Delta }\tau ^2\overline{\mathrm{}}_i\mathrm{}_i+\overline{\mathrm{}}_t\mathrm{}_t}$$
(A19)
is a three-dimensional lattice massless propagator.

The propagator can be expressed as a momentum integral:

$$G(r,t)=_\pi ^\pi \frac{dk_0}{2\pi }_\pi ^\pi \frac{dk_1}{2\pi }_\pi ^\pi \frac{dk_2}{2\pi }\frac{\mathrm{\Delta }\tau \text{e}^{ik.r+k_0t)}}{4\mathrm{sin}^2k_0/2+\mathrm{\Delta }\tau ^2\left(42\mathrm{cos}k_12\mathrm{cos}k_2\right)}$$
(A20)
For the isotropic lattice ($`\mathrm{\Delta }\tau =1`$), $`G(k=0)=0.2527`$; while in the limit $`\mathrm{\Delta }\tau \mathrm{}`$, $`G(k=0)`$ takes the value $`0.3214`$.
The resulting partition function may be written as

$$Z(J)=Z_{spinwave}Z_{ext}Z_{monopoles}$$
(A21)
where

$`Z_{spinwave}`$ $`=`$ $`{\displaystyle [d\chi ]\text{exp}\left[\frac{1}{\beta }\underset{r,i}{}\left(\mathrm{}_i\chi (r)\right)^2\right]}`$ (A22)
$`Z_{ext}`$ $`=`$ $`\text{exp}\left[{\displaystyle \frac{1}{2\beta }}{\displaystyle \underset{r,r^{},i}{}}J_i(r)G(rr^{})J_i(r^{})\right]`$ (A23)
$`Z_{monopole}`$ $`=`$ $`{\displaystyle \underset{r}{}}{\displaystyle \underset{m(r)=\mathrm{}}{\overset{\mathrm{}}{}}}\text{exp}\left[2\pi ^2\beta {\displaystyle \underset{r,r^{}}{}}m(r)G(rr^{})m(r^{})+2i\pi {\displaystyle \underset{r,r^{}}{}}G(rr^{})\eta (r)m(r)\right]`$ (A24)
The first factor in eq. (A20), the “spin-wave” piece, contributes the ordinary Coulomb part to the potential between the external charges. The rest of the eq. (A20) is the partition sum of the Coulomb gas of magnetic monopoles interacting with a stationary current loo via $`G(rr^{^{}})`$. This three-dimensional Coulomb gas is always in a plasma phase. The only thing that changes with $`\beta `$ is the density of the monopoles which decays exponentially to zero as $`\beta `$ increases. The monopole density generated by the external loop and coupling to the lattice monopoles m(r) is a dipole sheet along the surface spanning the contour C. Thus eq. (A20) represents a dipole sheet immersed in a monopole gas, which interacts directly with the dipole sheet. The gas becomes polarized, and the monopoles of opposite polarity accumulate on the two sides of the sheet, screening both the magnetic field away from the sheet and the monopole-monopole interaction, thus minimizing the energy of the system. This gives the correlation function a finite range, whose inverse is interpreted as a mass of the scalar field.
In the naive continuum limit, the Coulomb self-energy of the monopoles blows up and their density goes to zero exponentially. This causes the coefficient of the linear force law to go to zero and one regains the usual continuum theory, i.e., free electromagnetism in three dimensions; but if one renormalizes or rescales in the standard way so as to maintain the mass gap constant, then one obtains a confining theory of free massive bosons. The detailed analysis of Polyakov shows that the sheet cannot be screened completely and there is a finite energy left that binds the quarks through a potential which grows linearly with separation. Polyakov showed that for arbitrarily large finite $`\beta `$, theory has a mass gap and there appears a linear potential appears between static charges due to instantons in the lattice theory.
Splitting the exponent in eq. (A20) into pieces with $`r=r^{}`$ and $`rr^{}`$, the monopole partition sum can be written as

$`Z_{monopole}(J)`$ $`=`$ $`{\displaystyle \underset{r}{}}{\displaystyle }d\varphi (r)\text{exp}[{\displaystyle \frac{1}{8\pi ^2\beta }}{\displaystyle \underset{r}{}}\left(\mathrm{\Delta }_i\varphi (r)\right)^2{\displaystyle \underset{m=\mathrm{}}{\overset{\mathrm{}}{}}}2\pi ^2\beta G_0`$ (A26)
$`2i\pi {\displaystyle \underset{r}{}}\eta (r)m(r)+i{\displaystyle \underset{r}{}}\varphi (r)m(r)]`$
For large $`\beta `$, one can neglect all the monopoles with charge other than 0, $`\pm 1`$ in the sum, since the monopoles with a higher charge gets extra powers of $`\text{e}^\beta `$.  This argument becomes problematical when $`\mathrm{\Delta }\tau 0`$ . In this limit the monopoles of higher charges contribute to the action and hence are no longer suppressed. Applying the large-$`\beta `$ approximation to the monopole partition function, the above equation can be written, in the naive continuum limit, as

$$Z_{monopole}(J)D\varphi \text{exp}\left[\frac{g^2}{4\pi ^2}\left\{\frac{1}{2}(_i\varphi )^2M^2\mathrm{cos}\left(\varphi +2\pi \eta \right)\right\}d^3x\right]$$
(A27)
where the lattice photon mass $`M`$ is given by

$$M^2=\frac{8\pi ^2}{g^2a^3}\text{e}^{2\pi ^2G_0/g^2a}.$$
(A28)
At this point out analysis parallels the work of original work of Banks et al .
The the expectation value of the Wilson loop becomes

$`<W(C)>`$ $``$ $`\text{exp}[{\displaystyle \frac{1}{2\beta }}{\displaystyle \underset{r,r^{},i}{}}J_i(r)G(rr^{})J_i(r^{})]{\displaystyle \underset{r}{}}{\displaystyle }d\varphi (r)\text{exp}[{\displaystyle \frac{1}{8\pi ^2\beta }}{\displaystyle \underset{r}{}}\left(\mathrm{\Delta }_i\varphi (r)\right)^2`$ (A30)
$`2\text{e}^{2\pi ^2\beta G_0}{\displaystyle \underset{r}{}}(1\mathrm{cos}(\varphi (r)+2\pi \eta (r)))]`$
## B Large-$`\beta `$ prediction in Hamiltonian formulation

We now do the Hamiltonian study of U(1) model in (2+1) dimensions using weak coupling perturbation theory. The Hamiltonian formulation of compact U(1) gauge theory has been studied by Banks et al.  and Drell et al.  following the Kogut and Susskind approach .
Consider the anisotropic Wilson gauge action for U(1) model in (2+1) dimensions,

$$S=\beta \left[\frac{1}{\mathrm{\Delta }\tau }\underset{i}{}P_{it}+\mathrm{\Delta }\tau \underset{i<j}{}P_{ij}\right]$$
(B1)
where we made assumed the lattice symmetric in spatial directions with $`a_s=1`$ and $`a_t=\mathrm{\Delta }\tau `$.
The plaquette variable is given by

$$P_{ij}=1\mathrm{cos}\theta _{ij}^P$$

and $`\beta =1/g^2`$.

In the weak coupling approximation

$$\mathrm{cos}\theta _{ij}^P1\frac{1}{2}\theta _{ij}^{P}{}_{}{}^{2}$$
(B2)
The plaquette angles $`\theta ^P`$ are related to the link angles $`A_l`$ by

$$\theta ^P(\stackrel{~}{n}^{^{}},\stackrel{~}{\widehat{k}})=ϵ_{ijk}\left[A_l(\underset{¯}{\text{n}},\widehat{i})A_l(\stackrel{~}{n}+\widehat{j},\widehat{i})+A_l(\stackrel{~}{n}+\widehat{i},\widehat{j})A_l(\stackrel{~}{n},\widehat{j})\right]$$
(B3)
where $`(\stackrel{~}{n}^{^{}},\stackrel{~}{\widehat{k}})`$ denote a plaquette centered at $`\stackrel{~}{n}^{^{}}`$ and oriented in the $`\widehat{k}`$ while $`(\stackrel{~}{n}^{^{}},\widehat{i})`$ describes a link that starts at site $`\stackrel{~}{n}`$ in the direction $`\widehat{i}`$, with

$$\stackrel{~}{n}^{^{}}=\stackrel{~}{n}+(\frac{\widehat{i}+\widehat{j}}{2})$$

The Fourier transformations of the plaquette angles and the link angle give

$`\theta _{\stackrel{~}{k}}(\widehat{\mu })`$ $`=`$ $`{\displaystyle \frac{1}{\sqrt{N}}}{\displaystyle \underset{\stackrel{~}{n}^{^{}}}{}}\text{e}^{i\stackrel{~}{k}\stackrel{~}{n}^{^{}}}\theta _P(\stackrel{~}{n}^{^{}},\widehat{\mu })`$ (B4)
$`A_{\stackrel{~}{k}}q(\widehat{\mu })`$ $`=`$ $`{\displaystyle \frac{1}{\sqrt{N}}}{\displaystyle \underset{\stackrel{~}{n}}{}}\text{e}^{i\stackrel{~}{k}.\stackrel{~}{n}}A_l(\stackrel{~}{n},\widehat{\mu })`$ (B5)
where $`N=N_s^2N_t`$ is the total number of sites; then

$`A_{\stackrel{~}{k}}^{}(\widehat{\mu })`$ $`=`$ $`A_{\stackrel{~}{k}}(\widehat{\mu })`$ (B6)
$`\theta _{\stackrel{~}{k}}^{}(\widehat{\mu })`$ $`=`$ $`\theta _{\stackrel{~}{k}}(\widehat{\mu })`$ (B7)
Therefore eq. (B3) can be written as

$$\theta _{\stackrel{~}{k}}(\widehat{\mu })=\underset{\widehat{\nu }}{}G_{\stackrel{~}{k}}(\widehat{\mu },\widehat{\nu })A_{\stackrel{~}{k}}(\widehat{\nu })$$
(B8)
where

$$G_{\underset{¯}{\text{k}}}(\widehat{\mu },\widehat{\nu })=2i\underset{\nu ,\rho }{}ϵ_{\mu \nu \sigma }\text{e}^{ik_\nu /2}\mathrm{sin}(k_\sigma /2)$$
(B9)
Choosing the temporal gauge $`A_0=0`$, the time-like plaquette contributions, in the weak coupling approximation, become:

$`{\displaystyle \underset{\stackrel{~}{n}}{}}\theta _{i0}^2(\stackrel{~}{n})`$ $`=`$ $`{\displaystyle \underset{\stackrel{~}{k}}{}}\theta _{i0,\stackrel{~}{k}}^{}\theta _{i0,\stackrel{~}{k}}`$ (B10)
$`=`$ $`{\displaystyle \underset{\stackrel{~}{k}}{}}A_{\stackrel{~}{k}}^{}(\widehat{i})A_{\stackrel{~}{k}}(\widehat{i})\times 4\mathrm{sin}^2(k_0/2)`$ (B11)
while the space-like plaquette contribution becomes:

$`{\displaystyle \underset{\stackrel{~}{n}}{}}\theta _{12}^2(\stackrel{~}{n})`$ $`=`$ $`{\displaystyle \underset{\stackrel{~}{k}}{}}\theta _{12,\stackrel{~}{k}}^{}\theta _{12,\stackrel{~}{k}}`$ (B12)
$`=`$ $`{\displaystyle \underset{\stackrel{~}{k}}{}}A_{\stackrel{~}{k}}^{}R_{\stackrel{~}{k}}A_{\stackrel{~}{k}}`$ (B13)
where $`A_{\stackrel{~}{k}}`$ and $`R_{\stackrel{~}{k}}`$ are the vectors and matrices, $`(2\times 2)`$ in $`(\widehat{1},\widehat{2})`$ square, respectively, and $`R_{\stackrel{~}{k}}=G_{\stackrel{~}{k}}^{}G_{\stackrel{~}{k}}`$. Therefore

$$R_{\stackrel{~}{k}}=\left(\begin{array}{ccc}\mathrm{sin}^2(k_y/2)& \text{e}^{i(k_x/2k_y/2)}\mathrm{sin}k_x/2\mathrm{sin}k_y/2& \\ \text{e}^{i(k_y/2k_x/2)}\mathrm{sin}k_x/2\mathrm{sin}k_y/2& \mathrm{sin}^2(k_x/2)& \end{array}\right)$$
(B14)
The matrix $`R_{\stackrel{~}{k}}`$ is hermitian, and has eigenvalues:

$`\lambda _{\stackrel{~}{k}}(1)`$ $`=`$ $`4{\displaystyle \underset{i=1}{\overset{2}{}}}\mathrm{sin}^2(k_i/2)`$ (B15)
$`\lambda _{\stackrel{~}{k}}(1)`$ $`=`$ $`0`$ (B16)
We can define a unitary transformation

$$\omega _{\stackrel{~}{k}}=u_{\stackrel{~}{k}}A_{\stackrel{~}{k}}$$
(B17)
which diagonalizes the matrix $`R_{\stackrel{~}{k}}`$:

$$u_{\stackrel{~}{k}}=\frac{i}{sqrt\mathrm{sin}^2k_x/2+\mathrm{sin}^2k_y/2}\left(\begin{array}{ccc}\text{e}^{ik_x/2}\mathrm{sin}(k_y/2)& \text{e}^{ik_y/2)}\mathrm{sin}k_x/2& \\ \text{e}^{ik_x/2)}\mathrm{sin}k_x/2& \text{e}^{ik_y/2}\mathrm{sin}^2(k_x/2)& \end{array}\right)$$
(B18)
The components of $`\omega _{\stackrel{~}{k}}`$ are

$`\omega _{\stackrel{~}{k}}(1)`$ $`=`$ $`{\displaystyle \frac{i}{c_1}}\left[\text{e}^{ik_y/2)}\mathrm{sin}k_x/2A_{\underset{¯}{\text{k}}}(\widehat{2})\text{e}^{ik_x/2}\mathrm{sin}(k_y/2)A_{\underset{¯}{\text{k}}}(\widehat{2})\right]`$ (B19)
$`\omega _{\stackrel{~}{k}}(2)`$ $`=`$ $`{\displaystyle \frac{i}{c_1}}\left[\text{e}^{ik_x/2)}\mathrm{sin}k_x/2A_{\underset{¯}{\text{k}}}(\widehat{1})+\text{e}^{ik_y/2}\mathrm{sin}(k_y/2)A_{\underset{¯}{\text{k}}}(\widehat{2})\right]`$ (B20)
where $`c_1=\sqrt{\mathrm{sin}^2k_x/2+\mathrm{sin}^2k_y/2}`$.
The component of $`w_{\stackrel{~}{k}}`$ corresponding to the zero eigenvalue is

$$\stackrel{~}{\omega }_{\stackrel{~}{k}}=\frac{i}{c_1}\left[(1\text{e}^{ik_x})A_{\stackrel{~}{k}}(\widehat{x})+(1\text{e}^{ik_y})A_{\stackrel{~}{k}}(\widehat{y})\right].$$
(B21)
This vanishes at all orders and is an ”ignorable” coordinate.
Thus the action eq. (B1) can be written as

$`S`$ $`=`$ $`{\displaystyle \frac{\beta }{2}}{\displaystyle \underset{\stackrel{~}{k}}{}}[\mathrm{\Delta }\tau \left(\mathrm{}_{\stackrel{~}{k}}\omega _{\stackrel{~}{k}}^{}(1)\omega _{\stackrel{~}{k}}(1)\right]`$ (B23)
$`+{\displaystyle \frac{1}{\mathrm{\Delta }\tau }}\left[4\mathrm{sin}^2(k_0/2){\displaystyle \underset{\mu =1}{\overset{2}{}}}\omega _{\stackrel{~}{k}}^{}(\widehat{\mu })\omega _{\stackrel{~}{k}}(\widehat{\mu })\right]+\text{const}.`$
where

$$\mathrm{}_{\stackrel{~}{k}}=4\left[\mathrm{sin}^2k_1/2+\mathrm{sin}^2k_2/2\right]$$
(B24)
Therefore the partition function can be written as

$$Z=\underset{\stackrel{~}{k}}{}\{𝑑\omega \stackrel{~}{k}\}\text{e}^{S[\omega ]}$$
(B25)
where the action can be written in the following form:

$`s[\omega ]`$ $``$ $`{\displaystyle \frac{\beta }{2}}{\displaystyle \underset{\stackrel{~}{k}}{}}[\omega _{\stackrel{~}{k}}^{}(1)\omega _{\stackrel{~}{k}}(1)(\mathrm{\Delta }\tau \mathrm{}_{\stackrel{~}{k}}+{\displaystyle \frac{4}{\mathrm{\Delta }\tau }}\mathrm{sin}^2(k_0/2))`$ (B27)
$`+\omega _{\stackrel{~}{k}}^{}(2)\omega _{\stackrel{~}{k}}(2)\left({\displaystyle \frac{4}{\mathrm{\Delta }\tau }}\mathrm{sin}^2(k_0/2)\right)]`$
Now $`\omega _{\stackrel{~}{k}}`$ is the ignorable coordinate and will be ignored from here on.
The partition function becomes

$`Z`$ $``$ $`{\displaystyle \underset{\underset{¯}{\text{k}}}{}}\left[{\displaystyle \frac{2\pi }{\beta }}{\displaystyle \frac{1}{\left(\mathrm{\Delta }\tau \mathrm{}_{\underset{¯}{\text{k}}}+\frac{4}{\mathrm{\Delta }\tau }\mathrm{sin}^2(k_0/2)\right)}}\right]^{1/2}`$ (B28)
$`=`$ $`\text{exp}[{\displaystyle \frac{N}{2}}\text{ln}({\displaystyle \frac{2\pi }{\beta }}){\displaystyle \frac{1}{2}}{\displaystyle \underset{\underset{¯}{\text{k}}}{}}\text{ln}\left((\mathrm{\Delta }\tau \mathrm{}_{\underset{¯}{\text{k}}}+{\displaystyle \frac{4}{\mathrm{\Delta }\tau }}\mathrm{sin}^2(k_0/2))\right]`$ (B29)
This lead to the following free energy density:

$$f=\frac{1}{2N}\underset{\stackrel{~}{k}}{}\text{ln}\left[\frac{\beta }{2\pi }d(\stackrel{~}{k})\right]$$
(B30)
where

$$d(\stackrel{~}{k})=(\mathrm{\Delta }\tau \mathrm{}_{\stackrel{~}{k}}+\frac{4}{\mathrm{\Delta }\tau }\mathrm{sin}^2(k_0/2)$$
(B31)
### 1 Mean space-like plaquette

At weak coupling the cosine term in the space-like part may be expanded in the powers of $`\theta _P`$. Therefore

$`<\mathrm{cos}\theta _{ij}>`$ $`=`$ $`1{\displaystyle \frac{1}{2}}<\theta _{ij}^2>`$ (B32)
$`=`$ $`1{\displaystyle \frac{1}{2}}<{\displaystyle \frac{1}{N}}{\displaystyle \underset{\stackrel{~}{k}}{}}\mathrm{}_{\stackrel{~}{k}}\omega _{\stackrel{~}{k}}^{}(1)\omega _{\stackrel{~}{k}}(1)>`$ (B33)
where

$`<\omega _{\stackrel{~}{k}}^{}(1)\omega _{\stackrel{~}{k}}(1)>`$ $`=`$ $`{\displaystyle \frac{𝑑\omega _{\stackrel{~}{k}}\text{e}^{\frac{\beta }{2}d(\stackrel{~}{k})\omega _{\stackrel{~}{k}}^{}\omega _{\stackrel{~}{k}}}\omega _{\stackrel{~}{k}}^{}\omega _{\stackrel{~}{k}}}{𝑑\omega _{\stackrel{~}{k}}\text{e}^{\frac{\beta }{2}d(\stackrel{~}{k})\omega _{\stackrel{~}{k}}^{}\omega _{\stackrel{~}{k}}}}}`$ (B34)
Using the integrals

$`{\displaystyle _{\mathrm{}}^{\mathrm{}}}𝑑\omega \text{e}^{a\omega ^2}`$ $`=`$ $`\sqrt{{\displaystyle \frac{\pi }{2}}}`$ (B35)
$`{\displaystyle _{\mathrm{}}^{\mathrm{}}}𝑑\omega \omega ^2\text{e}^{a\omega ^2}`$ $`=`$ $`{\displaystyle \frac{1}{2a}}\sqrt{{\displaystyle \frac{\pi }{2}}}`$ (B36)
we get

$$<\omega _{\stackrel{~}{k}}^{}\omega _{\stackrel{~}{k}}>=\frac{1}{\beta d(\stackrel{~}{k}}$$
(B37)
Therefore

$$<\mathrm{cos}\theta _{ij}>=1\frac{1}{2\beta N}\underset{\stackrel{~}{k}}{}\frac{\mathrm{}_{\stackrel{~}{k}}}{d(\stackrel{~}{k})}$$
(B38)
Now consider the “cylindrical limit”, $`N_t\mathrm{}`$, $`\mathrm{}\tau 0`$: In this limit the sum term in the right hand side of the above equation becomes;

$`{\displaystyle \underset{\underset{¯}{\text{k}}}{}}{\displaystyle \frac{\mathrm{}_{\stackrel{~}{k}}}{d(\stackrel{~}{k})}}`$ $``$ $`{\displaystyle \underset{\stackrel{}{k}}{}}\mathrm{}(\stackrel{}{k}){\displaystyle \frac{N_t}{2\pi }}{\displaystyle _0^{2\pi }}𝑑k_0{\displaystyle \frac{1}{\mathrm{}\tau +\frac{4}{\mathrm{}\tau }\mathrm{sin}^2(k_0/2)}}`$ (B39)
$`=`$ $`{\displaystyle \underset{\stackrel{}{k}}{}}{\displaystyle \frac{N_t}{2\pi }}2{\displaystyle _0^{2\pi }}{\displaystyle \frac{dx}{\left(\mathrm{}\tau \mathrm{}(\stackrel{}{k})+\frac{x^2}{\mathrm{}\tau }\right)}}`$ (B40)
$`=`$ $`{\displaystyle \underset{\stackrel{}{k}}{}}{\displaystyle \frac{N_t}{\pi }}{\displaystyle \frac{\pi }{2\sqrt{\mathrm{}(\stackrel{}{k})}}}`$ (B41)
In the bulk limit, $`N\mathrm{}`$, this becomes;

$``$ $`{\displaystyle \frac{N_t2.(\frac{N_s}{2\pi })^2}{}}{\displaystyle _0^{2\pi }}𝑑k_x𝑑k_y\sqrt{4(\mathrm{sin}^2k_x/2+\mathrm{sin}^2k_y/2}`$ (B42)
$`=`$ $`{\displaystyle \frac{N_t}{2}}.2N_s^2C_1^{(2)}`$ (B43)
where

$`C_1^{(2)}`$ $`=`$ $`{\displaystyle \frac{1}{N_s^2}}{\displaystyle \underset{\stackrel{}{k}}{}}(1\gamma _\stackrel{}{k})^{1/2}`$ (B44)
$`=`$ $`{\displaystyle \frac{1}{2N_s^2}}{\displaystyle \underset{\stackrel{}{k}}{}}\mathrm{\Delta }^{1/2}(\stackrel{}{k})`$ (B45)
$`=`$ $`0.958`$ (B46)
is a standard lattice sum.
Therefore

$`<\mathrm{cos}\theta _{ij}>`$ $`=`$ $`1{\displaystyle \frac{1}{2\beta }}C_1^{(2)}`$ (B47)
$`=`$ $`1{\displaystyle \frac{g^2}{2}}C_1^{(2)}`$ (B48)
$`=`$ $`10.479g^2`$ (B49)
This result agrees with the estimates in the Hamiltonian limit  Now taking the bulk limit in isotropic case, we get

$$\underset{\stackrel{~}{k}}{}\frac{\mathrm{\Delta }(\stackrel{}{k})}{d(\stackrel{~}{k})}\frac{N}{(2\pi )^3}_0^{2\pi }𝑑k_x𝑑k_y𝑑k_0\frac{(\mathrm{sin}^2k_x/2+\mathrm{sin}^2k_y/2)}{(\mathrm{sin}^2k_x/2+\mathrm{sin}^2k_y/2+\mathrm{sin}^2k_0/2)}$$
(B50)
By symmetry, we can easily find that

$$\underset{\stackrel{~}{k}}{}\frac{\mathrm{\Delta }(\stackrel{}{k})}{d(\stackrel{~}{k})}=\frac{2N}{3}$$
(B51)
Therefore

$$<\mathrm{cos}\theta _{ij}>=1\frac{g^2}{3}$$
(B52)
### 2 Mean time-like plaquette

In weak coupling approximation the time-like plaquette can be written as

$`{\displaystyle \frac{1}{2}}{\displaystyle \underset{i=1}{\overset{2}{}}}<\mathrm{cos}\theta _{i0}>`$ $`=`$ $`1{\displaystyle \frac{1}{4}}{\displaystyle \underset{i=1}{\overset{2}{}}}<\theta _{i0}^2>`$ (B53)
$`=`$ $`1{\displaystyle \frac{1}{4}}<{\displaystyle \frac{1}{N}}{\displaystyle \underset{\stackrel{~}{k}}{}}4\mathrm{sin}^2k_0/2A_{\stackrel{~}{k}}^{}(\widehat{i})A_{\stackrel{~}{k}}(\widehat{i})>`$ (B54)
$`=`$ $`1{\displaystyle \frac{1}{4N}}<{\displaystyle \underset{\stackrel{~}{k}}{}}4\mathrm{sin}^2k_0/2(\omega _{\stackrel{~}{k}}^{}(\widehat{1})\omega _{\stackrel{~}{k}}(\widehat{1})+\omega _{\stackrel{~}{k}}^{}(\widehat{2})\omega _{\stackrel{~}{k}}(\widehat{2})>`$ (B55)
Now

$$<\omega _{\stackrel{~}{k}}^{}(\widehat{1})\omega _{\stackrel{~}{k}}(\widehat{1})+\omega _{\stackrel{~}{k}}^{}(\widehat{2})\omega _{\stackrel{~}{k}}(\widehat{2})>=\frac{1}{\beta d(\stackrel{~}{k})}+\frac{\mathrm{\Delta }\tau }{4\beta \mathrm{sin}^2k_0/2}$$
(B56)
$$\frac{1}{2}\underset{i=1}{\overset{2}{}}<\mathrm{cos}\theta _{i0}>=1\frac{g^2}{4N}\underset{\stackrel{~}{k}}{}4\mathrm{sin}^2k_0/2\left[\frac{1}{d(\stackrel{~}{k})}+\frac{\mathrm{\Delta }\tau }{4\mathrm{sin}^2k_0/2}\right]$$
(B57)
The term in the square brackets can be written as

$`{\displaystyle \underset{\stackrel{~}{k}}{}}4\mathrm{sin}^2k_0/2\left[{\displaystyle \frac{1}{d(\stackrel{~}{k})}}+{\displaystyle \frac{\mathrm{\Delta }\tau }{4\mathrm{sin}^2k_0/2}}\right]`$ $`=`$ $`{\displaystyle \frac{N}{(2\pi )^3}}{\displaystyle }{\displaystyle _0^{2\pi }}{\displaystyle }dk_xdk_ydk_0\mathrm{\Delta }\mathrm{sin}^2k_0[{\displaystyle \frac{1}{\mathrm{sin}^2k_0/2}}`$ (B59)
$`+{\displaystyle \frac{1}{\mathrm{sin}^2k_0/2+\mathrm{sin}^2k_0/2+\mathrm{\Delta }\tau ^2(\mathrm{sin}^2k_x/2+\mathrm{sin}^2k_y/2)}}]`$
$`=`$ $`\{\begin{array}{cc}0\hfill & (\mathrm{\Delta }\tau =0)\hfill \\ N(1+1/3)\hfill & (\mathrm{\Delta }\tau =1)\hfill \end{array}`$ (B62)
Thus the mean time-like plaquette gives;

$$\frac{1}{2}\underset{i=1}{\overset{2}{}}<\mathrm{cos}\theta _{i0}>=\{\begin{array}{cc}0\hfill & (\mathrm{\Delta }\tau =0)\hfill \\ 1\frac{g^2}{3}\hfill & (\mathrm{\Delta }\tau =1)\hfill \end{array}$$
(B63)
as is expected.

### 3 Expectation value of Wilson loop

The expectation value of the Wilson loop is given by

$$<W(R,T)>=\frac{1}{Z}\underset{\stackrel{~}{k}}{}𝑑\omega _{\stackrel{~}{k}}W_c\text{e}^{S[\omega ]}$$
(B64)
where the Wilson loop operator is given by

$$W_C=\frac{1}{2}\left[\underset{lC}{}\text{e}^{A_l}+\text{h.c}\right]$$
(B65)
Taking the separation R in the $`x`$ direction, for simplicity, the Fourier transformation of above equation) is

$$W_C=\text{exp}[\frac{i}{\sqrt{N}}\underset{\stackrel{~}{k}}{}\underset{\stackrel{~}{x},\mu C}{}\text{e}^{i\stackrel{~}{k}.\stackrel{~}{x}}A_{\stackrel{~}{k}}(\widehat{\mu }]+\text{h.c}$$
(B66)
We now consider the contribution from class of “zero modes” with $`k=(k_0,\stackrel{}{k}=0)`$. These modes are likely to be responsible for the dominant finite-size corrections , i.e.,

$$W_0(R,t)=\text{exp}\left[\frac{i}{\sqrt{N}}\underset{k_0}{}A_{k_0,\stackrel{}{0}}(\widehat{x})R\left(\text{e}^{ik_0x_{01}}\text{e}^{ik_0x_{02}}\right)\right]$$
(B67)
Recalling that $`A_0=0`$, the above equation becomes

$`W_0(R,t)`$ $`=`$ $`\text{exp}\left[{\displaystyle \frac{2}{\sqrt{N}}}{\displaystyle \underset{k_0}{}}A_{k_0,\stackrel{}{0}}(\widehat{x})R\text{e}^{ik_0(\frac{x_{01}+x_{02}}{2})}\mathrm{sin}k_0T/2\right]`$ (B68)
$`=`$ $`\text{exp}\left[{\displaystyle \frac{2}{\sqrt{N}}}{\displaystyle \underset{k_0}{}}\omega _{\stackrel{~}{k}}(\widehat{1})R\mathrm{sin}k_0T/2\right]`$ (B69)
where $`T=x_{02}x_{01}`$ and we have ignored the $`\text{e}^{ik_0(\frac{x_{01}+x_{02}}{2})}`$. Therefore

$$<W_0(R,T)>=\frac{1}{Z(0)}\underset{k_0}{}𝑑\omega _{k_0}(\widehat{1})\text{exp}\left[\frac{2}{\sqrt{N}}\omega _{k_0}R\mathrm{sin}(k_0T/2)\frac{\beta }{2\mathrm{\Delta }\tau }4\mathrm{sin}^2(k_0/2)\omega _{k_0}^{}\omega _{k_0}\right]$$
(B70)
Let

$$\frac{\beta }{2\mathrm{\Delta }\tau }4\mathrm{sin}^2(k_0/2)=c(k_0);$$

then the numerator of above equation can be written as

$$I_{k_0}=_{\mathrm{}}^{\mathrm{}}d\omega \text{exp}[c(k_0(\omega i\frac{R\mathrm{sin}(k_0T/2)}{\sqrt{N}c(k_0)})^2]\text{exp}[\frac{R^2\mathrm{sin}^2(k_0T/2)}{\sqrt{N}c(k_0)}]$$
(B71)
Therefore the expectation value of the Wilson loop is given by

$`<W_0(R,T)>`$ $`=`$ $`{\displaystyle \underset{k_0}{}}\text{exp}\left[{\displaystyle \frac{R^2\mathrm{sin}^2(k_0T/2)}{\sqrt{N}2\beta \mathrm{sin}^2(k_0)}}\right]`$ (B73)
$`\text{exp}\left[{\displaystyle \frac{R^2\mathrm{\Delta }\tau }{2N_s^2N_t\beta }}{\displaystyle \frac{N_t}{2\pi }}{\displaystyle _0^{2\pi }}𝑑k_0{\displaystyle \frac{\mathrm{sin}^2(k_0T/2}{\mathrm{sin}^2(k_0/2}}\right]`$
Now consider the cylindrical limit, $`N_t\mathrm{}`$, $`\mathrm{\Delta }\tau 0`$. In this limit

$$\frac{1}{2\pi }_0^{2\pi }𝑑k_0\frac{\mathrm{sin}^2(k_0T/2)}{\mathrm{sin}^2(k_0/2)}T$$
(B74)
Therefore

$`<W_0(R,T)>`$ $`=`$ $`\text{exp}\left[{\displaystyle \frac{g^2}{2}}({\displaystyle \frac{R}{N_s}})^2T\mathrm{\Delta }\tau \right]`$ (B76)
$`\text{exp}\left(\sigma \times \text{area}\right)`$
Where

$$\sigma _{eff}=\frac{g^2}{2}(\frac{R}{N_s^2})=\frac{g^2}{2}(\frac{R}{L^2})$$
(B77)
For R=L, i.e., for pair of Polyakov loops mapping around the lattice in the spatial direction, we get

$$\sigma _{eff}=\frac{g^2}{2L}$$
(B78)
This corresponds to the value calculated by Hamer and Zheng , $`\sigma =1/L`$, for a string wrapping around the entire lattice, if we multiply by an extra factor $`g^2/2`$ coming from renormalizing Hamiltonian. As is obvious that for smaller loops, the effective string tension is less.