Inhomogeneous phases in the Gross-Neveu model in 1+1 dimensions at finite number of flavors
Julian Lenz, Laurin Pannullo, Marc Wagner, Björn Wellegehausen, Andreas Wipf
IInhomogeneous phases in the
Gross-Neveu model in dimensions at finite number of flavors
Julian Lenz , Laurin Pannullo , Marc Wagner , Bj¨orn Wellegehausen ,Andreas Wipf Friedrich Schiller-Universit¨at Jena, Theoretisch Physikalisches Institut, Fr¨obelstieg 1,D-07743 Jena, Germany Goethe-Universit¨at Frankfurt am Main, Institut f¨ur Theoretische Physik,Max-von-Laue-Straße 1, D-60438 Frankfurt am Main, GermanyApril 1, 2020
Abstract
We explore the thermodynamics of the 1 + 1-dimensional Gross-Neveu (GN) model at finitenumber of fermion flavors N f , finite temperature and finite chemical potential using latticefield theory. In the limit N f → ∞ the model has been solved analytically in the continuum. Inthis limit three phases exist: a massive phase, in which a homogeneous chiral condensatebreaks chiral symmetry spontaneously, a massless symmetric phase with vanishing condensateand most interestingly an inhomogeneous phase with a condensate, which oscillates in thespatial direction. In the present work we use chiral lattice fermions (naive fermions and SLACfermions) to simulate the GN model with 2, 8 and 16 flavors. The results obtained with bothdiscretizations are in agreement. Similarly as for N f → ∞ we find three distinct regimes in thephase diagram, characterized by a qualitatively different behavior of the two-point function ofthe condensate field. For N f = 8 we map out the phase diagram in detail and obtain aninhomogeneous region smaller as in the limit N f → ∞ , where quantum fluctuations aresuppressed. We also comment on the existence or absence of Goldstone bosons related to thebreaking of translation invariance in 1 + 1 dimensions. a r X i v : . [ h e p - l a t ] A p r Introduction
The GN model describes Dirac fermions with N f flavors interacting via quartic interactions in1 + 1 dimensions. It was originally introduced as a toy model that shares several fundamentalfeatures with QCD [1]: it is renormalizable, asymptotically free, exhibits dynamical symmetrybreaking of the Z chiral symmetry, and has a large N f limit that behaves like the ’t Hooft large N c limit of QCD. The particle spectrum and thermodynamics of the theory in the N f → ∞ limitis known analytically. Similarly, the 1-flavor model is equivalent to the 1-flavor Thirring modelwhich can be solved analytically in the massless limit [2] (it has a vanishing β -function). But forintermediate numbers of flavors 1 < N f < ∞ there is – despite many analytical and numericalstudies – no complete understanding of the thermodynamics and particle spectrum.The GN model and related four-Fermi theories in 1 + 1 dimensions have been used in particlephysics, condensed matter physics and quantum information theory. For example, in condensedmatter physics the GN model describes the charge-soliton-conducting to metallic phase transitionin polyacetylene (CH) x as a function of a doping parameter [3]. It is equivalent to the Takayama-Lin-Liu-Maki model [4] which describes the electron-phonon interactions in CH in an effectivelow-energy continuum description, see Ref. [5]. The multi-flavor chiral GN model is relatedto the interacting Su-Schrieffer-Heeger model [6] which is used to investigate cold atoms in anoptical lattice. Four-Fermi models are intensively studied to better understand and classifysymmetry-protected topological phases of strongly interacting systems. For more details werefer to the nice summary in Ref. [7].Recently we have seen a renewed interest in the physics of the GN model at low temperature andhigh baryon density, because it is the region of the QCD phase diagram which is particularlychallenging for first-principles QCD approaches. Corresponding results are only available atasymptotically high densities, where the QCD coupling constant is small so that perturbationtheory can be applied, at vanishing density, where lattice QCD does not suffer from the signproblem, or for unphysically large quark masses, where effective theories exist that mitigate thesign problem. At moderate densities and realistic values of the quark masses, i.e. the regime,which is probed by heavy-ion experiments and which is relevant for supernovae and compactstars, neither approach can be applied. In this regime our current picture of the QCD phasediagram is, thus, mostly based on QCD-inspired models, e.g. the GN model, the Nambu-Jona-Lasinio (NJL) model or the quark-meson model. In early calculations within these models itwas assumed that the chiral condensate is homogeneous, i.e. constant with respect to the spatialcoordinate(s). However, allowing for spatially varying condenates it turned out that there existregions in the phase diagram, where inhomogeneous chiral phases are favored [8,9]. The majorityof the existing calculations have been performed in the limit N f → ∞ or, equivalently, the mean-field approximation (see Ref. [10] for a review and Refs. [11–18] for examples of recent work).Using lattice field theory and related numerical methods and considering N f → ∞ , the GNmodel has been explored in 1 + 1 and 2 + 1 dimensions, the chiral GN model in 1+1 dimensionsand the NJL model in 1 + 1 and 3 + 1 dimension [19–23]. However, a full lattice simulation andinvestigation of the phase diagram of any of these models at finite number of fermion flavors,where quantum fluctuations are taken into account, is still missing. The main goal of the presentwork is to make a step in this direction and to explore, whether such inhomogeneous phases alsoexist in the 1 + 1-dimensional GN model at finite N f .In this work we shall use naive fermions and SLAC fermions to study the multi-flavor GN1odel. These fermion discretizations are all chiral and no fine tuning is required to end up witha chirally symmetric continuum limit. But the theorem of Nielsen and Ninomyia [24] tells usthat we have to pay a price for using (strictly) chiral fermions. And indeed, with naive fermionswe can only simulate 4 , , . . . flavors. With SLAC fermions we can simulate 1 , , . . . flavors, butthe associated Dirac operator is non-local. It has been argued elsewhere that there is no problemwith SLAC fermions for lattice systems without local symmetries, see for example Ref. [25]. Weshall consider GN models with 2 , µ - T parameter space and carefully check for discretization and finite size effects. WithGinsparg-Wilson fermions this would be would be too time-consuming.This paper is organized as follows. In section 2 we summarize some known features of theGN model that are relevant for its thermodynamical properties. These include properties ofthe fermion determinant in the continuum and on the lattice, homogeneous and inhomogeneousphases in the N f → ∞ limit and some comments concerning the spontaneous symmetry breaking(SSB) of translation invariance. In section 3 we discuss different lattice discretizations, the scalesetting and some details of the simulations. Our numerical results are presented in section 4.The main focus concerns the behavior of the two-point function of the order parameter forchiral symmetry breaking and the resulting consequences for the phase diagram in the planespanned by the chemical potential µ and the temperature T . We shall see that the GN modelwith N f = 2, 8 and 16 flavors behaves qualitatively similar to the model in the N f → limit. Inparticular we localize three regions in µ - T parameter space, where the two-point function showsa qualitatively different dependence on the spatial separation. The model with N f = 2 is alsosimulated on rather large lattices with spatial extent N s up to 725 lattice points to carefullyinvestigate the long-range behavior of the correlator. In appendix A we discuss, why the latticeGN model with naive fermions may have an incorrect continuum limit, and how to modify theinteraction term to end up with an (almost) naive fermion discretization with correct continuumlimit. 2 Theoretical basics
The Gross-Neveu model (GN model) is a relativistic quantum field theory describing N f flavorsof Dirac fermions with a four fermion interaction. In this work we investigate this asymptoticallyfree model in 2 spacetime dimensions. The fermions are described by a field ψ = ( ψ , . . . , ψ N f ),the components of which are two-component Dirac spinors. Originally it has been studied inthe 1 /N f expansion. The action and partition function are S ψ = (cid:90) d x (cid:18) ¯ ψ i /∂ψ + g N f ( ¯ ψψ ) (cid:19) , Z = (cid:90) D ¯ ψ D ψ e − S ψ , (1)where the fermion bilinears contain sums over flavor indices, e.g. ¯ ψψ = (cid:80) i ¯ ψ i ψ i .To be able to perform the fermion integration one follows Hubbard and Stratonovich by in-troducing a fluctuating auxiliary scalar field σ to linearize the operator ¯ ψψ in the interactionterm, S σ = (cid:90) d x (cid:18) ¯ ψ i Dψ + N f g σ (cid:19) , Z = (cid:90) D ¯ ψ D ψ D σ e − S σ , (2)where D = /∂ + σ + µγ (3)is the Dirac operator. The four-Fermi term in (1) is recovered after eliminating σ by its equationof motion or equivalently by integrating over σ in the functional integral. In eqs. (2) and (3) wealso introduced a chemical potential µ to study the system at finite fermion density. Expectationvalues of operators O ( ψ, ¯ ψ, σ ) in the grand canonical ensemble are given by (cid:104)O(cid:105) = 1 Z (cid:90) D ¯ ψ D ψ D σ e − S σ O ( ψ, ¯ ψ, σ ) . (4)Note that the integration is over fermion fields, which are anti-periodic in the Euclidean timedirection, with period β = 1 /T , while the auxiliary scalar field is periodic.Integrating over the fermion fields leads to S eff = 12 g (cid:90) d x σ − log det D , Z = (cid:90) D σ e − N f S eff (5)with expectation values of operators O ( σ ) given by (cid:104)O(cid:105) = 1 Z (cid:90) D σ e − N f S eff O ( σ ) . (6)Of particular interest in the present work is the chiral condensate, which distinguishes the differ-ent phases of the GN model. Translation invariance of the integral over d σ x in the (well-defined)functional integral on the lattice implies a Ward-identity which states, that the condensate isproportional to the average auxiliary field, (cid:104) ¯ ψ ( x ) ψ ( x ) (cid:105) = i N f g (cid:104) σ ( x ) (cid:105) , x = ( x µ ) = (cid:18) tx (cid:19) . (7)3lso of interest is the Ward-identity relating the two-point function of the condensate to thetwo-point function of the auxiliary field, (cid:104) ( ¯ ψψ )( x )( ¯ ψψ )( y ) (cid:105) = N f g δ ( x − y ) − (cid:18) N f g (cid:19) (cid:104) σ ( x ) σ ( y ) (cid:105) . (8)In our analysis of the phase diagram at finite temperature and density, the two-point functionof the auxiliary field on the right hand side will play a crucial role (see section 4.3). In this section we study some relevant spectral properties of the Euclidean Dirac operator D with auxiliary field and chemical potential as defined in eq. (3).Clearly, in the continuum the free massless Dirac operator /∂ = γ µ ∂ µ and the partial derivatives ∂ µ have the following properties:a) the differential operator /∂ is anti-hermitean and anti-commutes with γ ∗ = i γ γ ,b) the partial derivatives ∂ µ are real, ∂ ∗ µ = ∂ µ .Later on we shall discretize Euclidean spacetime on a lattice such that ∂ µ turns into a differenceoperator. For most discretizations one does not retain the above properties without introducingdoublers – this is what the celebrated Nielsen-Ninomiya theorem tells us [24]. In the presentwork, however, we shall use chiral lattice fermions with the above properties, naive fermions(having doublers) and SLAC fermions. Now we investigate the spectral properties of the neitherhermitian nor anti-hermitian operator D with eigenvalue equation Dψ = ( /∂ + σ + µγ ) ψ = λψ (9)in the continuum or on the lattice under the assumption that a) and b) hold true. Charge conjugation:
In 2 Euclidean spacetime dimensions there exists a symmetric charge con-jugation matrix C with C − γ µ C = γ µT . Since the (Euclidean) γ -matrices are hermitian we have γ µ ∗ = γ µT , and property b) implies D ∗ = γ µ ∗ ∂ µ + σ + µγ ∗ = C − D C . (10)It follows that all non-real eigenvalues come in complex conjugated pairs ( λ, λ ∗ ) such that det D is real. Hence there is no sign problem for an even number of flavors, since then (det D ) N f isnon-negative. Chiral symmetry:
The four-Fermi term breaks the U A ( N f ) chiral symmetry of the kinetic term.But a discrete discrete Z chiral symmetry still remains under which ¯ ψψ and σ change theirsigns. Under this discrete chiral symmetry the Dirac operator is conjugated with γ ∗ = i γ γ : γ ∗ Dγ ∗ = − /∂ + σ − µγ . (11)Since (on a finite lattice) the number of eigenvalues of D is even, we conclude that the determi-nant is an even function of the auxiliary field,(det D )[ σ, µ ] = (det D )[ − σ, µ ] . (12)4 ermitean conjugation: The Dirac operator in eq. (9) is the sum of one anti-hermitean and twohermitean terms and D † [ σ, µ ] = − D [ − σ, − µ ] . (13)Since the determinant is real (and the number of eigenvalues is even) it follows that det D isinvariant under a simultaneous sign change of σ and µ ,(det D )[ σ, µ ] = (det D )[ − σ, − µ ] . (14)Together with (12) this leads to a determinant which is an even function of the chemical potential,(det D )[ σ, µ ] = (det D )[ σ, − µ ] . (15)Note that for most lattice fermions not all of the above properties hold true, for example forWilson fermions the Z chiral symmetry is explicitly broken. N f limit Before 2003 most field-theoreticians and particle physicists took for granted that in thermalequilibrium translation invariance is realized such that the chiral condensate (cid:104) ¯ ψψ (cid:105) is constant.Assuming translation invariance one can analytically determine the phase diagram of the GNmodel at finite temperature and fermion density in the large - N f limit [27]. But in the con-densed matter community it has been known for a while that the Peierls instability may triggera breaking of translation invariance. This explains, for example, the inhomogenous Fulde-Ferrell-Larkin-Ovchinnikov equilibrium state for ultracold fermions [28,29] (see [30] for a recent review).Subsequently it was shown that for N f → ∞ the relativistic GN model exhibits an inhomoge-neous condensate at low temperature and high density. In a series of interesting papers [8, 9] anexplicit expression for the condensate in terms of Jacobi elliptic functions has been derived. N f -limit For N f → ∞ the saddle point approximation to the functional integral with integrandexp( − N f S eff ) in eq. (5) becomes exact. This means that the condensate (cid:104) σ (cid:105) is identical to thefield σ which minimizes S eff . In particular, if we assume translation invariance then we mayminimize S eff [ σ ] on the set of constant fields. But for a constant σ the regularized action isproportional to the Euclidean spacetime volume, S eff = ( βL ) U eff , U eff = 12 g σ − βL log det Λ D . (16)We renormalize the theory such that the minimum of the effective potential at zero temperatureand zero chemical potential is at some σ >
0. This determines the bare coupling g Λ as functionof the dimensional parameter σ and the momentum-cutoff Λ,1 g = 12 π log (cid:18) σ (cid:19) . (17)5 .0 0.2 0.4 0.6 0.8 1.0 1.2/ T / homogeneouslybroken phase inhomogeneousphasesymmetric phase full phase diagram homomogenous condensateLifshitz point Figure 1: (left) The symmetric phase and the broken phase of the large- N f GN model assuming ahomogeneous condensate σ . The value of the condensate in units of σ is color-coded. For largetemperature or large chemical potential (the green region) the condensate vanishes. (right) Thecorrected phase diagram of the large- N f GN model with homogeneous and inhomogeneous con-densates. The first order line from the Lifshitz point to T = 0 (red dashed line) obtained, whenassuming a homogeneous condensate (see plot on the left side), turns into two second order lines(orange and green). In the region with large µ and low T the condensate is inhomogeneous.The renormalized potential has the simple form U eff = σ π (cid:16) log σ σ − (cid:17) − π (cid:90) ∞ d k k ε k (cid:18)
11 + e β ( ε k + µ ) + 11 + e β ( ε k − µ ) (cid:19) (18)with one-particle energies (cid:15) k = √ k + σ . In accordance with our previous discussion it is aneven function of the auxiliary field and of the chemical potential.The minimizing field as function of the temperature and chemical potential is depicted in Figure1, left. Throughout the present work we use σ to set the scale and thus measure the chemicalpotential, temperature and condensate field in units of σ . The phase diagram shows a symmetricphase with vanishing condensate and a broken phase with homogeneous condensate.The system undergoes a phase transition from the symmetric phase at high temperature or largechemical potential to the broken phase at low T and small µ [27, 31]. At vanishing chemicalpotential the transition happens at T c = e γ /π ≈ . µ , T ) ≈ (0 . , . µ c = 1 / √ ≈ . µ c and T = 0 the condensatejumps from 1 to 0. N f -limit At low temperature and large chemical potential the minimum of S eff does not correspond to ahomogeneous but to a spatially inhomogeneous condensate (cid:104) σ ( x ) (cid:105) . For a time-independent butspatially varying auxiliary field σ the eigenfunctions of D have the form ψ nm ( x ) = e i ω n t ψ m ( x )with Matsubara frequency ω n . Summing over these frequencies in log det D one arrives at the6enormalized effective action S eff [ σ ] = βL π σ (cid:16) log σ σ − (cid:17) + β (cid:16) (cid:88) n : ε m < ε m − (cid:88) m :¯ ε m < ¯ ε m (cid:17) − (cid:88) m : ε m > (cid:16) log (cid:0) e − β ( ε m + µ ) (cid:1) + log (cid:0) e − β ( ε m − µ ) (cid:1)(cid:17) , (19)where the same renormalization prescription as in the homogeneous case has been adopted. The ε m are the real eigenvalues of the hermitean Dirac-Hamiltonian h σ ψ m = ε m ψ m , h σ = γ γ ∂ x + γ σ (20)which appears in the decomposition γ D = ∂ + µ + h σ . (21)The ¯ ε m are the eigenvalues of the Dirac Hamiltonian with constant auxiliary field ¯ σ given by¯ σ = 1 L (cid:90) d x σ ( x ) . (22)The two sums over the negative one-particle energies in the first line of (19) are easily identifiedas difference of two divergent vacuum energies: one for the prescribed auxiliary field σ ( x ) andthe other for the constant reference field σ defined in (22). A heat kernel regularization revealsthat the first line in (19) is UV-finite if and only if the reference field is chosen as in (22). The T - and µ - dependent traces in the second line in (19) are manifestly UV-finite and represent thefinite temperature and density corrections.In the large- N f limit only auxiliary fields which minimize S eff [ σ ] contribute to the functionalintegral in (2). With the Hellman-Feynman formula for the expectation values ε m = (cid:104) ψ m | h σ | ψ m (cid:105) the variational derivative of S eff with respect to σ can be calculated and one ends up with theGap equation12 π σ ( x ) log ¯ σ σ + (cid:88) m : ε m < ψ † m ( x ) γ ψ m ( x ) − (cid:88) m :¯ ε m < ¯ ψ † m ( x ) γ ¯ ψ m ( x )+ (cid:88) m : ε m > (cid:18)
11 + e β ( ε m + µ ) + 11 + e β ( ε m − µ ) (cid:19) ψ † m ( x ) γ ψ m ( x ) = 0 . (23)This renormalized self-consistency equation is a complicated functional equation, whose solu-tions have been investigated at various times in the literature. Most derivations given previouslyderived the regularized gap equation from the regularized trace of the Green-function with barecoupling constant and cutoff parameter [32–35]. Here the point of departure is the renormal-ized effective action (19) with physical scale parameter σ and only finite quantities enter thederivation of the gap equation.To summarize, to calculate the chiral condensate at finite temperature and finite density inthe large- N f limit one must solve the spectral problem for the σ -dependent Dirac Hamiltonian(20) and find a self-consistent solution σ ( x ) of the gap equation (23). At zero temperature andfermion density Dashen et al. indeed could solve the coupled system for the modes ψ n ( x ) and7he scalar field σ ( x ) by using powerful inverse scattering methods [32]. They observed that ascalar field could only solve the gap equation if the solutions of the Dirac equation ψ n are notreflected. Their space-dependent solutions describe n -particle bound states with filled Dirac seaand masses m B = 2 σ π (cid:18) sin θθ (cid:19) , θ = nπ N f , n = 2 , . . . , N f − . (24)Self-consistent solutions at finite temperature and fermion density have been constructed byThies et al. [8, 9] by some (nonlinear) superposition of kink-antikink solutions. They succeededto construct periodic solutions σ ( x ) with associated Bloch waves ψ n ( x ) of the coupled system(20) and (23) in a certain region of the ( T, µ ) phase diagram. The Bloch waves (they are solutionsof the Lam´e equation) and the scalar field σ ( x ) are given in terms of Jacobi’s elliptic functions,see Ref. [9]. The associated Dirac-Hamiltonian has one gap in the spectrum and the periodicand anti-periodic states at the band-edges are given by particular simple Jacobi functions. Thus,the property that h σ shows no reflection for baryon excitations above the vacuum is replaced bythe property of having exactly one band gap in the spectrum if the system has high density.In the large- N f limit where the saddle point approximation to the functional integral (5) be-comes exact the inhomogeneous condensate (cid:104) σ ( x ) (cid:105) minimizes the effective action (19) and thusis given by the solution of the gap equation, i.e. by a Jacobi elliptic function. For points inthe phase diagram where the inhomogeneous solution has a lower effective action as any homo-geneous solution the system is in a inhomogeneous phase. The correct phase diagram in thelarge- N f limit is depicted in Figure 1, right. Note that the metastable phases and first ordertransition line (to guide the eye this line is kept as dashed line) disappear and are replaced bytwo second order transition lines. At low temperature and small chemical potential there is ahomogeneous phase with broken chiral symmetry, at sufficiently high temperature we are in thehomogeneous symmetric phase and at low temperature and large chemical potential we are inthe inhomogeneous phase with an oscillating chiral condensate.The wave length and amplitude of the condensate in the inhomogeneous phase are determinedby the chemical potential or equivalently by the Fermi-momentum and by the temperature. Ifone moves within the inhomogeneous phase towards the symmetric phase, the amplitude of thecondensate vanishes. If one moves towards the homogeneously broken phase, then the wavelength of the condensate increases. In this work we mainly address the question whether thereexists an inhomogeneous phase for a finite number of flavors N f or whether such a phase is anartifact of the large- N f limit. A well-known theorem by N. D. Mermin and H. Wagner in statistical mechanics states thata continuous symmetry cannot be spontaneously broken at finite temperature in 1- and 2-dimensional statistical systems with short range interaction [36]. A similar theorem has beenproven by S. Coleman for relativistic quantum field theory in d ≤ G shouldbe broken to a subgroup H . Then there exists one massless scalar particle for each broken sym-metry (or broken generator) such that there are n BG = dim( G/H ) massless Goldstone bosons.In non-relativistic systems and for spacetime symmetries the situation is more intricate, sincesometimes NGBs have unusual dispersion relations or they are even redundant. • For example, in a ferromagnet and antiferromagnet we may have the same spontaneoussymmetry breaking O(3) → O(2), but in the first case we have only one NGB (the magnon)whereas in the second case there are two NGB. Since quantum field theories at highdensities may be described by quasi-excitations with non-relativistic dispersion relationsa similar reduction of NGB may happen in the high-density GN model. • In addition, for a breaking of spacetime symmetries the simple counting rule does notapply. For example, crystals have phonons for spontaneously broken translations but nogapless excitations for equally spontaneously broken rotations. Again a reduction of thenumber of NGB may happen if we are dealing with spacetime symmetries instead of innersymmetries. • Finally, it may happen that the NGB completely decouple from the rest of the system.Then one may evade the conclusion of Coleman’s theorem about the non-existence ofNGB in 2 spacetime dimensions. This seems to happen in the large N f -limit of the GNmodel [40], where translation invariance is definitely broken for high fermion density.In 1976, Nielsen and Chadha [41] presented a general counting rule of NGBs valid either with orwithout relativistic invariance. They divided the modes into two classes, based on the behaviorof their dispersion relations for small | k | : ε k ∝ (cid:40) | k | n +1 n ≥ , type I | k | n n > , type II . (25)Relativistic modes are of type I and non-relativistic modes are of type II. By examining analyticproperties of correlation functions they showed that n NGB ≤ n BG ≤ n I + 2 n II , (26)where n NGB is the total number of NGBs and n I and n II are the number of type I and typeII NGBs. The number of broken generators n BG = dim( G/H ) agrees with the number of flatdirections of fluctuations of the order parameter.In passing we note that there exists a related, but in general slightly different division of NGBsinto type A and B [42, 43]. It is an algebraic classification based on the Lie algebra of symmetrygenerators. To each pair of non-commuting symmetry generators ˆ Q i , ˆ Q j (the conserved chargesbelonging to the symmetry group G ) a NGB of type B is associated. The NGBs of type A tendto be linearly dispersive for small | k | . There is a simple counting for these NGBs: n A = n BG − rank ρ, n B = 12 rank ρ such that n NGB = n BG −
12 rank ρ , (27)where ρ is the Watanabe-Brauner matrix build from the conserved charges related to the sym-metry group G , ρ ij ∝ (cid:104) [ ˆ Q i , ˆ Q j ] (cid:105) . For more details we refer to Refs. [44, 45].9trictly speaking the above results hold for internal symmetries only. But it is believed thatthe NGBs originating from a spontaneously broken translation symmetry can be treated inessentially the same way as those associated with internal symmetries [44]. This leaves us withthe following scenarios for the finite-temperature GN model at high density: • Only the (abelian) spatial translation symmetry is broken such that the Watanabe-Braunermatrix ρ vanishes. Then the results (27) would imply that there is just one type A NGB.If this would be – as expected – a NGB of type I with relativistic dispersion relation thenwe would be confronted with infrared divergences. The way out could be that it is not oftype I but of type II or that it fully decouples from the system. • Alternatively, if the above results do not apply to the breaking of translation invarianceat high density systems, then we may as well find no NGB or a NGB of type II with non-relativistic dispersion relation ε k ∼ | k | . Its correlation function is not infrared divergentand the problem with the spontaneous breaking would go away.Since the inhomogeneous condensate appears (at least in the large- N f limit) at high density, anon-relativistic dispersion relation seems to be more likely than a relativistic one. Unfortunately,with the available ensembles on lattices of spatial extent up to N s = 725 lattice sites we cannotreliably measure the dispersion relation of the NGB (if it exists) in the model with N f = 2flavors in the infrared and thus cannot decide whether the NGB has a non-relativistic or arelativistic dispersion relation. It may even be that for finite N f there is no SSB of translationinvariance in the strict sense and that the model behaves like a simple atomic liquid, for exampleas liquid argon. Indeed, the correlator of the condensate on large lattices as presented in section4 resembles the radial pair correlation function in an atomic fluid, see the reviews [46, 47]. Ina forthcoming accompanying publication we will further substantiate, by studying the baryonnumber as function of the chemical potential, that the GN model at high density is either acrystal or an extremely viscous fluid. 10 Lattice field theory techniques
The number of lattice sites in temporal and spatial direction are denoted by N t and N s , respec-tively. Consequently, the temperature is given by T = 1 /N t a and the extent of the periodicspatial direction by L = N s a , where a is the lattice spacing.In the following we consider spacetime averages of observables O [ σ ] for given field configurations σ ( x ) with x = ( t, x ), O = 1 N t N s (cid:88) x O [ σ ] . (28)We also compute ensemble averages, (cid:10) O (cid:11) = 1 N conf (cid:88) σ O [ σ ] ≈ Z (cid:90) Dσ e − S eff [ σ ] O [ σ ] , (29)where the sum over σ extends over N conf field configurations σ ( x ) generated by Monte Carlosampling according to e − S eff [ σ ] . The “ ≈ ” sign in eq. (29) becomes the identity “=” in the limit N conf → ∞ .Moreover, we use the discrete Fourier transform either with respect to spacetime˜ F ( k ) = F ( F )( k ) = 1 √ N t N s (cid:88) t,x e − i k · x F ( x ) , k = (cid:16) ωk (cid:17) (30)or space ˜ F ( k ) = F x ( F )( k ) = 1 √ N s (cid:88) x e − i kx F ( x ) , (31)with time and space restricted to the lattice sites, i.e. x = a m , m = 0 , , . . . , N t − , m = 0 , , . . . , N s − . (32)The corresponding momenta are k ≡ ω ∈ (cid:26) πaN t n (cid:27) , k ≡ k ∈ (cid:26) πaN s n (cid:27) . (33)For fermions we impose antiperiodic boundary conditions (BC) in t -direction such that theinteger-spaced n are half-integer valued. For bosons we impose periodic BC in t -directionsuch that the integer-spaced n are integer valued. Both bosons and fermions are periodic in x -direction such that the integer-spaced n are integer valued. We use two different lattice discretizations of fermions, naive fermions and SLAC fermions. Bothdiscretizations have certain advantages but come also with subtleties, which are discussed in thefollowing.In section 4 we present numerical results for both discretizations and find agreement, which weconsider an important cross check. Another comparison of the two discretizations can be foundin section 3.2.3, where we show the N f → ∞ phase diagram with the restriction to homogeneouscondensate σ . 11 .2.1 Naive fermions The naive discretization is at first glance the most straightforward lattice discretization offermions (see e.g. the textbooks [48–50]). In contrast to other common fermion discretizations,e.g. Wilson fermions, the free massless naive fermion action is chirally symmetric, which is anessential and necessary property in the context of this work. Naive fermions, however, lead tofermion doubling according to the Nielsen-Nynomia theorem [24]. Thus, for most applications,e.g. QCD with 2, 2+1 or 2+1+1 quark flavors, naive fermions are not appropriate. In our caseof 1 + 1 spacetime dimensions the number of fermion flavors is restricted to multiples of 4. Thisis not a severe limitation, since we are not interested in a particular number of flavors <
4, butmostly in simulating finite numbers of flavors, e.g. N f = 8 or N f = 16.Besides fermion doubling there are, however, further pitfalls, which might lead to a continuumlimit different from the theory of interest. In the context of the GN model, this was first observedand discussed in Ref. [51]. In appendix A we reproduce the arguments of Ref. [51] and we derivea modification of the straightforward naive discretization of the GN model, which has the correctcontinuum limit. This modified naive lattice action is S GN = N f / (cid:88) i =1 S free [ χ i , ¯ χ i ] + i √ N t N s N f / (cid:88) i =1 (cid:88) x , y ¯ χ i ( x ) W ( x − y ) σ ( y ) χ i ( x ) + N f g (cid:88) x σ ( x ) (34)with the well-known action for naive free fermions coupled to a chemical potential S free [ χ, ¯ χ ] = (cid:88) x (cid:88) ν =0 , i2 a ¯ χ ( x ) γ ν (cid:16) e aµ δ ν χ ( x + a e ν ) − e − aµ δ ν χ ( x − a e ν ) (cid:17) (35)and the auxiliary field summed over neighbors with separation a and √ a √ N t N s (cid:88) y W ( x − y ) σ ( y ) = 14 σ ( x ) + 18 (cid:88) y : | y − x | = a σ ( y ) + 116 (cid:88) y : | y − x | = √ a σ ( y ) , (36)where N f is a multiple of 4 and N t , N s are even integers such that all doublers obey the same BC(for details see appendix A). For all computations with naive fermions presented in the followingsections we use the action (34). For SLAC-fermions the non-local derivatives in the Dirac operator are easily characterized inmomentum space [54], F ( ∂ slac µ ψ )( k ) = i k µ F ( ψ )( k ) (37)with the Fourier transform F as defined in eq. (30). We choose the discrete momenta k µ symmetric to the origin to end up with a real and antisymmetric matrix ∂ µ [55]. This meansthat in spatial direction (with periodic BC) the lattice has an odd number N s of lattice sites,whereas in temporal direction (with antiperiodic BC) the lattice has an even number N t of At an early stage of this work we used the straightforward naive discretization [52, 53]. k = 2 πN t a n , n = N t − , N t − , . . . , − N t , N t even k = 2 πN s a n , n = N s − , N s − , . . . , − N s , N s odd . (38)Both the naive and the SLAC derivative define chiral fermions, for which i /∂ is hermitean andanticommutes with γ ∗ = i γ γ . In contrast to naive fermions, however, there are no doublers forSLAC-fermions. Thus they can be used to study any positive integer number of fermion flavors. We point out that the non-local SLAC-derivative must not be used in gauge theories, where theedge of the Brillouin zone (where k µ jumps) is probed in the functional integral, which leadsto a clash with Euclidean Lorentz invariance [56]. But SLAC-fermions have been successfullyapplied to non-gauge theories, for example scalar field theories [57], non-linear sigma models [58],supersymmetric Yukawa models [57] and more recently to interacting fermion systems [25].For SLAC-fermions the chemical potential µ is introduced as in the continuum theory,¯ ψ ( x ) γ (cid:16) ∂ + µ (cid:17) ψ ( x ) → ¯ ψ ( x ) γ (cid:16) ∂ slac0 + µ (cid:17) ψ ( x ) . (39)Note that the chemical potential µ enters linearly and not via exp( ± aµ ) multiplying a hoppingterm as e.g. for naive fermions (see e.g. eq. (35)) For some observables, for example the fermiondensity, this introduces an additional term ∝ µ in the continuum limit, which can be easilysubtracted, since (in 2 spacetime dimensions) the term is finite and can be calculated analytically.We emphasize that this term is not a lattice artifact – it also exists in continuum theories whenappropriate care is taken in manipulating divergent integrals [59]. After the subtraction isperformed, results obtained with SLAC-fermions converge much faster to the continuum limitthan for other fermion discretizations. A similar observation has been made when using a linear µ for naive fermions in 4 spacetime dimensions [59]. All observables considered in the presentwork need no such subtraction. N f → ∞ and homogeneous con-densate In Figure 2 we show the N f → ∞ phase diagram with the restriction to a homogeneous con-densate σ for naive and for SLAC fermions. The corresponding computations are straight-forward and computationally rather cheap, when using techniques similar to those discussedin Refs. [19–21]. For both discretizations we performed computations for several significantlydifferent values of the lattice spacing, a ≈ . /σ and a ≈ . /σ (naive and SLAC) and a ≈ . /σ (only naive), but similar spatial extent L . When decreasing the lattice spacing,the results obtained with each of the two discretizations approach the continuum result fromRef. [27]. Note, however, that discretization errors for SLAC fermions are almost negligible, i.e.significantly smaller than discretization errors for naive fermions. Because of the sign problem our numerical simulations are restricted to even N f . These findings will be published elsewhere. .0 0.1 0.2 0.3 0.4 0.5 0.6 0.7/ T / naive, a = 0.41/ , N s = 64SLAC, a = 0.41/ , N s = 63naive, a = 0.2/ , N s = 128SLAC, a = 0.2/ , N s = 127naive, a = 0.1/ , N s = 256continuum result Figure 2: N f → ∞ phase diagram with the restriction to homogeneous condensate σ for naiveand for SLAC fermions (three different lattice spacings, a ≈ . /σ , . /σ , . /σ , similarspatial extent L ). The solid grey line represents the continuum result from Ref. [27]. We use a standard RHMC (Rational Hybrid Monte Carlo) algorithm [60] to perform numericalsimulations. In detail we use the implementation described in Ref. [61], which was also used inRefs. [25, 62, 63].
We assume that at chemical potential µ = 0 and temperature T = 0 the system is in a homoge-neously broken phase and use the (positive) expectation value σ = lim L →∞ (cid:104) σ (cid:105) = lim L →∞ (cid:104)| σ |(cid:105) (cid:12)(cid:12)(cid:12) µ =0 ,T =0 (40)to set the scale. In other words, we express all dimensionful quantities in units of σ , e.g. weuse for the chemical potential µ/σ , for the temperature T /σ , etc. Setting the scale via σ wasalso done in previous analytical and numerical studies of the phase diagram of the GN modelin the N f → ∞ limit (see e.g. [8, 9, 19–21]), i.e. expressing dimensionful quantities in units of σ allows a straightforward comparison of our results at finite N f to existing N f → ∞ results.The determination of σ in lattice units is technically straightforward. When increasing thenumber of lattice sites in temporal direction N t as well as in spatial direction N s at fixedcoupling g , the ensemble average (cid:104)| σ |(cid:105)| µ =0 ,T quickly approaches the constant σ . Thus, inpractice, one just has to compute (cid:104)| σ |(cid:105)| µ =0 ,T on a lattice with sufficiently large N t and N s ,where (cid:104)| σ |(cid:105)| µ =0 ,T ≈ σ . This is illustrated in Figure 3 for N f = 8, SLAC fermions and twodifferent g .As e.g. in lattice simulations of 4-dimensional Yang-Mills theory or QCD, the lattice spacing a is a function of the dimensionless coupling g and can be set by choosing appropriate valuesfor g . This is reflected by the two plateau values at small 1 /N t in Figure 3 representing aσ (the lattice spacing in units of σ ), which correspond to g = 0 .
192 (larger lattice spacing) and g = 0 .
161 (smaller lattice spacing). 14 .00 0.05 0.10 0.15 0.20 0.251/ N t a ||| = , T g N s = 63 g N s = 127 g N s = 63 g N s = 127 Figure 3: a (cid:104)| σ |(cid:105)| µ =0 ,T as a function of 1 /N t = T a for µ = 0, N f = 8, SLAC fermions, twodifferent g and two different N s = L/a . The plateau values at small 1 /N t correspond to aσ . To explore the µ - T phase diagram of the 1 + 1-dimensional GN model and its dependenceon the number of fermion flavors N f and to exclude sizable lattice discretization and finitevolume corrections, we generated a large number of ensembles of field configurations σ ( x ).These ensembles are listed in Table 1.For given coupling g , i.e. for fixed lattice spacing a , we vary the temperature T = 1 /N t a bychanging N t , the number of lattice sites in temporal direction. Thus, at fixed g the temperature T can only be changed in discrete steps. The chemical potential µ , on the other hand, is notrestricted in such a way and can be set to any value.The majority of simulations were carried out for N f = 8: • We simulated at several different spatial extents with 31 ≤ N s ≤
128 corresponding to L = N s a to check for finite volume corrections. • We simulated at several different values of the coupling g corresponding to four differentlattice spacings a ≈ . /σ , . /σ , . /σ , . /σ (the lattice spacing is listed in unitsdefined by σ in the column “ aσ ” of Table 1). • We simulated at many different values of the chemical potential, to explore the phasediagram. • We carried out a sizable amount of these simulations using both fermion discretizations, i.e.SLAC fermions and naive fermions, to cross-check our results (the corresponding couplingconstants g have been tuned in such a way, that the simulated lattice spacings are almostidentical).Simulations at N f = 2 and N f = 16 were done with SLAC fermions, but not with naive fermions.For each ensemble between 300 and 10 000 configurations were generated.15 f N s = L/a N t = 1 /T a µ/σ g aσ SLAC, thermodynamics , , . . . , , , , , . . . , ,
80 0 . , . . . , . . . , , ,
127 4 , , . . . , , , , , . . . , ,
80 0 . , . . . , . . . . . . . , , . . . , , , , , . . . ,
64 0 1 . . Naive, thermodynamics , , . . . , , , , . . . ,
64 0 . , . . . , . . . . , . . . , . . . . , . . . SLAC, long-range behavior of the correlation function , , , , , ,
725 80 0 . . . Numerical results
The majority of results shown in the following (section 4.2 to section 4.3) correspond to N f =8, the minimal number of flavors where computations are possible for both naive and SLACfermions. Results for N f = 2 and N f = 16 are presented in section 4.4 and section 4.5. In the 1 + 1-dimensional GN model in the limit N f → ∞ there are three phases, a symmetricphase, a homogeneously broken phase and an inhomogeneous phase (see the discussion in sec-tion 2.3). The structure of the field configurations σ ( x ) generated during our simulations bythe HMC algorithm suggest that there is a similar phase structure at finite N f . At large T thefield σ ( x ) mostly fluctuates around zero, while at small T and small µ either σ ( x ) ≈ + σ or σ ( x ) ≈ − σ . Most interestingly, however, at small T and large µ the field σ ( x ) exhibits spatialperiodic oscillations similar to a cos-function, which might signal an inhomogeneous phase. An example of a typical field configuration at ( µ/σ , T /σ ) ≈ (0 . , . σ ( x ) generated by the Monte Carlo algorithm arecrudely described by the following model: • Inside a symmetric phase σ ( x ) = εη ( x ) . (41) • Inside a homogeneously broken phase σ ( x ) = ± ζσ + εη ( x ) . (42) • Inside an inhomogeneous phase σ ( x ) = A cos (cid:18) π ( x + δx ) λ (cid:19) + εη ( x ) . (43) ε ≥ ζ > A ≥ µ and T . η ( x ) are independentcontinuous random variables with Gaussian probability distributions p ( η ( x )) ∝ exp( − η ( x ) / λ = L/ ( q + δq ) is the wavelength of σ in an inho-mogeneous phase, where q ≥ δq is an integer-valued discreterandom variable with Gaussian probabilities p ( δq ) ∝ exp( − δq / q ). ∆ q (cid:28) q , the width ofthe Gaussian, is another real parameter. δx ∈ [0 , L ) is also a random variable, where it is apriori not clear, what kind of distribution to expect. The distribution could depend on thedetails of the HMC algorithm, and whether translation symmetry is spontaneously broken ornot. Note, however, that the observables we are studying are constructed in such a way thatthey are independent of this distribution. To summarize, the model defined by eqs. (41) to (43)describes field configurations σ ( x ), which fluctuate around 0 in a symmetric phase, around ± ζσ The existence of kink-antikink structures in simulations of the 1 + 1-dimensional GN model with N f = 12 wasobserved already many years ago [64].
20 40 60 80 100 120 x / a t / a -2.0-1.0-0.50.00.51.02.0 / Figure 4: A typical field configuration σ ( x ) /σ generated by the HMC algorithm at large µ/σ ≈ .
450 and small
T /σ ≈ . N f = 8, SLACfermions, a ≈ . /σ , N s = 127). The clearly visible vertical stripes indicate six oscillationsin spatial direction.in a homogeneously broken phase and around a cos-function with varying wavelength λ in aninhomogeneous phase.With this model in mind, which is based on existing results in the N f → ∞ limit [8, 9], wedesigned several observables, which are able to distinguish the three phases. Note that the solepurpose of this model is to provide some guidelines for the construction of observables and todevelop expectations, in which way they characterize the three phases. The model is not usedelsewhere in this work, in particular not for the analysis of our numerical results. σ ( x ) A rather simple observable is Σ = (cid:104) σ (cid:105) σ , (44)the normalized ensemble average of the squared spacetime average of σ ( x ). It is not suited tocharacterize an inhomogeneous phase, but still useful to distinguish a homogeneously brokenphase from a symmetric phase and an inhomogeneous phase. Note that our numerical resultsdo not allow to decide, whether these regions in the µ - T plane are phases in a strict thermody-namical sense or rather regimes, which strongly resemble phases. In any case, throughout thispaper we denote these regions as “phases”.Within the model defined in section 4.1 one finds • inside a symmetric phase Σ = ε /N t N s σ . (45) • inside a homogeneously broken phaseΣ = ζ + ε /N t N s σ . (46)18 .0 0.2 0.4 0.6 0.8 1.0/ T / naive, a , N s = 64 large- N f boundary 0.00.20.40.60.81.0 T / SLAC, a , N s = 63 large- N f boundary 0.00.20.40.60.81.0 Figure 5: Σ in the µ - T plane for naive fermions (left plot) and SLAC fermions (right plot)( N f = 8). The red regions correspond to a homogeneously broken phase and the green regionsto a symmetric and/or an inhomogenous phase. The gray lines represent the N f → ∞ phaseboundary of the homogeneously broken phase [8, 9]. • inside an inhomogeneous phase Σ = ε /N t N s σ . (47)Thus, we expect Σ ≈ ≈
1, inside a homogeneously broken phase.In Figure 5 we show Σ in the µ - T plane for naive fermions and SLAC fermions. The red regionsclearly indicate a homogeneously broken phase, while the green regions represents a symmetricand/or an inhomogenous phase. The two plots are very similar. The main reason for the smalldiscrepancies are lattice discretization errors, which are expected to be significantly larger fornaive fermions than for SLAC fermions (see section 3.2.3). To ease comparison with N f → ∞ results, we included the corresponding phase boundary of the homogeneously broken phase fromRefs. [8, 9]. It is obvious that the homogeneously broken phase at finite N f = 8 is of similarshape, but of smaller size than its analog at N f → ∞ . Such a reduction in size is expected,because at finite N f there are fluctuations in σ ( x ), which increase disorder and, thus, favor asymmetric phase. Note that at small T the boundary between the red and the green regionstarts to deviate significantly from the N f → ∞ boundary and turns towards ( µ, T ) = (0 , N f → ∞ lattice study of the GN model [19]. σ ( x ) In the limit N f → ∞ in the inhomogeneous phase, σ ( x ) is a periodic function of the spatialcoordinate x . It has a kink-antikink structure with large wavelength close to the boundary tothe homogeneously broken phase and is sin-like with smaller wavelength for larger µ . We expecta similar behavior also at finite N f (see also eq. (43)).Since the action S eff is invariant under spatial translations, field configurations, which are spa-tially shifted relative to each other, i.e. σ ( t, x ) and σ ( t, x + δx ), contribute with the same weight19 − S eff to the partition function and, thus, should be generated with the same probability bythe HMC algorithm (see also the discussion on the distribution of δx in section 4.1). Conse-quently, simple observables like (cid:104) σ ( x ) (cid:105) are not suited to detect an inhomogeneous phase in afinite system, because destructive interference should lead to (cid:104) σ ( x ) (cid:105) = 0 in a finite system, evenin cases, where all field configurations exhibit spatial oscillations with the same wavelength. Anobservable, which does not suffer from destructive interference and is able to exhibit informationabout possibly present inhomogeneous structures, is the spatial correlation function of σ ( x ), i.e. C ( x ) = (cid:104) σ ( t , x ) σ ( t , (cid:105) = 1 N t N s (cid:88) t,y (cid:10) σ ( t, y + x ) σ ( t, y ) (cid:11) . (48)The equality holds in thermal equilibrium and a finite box of length L with periodic boundaryconditions since (cid:104) σ ( t, x + y ) σ ( t, y ) (cid:105) does neither depend on t nor on y . Actually, our HMC algo-rithm is able to sample all field configurations and, thus, to produce y -independent expectationvalues (cid:104) σ ( t, x + y ) σ ( t, y ) (cid:105) . We use the sum over t and y in Eq. (48) to decrease statistical errorsin the Monte Carlo average.The correlator C ( x ) is our main observable to detect and to distinguish the three expectedphases, in particular an inhomogeneous phase. In contrast to the typical exponential decay ofcorrelation functions, C ( x ) is expected to oscillate in an inhomogeneous phase, i.e. C ( x ) shouldbe positive, if x/λ is close to an integer, and negative, if x/λ is close to a half-integer, where λ denotes the wavelength of the spatial periodic structure of σ ( x ). Such oscillations are alsofound in the N f → ∞ limit [65]. The expectation is also supported by analytical calculationswithin our model defined in section 4.1, where • inside a symmetric phase C ( x ) = ε δ x, (49) • inside a homogeneously broken phase C ( x ) = ( ζ σ ) + ε δ x, (50) • inside an inhomogeneous phase C ( x ) ≈ A ϑ ( x/L, i/ π ∆ q ) ϑ (0 , i/ π ∆ q ) cos (cid:18) πqxL (cid:19) + ε δ x, (51)with the Jacobi ϑ function ϑ ( z, τ ) = 1 + 2 ∞ (cid:88) n =1 e iπn τ cos(2 πnz ) . (52)The cos-term in eq. (51) leads to oscillations with wave length L/q , while the factor includingthe ϑ function causes a damping of these oscillations for increasing separations x . This dampingis due to the random fluctuations of the wave number q + δq entering eq. (43) via λ , which cause An alternative would be to break translation invariance explicitly, for example by imposing Dirichlet boundaryconditions on σ ( x ). x . The damping is strong for large fluctuations, i.e. large ∆ q ,and not present in the limit ∆ q →
0. Note that the result for C ( x ) inside an inhomogeneousphase, i.e. the right hand side of eq. (51), is independent of the distribution of the randomvariable δx introduced in section 4.1.Of similar interest as C ( x ) is its Fourier transform˜ C ( k ) = F x ( C )( k ) (53)(see eq. (31)). The expected behavior is the following: • inside a symmetric phase ˜ C ( k ) is rather small and smooth without any pronounced peak. • inside a homogeneously broken phase ˜ C ( k ) has a pronounced peak at k = 0 and is rathersmall and smooth at k (cid:54) = 0. • inside an inhomogeneous phase ˜ C ( k ) has pronounced peaks at k = ± q , where q is relatedto the wavelength of the spatial oscillations of C ( x ) via λ = L/q .This expectation is in agreement with results obtained within our model from section 4.1, where • inside a symmetric phase ˜ C ( k ) = 1 √ N s ε , (54) • inside a homogeneously broken phase˜ C ( k ) = (cid:112) N s ( ζσ ) δ k, + 1 √ N s ε , (55) • inside an inhomogeneous phase˜ C ( k ) ≈ √ N s A ϑ (0 , i/ π ∆ q ) (cid:18) exp (cid:18) − ( k − q ) q (cid:19) + exp (cid:18) − ( k + q ) q (cid:19)(cid:19) + 1 √ N s ε . (56)Exemplary results for C ( x ) and ˜ C ( k ) are shown in Figure 6 inside the symmetric phase(( µ/σ , T /σ ) ≈ (0 , . µ/σ , T /σ ) ≈ (0 , . µ/σ , T /σ ) ≈ (0 . , . µ/σ , T /σ ) ≈ (0 . , . a ≈ . /σ and a ≈ . /σ . Those corresponding to the finer lattice spacing are closer to the SLACresults, where a ≈ . /σ . We interpret this as indication that both discretizations agreein the continuum limit. Moreover, on a qualitative level there is agreement with our crudeexpectations summarized by eqs. (49) to (51) and eqs. (54) to (56), when the parameters inthese equations are chosen appropriately. In particular the plots in the lower half of Figure 6clearly indicate the existence of an inhomogeneous phase. C ( x ) exhibits cos-like oscillations withdecreasing wavelength λ for increasing µ , as observed in the N f → ∞ limit. This is also reflected21 ymmetric phase x C / ( / , T / ) (0.000, 0.993) k / C / N s ( / , T / ) (0.000, 0.993) SLAC, a , N s = 63naive, a , N s = 64 homogeneously broken phase x C / ( / , T / ) (0.000, 0.083) k / C / N s ( / , T / ) (0.000, 0.083) SLAC, a , N s = 63naive, a , N s = 64 inhomogeneous phase x C / ( / , T / ) (0.700, 0.083) k / C / N s ( / , T / ) (0.700, 0.083) SLAC, a , N s = 63naive, a , N s = 64naive, a , N s = 128model0 1 2 3 4 5 6 7 8 x C / ( / , T / ) (0.900, 0.083) k / C / N s ( / , T / ) (0.900, 0.083) SLAC, a , N s = 63naive, a , N s = 64naive, a , N s = 128model Figure 6: The spatial correlation function C ( x ) /σ (left column, where the x − range extends tohalf the box length) and its Fourier transform ˜ C ( k ) /σ N s (right column) for SLAC fermions (bluedots) and naive fermions (orange and green dots) in the symmetric, the homogeneously brokenand the inhomogeneous phase ( N f = 8) together with model expectations for the inhomogeneousphase (eqs. (51) and (56); gray curves). 22 .0 0.2 0.4 0.6 0.8 1.0/ T / naive, a , N s = 64 large- N f boundary 0.20.00.20.40.60.81.0 C m i n / T / SLAC, a , N s = 63 large- N f boundary 0.20.00.20.40.60.81.0 C m i n / Figure 7: C min /σ in the µ - T plane for naive fermions (left plot) and SLAC fermions (right plot)( N f = 8). The red regions correspond to a homogeneously broken phase, the green regions to asymmetric phase and the blue regions to an inhomogenous phase. The gray lines represent the N f → ∞ phase boundaries [8, 9].by the symmetric pair of peaks of ˜ C ( k ) at the corresponding wave numbers q = L/λ . The graycurves in these plots represent the model expectations (eqs. (51) and (56)) with parameters A , q , ∆ q and (cid:15) determined by fits to the lattice results for C ( x ) and ˜ C ( k ).A straightforward calculation leads to˜ C ( k ) = (cid:28) N t √ N s (cid:88) t | ˜ σ ( t, k ) | (cid:29) , (57)where ˜ σ ( t, k ) = F x ( σ )( t, k ) (58)is the Fourier transform of σ ( x ) with respect to the spatial coordinate (see eq. (31)). Theabsolute values | ˜ σ ( t, k ) | are invariant under spatial translations x → x + δx , because F x ( σ ( t, x + δx ))( t, k ) = e − i kδx F x ( σ ( t, x ))( t, k ) . (59)This shows again that both C ( x ) and ˜ C ( k ) do not suffer from destructive interference, as alreadydiscussed at the beginning of this subsection. Moreover, eq. (57) shows in an explicit way thatthe Fourier transformed correlation function ˜ C ( k ) also provides information about the absolutevalues of the Fourier coefficients of the field σ ( x ). In particular the peaks in ˜ C ( k ) at non-vanishing k in the plots in the lower half of Figure 6 indicate that inside an inhomogeneousphase strong oscillations with the same wavelength are present in the majority of the generatedfield configurations.From Figure 6 one can see C min = min x C ( x ) (cid:29) ≈ (cid:28) . (60)Thus, the minimum of the correlation function C ( x ) is suited to plot a crude phase diagram asshown in Figure 7 both for naive and for SLAC fermions.23 .0 0.2 0.4 0.6 0.8 1.00.00.20.40.60.81.0 T / a = 0.4100/ a = 0.250/ a = 0.195/ N s = N s = N s = N s = C m i n / Figure 8: C min /σ in the µ - T plane for three different values of the lattice spacing a (thecolumns) and four different numbers of lattice sites in spatial direction N s (the rows) ( N f = 8,SLAC fermions). 24 omogeneously broken ↔ symmetric (blue points)inhomogeneous ↔ symmetric (orange points) T / C m i n / = 0/ homogeneously broken ↔ inhomogeneous C m i n / T / Figure 9: C min /σ as function of T for µ = 0 and for µ/σ ≈ .
600 (left plot) and C min /σ asfunction of µ for T /σ ≈ .
102 (right plot) ( N f = 8, SLAC fermions, a ≈ . /σ , N s = 63).The red region indicates a homogeneously broken phase, the green region a symmetric phase andthe blue region an inhomogeneous phase. As before, results obtained with these two differentfermion discretizations are in fair agreement. Moreover, the phase diagram is qualitativelysimilar to the N f → ∞ phase diagram from Refs. [8, 9], whose phase boundaries are also shownin Figure 7. The homogeneously broken phase and the inhomogeneous phase are, however,somewhat smaller for finite N f than for N f → ∞ , presumably because quantum fluctuationsat finite N f increase disorder and, thus, favor a symmetric phase. Note, however, that C min is not the expectation value of a product of local operators as for example C ( x ), which is thetwo-point function of the order parameter. In general one must be cautious, when using non-local quantities like C min , since they can fake non-existing phase transitions [66, 67]. But in thepresent case transition lines have been localized by C min as well as the correlator C ( x ) of thelocal field σ .We also checked the stability of the phase diagram with respect to variations of the latticespacing and the spatial volume. To this end we performed simulations using SLAC fermions atthree different values of the lattice spacing, a ≈ . /σ , . /σ , . /σ (the columns inFigure 8), and for four different numbers of lattice sites in spatial direction, N s = 31 , , , N f → ∞ the phase boundaries between the three phases are of second order. Inthe following we present selected results for finite N f and discuss, whether there are also phasetransitions or rather crossovers. • For the transition between the symmetric and the homogeneously broken phase we com- As already pointed out on page 18, we call the three regions in µ - T plane “phases”, although they may notbe phases according to the Ehrenfest classification. C min as function of the temperature T for vanishing chemical potential µ = 0 (seeFigure 9, left plot, blue points). There is a rapid decrease of C min at around T /σ = 0 . N f → ∞ result and, thus, might indicate thatthere is also a second order phase transition at finite N f . • For the transition between the homogeneously broken and the inhomogeneous phase wecomputed C min as function of the chemical potential µ for rather low temperature T /σ ≈ .
102 (see Figure 9, right plot). We observe a rapid decrease of C min at around µ/σ = 0 . N f → ∞ case. • As can also be seen from the phase diagram in Figure 7, the transition between thesymmetric and the inhomogeneous phase is somewhat washed-out. This is also reflected byFigure 9, left plot, where the orange points represent C min as function of the temperature T for chemical potential µ/σ ≈ . As has been pointed out earlier, these conclusions may not be fully coherent since C min is anon-local quantity. C ( x ) We have argued in section 2.4 that SSB of a continuous (spacetime) symmetry in 1+1 dimensionsis a delicate issue. To decide, whether there are NGB and if so, what type of NGB, we investigatethe long-range correlations of the GN model with N f = 2 flavors. The question is, whether thelong-range order (which in the present context means that C ( x ) oscillates with a constant non-zero amplitude for arbitrarily large | x | ), which is necessary to form a crystal at N f = ∞ , remainsat finite N f long-range or becomes almost long-range `a la Berezinskii, Kosterlitz and Thouless(BKT) [68, 69] (in which case the amplitude decreases with distance like an inverse power).Indeed, by studying the long-range behavior of the SU( N f ) Thirring model (this is a four-Fermitheory with current-current interaction) in 1 + 1 dimensions, Witten argued that for finite N f the correlations have almost long-range order, (cid:104) ¯ ψψ ( x ) ¯ ψψ (0) (cid:105) ∼ | x | /N f , | x | → ∞ , (61)such that only for N f → ∞ the continuous chiral symmetry is broken, i.e. (cid:104) ¯ ψψ (cid:105) (cid:54) = 0 [70]. Thisway the system circumvents the no-go theorems for SSB of continuous inner symmetries. In whatfollows we try to answer the question whether a similar mechanism is at work for translationsymmetry in the GN model.To detect SSB of translation invariance directly we could break translation symmetry explicitlyin a finite box with periodic BC, for example by adding a term ε ( x ) σ ( t, x ) to the Lagrangian,perform the infinite volume limit and finally remove the source ε ( x ). Assuming clustering inthermal equilibrium andlim ε → lim L →∞ (cid:104) σ ( t, x ) σ ( t, (cid:105) ε = lim L →∞ (cid:104) σ ( t, x ) σ ( t, (cid:105) ε =0 = C ( x ) (62) A more detailed investigation of the long-range behavior of C ( x ) presented in section 4.4 points towards aphase transition.
26e conclude that C ( x ) is for large | x | proportional to the condensate (calculated with first L → ∞ and afterwards an adapted ε → C ( x ) = (cid:10) σ ( t, x ) σ ( t, (cid:11) → (cid:10) σ ( t, x ) (cid:11)(cid:10) σ ( t, (cid:11) for | x | → ∞ . (63)In the inhomogeneous phase we can write C ( x ) = A ( x ) C periodic ( x ) , (64)where A ( x ) is the non-increasing amplitude function, while C periodic ( x ) represents the periodicoscillations. If the system forms a crystal, the amplitude function A ( x ) should approach a non-zero constant for sufficiently large separations | x | . In case the system has almost long-rangeorder `a la BKT, the amplitude function decreases with | x | as an inverse power. To distinguishthe two scenarios we study C ( x ) for small N f = 2, to detect a possible deviation of A ( x ) fromthe asymptotically constant behavior in the large- N f limit (for small N f quantum fluctuationsmight be strong enough to change long-range into almost long-range order as it happens in thechirally invariant SU( N f ) Thirring model for finite N f [70]). Thus, we expect one of the followingamplitude functions:1. In a BKT-like phase without SSB we expect the amplitude function A ( x ) to have thefollowing behavior for large | x | : A ( x ) = A BKT ( x ) ∼ α | x | β + . . . (65)The dots indicate sub-leading terms and terms arising from the finite spatial extent of thesystem.2. If there is SSB of translation invariance, C ( x ) oscillates with constant non-zero amplitudeat large | x | , where the short-ranged contributions from excited states are suppressed. Theamplitude function would then be A ( x ) = A SSB ( x ) ∼ γ + αe − m | x | + . . . or A ( x ) = A SSB (cid:48) ( x ) ∼ γ + α | x | β + . . . (66)depending on whether the excitations over the oscillating condensate are massive or mass-less. If the NGB decouple from the system, the amplitude function approaches the constantamplitude γ (cid:54) = 0 exponentially fast. If not, A ( x ) will approach γ (cid:54) = 0 with an inverse powerof | x | .Before discussing the long-range behavior of C ( x ) we present the phase diagram from C min for N f = 2 in Figure 10. Again we recognize the same phases as in the large- N f limit: ahomogeneously broken phase which (as expected) is significantly smaller than for N f = 8 and N f → ∞ , a symmetric phase for sufficiently large temperature and a region, where C min is clearlynegative. One unexpected and striking feature of the latter phase is that the temperature rangewith negative C min grows with increasing µ , which is qualitatively different from the situationat large N f . This already happens for N f = 8 in a small region of parameter space (see Figure 8)but is so pronounced at N f = 2 that up to µ/σ = 1 . µ axis.27 .0 0.2 0.4 0.6 0.8 1.0/ T / large- N f boundary 0.20.00.20.40.60.81.0 C m i n / Figure 10: C min /σ in the µ - T plane for N f = 2. The red region corresponds to a homogeneouslybroken phase, the green region to a symmetric phase and the blue region to an inhomogenousphase. The grey lines represent the N f → ∞ phase boundaries [8, 9].Now we investigate the long-range behavior of the amplitude function A ( x ) defined in eq. (64)for N f = 2 on lattices with a rather large number of sites in spatial direction. Figure 11 showsthe correlator C ( x ) and its Fourier transform ˜ C ( k ) for ( µ/σ , T /σ ) = (0 . , . N s up to 725 corresponding to spatial extents up to L ≈ . /σ . We observe 32 periods ofstatistically significant oscillations in C ( x ) over the whole range of separations and a pronouncedpeak of ˜ C ( k ) at the corresponding wave number. The position of the peak is essentially the samefor all L demonstrating once more that the wave length is independent of the spatial extent.The amplitude function A ( x ) is extracted from the peaks of the correlation function C ( x ) (whichwe identified using the scipy.signal.find peaks method [71] with prominence=0.01 ) for all N s , as exemplified for N s = 725 in the left plot of Figure 12. The peaks of C ( x ) for various N s are depicted in the right plot of Figure 12. There is a rapid drop for small separations x that flattens out for asymptotically large x . We performed χ -minimizing fits of symmetrized x C /
40 60 80 1000.0250.0000.025 0 1 2 3 4 5 6 7 8 k / C N s / N s = 65 N s = 125 N s = 185 N s = 255 N s = 525 N s = 7250.65 0.70 0.7514161820 Figure 11: The spatial correlation function C ( x ) /σ (left) and its Fourier transform ˜ C ( k ) √ N s /σ (right) at ( µ/σ , T /σ ) = (0 . , . N f = 2 and N s = 65 , , , , ,
20 40 60 80 100 120 140 x C / dataPeaks SSB , x min = 15SSB', x min = 25BKT , x min = 3780 100 120 1400.010.000.01 0 20 40 60 80 100 120 140 x A N s = 65 N s = 125 N s = 185 N s = 255 N s = 525 N s = 725 Figure 12: The left plot shows the correlation function C ( x ) /σ at ( µ/σ , T /σ ) = (0 . , . N s = 725, the extracted peaks and fits with the functions defined in eqs. (65) and (66)(the parameters obtained by these fits are shown in Table 2). For clarity interpolating linesare shown instead of the data points. The right plot shows the extracted peaks for N s =65 , , , , , x min α β resp. m γ χ SSB 15 0 . ± .
026 0 . ± . . ± . . . ± . · . ± .
50 0 . ± . . . ± .
082 0 . ± .
13 - 0 . A ( x ) in eqs. (65) and (66) obtained by fits tothe extracted peaks for N s = 725.versions of the expectations for the amplitude function A ( x ) (eqs. (65) and (66)) to the extractedpeaks for N s = 725 (see Figure 12, left plot). We only used data points with x ≥ x min in thefitting procedure, because all three fit functions are expected to model the amplitude functionfor large | x | . Since χ as a function of x min is almost constant for large | x | (see Figure 13,upper plot), we chose (independently for each of the three fits) the minimal value x min , where χ is consistent with the asymptotic constant, i.e. that value, where the constant behaviorsets in. In the lower plot of Figure 13 we show the extracted parameters γ appearing in thetwo fit functions in eq. (65) as functions of x min . It is reassuring that both parameters γ areessentially independent of x min as long as the corresponding χ is small. Since there are massiveexcitations in the SSB model and massless excitations in the SSB’ and BKT models, one shouldnot expect to find the same x min for the three models, but a smaller x min for the SSB model,where the excitations are short-range. This expectation is confirmed by our numerical results(see column “ x min ” in Table 2).The fit results for the parameters of the amplitude functions A ( x ) from eqs. (65) and (66) for N s = 725 are collected in Table 2. There are several points to note: • The SSB model admits the most stable fits with resulting parameters almost independentof x min and the initial values used in the fitting algorithm. γ is cleary different from zero. • Fits for the SSB’ model are less stable, which is reflected by the large uncertainties obtained29 .20.40.6 r e d SSBSSB'BKT10 15 20 25 30 35 40 45 x min 0 Figure 13: Top: χ as a function of x min for all three models ( N s = 725). Bottom: Extractedparameters γ for the SSB and SSB’ model as functions of x min .for α and β . γ , however, can be determined in a reliable way and the result is againdifferent from zero. Moreover, it is in excellent agreement with the corresponding resultfor the SSB model. It is also interesting to note that the SSB’ model, which differs fromthe BKT model by the additive constant γ , leads to significantly smaller χ (and γ (cid:54) = 0)than the BKT model, when the same x min is used. • The BKT model is only able to describe the amplitude function for rather large | x | , i.e.the corresponding x min is significantly larger than for the SSB model and SSB’ model.The fit result for the exponent is β = 0 . ± .
13, which is similar to the correspondinganalytically known value β = 1 /N f = 1 / N f ) Thirring model.To summarize, it seems that the excitations are probably massive, because the SSB model isable to describe the extracted peaks for significantly smaller separations x , than it is possiblewith the SSB’ model or the BKT model. When allowing for a constant γ (as in the SSB modeland the SSB’ model), the fits lead to stable and clearly non-zero results. Thus, our current datais best described by the SSB scenario, where the NGB completely decouples from the system.However, this does not imply that this is the physical situation in the thermodynamic limit,since we saw that even N s = 725 lattice sites in spatial direction are not sufficient, to rule outthe 1 / | x | β almost long-range behavior of the BKT scenario. In other words, a constant behavior ∝ γ and an inverse power ∝ / | x | β are very similar for large | x | , in particular in a periodicspatial volume, such that even larger lattices or higher accuracy is needed to clearly distinguishbetween these scenarios.Despite the fact that we could not fully reveal the nature of the inhomogeneous phase, we canstill argue that there exists a phase transition between the inhomogeneous low-temperaturephase and the symmetric high temperature phase. This can be seen from Figure 14, wherewe show the spatial correlation function C ( x ) together with cosh-fits for three different ( µ, T )in the symmetric phase. It is evident that the cosh-functions perfectly fit the data points,which indicates that in the symmetric phase the excitations are massive. In contrast to that,at low temperature in the inhomogeneous phase C ( x ) is long-range or almost long-range (seeFigure 11). Since C ( x ) behaves qualitatively different in the inhomogeneous phase and in thesymmetric phase, i.e. long-range or almost long-range versus exponentially decaying, we expecta phase transition, either a symmetry-restoring transition or a BKT-like transition from thelow-temperature inhomogeneous phase to the high-temperature symmetric phase.30 x C / datafit T / = 0.610, / = 0.00 T / = 0.305, / = 0.10 T / = 0.305, / = 0.30 Figure 14: The spatial correlation function C ( x ) /σ for three different ( µ, T ) in the symmetricphase and fits with A cosh( m ( x − L/ A cosh( m ( x − L/ ≤ xσ ≤
16 because of large systematic errors due to autocorrelations. T / N f = N f = 2 N f = 8 N f = 16 Figure 15: Σ as a function of T for µ = 0 and N f = 2 , , , ∞ (SLAC fermions, a ≈ . /σ , N s = 63). N f → ∞ results with computations at finite N f In sections 4.2 to 4.4 we have presented results at N f = 2 and N f = 8, which are similar to theanalytically known N f → ∞ results [8, 9], e.g. for the phase diagram. To check and to confirmthat results at finite N f approach for increasing N f the N f → ∞ results, we also performedsimulations at N f = 16. An exemplary plot is shown in Figure 15, where Σ is shown as afunction of the temperature T for vanishing chemical potential µ = 0 and N f = 2 , , , ∞ .While results for N f = 2 agree with the N f → ∞ result only for rather small T /σ , there isagreement also for larger T /σ , when N f is increased, indicating that one can approach theanalytically known N f → ∞ results with computations at finite N f .31 Conclusions
In the present work we could localize three regimes in the space of thermodynamic controlparameters T and µ , in which the two-point function of the order parameter shows qualitativelydifferent behaviors. We spotted a homogeneously broken phase, a symmetric phase and a regionwith oscillating correlation function C ( x ) defined in eq. (48). The results of our Monte Carlosimulations with two different types of chiral fermions for systems with 2, 8 and 16 flavorson lattices with sizes up to N t = 80 and N s = 725 were presented, analyzed and discussedin the main body of the text. Although we could not answer the question, whether in GNmodels with a finite number of flavors translation invariance is spontaneously broken at low T and large µ , or whether the system is in a Berezinskii-Kosterlitz-Thouless like phase, weclearly spotted a low-temperature and high-density region, where the model exhibits oscillatingspatial correlators. The wave-length of the spatial oscillation is determined by the chemicalpotential and temperature and not by the lattice spacing and the spatial extent of the lattice,i.e. spatial oscillations are neither a lattice artifact nor a finite size effect. We argued thatthere is a transition between the inhomogeneous phase and the symmetric phase which couldbe an infinite order transition (according to the Ehrenfest classification). In an accompanyingpaper we shall demonstrate that the ratio of the system size and the dominant wave lengthof the condensate oscillations is equal to the number of baryons in the systems. This furthersubstantiates the physical picture that the GN model in equilibrium at low temperature andlarge fermion density either forms a crystal of baryons or a viscous fluid of baryons. In thiswork we also showed that the amplitude of oscillations stays constant or decays very slowly assuggested by a related result [70]. The first behavior is expected for a baryonic crystal the secondbehavior for a viscous baryonic fluid. To better understand, how the long range behavior atlow temperature and high density does not clash with the absence of Nambu-Goldstone bosonsin 1 + 1 dimensions, needs further high-precision results on the two-point function of the orderparameter. If the dispersion relation is non-relativistic or if the massless modes fully decouplefrom the system then there should be no problem.Independent of whether the oscillating correlator C ( x ) points to a baryonic crystal or a baryonicliquid for finite N f we have seen that mean-field/large- N f approximations may keep more infor-mation on the physics at finite N f than one might expect. This is reassuring since in particlephysics and even more so in solid state physics we often rely on mean-field type approximations.An important question is, of course, whether our results have any relevance for QCD at finitebaryon density. On the one hand, we established that the interpretation as baryonic matter isnot spoiled by taking quantum fluctuations into account. On the other hand, although recent numerical investigations of four-Fermi theories in 2 + 1 dimensions and for N f → ∞ spottedinhomogenous condensates [22], the spatial modulation is related to the cutoff scale and seems todisappear in the continuum limit [23]. Clearly, if this happens then we cannot expect a breakingof translation invariance for a finite number of flavors. Thus, extending our lattice studies tohigher dimensions is of relevance for QCD. Simulations of interacting Fermi theories in 2 + 1dimensions are under way and we hope to report on our findings soon.32 Lattice discretization of the GN model with naive fermions
A.1 Free naive fermions
The action of free naive fermions with chemical potential µ is given in eq. (35). The Fourierrepresentations of the fermion fields are χ ( x ) = 1 √ N t N s (cid:88) k e i k · x ˜ χ ( k ) , ¯ χ ( x ) = 1 √ N t N s (cid:88) k e − i k · x ˜¯ χ ( k ) , (67)where the discrete 2-momenta k = ( k , k ) are in the first Brillouin zone, − π ≤ k µ a ≤ π , andare chosen such that BC in t direction are antiperiodic and in x direction are periodic (seesection 3.1, in particular eq. (30)). Inserting these Fourier representations into eq. (35) leads to S free [ ˜ χ, ˜¯ χ ] = − (cid:88) k ˜¯ χ ( k ) (cid:0) γ ˚ k + γ ˚ k (cid:1) ˜ χ ( k ) , (68)where we abbreviated˚ k = cosh( µa ) sin( k a ) a − i sinh( µa ) cos( k a ) a , ˚ k = sin( k a ) a . (69)In the limit a → k and k can be restricted to the “soft modes”, where both | ˚ k a | (cid:28) | ˚ k a | (cid:28)
1. There are four regions of soft modes in the first Brillouin zone, andthey are denoted by R uv with u, v ∈ { , } . The momenta of the soft modes in region R uv arein the neighborhood of the four momenta k uv = πa (cid:16) uv (cid:17) , u, v ∈ { , } , (70)at which (for µ = 0) the lattice momenta ˚ k and ˚ k vanish. For k ∈ R uv we have˚ k = ( − u k − i µ + O ( a ) , ˚ k = ( − v k + O ( a ) . (71)Now we define the soft modes in the four regions according to ˜ χ uv ( k ) = ˜ χ ( k uv + k ) (and analogousfor ¯˜ χ ) with small | k µ a | . Neglecting the O ( a ) corrections in (71) we can approximate the freelattice action (68) by S free [ ˜ χ, ˜¯ χ ] ≈ − (cid:88) u,v (cid:88) k ∈R uv ˜¯ χ uv ( k ) (cid:16) γ (( − u k − i µ ) + γ ( − v k (cid:17) ˜ χ uv ( k ) . (72)This short calculation exhibits the well-known fermion flavor doubling for each spacetime di-mension. It also shows that both N t and N s must be even to obey anti-periodic boundaryconditions in t direction and periodic boundary conditions in x direction for each of the fourfermion flavors.It is important to note that the action (72) differs in a couple of minus signs in front of the γ matrices for flavors ( u, v ) (cid:54) = (0 ,
0) from the corresponding continuum expression for four freefermion flavors. These minus signs can be eliminated by changing field coordinates via˜ χ uv = ( γ ) u ( γ ) v ˜ ψ uv , ¯˜ χ uv = ¯˜ ψ uv ( γ ) v ( γ ) u , u, v ∈ { , } . (73)33hen eq. (72) becomes S free [ ˜ ψ, ˜¯ ψ ] ≈ − (cid:88) u,v (cid:88) k ˜¯ ψ uv ( k ) (cid:0) γ ( k − i µ ) + γ k (cid:1) ˜ ψ uv ( k ) . (74)This shows that the lattice action (35) corresponds in the continuum limit to four masslessnon-interacting fermion flavors. A.2 Naive fermions and the GN model
Discretizing the GN model (2) with N f = 4 flavors in a straightforward way using the fields χ and ¯ χ via S σ [ χ, ¯ χ, σ ] = S free [ χ, ¯ χ ] + i (cid:88) x ¯ χ ( x ) σ ( x ) χ ( x ) + N f g (cid:88) x σ ( x ) , N f = 4 (75)actually results in a theory different from the GN model. To show this, we insert again theFourier representations of the fermionic fields (67) as well as of the real scalar field σ , σ ( x ) = 1 √ N t N s (cid:88) k e i k · x ˜ σ ( k ) , (76)where ˜ σ ( − k ) = ˜ σ ∗ ( k ) and the discrete momenta k are chosen such that σ is periodic in x and x direction. The action (75) becomes S σ [ ˜ χ, ˜¯ χ, ˜ σ ] = S free [ ˜ χ, ˜¯ χ ] + i √ N t N s (cid:88) k (cid:88) k (cid:48) ˜¯ χ ( k )˜ σ ( k − k (cid:48) ) ˜ χ ( k (cid:48) ) + N f g (cid:88) k (cid:12)(cid:12) ˜ σ ( k ) (cid:12)(cid:12) , N f = 4 . (77)In the limit a → σ and, thus, no corresponding suppression of σ modes. Consequently, for a → √ N t N s (cid:88) k , k (cid:48) (cid:88) u,v,u (cid:48) ,v (cid:48) ¯˜ χ uv ( k )˜ σ uv,u (cid:48) v (cid:48) ( k − k (cid:48) ) ˜ χ u (cid:48) v (cid:48) ( k (cid:48) ) , (78)with symmetric kernel in momentum space˜ σ uv,u (cid:48) v (cid:48) ( k ) = ˜ σ u (cid:48) v (cid:48) ,uv ( k ) = ˜ σ ( k uv − k u (cid:48) v (cid:48) + k ) . (79)In terms of the usual field coordinates ˜ ψ uv , related to ˜ χ uv via eq. (73), the interaction term isi √ N t N s (cid:88) k , k (cid:48) (cid:88) u,v,u (cid:48) v (cid:48) ¯˜ ψ uv ( k )( γ ) v ( γ ) u σ uv,u (cid:48) v (cid:48) ( k − k (cid:48) )( γ ) u (cid:48) ( γ ) v (cid:48) ˜ ψ u (cid:48) v (cid:48) ( k (cid:48) ) . (80)Now it is obvious that the action (75) is not a discretization of the GN model with N f = 4fermion flavors. While the four terms with uv = u (cid:48) v (cid:48) in eq. (80) represent the correct GNinteraction for the four fermion flavors,i √ N t N s (cid:88) u,v (cid:88) k , k (cid:48) ˜¯ ψ uv ( k )˜ σ ( k − k (cid:48) ) ˜ ψ u (cid:48) v (cid:48) ( k (cid:48) ) , (81)34 .0 0.1 0.2 0.3 0.4 0.5 0.6 T / SLAC, N s = 63naive, N s = 64incorrect naive, N s = 64 Figure 16: Σ as a function of T at µ = 0 for SLAC fermions, naive fermions and the straight-forward, but incorrect naive discretization (75) ( N f = 8, a ≈ . /σ ).there are twelve additional terms not present in the GN model, where the field σ couples twodifferent fermion flavors,i √ N t N s (cid:88) k , k (cid:48) ˜¯ ψ ( k ) γ ˜ σ ( π/a + k − k (cid:48) , k − k (cid:48) ) ˜ ψ ( k (cid:48) ) + eleven more terms , (82)as was already pointed out in Ref. [51]. Including these twelve terms in a numerical simulation,i.e. using the action (75), corresponds to studying a different theory and leads to results signif-icantly different from those obtained with a correct discretization of the GN model (examplesare shown at the end of this section).Now we derive a proper lattice discretization of the GN model. To this end, we note that onlythe soft fermion modes contribute in the limit a → σ uv,uv , while the twelve spurious interaction terms are proportional to ˜ σ uv,u (cid:48) v (cid:48) with uv (cid:54) = u (cid:48) v (cid:48) . Thus, one can eliminate the spurious terms by replacing ˜ σ ( k ) in the interactionterm in eq. (77) by ˜ W ( k )˜ σ ( k ). Here ˜ W is a weight-function with • ˜ W ( k ) → k ≈ k = (0 ,
0) (i.e. in region R ), • ˜ W ( k ) → k ≈ k uv with ( u, v ) (cid:54) = (0 ,
0) (i.e. in the other regions R uv ).A simple choice, which we use for our numerical simulations, is˜ W ( k ) = 14 (cid:89) µ =0 (cid:16) ak µ ) (cid:17) . (83)Expressing this modified action in terms of χ ( x ), ¯ χ ( x ) and σ ( x ) is straightforward, S GN [ χ, ¯ χ, σ ] = S free [ χ, ¯ χ ] + i √ N t N s (cid:88) x , x (cid:48) ¯ χ ( x ) W ( x − x (cid:48) ) σ ( x (cid:48) ) χ ( x ) + N f g (cid:88) x σ ( x ) , N f = 4 , (84)35here W is the Fourier transform of ˜ W and given in eq. (36). All calculations and argumentspresented in this section apply to N f = 8 , , , . . . flavors in a trivial way. The generalizationof eq. (84) is eq. (34).In Figure 16 we show numerical evidence that using the straightforward naive discretization ofthe GN model (75) leads to incorrect results, i.e. results not corresponding to the GN model.We plot Σ as a function of the temperature T for chemical potential µ = 0. The blue andorange curves correspond to the SLAC discretization (see section 3.2.2) and the correct naivediscretization (34) (or equivalently (84)). These curves are rather similar and get closer, whendecreasing the lattice spacing, indicating that they have the same continuum limit. The greencurve, on the other hand, corresponding to the straightforward naive discretization (75) is quitedifferent and does not approach the blue and orange curves, when decreasing the lattice spacing.We obtained similar results also for non-vanishing chemical potenial.36 cknowledgements We acknowledge useful discussions with Michael Buballa, Holger Gies, Felix Karbstein, AdrianK¨onigstein, Maria Paola Lombardo, Dirk Rischke, Alessandro Sciarra, Lorenz von Smekal, StefenTheisen, Michael Thies, Marc Winstel and Ulli Wolff on various aspects of fermion theories andspacetime symmetries. Special thanks go to Philippe de Forcrand for his constructive remarksand to Martin Ammon who shared his knowledge on SSB and Goldstone bosons and for hissteady encouragement in the past 15 months.J.J.L. and A.W. have been supported by the Deutsche Forschungsgemeinschaft (DFG) underGrant No. 406116891 within the Research Training Group RTG 2522/1. L.P. and M.W. ac-knowledge support by the Deutsche Forschungsgemeinschaft (DFG, German Research Founda-tion) through the CRC-TR 211 “Strong-interaction matter under extreme conditions” – projectnumber 315477589 – TRR 211. M.W. acknowledges support by the Heisenberg Programme ofthe Deutsche Forschungsgemeinschaft (DFG, German Research Foundation) – project number399217702.This work was supported in part by the Helmholtz International Center for FAIR within theframework of the LOEWE program launched by the State of Hesse. Calculations on theGOETHE-HLR and on the on the FUCHS-CSC high-performance computers of the Frank-furt University were conducted for this research. We would like to thank HPC-Hessen, fundedby the State Ministry of Higher Education, Research and the Arts, for programming advice.We used the python programming language, most notably pandas [72] and numpy [73] for dataanalysis and matplotlib [74] for plotting but also further algorithms from the SciPy ecosystem[71]. 37 eferences [1] D. J. Gross and A. Neveu,
Dynamical Symmetry Breaking in Asymptotically Free FieldTheories , Phys. Rev.
D10 (1974) 3235.[2] I. Sachs and A. Wipf,
Generalized Thirring models , Annals Phys. (1996) 380, arXiv:hep-th/9508142.[3] A. Chodos and H. Minakata,
The Gross-Neveu model as an effective theory for polyacetylene ,Phys. Lett.
A191 (1994) 39.[4] H. Takayama, Y. R. Lin-Liu and K. Maki,
Continuum model for solitons in polyacetylene ,Phys. Rev.
B21 (1980) 2388.[5] H. Caldas,
An Effective Field Theory Model for One-Dimensional CH Chains: Effectsat Finite Chemical Potential, Temperature and External Zeeman Magnetic Field , J. Stat.Mech. (2011), no. 10 P10005, arXiv:1106.0948.[6] Y. Kuno,
Phase structure of the interacting Su-Schrieffer-Heeger model and the relationshipwith the Gross-Neveu model on lattice , Phys. Rev.
B99 (2019) 064105, arXiv:1811.01487.[7] A. Bermudez, E. Tirrito, M. Rizzi, M. Lewenstein and S. Hands,
Gross-Neveu-Wilsonmodel and correlated symmetry-protected topological phases , Annals Phys. (2018) 149,arXiv:1807.03202.[8] M. Thies and K. Urlichs,
Revised phase diagram of the Gross-Neveu model , Phys. Rev.
D67 (2003) 125015, arXiv:hep-th/0302092.[9] O. Schnetz, M. Thies and K. Urlichs,
Phase diagram of the Gross-Neveu model: Exact re-sults and condensed matter precursors , Annals Phys. (2004) 425, arXiv:hep-th/0402014.[10] M. Buballa and S. Carignano,
Inhomogeneous chiral condensates , Prog. Part. Nucl. Phys. (2015) 39, arXiv:1406.1367.[11] G. Basar, G. V. Dunne and M. Thies, Inhomogeneous Condensates in the Thermodynamicsof the Chiral NJL(2) model , Phys. Rev.
D79 (2009) 105012, arXiv:0903.1868.[12] D. Nickel,
Inhomogeneous phases in the Nambu-Jona-Lasino and quark-meson model , Phys.Rev.
D80 (2009) 074025, arXiv:0906.5295.[13] S. Carignano, D. Nickel and M. Buballa,
Influence of vector interaction and Polyakovloop dynamics on inhomogeneous chiral symmetry breaking phases , Phys. Rev.
D82 (2010)054009, arXiv:1007.1397.[14] S. Carignano and M. Buballa,
Two-dimensional chiral crystals in the NJL model , Phys.Rev.
D86 (2012) 074018, arXiv:1203.5343.[15] A. Heinz, F. Giacosa and D. H. Rischke,
Chiral density wave in nuclear matter , Nucl. Phys.
A933 (2015) 34, arXiv:1312.3244.[16] J. Braun, F. Karbstein, S. Rechenberger and D. Roscher,
Crystalline ground states inPolyakov-loop extended Nambu-Jona-Lasinio models , Phys. Rev.
D93 (2016), no. 1 014032,arXiv:1510.04012. 3817] M. Buballa and S. Carignano,
Inhomogeneous chiral phases away from the chiral limit ,Phys. Lett.
B791 (2019) 361, arXiv:1809.10066.[18] S. Carignano and M. Buballa,
Inhomogeneous chiral condensates in three-flavor quark mat-ter , Phys. Rev.
D101 (2020), no. 1 014026, arXiv:1910.03604.[19] P. de Forcrand and U. Wenger,
New baryon matter in the lattice Gross-Neveu model , PoS
LAT2006 (2006) 152, arXiv:hep-lat/0610117.[20] M. Wagner,
Fermions in the pseudoparticle approach , Phys. Rev.
D76 (2007) 076002,arXiv:0704.3023.[21] A. Heinz, F. Giacosa, M. Wagner and D. H. Rischke,
Inhomogeneous condensation in effec-tive models for QCD using the finite-mode approach , Phys. Rev.
D93 (2016), no. 1 014007,arXiv:1508.06057.[22] M. Winstel, J. Stoll and M. Wagner,
Lattice investigation of an inhomogeneous phase ofthe 2+1-dimensional Gross-Neveu model in the limit of infinitely many flavors (2019),arXiv:1909.00064.[23] R. Narayanan,
Phase diagram of the large N Gross-Neveu model in a finite periodic box (2020), arXiv:2001.09200.[24] H. B. Nielsen and M. Ninomiya,
No Go Theorem for Regularizing Chiral Fermions , Phys.Lett. (1981) 219.[25] B. H. Wellegehausen, D. Schmidt and A. Wipf,
Critical flavor number of the Thirring modelin three dimensions , Phys. Rev.
D96 (2017), no. 9 094504, arXiv:1708.01160.[26] S. Aoki and K. Higashijima,
The Recovery of the Chiral Symmetry in Lattice Gross-NeveuModel , Prog. Theor. Phys. (1986) 521.[27] U. Wolff, The phase diagram of the inifinite N Gross-Neveu model at finite temperature andchemical potential , Phys. Lett. (1985) 303.[28] P. Fulde and R. A. Ferrell,
Superconductivity in a Strong Spin-Exchange Field , Phys. Rev. (1964) A550.[29] A. I. Larkin and Y. N. Ovchinnikov,
Nonuniform state of superconductors , Zh. Eksp. Teor.Fiz. (1964) 1136.[30] J. J. Kinnunen, J. E. Baarsma, J.-P. Martikainen and P. T¨orm¨a, The Fulde-Ferrell-Larkin-Ovchinnikov state for ultracold fermions in lattice and harmonic potentials: a review , Rept.Prog. Phys. (2018), no. 4 046401, arXiv:1706.07076.[31] A. Barducci, R. Casalbuoni, M. Modugno, G. Pettini and R. Gatto, Thermodynamics ofthe massive Gross-Neveu model , Phys. Rev.
D51 (1995) 3042, arXiv:hep-th/9406117.[32] R. F. Dashen, B. Hasslacher and A. Neveu,
Semiclassical Bound States in an AsymptoticallyFree Theory , Phys. Rev.
D12 (1975) 2443.[33] R. Pausch, M. Thies and V. L. Dolman,
Solving the Gross-Neveu model with relativisticmany body methods , Z. Phys.
A338 (1991) 441.3934] J. Feinberg,
All about the static fermion bags in the Gross-Neveu model , Annals Phys. (2004) 166, arXiv:hep-th/0305240.[35] G. Basar and G. V. Dunne,
Self-consistent crystalline condensate in chiral Gross-Neveuand Bogoliubov-de Gennes systems , Phys. Rev. Lett. (2008) 200404, arXiv:0803.1501.[36] N. D. Mermin and H. Wagner,
Absence of ferromagnetism or antiferromagnetism in one-dimensional or two-dimensional isotropic Heisenberg models , Phys. Rev. Lett. (1966)1133.[37] S. R. Coleman, There are no Goldstone bosons in two-dimensions , Commun. Math. Phys. (1973) 259.[38] J. Goldstone, Field Theories with Superconductor Solutions , Nuovo Cim. (1961) 154.[39] Y. Nambu, Axial vector current conservation in weak interactions , Phys. Rev. Lett. (1960)380.[40] S.-S. Shei, Semiclassical Bound States in a Model with Chiral Symmetry , Phys. Rev.
D14 (1976) 535.[41] H. B. Nielsen and S. Chadha,
On How to Count Goldstone Bosons , Nucl. Phys.
B105 (1976) 445.[42] H. Watanabe and T. Brauner,
On the number of Nambu-Goldstone bosons and its relationto charge densities , Phys. Rev.
D84 (2011) 125013, arXiv:1109.6327.[43] H. Watanabe and H. Murayama,
Unified Description of Nambu-Goldstone Bosons withoutLorentz Invariance , Phys. Rev. Lett. (2012) 251602, arXiv:1203.0609.[44] H. Watanabe,
Counting Rules of Nambu-Goldstone Modes , Ann. Rev. Condensed MatterPhys. (2020) 169, arXiv:1904.00569.[45] M. Ammon, M. Baggioli and A. Jimnez-Alba, A Unified Description of Translational Sym-metry Breaking in Holography , JHEP (2019) 124, arXiv:1904.05785.[46] J. A. Barker and D. Henderson, What is ’liquid’ ? Understanding the states of matter , Rev.Mod. Phys. (1976) 587.[47] D. Chandler, Introduction to modern statistical mechanics , Oxford University Press, NewYork (1987).[48] H. J. Rothe,
Lattice Gauge Theories , World Scientific, 4th edn. (2012).[49] J. Smit,
Introduction to quantum fields on a lattice: A robust mate , Cambridge Lect. NotesPhys. (2002) 1.[50] A. Wipf, Statistical Approach to Quantum Field Theory , Lecture Notes in Physics (2013).[51] Y. Cohen, S. Elitzur and E. Rabinovici,
A Monte Carlo study of the Gross-Neveu model ,Nucl. Phys.
B220 (1983) 102. 4052] L. Pannullo, J. Lenz, M. Wagner, B. Wellegehausen and A. Wipf,
Inhomogeneous phasesin the 1+1 dimensional Gross-Neveu model at finite number of fermion flavors , Acta Phys.Polon. Supp. (2020) 127, arXiv:1902.11066.[53] L. Pannullo, J. Lenz, M. Wagner, B. Wellegehausen and A. Wipf, Lattice investigation ofthe phase diagram of the 1+1 dimensional Gross-Neveu model at finite number of fermionflavors , in (2019) arXiv:1909.11513.[54] S. D. Drell, M. Weinstein and S. Yankielowicz,
Variational Approach to Strong CouplingField Theory. 1. Phi**4 Theory , Phys. Rev.
D14 (1976) 487.[55] G. Bergner, T. Kaestner, S. Uhlmann and A. Wipf,
Low-dimensional Supersymmetric Lat-tice Models , Annals Phys. (2008) 946, arXiv:0705.2212.[56] L. H. Karsten and J. Smit,
The Vacuum Polarization With SLAC Lattice Fermions , Phys.Lett. (1979) 100.[57] C. Wozar and A. Wipf,
Supersymmetry Breaking in Low Dimensional Models , Annals Phys. (2012) 774, arXiv:1107.3324.[58] R. Flore, D. Korner, A. Wipf and C. Wozar,
Supersymmetric Nonlinear O(3) Sigma Modelon the Lattice , JHEP (2012) 159, arXiv:1207.6947.[59] R. V. Gavai and S. Sharma, Divergences in the quark number susceptibility: The origin anda cure , Phys. Lett.
B749 (2015) 8, arXiv:1406.0474.[60] M. A. Clark and A. D. Kennedy,
The RHMC algorithm for two flavors of dynamical stag-gered fermions , Nucl. Phys. Proc. Suppl. (2004) 850, arXiv:hep-lat/0309084.[61] B. H. Wellegehausen,
Phase diagrams of exceptional and supersymmetric lattice gauge the-ories , Ph.D. thesis, University Jena (2012).[62] J. J. Lenz, B. H. Wellegehausen and A. Wipf,
Absence of chiral symmetry break-ing in Thirring models in 1+2 dimensions , Phys. Rev.
D100 (2019), no. 5 054501,arXiv:1905.00137.[63] D. August, M. Steinhauser, B. H. Wellegehausen and A. Wipf,
Mass spectrum of -dimensional N = (2 , super Yang-Mills theory on the lattice , JHEP (2019) 099,arXiv:1802.07797.[64] F. Karsch, J. B. Kogut and H. W. Wyld, The Gross-Neveu Model at Finite Temperatureand Density , Nucl. Phys.
B280 (1987) 289.[65] M. Thies (2019), unpublished notes.[66] J. Kertesz,
Existence of weak singularities when going around the liquid-gas critical point ,Physica A (1989) 58.[67] S. Wenzel, E. Bittner, W. Janke, A. M. J. Schakel and A. Schiller, Kertesz line in thethree-dimensional compact U(1) lattice Higgs model , Phys. Rev. Lett. (2005) 051601,arXiv:cond-mat/0503599. 4168] V. L. Berezinskii, Destruction of Long-range Order in One-dimensional and Two-dimensional Systems having a Continuous Symmetry Group I. Classical Systems , Sov. Phys.JETP (1971) 493.[69] J. M. Kosterlitz and D. J. Thouless, Ordering, metastability and phase transitions in two-dimensional systems , J. Phys. C6 (1973) 1181.[70] E. Witten, Chiral Symmetry, the 1/n Expansion, and the SU(N) Thirring Model , Nucl.Phys.
B145 (1978) 110.[71] SciPy Contributors,
SciPy 1.0–Fundamental Algorithms for Scientific Computing inPython , Nature Meth. (2020), arXiv:1907.10121.[72] W. McKinney,
Data Structures for Statistical Computing in Python , in
Proceedings of the9th Python in Science Conference (2010) pages 51 – 56.[73] S. van der Walt, S. C. Colbert and G. Varoquaux,
The NumPy Array: A Structure forEfficient Numerical Computation , Comput. Sci. Eng. (2011), no. 2 22, arXiv:1102.1523.[74] J. D. Hunter, Matplotlib: A 2D Graphics Environment , Comput. Sci. Eng.9