Analyzing the Higgs-confinement transition with non-local operators on the lattice

In gauge theories coupled with fundamental matter, the Higgs phase and the confinement phase cannot be distinguished by any local order parameter because they have the same symmetry in the Ginzburg-Landau sense [1, 2, 3]. When the theory has topological order, the two phases are distinguishable by non-local order parameters and separated by a first-order phase transition. Otherwise, the two phases are continuously connected. A prominent example of the Higgs-confinement continuity is quantum chromodynamics (QCD) at nonzero density. Quarks form color singlet hadrons in the vacuum while quarks form diquark Cooper pairs with color-flavor locking in high density limit. The hadronic phase and the color-flavor locked phase have the same symmetry breaking pattern, and thus are expected to be smoothly connected without phase boundary. This is referred to as the quark-hadron continuity [4]. The quark-hadron continuity might be relevant for phenomenology of dense QCD matter, such as neutron star physics [5].

Non-local operators, such as loop and surface operators, play a key role in the Higgs-confinement phase transition. Well-known examples of the non-local operators are the Polyakov and ’t Hooft loops. These loops are (dis)order parameters of electric and magnetic confinement, respectively. In general, there is more than one order parameter for one phase transition. Among several candidates, a topological vortex [6, 7, 8, 9] and its Aharonov-Bohm phase [10, 11] were intensely discussed in the context of the quark-hadron continuity. Also, it was recently suggested that a correlation function on a vortex can distinguish the Higgs and confinement phases even if the transition is a smooth crossover [12]. Numerical simulation of lattice gauge theory is a powerful way to compute these non-local operators nonperturbatively. Some of them were already established and others were not. We need to deepen our understanding of their lattice implementation.

In this paper, we study the lattice implementation of a few non-local operators for analyzing the Higgs-confinement phase transition. We aim to showcase numerical analysis and to discuss practical difficulty etc. For this purpose, we should choose a simple gauge-Higgs model as a testing ground. We adopt the charge-$q$ Abelian Higgs model introduced by Fradkin and Shenker long ago [1] because its phase structure is well understood, thanks to many earlier works [13, 14, 15, 16, 17, 18, 19, 20]. Once we understand the implementation in this simple model, we can proceed to the next stage in forthcoming works; e.g., the SU(3) gauge-Higgs model, which is the Ginzburg-Landau theory of dense QCD.

The rest of the paper is organized as follows. We begin with a brief introduction of the Abelian Higgs model and the dual lattice in Sec. II. We discuss three non-local operators: the Polyakov loop in Sec. III, the ’t Hooft loop in Sec. IV, and the Aharonov-Bohm phase in Sec. V. All the equations are written in the lattice unit throughout the paper.

We introduce the lattice formulation of the charge-$q$ Abelian Higgs model in $3+1$ dimensions. We employ the compact formulation with the unit-radius Higgs field. In the compact formulation, there are only two kinds of field variables: the Higgs field angle $\theta(x)$ and the Abelian gauge field $A_{\mu}(x)$. The action takes the form

where $\hat{\mu}$ is the unit lattice vector in $\mu$ direction and $F_{\mu\nu}(x)$ is the discretized field strength

This action is invariant under the local U(1) gauge transformation

by an arbitrary function $\lambda(x)$ and the global $Z_{q}$ transformation

The global $Z_{q}$ symmetry exists only for multiple charge $q\geq 2$. We set $q=2$ in our numerical simulation.

The model has two parameters: the gauge coupling constant $\beta$ and the hopping parameter $\kappa$. These parameters control the phase of matter. The small $\beta$ regime is the confinement phase. The fields are strongly coupled and electrically charged particles are confined. The large $\beta$ and large $\kappa$ regime is the Higgs phase. Since the hopping term is dominant, $\theta(x)$ is frozen and Higgs condensation is induced. The large $\beta$ and small $\kappa$ regime is the Coulomb phase. The vacuum is the weakly-coupled perturbative one. For $q\geq 2$, these three regimes are separated by the first-order phase transitions. In particular, the Higgs and confinement regimes are distinguishable by the global $Z_{q}$ symmetry. This is in stark contrast to the $q=1$ model, where the global $Z_{q}$ symmetry is absent and the confinement and Higgs regimes are smoothly connected without phase transition. The phase diagram of the charge-$q$ Abelian Higgs model was investigated in several previous lattice calculations. The calculation was first done with $Z_{N}$ approximation [13], and then with exact U(1) for the specific heat [14], hysteresis [15], and magnetic monopole flow [16]. The location of phase boundaries was recently obtained with large lattice sizes [17]. The $(2+1)$-dimensional model was also studied for $q=1$ and $2$ [18].

We also introduce the dual lattice for later convenience. The above model is defined on a four-dimensional hypercubic lattice with the volume $V=N_{x}N_{y}N_{z}N_{\tau}$. The dual lattice is also a four-dimensional hypercubic lattice with the volume $V$. The dual sites are defined by the geometric centers of hypercubes of the original lattice. A $n$-dimensional manifold on the original lattice has one-to-one correspondence to a $(4-n)$-dimensional manifold on the dual lattice. We use a star symbol $*$ for describing dual manifolds. For example, an one-dimensional link $\mathcal{L}$ penetrates a three-dimensional dual cube $\mathcal{V}^{*}$, a two-dimensional face $\mathcal{S}$ intersects a two-dimensional dual face $\mathcal{S}^{*}$, and a three-dimensional cube $\mathcal{V}$ is penetrated by an one-dimensional dual link $\mathcal{L}^{*}$. Magnetic objects, such as magnetic monopoles and vortices, are not elemental variables in the path integral but defined on the dual manifolds. The operators to create the magnetic objects are defined so as to act on the corresponding manifolds on the original lattice. The explicit operator constructions are demonstrated in Secs. IV and V. Note that the Abelian Higgs model can be rewritten and simulated by using dual variables [21, 22, 23], but we do not adopt such dual formulation in this work.

The Polyakov loop is the most established order parameter of confinement-deconfinement phase transition [24]. Since the Polyakov loop is zero in the confinement phase and nonzero in the deconfinement phase, it is sometimes called a disorder parameter. The Polyakov loop was calculated in many studies of lattice gauge theory, mainly in SU($N$) gauge theory. As for the charge-$q$ Abelian Higgs model, there is only one calculation and brief description [17]. Here we provide more detailed analysis.

We define the averaged Polyakov loop

from the Polyakov loop at the spatial coordinate $\bm{r}$

The expectation value

and $\langle p(\bm{r})\rangle$ are order parameters of spontaneous $Z_{q}$ symmetry breaking. In numerical simulation, however, they are zero regardless of the phase because symmetry is not spontaneously broken in finite volume. The difference of the phase is seen in probability distribution. Figure 1 is scattering plots of the averaged Polyakov loop for each configuration. The number of gauge configuration is $N_{\mathrm{conf}}=100$. The left panel is the case of $\kappa=0$ and $\beta=0.5$, which corresponds to the confinement phase of pure U(1) gauge theory. The distribution is very close to zero. For large $\beta$, the system goes to the Coulomb phase. The distribution is away from zero and spreads into a circle, as shown in the middle panel. The circular distribution reflects the global U(1) symmetry in the pure gauge theory. The right panel corresponds to the Higgs phase. For $\kappa>0$, the hopping term in Eq. (6) breaks the global symmetry from U(1) to $Z_{2}$. When $\kappa$ is very large, the hopping term restricts link variables to $\pm 1$ and the model is reduced to $Z_{2}$ gauge theory. The distribution is concentrated in two minima on the real axis. The change of the distribution from U(1) to $Z_{2}$ characterizes Coulomb-Higgs phase transition.

For determining transition points and drawing a phase diagram, we use $\expectationvalue{\absolutevalue{P}}$ instead of $\expectationvalue{P}$. (The non-averaged one $\expectationvalue{p(\bm{r})}$ cannot be used for this purpose because of $|p(\bm{r})|=1$ in U(1) gauge theory. This is the reason why we take the spatial average in Eq. (6).) The results are shown in Fig. 2. The left panel shows the development of the Polyakov loop against $\beta$ at fixed $\kappa$. While $\expectationvalue{P}$ is always zero, $\expectationvalue{\absolutevalue{P}}$ is nonzero and changes from small to large values. There is a sharp jump between $\beta=0.8$ and $\beta=0.9$. This is the confinement-deconfinement phase transition and its critical coupling constant is $0.8<\beta_{c}<0.9$ at $\kappa=0.8$. The right panel is the phase diagram in $\beta$-$\kappa$ plane, where the values of $\expectationvalue{\absolutevalue{P}}$ are highlighted in different colors. The phase boundaries of the three phases are clearly seen. At $\kappa=0$, the system is pure U(1) gauge theory and the critical coupling constant is $\beta_{c}\simeq 1.0$, which is consistent with pure gauge simulation [25, 26, 27]. The critical coupling constant $\beta_{c}$ gradually decreases as $\kappa$ increases. The deconfinement phase is further decomposed into the Higgs and Coulomb regimes. They can be distinguished by the values of the Polyakov loop; the Higgs regime has larger values and the Coulomb regime has smaller values. This can be understood from the different distribution in Fig. 1. The Higgs-Coulomb phase transition occurs at $\kappa_{c}=0.4$-$0.5$. All the results are consistent with previous studies [14, 15, 16, 17].

The ’t Hooft loop is dual to the Polyakov loop [28]. (In this paper, we consider the ’t Hooft loop encircling imaginary time direction.) The ’t Hooft loop is the worldline of a magnetic monopole while the Polyakov loop is the worldline of an electric charge. Since the ’t Hooft loop is nonzero in the confinement phase and zero in the Higgs phase, it is an order parameter of Higgs-confinement phase transition.

The lattice implementation of the ’t Hooft loop has been known since long ago [29, 30, 31] and applied to SU($N$) gauge theory [32, 33, 34, 35, 36, 37, 38, 39, 40, 41]. The implementation generates the worldline of a magnetic monopole in $2+1$ dimensions. Here we apply it to the Abelian gauge theory in $3+1$ dimensions. In $(3+1)$-dimensional case, the implementation generates the worldsheet of a magnetic string, not the worldline of a single monopole, because of extra spatial direction, i.e., $z$ direction. Let us consider the worldsheet on a two-dimensional dual surface $\mathcal{S}^{*}$ in $z$-$\tau$ plane. The surface $\mathcal{S}^{*}$ has unit length in $z$ direction and wraps around $\tau$ direction. The corresponding magnetic string has unit length and its the endpoints are identified as a monopole and an antimonopole, as shown in the left panel of Fig. 3. We define the surface operator $T$ to create the worldsheet. The operation of $T$ multiplies a nontrivial $Z_{q}$ element to all the plaquette variables intersecting $\mathcal{S}^{*}$. The expectation value of $T$ is given by

with the modified action

The surface $\mathcal{S}$ is defined as a set of the $x$-$y$ plaquettes intersecting $\mathcal{S}^{*}$. The expectation value is positive real. It is interpreted as $\exp(-F_{m})$, where $F_{m}$ is the free energy of a unit-length magnetic string, i.e., a monopole-antimonopole pair.

As shown in the right panel of Fig. 3, the simulation results clearly exhibit a jump of the first-order phase transition in spite of finite volume. The free energy of the monopole-antimonopole pair (two monopoles for $q=2$) is finite in the confinement phase and nearly infinite in the Higgs phase. The transition point $\beta_{c}=0.8$-$0.9$ is consistent with that estimated from the Polyakov loop in the previous section. This indicate that the above formulation works well for analyzing the Higgs-confinement phase transition. Strictly speaking, however, the expectation value is not an order parameter of any symmetry. In $3+1$ dimensions, it is not the free energy of a single monopole but a monopole-antimonopole pair. The monopole-antimonopole pair has the same symmetry as the vacuum and its free energy will be large but finite even if a monopole is confined. The expectation value seems nearly zero but is not exactly zero in Fig. 3.

In practice, we computed Eq. (9) by the Monte Carlo sampling with sequential reweighting [38]. The computational cost of the reweighting method rapidly increases as the area of $\mathcal{S}$ increases. The area of $\mathcal{S}$ is proportional to $N_{\tau}$ and the length of the magnetic string in $z$ direction. Above we considered the magnetic string with unit length because the computational cost is small. The magnetic string wrapping around $z$ direction, i.e., without endpoints, can be considered and its free energy can be computed. This will show similar tendency although the cost is larger than ours.

The difference between the Higgs and confinement regimes also emerges in nontrivial correlation between non-local operators. It was recently conjectured that the two regimes are distinguished by the Aharonov–Bohm phase around a superfluid vortex [10]. This conjecture was analytically examined in the strong coupling and deep Higgs limits [11]. Lattice simulation will be necessary for full quantum examination. Here we formulate the lattice simulation of the Aharonov–Bohm phase and test it in the Abelian Higgs model. Note that the vortex in the Abelian Higgs model is not a global vortex of superfluid but a gauged vortex of superconductor. Although physics is different, the model can be used for testing formulation. The formulation is straightforwardly applicable to lattice gauge theories with a superfluid vortex, such as the non-Abelian Higgs model [42].

A vortex of the Higgs field is defined by nonzero circulation of $\theta(x)$. A vortex worldsheet can be generated by a cuboid operator on the dual lattice [43]. The setup is schematically depicted in the left panel of Fig. 4. Boundary conditions are open in $x$ and $y$ directions and periodic in $z$ and $\tau$ directions. Let $\mathcal{V}^{*}$ denote the three-dimensional volume which lies from the center of the box to the boundary, and $\mathcal{L}$ denote the links penetrating $\mathcal{V}^{*}$. The cuboid operator on $\mathcal{V}^{*}$ is defined by multiplying the factor $\exp(i\pi)=-1$ to the hopping terms in $\mathcal{L}$. This results in the $\pi$ circulation around the center of the box. The vortex exists at the center of the box and its worldsheet wraps around $z$ and $\tau$ directions. We define the charge-$q$ spatial Wilson loop

encircling the vortex, as shown in the figure. The path $\mathcal{C}$ is a square on the $x$-$y$ plane and its center is located at the vortex position. We computed the expectation value

with

Since the spatial Wilson loop is the trajectory of a charged particle, it picks up the Aharonov–Bohm phase when it encircles the vortex. We took the average of the translationally invariant $N_{z}N_{\tau}$ loops in Eq. (12).

The right panel of Fig. 4 shows the real parts of the $3\times 3$ and $5\times 5$ Wilson loops. The imaginary parts are always zero in case of $q=2$. The results of the $3\times 3$ loop are clear. The expectation value is negative in the Higgs phase. This corresponds to the non-trivial Aharonov–Bohm phase $\pi$ under coherent Higgs condensation. In the confinement regime, the expectation value is zero within statistical error and the Aharonov–Bohm phase is not calculable. The $5\times 5$ loop is qualitatively the same but not so clear due to statistical noise. In general, the expectation value of the Wilson loop is an exponentially damping function of loop size; the area law in the confinement regime and the perimeter law in the deconfinement regime. The signal to noise ratio drastically gets worse for larger loop. If one wants to take large loop limit, it will be computationally demanding. Although the asymptotic value of the Aharonov–Bohm phase is considered to be an order parameter, it will be hard to use it in practice.

For looking more closely at the Aharonov–Bohm phase, let us quantify the angle of the Wilson loop. The expectation value $\langle{\rm arg}W\rangle$ is however ill-defined in quantum theory because it depends on the domain of argument. For example, when ${\rm arg}W$ fluctuates around $\pi$, $\langle{\rm arg}W\rangle\simeq\pi$ if the argument is defined in $(0,2\pi]$, but $\langle{\rm arg}W\rangle\simeq 0$ if the argument is defined in $(-\pi,\pi]$. Instead, we here introduce probability distribution of ${\rm arg}W$. Because of $|W|=1$ in U(1) gauge theory, probability distribution of ${\rm arg}W$ uniquely determines the expectation value. We generated $N_{\rm conf}=10^{3}$ gauge configurations, computed $N_{z}N_{\tau}=10^{2}$ translationally invariant Wilson loops, and made one histogram, i.e., how many data exist in each domain of ${\rm arg}W$, from $10^{3}\times 10^{2}$ data. We repeated these steps for 10 independent sets of gauge configurations to evaluate statistical average and error. Thus, totally, we made a histogram from the $10^{3}\times 10^{2}\times 10=10^{6}$ data. As seen from the histograms in Fig. 5, the probability distribution changes between the Higgs and confinement regimes. In the confinement regime, the phase is random because of strong fluctuation of gauge fields. In the Higgs regime, the probability distribution has a peak at $\pi$, but the peak is tiny, only about 0.1 %. This means that quantum fluctuation is still strong and semi-classical approximation is invalid. The probability distribution will go to delta function in classical limit, such as the deep Higgs limit $\kappa\to\infty$.