Baryons in the Gross-Neveu model in 1+1 dimensions at finite number of flavors

In a recent work we studied the phase structure of the Gross-Neveu (GN) model in $1+1$ dimensions at finite number of fermion flavors $N_\mathrm{f} = 2, 8, 16$, finite temperature and finite chemical potential using lattice field theory. Most importantly, we found an inhomogeneous phase at low temperature and large chemical potential, quite similar to the analytically solvable $N_\mathrm{f} \rightarrow \infty$ limit. In the present work we continue our lattice field theory investigation of the finite-$N_\mathrm{f}$ GN model by studying the formation of baryons, their spatial distribution and their relation to the chiral condensate. As a preparatory step we also discuss a linear coupling of lattice fermions to the chemical potential.


I. INTRODUCTION
In recent years experiments provided many interesting insights concerning strongly interacting matter at high density (see e.g. Ref. [2] for a comprehensive review). On the theoretical side our present understanding of the QCD phase diagram at non-zero chemical potential µ is to a large extent based on conjectures relying on physical intuition, on model calculations and on effective low energy descriptions [3,4], while reliable ab-initio results are still missing, mostly due to the infamous sign problem in lattice-QCD. Even though there are a number of interesting approaches, which led to considerable progress to mitigate or solve the sign problem, such as using complex Langevin algorithms [5][6][7][8] or thimble methods [9][10][11], finding more suitable variables [12][13][14] or refining the density of states approach [15][16][17][18][19], a better understanding of lattice-QCD at finite baryon density is certainly an urgent problem. Urgent, for example, since our colleagues from gravitational wave astronomy and astrophysics are in need of more reliable equations of state of strongly interacting matter at baryon density n B up to several times the nuclear density n 0 ≈ 0.17 fm −3 .
It has been conjectured that in QCD at low temperature and large baryon density there is an inhomogeneous crystalline phase. This conjecture is based on mean field calculations in various effective four Fermi theories indicating the existence of such an inhomogeneous phase [20][21][22][23][24]. The underlying mean-field (or Hartree-Fock like) approximation becomes exact in the limit of an infinite * julian.johannes.lenz@uni-jena.de † bjoern.wellegehausen@uni-jena.de ‡ wipf@tpi.uni-jena.de § pannullo@itp.uni-frankfurt.de ¶ mwagner@itp.uni-frankfurt.de number of fermion flavors N f , since in this limit quantum fluctuations are negligible. Mean field approximations are also common in condensed matter physics. For example, for the GN model considered in the present work, the mean field phase diagram with homogeneous and inhomogeneous phases has been known in the condensed matter community [25,26] long before it has been rediscovered in particle physics [20,27].
More recently, interesting models implementing the breaking of translational invariance -for example by charge density waves, dynamical defects or by magnetic fields -have been proposed and studied within the holographic framework [28,29].
At present it is largely unknown, whether crystalline phases exist in effective four Fermi theories at finite number of fermion flavors, or whether quantum fluctuations lead to a qualitatively different phase structure. In a recent work [1] we performed lattice field theory simulations of the GN model in 1 + 1 dimensions with N f = 2, 8, 16 and found clear evidence for the existence of an inhomgeneous phase, qualitatively similar to that in the limit N f → ∞. In the present work we continue our investigation of the GN model in 1 + 1 dimensions at a finite number of fermion flavors and focus on baryonic excitations at low temperature and large chemical potential. We investigate their spatial distribution as well as their relation to the chiral condensate

II. CHEMICAL POTENTIAL FOR LATTICE FERMIONS
The continuum Lagrangian density of the (Euclidean) GN model with vanishing bare mass is given by where µ denotes the chemical potential for the conserved baryon number. To be able to perform the fermion integration one follows Hubbard and Stratonovich by introducing a fluctuating auxiliary scalar field σ to linearize the operatorψψ in the interaction term The four Fermi term in Eq. (1) is recovered after eliminating σ by its equation of motion or equivalently by integrating over σ in the functional integral. Translation invariance of dσ x in the (well-defined) functional integral dσ x for the lattice model implies the Ward identity Keeping as many global symmetries of the continuum model as possible in a discretization can be crucial to obtain a lattice model with the correct continuum limit. Thus we shall discretize the operator D using the chiral and doubler-free SLAC derivative [30,31]. While nonlocal SLAC fermions must not be used to discretize a field theory with local gauge symmetries [32][33][34], they have been used successfully in various scalar-field theories and fermionic theories with global symmetries only [35][36][37][38][39][40]. In addition to using SLAC fermions to simulate the GN model at finite N f , we have cross-checked our results with a discretization based on naive fermions. This fermion species is chiral as well but describes 2 d doublers in d dimensions. More details can be found in Ref. [1].
Besides our preceding paper [1], we are not aware of any work, in which SLAC fermions have been used to study fermion systems at finite density. Thus, we begin with comparing the thermodynamics of a gas of free massive fermions in a spatial box of size L in the continuum and on the lattice with different fermion discretiztations. A straightforward calculation yields the grand partition function at inverse temperature β and chemical potential µ in the continuum with single-particle energies E 2 k = k 2 + m 2 depending on the spatial wave number k = 2πn/L, n ∈ Z and mass m. The corresponding baryon density is Note that the sum over all Matsubara frequencies has already been performed, i.e. the continuum limit in (imaginary) time direction is already implied, while truncating the sum over k is conceptually similar to a finite lattice spacing in spatial direction. Since the SLAC derivative discretizes the continuum dispersion relation up to the maximal momentum given by the inverse lattice spacing, the finite (truncated) sum is also the result for free, massive SLAC fermions discretized in spatial direction only.
This is not what is implemented in lattice Monte Carlo simulations where also the (imaginary) time is discretized. Asymptotically, the error in truncating the Matsubara sum after N t terms (with Matsubara frequencies symmetrically about the origin) is ∼ β (E k ± µ) /π 2 N t . Letting the lattice constant in time direction β/N t → 0 at fixed k, we can neglect this error. However, E k will eventually become large in the sum over k and higher order corrections will contribute, if the temporal "cutoff" N t /β is not sent to infinity before taking the limit L/N s → 0. In the particular case of a uniform continuum limit N t = N s → ∞, we pick up the following correction terms in 1 + 1 dimensions: As argued above, the expressions on the left correspond to the continuum limit of free massive SLAC fermions, 1 for which the chemical potential enters the Lagrangian linearly via iµψγ 0 ψ as it does in the continuum, see e.g. Eq.
(1). We conclude that introducing the chemical potential linearly as in the continuum theory yields the correct partition function up to a (µ-dependent) constant and can thus be used in Monte Carlo simulations. In fact, since SLAC fermions couple fermion fields at well-separated lattice sites it would be difficult to introduce an exponentially coupled µ for all hopping terms in the Lagrangian.
In appendix A we show Eq. (7), i.e. that the baryon density computed from a lattice action with linearly coupled µ has to be corrected by the constant +µ/2π. This is actually not a defect of SLAC fermions but just expresses the fact that conventionally one first performs the continuum limit in time direction and afterwards in the spatial directions to arrive at the well-known expression for the thermodynamic potentials at finite temperature and density. A comparison of free massive baryon densities for various commonly used lattice discretizations is depicted in Fig. 1. The (properly corrected) SLAC result is almost indistinguishable from the continuum result. In passing, we note the following: • A similar analysis for naive fermions reveals that introducing a chemical potential as additive linear term requires the same correction as for SLAC fermions (when β/N t and L/N s approach zero simultaneously) and thus could also be used. the additive linear coupling according to Eq. (7) has been used, for Wilson fermions and for staggered fermions only exponential coupling and for naive fermions we distinguish between conventional exponential coupling (exp) and additive linear coupling (lin). Note that staggered fermions and naive fermions with exponential coupling lead to identical results.
• For d > 2 space-time dimensions the known correction term diverges in the continuum limit. In 4 spacetime dimensions it diverges quadratically [41]. A practical renormalization scheme on the lattice could then be to determine the constant c via (c is finite at finite lattice spacing a) and subtracting cµ 2 from the partition function. It has been demonstrated in Ref. [42] that this method of divergence removal works in (quenched) QCD. This is different from earlier attempts to eliminate the divergencies by suitably modifying the lattice action [43][44][45].
• The standard exponential coupling of the chemical potential for naive fermions has significantly larger discretization errors than the linear coupling (if corrected properly). In particular, a linearly coupled chemical potential yields quite accurately the position of the first step, i.e. the fermion mass.
• In appendix A we explicitly calculate the correction term for N t ≈ N s 1 for non-interacting fermions. We have observed that this correction is insensitive to the interaction and hence it is sufficient to subtract the same term in the GN model [1]. It seems that this statement holds true in higher dimensions as well, where the correction terms are UV-divergent. This has been observed in numerical simulations [42], but a general proof of this interesting observation in QCD and interacting GN models appears to be still missing.
We emphasize once more, that the µ-dependent correction terms are not lattice artifacts -they may also appear in continuum theories, depending on how divergent integrals are treated (cf. the detailed discussion below Eq. (7)). We have crosschecked results obtained within this project using naive fermions with a conventional exponentially coupled chemical potential (see Ref. [1]) as this is the established method of introducing a chemical potential in lattice field theory. Since the exponentially coupled µ couples to the exactly conserved charge on the lattice, no correction terms are needed [43].

III. BARYONIC MATTER
In a recent paper [1] we studied the phase structure of the GN model in 1 + 1 dimensions for N f = 8 using lattice field theory and SLAC fermions with exact chiral symmetry. The resulting phase diagram is shown in Fig. 2 in units of As expected, we identified a homogeneously broken phase with non-zero constant chiral condensate at small chemical potential and low temperature and a symmetric phase at high temperature. 2 Most interestingly, however, we also found a phase, where the spatial correlator is an oscillating function (see e.g. Fig. 4b, bottom). As a simple observable to distinguish the three phases we used where 0 inside the homogeneously broken phase ≈ 0 inside the symmetric phase < 0 inside the inhomogeneous phase (13) (see Fig. 2). The GN model can be solved analytically in the semiclassical approximation or, equivalently, in the limit N f → ∞ (see e.g. Refs. [20,27,46,47]) and it is known that extrema of the effective action which one obtains by using L σ from Eq. (2) and integrating over the fermions in the partition function, are not only given by σ = const. For example in Ref. [20] it was shown that at large chemical potential and small temperature a spatially oscillating function σ(x) minimizes the free energy. For each cycle of the oscillation N f fermions or antifermions, which can be interpreted as baryons, are located in the region of minimal σ 2 , i.e. where the sign of σ changes. This implies breaking of translational symmetry and a crystal of baryons (as shown in Fig. 3).
In the present work we investigate, whether traces of such a baryonic crystal are also present in the GN model with a finite number of fermion flavors. We use the same lattice setup as in our preceding paper [1]. For all plots shown in the following we performed computations with N f = 8 flavors of SLAC fermions, lattice spacing a ≈ 0.410/σ 0 and N s = 63 lattice sites in spatial direction, corresponding to a periodic spatial direction of extent L = N s a ≈ 25.8/σ 0 . From the extensive set of simulations we carried out in Ref. [1] for different a and L, we expect that both lattice discretization errors and finite volume corrections are small. Note that in Ref. [1] we also performed computations with N f = 8 flavors of naive fermions, to check and to confirm our numerical results.  54) and (80)).

A. Correlation of the baryon density and the condensate
We start by investigating the location of the fermions relative to the spatially oscillating condensate inside the inhomogeneous phase. It is important to note that the effective action (14) is invariant under spatial translations. Therefore, field configurations, which are spatially shifted relative to each other, i.e. σ(t, x) and σ(t, x + δx), contribute with the same weight e −S eff to the partition function and, thus, should be generated with the same probability by the HMC algorithm. Consequently, simple observables like σ(x) or n B (x) , where are not suited to detect an inhomogeneous condensate or baryon density in a lattice simulation, because destructive interference should lead to σ(x) = 0 and n B (x) = const, even in cases, where all field configurations exhibit spatial oscillations with the same wavelength. An observable, which does not suffer from destructive interference and is able to exhibit information about possibly present inhomogeneous structures, is the spatial correlation function of σ(x), as defined in Eq. (11) (for a more detailed discussion see section 4.3 of Ref. [1]). Similarly, the spatial correlation function of the baryon density and the squared condensate, can provide insights on the location of the fermions relative to the extrema of the condensate inside an inhomogeneous phase. For N f → ∞ the maxima of σ 2 (x) coincide with the minima of n B (x) and vice versa, as one can read off from Fig. 3. Thus, for N f → ∞ the correlator C n B σ 2 (x) has minima at nλ/2 and maxima at    16) and (11). The vertical gray lines indicate the roots of C(x).
(n + 1/2)λ/2, where λ is the wavelength of both σ 2 (x) and n B (x) (see Fig. 4a, where C n B σ 2 (x) and C(x) are shown for (µ/σ 0 , T /σ 0 ) = (0.700, 0.038)). Our corresponding lattice results for N f = 8 at the same chemical potential and temperature exhibit an almost identical behavior (see Fig. 4b). We interpret this as clear signal that baryons are centered at the roots of the condensate σ, where σ is a periodically oscillating function (we have investigated the latter in detail in our preceding work [1]). Thus separations between neighboring baryons should all be similar, which is reminiscent to the baryonic crystal found at N f → ∞.

B. Baryon number and its relation to the condensate
In this section we study the baryon number with the baryon density n B as defined in Eq. (15) and investigate its relation to the average number of cycles of the oscillating condensate σ. Inside the finite periodic lattice with extent L we have defined and computed By k max we denote the dominant momentum of c(x) (see Eq. (11)), i.e. that k, which maximizes the absolute value of the Fourier transformc(k). The central result of this subsection is that B and ν max are almost identical on each field configuration, i.e. even when omitting the average . . . over all generated field configurations in the definitions (17) and (18). In particular at small T there is almost perfect agreement. This is illustrated in Fig. 5 In the following we elaborate on this further by presenting and discussing results obtained at two different values for the temperature, T /σ 0 ≈ 0.076 and T /σ 0 ≈ 0.038.
with σ 0 and σ as defined in Eq. (9) and Eq. (10)   transition between the homogeneously broken and the inhomogeneous phase. While B ≈ 0 for small µ, B suddenly starts to increase at µ/σ 0 ≈ 0.51. At roughly the same µ value Σ 2 rapidly drops from around 1 to 0. Note that B is quite similar to n B | N f ,L→∞ L (the green solid line in Fig. 6), where n B | N f ,L→∞ is the analytical infinite volume result for the N f → ∞ baryon density according to Ref. [20]. However, the phase transition for N f = 8 takes place at smaller µ/σ 0 ≈ 0.51 compared to µ/σ 0 ≈ 2/π ≈ 0.64 for N f → ∞, where (at T = 0) a baryon corresponds to a kink-antikink field configuration σ with energy 2σ 0 /π per fermion (see e.g. Refs. [46,49]). This can also be seen in the phase diagram shown in Fig. 2 and implies that the baryon mass (per fermion flavor) is somewhat smaller compared to N f → ∞. Note that at finite N f a smaller homogeneously broken phase is expected, because of fluctuations in σ, which increase disorder (see the detailed discussion in our preceding work [1]). The careful reader might already note the remnants of the stair-like low temperature behavior discussed later on in section III B 2. Moreover, at large µ the N f = 8 result and the N f → ∞ result for the baryon number B are very similar, where the latter approaches the corresponding result for free fermions, as noted in Ref. [27].
In Fig. 7 we show B and ν max as functions of the chemical potential µ. The two curves are quite similar, which is a strong indication that the finite-N f theory is qualitatively well-described by the semi-classical N f → ∞ picture. Note that ν max is slightly below B for µ/σ 0 > ∼ 0.51. The reason can be seen in Fig. 8  stable quantity, while ν max exhibits sizable fluctuations on around 25% of the generated field configurations, mostly fluctuations towards small values, significantly below the median of ν max . Such fluctuations are more common for larger values of µ. This is expected from the Fourier transformed spatial correlation function, which we investigated in Ref. [1] in detail, and also reflected by the corresponding histogram for ν max shown in the right part of Fig. 8.
To summarize, Fig. 7 indicates that the baryon number B is quite similar to the number of cycles of the spatial oscillation of the condensate σ, as for N f → ∞. This observation supports our conclusions above that, also at finite N f , baryons are the relevant excitations of the GN model and that their number is closely related to the shape of the condensate σ.  Autocorrelatios at T /σ 0 ≈ 0.038 turned out to be rather large in our simulations, in particular near the boundary of the homogeneously broken phase and the inhomogeneous phase, in the region 0.40 < ∼ µ/σ 0 < ∼ 0.65. This is illustrated in Fig. 9, where we compare Monte Carlo histories of B at µ/σ 0 = 0.65 for a cold start (each field variable σ(t, x) = 1) and a hot start (each field variable drawn randomly from a Gaussian distribution with mean 0). After around 1500 Monte Carlo sweeps the two Monte Carlo histories eventually converge and the simulations seem to have thermalized. Nevertheless the autocorrelation time is quite large, of the order of the average number of Monte Carlo sweeps needed to create or annihilate a baryon, i.e. > ∼ 500. This is sizable compared to the typical number of Monte Carlo sweeps, between 2000 and 10000, we are able to carry out for each simulation with our available HPC resources. Thus, for µ/σ 0 ≈ 0.65 the errors we show for our results might be somewhat underestimated.
For larger µ/σ 0 , i.e. farther away from the phase boundary, autocorrelation times become smaller. For example, in Fig. 5 we show the Monte Carlo histories of B and of ν max at µ/σ 0 = 1.10 after thermalization. For almost all field configurations B ≈ ν max and their values are either close to 8 or close to 9. Even though we only show 150 Monte Carlo sweeps, there are many transtions between B ≈ ν max = 8 and B ≈ ν max = 9, indicating that the HMC algorithm is able to frequently increase or decrease the number of cycles of the spatial oscillation of the condensate σ. In Fig. 10 we show B and ν max as functions of µ in the range 0.80 ≤ µ/σ 0 ≤ 1.20. Both quantities exhibit a very clear stair-like behavior with the steps corresponding to integers. The steps of ν max (which are more pronounced at low temperature due to less thermal fluctuations) are explained by the semiclassical commensurability constraint, namely that the wave length of the periodic condensate must divide the box length L. In the limit N f → ∞ the baryon number is equal to the cycles of the condensate. Thus, a second baryon appears at chemical potential µ 2 > µ c ≈ (2/π)σ 0 , a third baryon at µ 3 > µ 2 , etc. Hence, with increasing µ the mean separations of baryons decreases which leads to more interaction. That the semi-classical picture explains the simulations so well is a further indication that the GN model at N f = 8 is quite similar to the GN model for N f → ∞, i.e. quantum fluctuations at N f = 8 seem to be rather weak. We note that the stair-like behavior is a consequence of the finite spatial extent L. To support this we minimized the SLAC regularized GN action with a specific ansatz for the chiral condensate, in the variables A and q. This ansatz is a reasonable approximation for the considered values of µ as the analytically known chiral condensate in the N f , L → ∞ case rapidly reduces to a cos-shape for increasing chemical potential. This enables us to calculate B for N f → ∞ and finite L. The result for a ≈ 0.410/σ 0 and N s = 63 corresponding to L = N s a ≈ 25.8/σ 0 (i.e. the same lattice spacing and extent used in our simulations) exhibits clear steps as shown in Fig. 10. The steps disappear for L → ∞ as shown by n B | N f ,L→∞ L, which is also plotted in Fig. 10. For µ/σ 0 < ∼ 0.65 autocorrelation times are extremely large and we were not able to reach thermal equilibrium in our simulations. This is shown in Fig. 11, where we present results for B obtained from cold starts and from hot starts. For 0.4 < ∼ µ/σ 0 < ∼ 0.65 the cold and the hot curves differ and the largest discrepancy is observed close to the phase transition around µ/σ 0 ≈ 0.51. We expect the true result for B(µ) to be somewhere between the two curves, i.e. the cold and the hot results represent lower and upper bounds for B(µ). This is also supported by our simulation results for B(µ) at T /σ 0 ≈ 0.076 (see section III B 1), which is bounded by the cold and the hot curves obtained at T /σ 0 ≈ 0.038 in this range. Fig. 11 as well as the exceedingly long autocorrelation times remind us of hysteresis effects near a first order transition. From the good agreement with the semi-classical picture, one could conjecture that in a finite volume the probability distribution e −S eff /Z has two peaks (due to the commensurability constraint) which leads to the observed hysteresis effects. The problem will probably go away for very large volumes, but for high-precision simulations on finite lattices near the transition improved algorithms are needed, which support the creation and annihilation of extended baryons.
We remark that consequences of the large autocorrelations are also visible in the phase diagram shown in Fig. 2. In the problematic region, i.e. for 0.40 < ∼ µ/σ 0 < ∼ 0.65 and T /σ 0 < ∼ 0.05, the boundary between the homogeneously broken phase and the inhomogeneous phase suddenly turns towards the origin, which amounts to an inhomogeneous phase larger than expected and qualitatively different from the N f → ∞ boundary. Knowing that all simulations for this phase diagram were started with hot field configurations, this behavior can now be understood as a thermalization problem (see Fig. 11, where the result for B(µ) obtained with a hot start incorrectly indicates the phase boundary at a rather small value for µ.

IV. CONCLUSIONS
In this work we investigated the distribution of the baryon density n B in the N f = 8 GN model enclosed in a finite box of size L at finite chemical potential µ. The simulations were performed with chiral SLAC fermions. We compared with recent results on the spatial inhomogeneities in the 1 + 1-dimensional GN model [1], which is interpreted as modulation of baryonic matter density as observed in the N f → ∞ limit of the model [20] and which is well-known in solid state physics [25,26]. Since translation symmetry is inherent to Monte Carlo simulations on finite lattices (the phase of a quasi-periodic configuration is a collective parameter) we could not measure the baryon density n B (x) directly. Instead we found a strong correlation between the dominant wave number of the spatial inhomogeneities and the baryon number. This is clear evidence for a region in the phase diagram corresponding to a regime of modulated baryonic matter. Via this detour we explicitly circumvented the question about the breaking of translation symmetry. The delicate question whether we found a rigid baryon crystal (as seen in the large-N f limit) or rather a baryonic liquid [50], where the baryons have a preferred separation locally, but are disordered on large scales, has been addressed in Ref. [1] and needs further investigations with improved algorithms.
Our results shed further light on what can happen in quantum field theories at large fermion densities. In particular, it shows that mean-field and large-N f approximations may contain more (hidden) information on the physics at finite N f than one would expect. This is reassuring, since in particle physics and even more so in solid state physics we often rely on these approximations. Our results may also be of relevance in condensed matter systems, e.g. for large, almost 1-dimensional polymers [51].
However, it is still unclear if the above results are relevant for QCD. On the one hand, we established that the interpretation as baryonic matter is not spoiled by taking quantum fluctuations into account. On the other hand, although recent numerical lattice-studies of the GN model in 2 + 1 dimensions and for N f → ∞ spotted inhomogenous condensates, the spatial modulation is related to the cutoff scale and disappears in the continuum limit [52,53]. Clearly, if this happens in the limit N f → ∞, then we cannot expect a breaking of translation invariance for a finite number of flavours. Thus, extending our numerical studies and simulations to fermion systems in higher dimensions is an important task. Interesting candidates are for example the Nambu-Jona-Lasinio model or the quark-meson model in 3 + 1 dimensions.
We thank Philippe de Forcrand for bringing reference [50] to our attention.
We also thank Marc Winstel for providing the code base which was used to obtain the finite Volume N f → ∞ lattice results.
Furthermore we thank Michael Thies for providing us the analytic expression for the spatial correlator in the N f , L → ∞ limit.
J  In this appendix we calculate the correction term, which must be added to the baryon density, when using a linearly coupled chemical potential and removing the lattice cutoff as it is typically done in lattice field theory.
We consider non-interacting fermions linearly coupled to µ enclosed in a (d − 1)-dimensional spatial box of linear extent L subject to periodic boundary conditions. The allowed wave numbers of the corresponding fermion field at finite temperature T = 1/β are with Matsubara frequencies The corresponding eigenvalues of the free Dirac operator with chemical potential are such that Because of spin, the degeneracy of the eigenvalues is denotes the greatest integer less than or equal to x. The logarithm of the grand partition function Z(µ) divided by βV is the pressure p. The µ-derivative of p is the baryon density, with The sum defining v(E k ) is convergent for any fixed k . But when we sum over the spatial momenta, E k (which is the positive square root of E 2 k ) becomes arbitrarily large and we will show that removing the cutoffs in frequency space and momentum space does not commute.
To this end, we regularize the sum over the Matsubara frequencies (A6) by only admitting N t frequencies. It is convenient to choose these frequencies symmetric to the origin (which is only possible for even N t ), i.e. we restrict ω n according to To calculate the truncated series over the ω n , denoted by v(E k , N t ), we combine terms with ±ω n and obtain v(E k , N t ) = (A8) where (E k → −E k ) represents the previous term with opposite sign of the energy. v(E k , N t ) is finite for N t → ∞ and can be calculated, Inserting (A9) into (A5) yields the baryon density An integration with respect to the chemical potential gives the pressure of the Fermi gas. The integration constant is the divergent contribution of the quantum fluctuations at zero temperature and zero chemical potential, To summarize: if we first perform the continuum limit in Euclidean time direction, which means sum over all n ∈ Z in Eq. (A6), then we get the sum E k at µ = T = 0 plus the finite sum known from quantum statistics.
Since the sum over the Matsubara frequencies is only conditionally convergent and the sum over the k is divergent, we get a different result, when we remove the cutoff in the Matsubara frequencies together with the cutoff in the spatial momenta, as one does in lattice field theory. To show that, we consider the difference between the series (A9) and the finite sum (A8), where we introduced the function In the following we focus on free fermions in 1 + 1 dimensions, where only the first term gives a finite contribution to the error. The other terms are suppressed by inverse powers of N t . Thus, keeping the relevant term we arrive at For large spatial momenta we have µ E k and To study the error for the baryon density ∆n B = 1 βL k ∆v(E k , N t ) (A18) (in 2 spacetime dimensions C = 1) as a function of the temporal and spatial cutoffs, we cut off the spatial momenta as (A19) For convenience we choose the spatial momenta symmetric to the origin. Inserting Eq. (A17) into the regularized sum (A18) we can calculate the leading term with the help of (A21) If we regularize the system on a lattice with the same lattice constant a in temporal and spatial direction then the argument of atan is equal to 1 and we obtain To summarize: If we use a linearly coupled chemical potential for free fermions on the lattice, then we must subtract from the resulting baryon density n B a term linear in µ, in order to recover the result obtained in a conventional continuum calculation. In spacetime dimensions d > 2 the correction term actually is UV-divergent, since an analogous estimate reveals that the suppression by inverse powers of κ ∝ N t does not anymore balance the sum over the spatial momenta. Finally, we note that higher order terms in this correction (vanishing for N t → ∞) are now straightforward to compute. For future reference, we just show the O(a) correction on a hypercubic lattice with lattice spacing a (β = aN t , L = aN s ): In higher orders in a, also terms ∼ µ 2n+1 , n ∈ N, appear.