Phase diagram of 1+1D Abelian-Higgs model and its critical point

We determine the phase diagram of the Abelian-Higgs model in one spatial dimension and time (1+1D) on a lattice. We identify a line of first order phase transitions separating the Higgs region from the confined one. This line terminates in a quantum critical point above which the two regions are connected by a smooth crossover. We analyze the critical point and find compelling evidences for its description as the product of two non-interacting systems, a massless free fermion and a massless free boson. However, we find also some surprizing results that cannot be explained by our simple picture, suggesting this newly discovered critical point to be an unusual one.

Introduction. Gauge theories in 1+1 dimensions (1D in space and time) are ideal playgrounds to characterize the effects of strong-coupling between matter and gauge fields. Many of the non-perturbative aspects of 3+1 dimensional gauge theories relevant to our understanding of particle physics, such as quark confinement and chiral symmetry breaking, have a 1+1D analogue. Furthermore, 1+1D field theories can often be treated analytically [1,2] providing important insights to the physics of 1+1D systems relevant also to condensed matter physics.
In AHM 2 , a weak matter-field coupling limit [13] suggests that the phase diagram is shared by two phases characteristic of the Higgs mechanism [8][9][10][11], a superfluid phase 1 with the quasi-condensation of bosons (Higgs phase) and a Mott-insulating phase 2 with strong interactions. However, non-perturbative calculations show that the phenomenology in the phase 1 is the same as in the phase 2, and bosons are always tightly confined [13] (for a recent discussion see [14,15]).
One can certify the presence of a phase transition in d + 1 dimension for any d > 2 [16,17], but due to the (boring) expectation of a single phase in the continuum, the phase diagram of AHM 2 has never been computed on the lattice in 1+1D. This work aims at filling this gap.
Our work is strongly motivated by the current prospects of simulating lattice gauge theories using cold atomic setups [18][19][20][21][22][23]. Since bosons are easier to cool down than fermions, experiments along the lines proposed in [24][25][26][27][28][29][30] should soon explore the phase diagram of AHM 2 . Physicists have been working hard to measure the Higgs mode in experiments with cold atoms for a longtime, as reviewed in [31]. In 2+1D, an explicit particlehole symmetry protects the decay of the Higgs mode into Goldstone modes allowing a proper measurement of its mass. These conditions are only met at the tip of the lobe of the Mott insulator to superfluid transition in Bose-Hubbard systems [32][33][34] which is, unfortunately, in 1+1D of the Berezinskii-Kosterlitz-Thouless type (see e.g., [35][36][37][38]) and is not particle-hole symmetric. This observation seems to strengthen the picture emerging from the presence of a single phase in AHM 2 and seems to suggest that a proper Higgs mode does not exist in 1+1D.
The results we present indicate a different picture, still characterized by a single phase, but with a reach land-scape of transitions.
By performing matrix product states (MPS) [39,40] simulations of the Hamiltonian version of AHM 2 , we confirm the presence of a single phase for all the values of the mass of the bosons µ 2 /q 2 and their interaction strength λ/q 2 (in unit of the bosonic charge q) in agreement with the field theoretical analysis. However, unexpectedly, for small λ/q 2 we find a line of first order quantum phase transitions (FOQPT) between the "Higgs" and the "confined" regions. This line ends, for a finite value of λ/q 2 , in a critical second order quantum phase transition (SO-QPT), above which the two regions are continuously connected through a smooth crossover. Close to FOQPT line the two regions are well separated enabling identification of a "Higgs" mode and its analysis in the continuum limit (for a sufficiently small λ/q 2 ). One can indeed take a different continuum limit to the standard one by approaching, from the Higgs region, the newly discovered SOQPT.
We precisely identify the position and the nature of the new critical point. By assuming Lorentz invariance at the critical point and then using the machinery of conformal field theories (CFT), we can understand the critical point as the direct sum of two non-interacting fields: a free fermionic field describing the Higgs mode and a free bosonic field, a collective mode of the Goldstone modes and the gauge field. Still the a complete characterization of this critical point remains an outstanding challenge as some results do not fit the above picture.
The model. Following [3], we discretize AHM 2 on a finite 1D lattice with L sites (with spacing a) (see [41] for details) arriving at the Hamiltonian (with open boundary conditions): with x = 1/a 2 q 2 . The matter fields {φ j ,φ † j ,Π j ,Π † j } operators act in Hilbert space at sites j, while the gauge-field {L j ,Û j ,Û † j } objects act in Hilbert space defined on the bond linking sites j and j + 1. The operators fulfill the standard commutation relations [φ j , Here we fix x = 2 and characterize the phase diagram on the lattice.
We can define creation and annihilation operators for particles 'a' and anti-particles 'b' asâ j andb j fulfilling . We use the density matrix renormalization group (DMRG) algorithm [39,40,[43][44][45][46] to find the ground state of the Hamiltonian (1). Specifically, we employ a strictly single-site variant of DMRG with subspace expansion [47]. For numerics we limit the occupations of bosonic modes to at most n a 0 = n b 0 = 10 [48].
In absence of external charges, the local U(1) symmetry implies the Gauss lawĜ j = 0, ∀j, where the generators are [3]Ĝ j =L j −L j−1 −Q j , andQ j =â † jâ j −b † jb j encodes the density of dynamical charges. Using the Gauss law, we can integrate-out the gauge-fields in a chain with open-boundary conditions in favor of a long-range potential for the matter fields [4].
Confined and Higgs, two shades of the same phase. The long-range interactions among bosons destroy the phases of the standard 1+1D Bose Hubbard model -see the phase diagram in Fig. 1. For λ/q 2 ≥ 0, the system is in the confined region as far as µ 2 /q 2 ≤ 0. In this region the model has a finite mass gap, and the elementary excitations are mesons, bound pairs of particle-antiparticle. The gauge bosons are in the lowest eigenstate ofL j , so For µ 2 /q 2 0, the system enters the gapped Higgs region where the variance σ 2 (L j ) becomes large. The effective gauge-field mediated tunneling amplitude O tunn = 1 2L j φ † j+1Û jφj + h.c. increases, so that it can distinguish the Higgs region from the confined one.
We can also characterize the two regions by considering the behavior of the entanglement entropy of a block made of l constituents starting from the boundary, defined as where ρ l = Tr l+1,l+2,...,L |ψ ψ| is the reduced density matrix. In the Higgs region, the entanglement entropy is systematically larger than in the confined region. However, since both phases are gapped, the entropy follows area-law scaling with respect to the bipartition size (see [41]). For sufficiently small λ/q 2 , the two regions are separated by a FOQPT line characterized by discontinuous jumps both in the tunneling amplitude O tunn and in S L/2 measured across the central bond [41]. This line terminates at a critical SOQPT at a finite value of λ c /q 2 and µ 2 c /q 2 , identified by the red cross in Fig. 1. We discuss precise location and characterization of this critical point below. Above the SOQPT, the two regions are smoothly connected as revealed by smooth changes in all the physical quantities while moving from one region to the another.
Nature of the critical point. We can precisely locate and characterize the critical point assuming its Lorentz invariance, that implies an applicability of a CFT at low energies. In a CFT, the finite-size scaling of the entanglement entropy of a block of first l consecutive sites in a chain with open boundary conditions and length L is where c is the central charge of the corresponding CFT, b is a non-universal constant and the chord length W is a function of both L and l: W (l, L) = ln 2L π sin(πl/L) [49][50][51]. We pin-point the SOQPT by adapting the idea of the phenomenological renormalization group [52] to the scaling of the entropy in (3) as explained in [53,54]. At the critical point, the c value should be independent of the system's size. For each L, we obtain c fit by fitting (3) to our numerical data for S(l, L). The extracted values in the (µ 2 /q 2 , λ/q 2 )-plane depend on L and become independent of the system size only at the critical point. Our data suggest that the L dependent central charges collapse to a single value at (µ 2 /q 2 = 0.447(1), λ/q 2 = 0.0565(1)) (see Fig. 2). The central charge at the critical point (µ 2 c /q 2 = 0.447, λ c /q 2 = 0.0565) is found to be c = 1.49(1) [55].
The value of the central charge mentioned above suggests that we are not dealing with a minimal model. However, we want to argue here that we are in the presence of the direct sum of two different minimal models, each contributing to a piece of the total central charge, a c f = 1/2 for a free Majorana fermion and a c b = 1 for free boson. This scenario is strongly motivated by the standard Higgs mechanism. The complex Higgs field separates into its amplitude and its phase. The amplitude mode is effectively described by a real λφ 4 theory that undergoes the standard Ising phase transition (the c = 1/2 part). The phase, on the other hand, provides the longitudinal degree of freedom to the photon field. The latter becomes massless at the transition and provides the c = 1 free bosonic part. The value of c = 1.5 furthermore suggests, based on the c theorem, that the two parts should be non-interacting [56].
In order to confirm this scenario we compute the entanglement spectrum that is also known to encode the central charge of the theory [57]. The entanglement spectrum, denoted by ε s , is the spectrum of the entanglement Hamiltonian H l = − log(ρ l ). By assuming a factorized ground state, we should observe that the smallest eigenvalue of H l , ε 0 diverge logarithmically. In particular we should see that [57] By fitting our numerical data to (4), we observe a perfect collapse on the functional form predicted by CFT, but the numerical result for of c eff = 1.20 (1) is not compatible with 1.5 (see Fig. 3(a)). This disagreement between the scaling of the entanglement entropy and that for the entanglement ground state already suggests that the critical point is unusual and exotic in nature.
We thus turn to analyze the operator content of the model by studying the correlation functions of local operators. We should be able to identify a set of primary operators, by studying the large distance two point correlations function that should decay algebraically as The presence of gauge symmetry, however, strongly reduces the set of operators we can consider. Most of the candidates that should couple to primary operators are either trivial (due to the lowdimensionality of the system) or vanishing since they are not gauge invariant. The only non-vanishing operators are indeed Wilson lines terminating on a boson-antiboson pair, and electric field correlations. We also have access to local operators such asφ †φ andΠ †Π that couple both to the real part and the phase of the field. By assuming we are dealing with a CFT, we can use the conformal map that maps the profile of local operators to two points correlation functions on the full plane (see e.g. [58]).
At first we analyze the behavior of L 2 l as function of the chord coordinate W . The numerical results show that L 2 l diverges linearly as a function of W , unveiling thatL 2 behaves as a free-bosonic field. Furthermore, the slope of such linear scaling is found to be 1.20(4)/12, according to the CFT prediction (Eq. (5)). It couples both to the identity operator and one primary with scaling dimension ∆ that comes out to be ∆ = 0.51(2) from the fit, with a being 0.5474 (2) .
matching that of the entanglement ground state energy. Turning to analyzing the profile of Π † lΠ l as a function of W we find that where a and b encode the overlap of the above expectation value with the identity operator and one of the primaries. The numerical data ( Fig. 3(b)) suggest that ∆ 0.5, the conjugate operator to the one that would match to the derivative of the electric field.
Unfortunately, we do not find any operator that couples to the primary of the Ising part of the CFT. Summarizing, our data seem to confirm that we have one part of the system that behaves as free boson and suggest that ∂ xL 2 should have a large overlap with the primary operator (the derivative of the free bosonic field), whileΠ †Π should have a strong overlap with the conjugate primary operator. Now moving away from the critical point, we can use the standard scaling hypothesis to extract the exponent ν from the collapse of the fitted central charge as where f (.) is a continuous function and ν is the corresponding critical exponent. Performing the data collapse according to (6) in the neighborhood of the critical point µ 2 c /λ c (see Fig. 4) we find the critical exponent to be ν = 1/2 ± 0.02 that matches the value observed in the transition from polarized to critical phase in the XX model in a magnetic field [59][60][61]. The same transition can be understood in terms of free bosons that pass from their Fock vacuum to the superfluid regime as the chemical potential exceeds the width of the first band. In our case, the strange thing is that there is no superfluid regime, but just a single critical point where the gauge-boson condense, while away from the critical point our system passes from vacuum to a Mott insulator phase. Now using the standard scaling hypothesis, once we have figured out that ν = 1/2 we can deduce that Π †Π (µ 2 /λ − µ 2 c /λ c ), meaning that β = 1. It is worth pointing out that ν = 1/2, seems inconsistent with a CFT, where we would expect, that d − 1/ν = ∆ with ∆ the thermal critical exponent and d being 2 for the 1+1D quantum system, while we find d − 1/ν = 0. However, all our results so far have been obtained by assuming a full conformal invariance in mapping the correlation function of our finite system to the ones of an infinite plane by means of a conformal transformations.
The appearance of ν = 1/2, together with the failure to identify a local operator that couples to the primary field of the Ising part of the CFT, contrasts with the factorization on the critical point. However, by repeating a similar analysis in a Z 3 gauge theory coupled to bosonic matter we find a c 0.8 + 0.5 [62]. As a result we still believe that the factorization hypothesis is correct, but it requires a further analysis to be appropriately confirmed. In particular, it would be interesting to analyze the system under periodic boundary condition, which, however, is impractical at current computational capabilities using MPS ansatz but may become a possibility using next-generation tensor network algorithms.
Discussion and Conclusions. We have analyzed the phase diagram of AHM 2 on a discrete lattice in 1+1D. We have found two distinct regions, the confined and the Higgs regions that are separated by line of FOQPT that terminates at a SOQPT. Beyond the SOQPT the two regions are connected by a smooth crossover. The presence of a SOQPT allows one to construct an unorthodox continuum limit of the theory that should be described by free fermions and free bosons that do not interact.
This would result in a CFT with central charge c = 3/2, compatible with our numerical result and would have a compelling interpretation in terms of the standard Higgs mechanism -the real part of the complex field undergoes an Ising transition (the c = 1/2 part), while the phase of it provides the transverse degree of freedom to the photon that becomes dynamical and massless (the c = 1 part).
However, further numerical analyses unveil surprizing pieces of the puzzle that do not fit our interpretation.We did not find a local operator that couples to the c = 1/2 part of the CFT. The scaling of the entanglement ground state should follow a similar law to the one of the entanglement entropy. The numerical value of the central charge that we extract from it is c = 1.20(1). We also obtain ν = 1/2 analyzing the collapse of the data for the entanglement entropy close to the critical point, in contrast to the expected ν = 1.
Are we actually observing a Lorentz invariant critical point where the Higgs and photon mode factorize? We believe this is the case as also supported by the presence of a linear dispersion relation witnessed by the non-zero "sound velocity" extracted from a finite-size scaling analysis of the ground state energy [41]. Still our study leaves some questions unanswered. We strongly believe that our paper will open the debate, and that together with the broader scientific community we will soon have a final picture of the mechanism behind this newly observed critical point.
We are grateful to Marcello Dalmonte for the valuable suggestion regarding the analysis of the sound velocity.
After fixing the temporal gauge A t (x, t) = 0, we get the quantum 1+1D Hamiltonian (AHM 2 ) in the continuum aŝ whereÊ x (x),Π(x), andΠ † (x) are the canonical conjugate operators corresponding toÂ x (x),φ(x), andφ † (x) respectively, satisfying the canonical commutation relations: Note: The quantization of the |φ| 4 term may be done in various manners. Here we do not enforce any normal ordering but rather take Following [3], we can straightforwardly discretize the above Hamiltonian on a 1D spatial lattice with the lattice spacing a, and we ultimately arrive at the Hamiltonian (1) in the main text. regions are separated by a first order quantum phase transition (FOQPT) line. To clearly demonstrate the existence and features of the FOQPT, in Figs. 5 we present different observables across the transition as a function of µ 2 /q 2 for fixed λ/q 2 < λ c /q 2 . Specifically, we depict the behaviors of O tunn = 1 2L j φ † j+1Û jφj + h.c. , the energy density and its derivative with respect to the system parameter µ 2 /q 2 respectively in Figs. 5(a)-(c) for different system-sizes L ∈ [20,120]. The average tunneling amplitude O tunn clearly shows discontinuous jumps as we vary µ 2 /q 2 across the phase transition for sufficiently large system-sizes. On other hand, the ground state energy shows non-analytic kinks in its profile at the transition point, such that its derivative (with respect to the parameter µ 2 /q 2 ) manifests a discontinuity. To observe the discontinuous jumps or kinks in the profile of O tunn or of the energy, we need sufficiently large system-size as FOQPTs are associated with large (but finite) correlation lengths, as a result smaller system sizes cannot properly resolve the FOQPT. In Fig. 5(d), we also show the scaling of the entanglement entropy at the FOQPT that features the expected area-law behavior for bipartitions larger than the correlation length.

Entanglement scaling in the gapped phases
As mentioned in the main text, apart from the critical point, the Hamiltonian in every point of the phase diagram is gapped, therefore the entanglement entropy follows asymptotically, for large enough blocks, the standard area-law scaling. To show this, we plot the the entanglement entropy with respect to the cord length W in the confined and the Higgs regions in Fig. 6. In both cases, the entanglement entropy becomes flat with increasing values of W , showing that for large enough blocks it saturates to a value that does not depend on the block size (area-law).

Evidence for Lorentz invariance of the critical point
In the main text, we have analyzed the critical point (µ 2 c /q 2 = 0.447, λ c /q 2 = 0.0565) by means of predictions from conformal field theory (CFT), assuming Lorentz invariance. In order to verify this assumption, we calculate the "sound velocity" at the critical point from the scaling of ground state energy. with the system size L as [63][64][65] where ∞ 0 is the ground state energy density in the bulk and ∞ 1 is the surface free energy in the thermodynamic limit, c is the central charge of the corresponding CFT, and v s is the sound velocity. For a Lorentz invariant system v s must be non-zero, while it vanishes for Lorentz non-invariant critical points with a quadratic dispersion.
To estimate v s , we perform the finite-size scaling of the ground state energy according to Eq. (11) at the critical point (for similar analysis, see [66][67][68][69]). Such a finitesize scaling (Fig. 7) yields cv s = 6.02 ± 0.08 (v s 4 for c = 1.5), confirming the Lorentz invariance of the critical point analyzed in the main text.

Details about numerical simulations
The results reported in the article have been obtained using the density matrix renormalization group (DMRG) method [39,40,[43][44][45][46] with the matrix product state (MPS) ansatz [39,40]. Specifically, we use a strictly single-site variant of DMRG (DMRG3S) with the subspace expansion [47]. In our calculations the gauge fields are integrated out, and thus we do not need to use gaugeinvariant tensor network [70][71][72][73] but we can use standard globally symmetric U(1) MPS [74,75]. This residual global symmetry corresponds to the conservation of the total dynamical charge jQ j .
The results from DMRG simulations, that we report, have been performed with maximum MPS bond dimension of χ = 600. To confirm the convergence of the DMRG sweeps, we continue the DMRG3S iterations until the energy difference in subsequent sweeps falls below 10 −13 . On the other hand, Fig. 8 shows the convergence of different quantities with respect to the bond dimension χ ≤ 600 for systems of sizes L ≤ 120. For example, the energy density converges close to the machine precision within χ ≤ 500 in the gapped regions -confined and Higgs. On the other hand, as expected, the convergence is slower at the critical point due to the diverging correlation length. However, as shown in Fig. 8(b), the precision we attain at the critical point for χ = 500, 600 is sufficient to perform the precise scaling analysis that we report in the main text. Since, our system Hamiltonian does not conserve particle numbers for the individual bosonic species 'a' and 'b', but conserves their number difference, we verify the convergence with respect to the maximum bosonic occupancy n 0 for both the species. Specifically, we vary the bosonic cutoff n 0 in the range [4,15], and check for convergence of different observables (see Fig. 9). Clearly, as seen in Fig. 9(a), even a small n 0 = 4 is sufficient to faithfully capture the confined phase. However, this is not the case at the critical point or in the Higgs region. This is due to the fact that semiclassically the field φ attains non-zero expectation values in the Higgs region. It is to be noted that in our calculations φ is trivially zero since it violates the global symmetry, while φ †φ , and thereby individual n a and n b , attain large expectation values in the Higgs region. Moreover, Fig. 9(b) shows that the smaller values of the bosonic cutoff n 0 are insufficient to capture the critical entropy scaling with respect to the cord length W at the critical point. On the other hand, the entropy profile for n 0 = 10 is essentially identical to the one with n 0 = 15, capturing the proper entropy scal- ing with respect to the cord length. Therefore, Fig. 9 clearly demonstrates that the results obtained with the bosonic cutoff n 0 = 10, provide accurate enough results, and this is the reason we have used this cutoff in our analysis.
We do not exceed n 0 = 10 since the local Hilbert space dimension is d = (n 0 + 1) 2 = 121, which is very large, and the complexity of the DMRG3S algorithm scales as ∼ Ldχ 3 . For this reasons we are forced to consider moderate bond dimensions (χ ≤ 600) that also limit the range of system sizes (L ≤ 120) that we can analyze accurately.
On the other hand, Fig. 8(b) shows that the entanglement ground state energy ε 0 converges faster with respect to the bond dimension than the convergence in the entanglement entropy S (that requires the convergence of the full entanglement spectrum, not only of its ground state). Therefore, for this specific quantity, and just for a comparison, we extend our finite size scaling analysis to the range ε 0 with L ∈ [40,160] keeping the bond dimension fixed at χ = 600. In Fig. 10 we show the same scaling first in the whole range of L ∈ [40, 160] and then for longer chains L ∈ [100, 160]. One could indeed think that the discrepancy between the central charge extracted from the entanglement entropy and the one extracted from the scaling of the first eigenvalue is ultimately a finite-size effect and that this would disappear in the thermodynamics limit.
Unfortunately, in both the cases, the slopes are pretty robust to the choice of system sizes that we use in our fit, and remain fixed to 1.2/12. Even using only the larger systems does not seem to induce a systematic deviation towards the 1.5/12 extracted from the scaling of the entanglement entropy reported in the main text. As a result, the level of accuracy, that we are able to achieve with our computational resources, is extremely good and under control. These extra analysis seems to confirm the discrepancy between the value of the central charge extracted from the scaling of the entanglement ground state with respect to the one extracted from the scaling of the entanglement entropy.