Entanglement entropies of an interval in the free Schr\"odinger field theory at finite density

We study the entanglement entropies of an interval on the infinite line in the free fermionic spinless Schr\"odinger field theory at finite density and zero temperature, which is a non-relativistic model with Lifshitz exponent $z=2$. We prove that the entanglement entropies are finite functions of one dimensionless parameter proportional to the area of a rectangular region in the phase space determined by the Fermi momentum and the length of the interval. The entanglement entropy is a monotonically increasing function. By employing the properties of the prolate spheroidal wave functions of order zero or the asymptotic expansions of the tau function of the sine kernel, we find analytic expressions for the expansions of the entanglement entropies in the asymptotic regimes of small and large area of the rectangular region in the phase space. These expansions lead to prove that the analogue of the relativistic entropic $C$ function is not monotonous. Extending our analyses to a class of free fermionic Lifshitz models labelled by their integer dynamical exponent $z$, we find that the parity of this exponent determines the properties of the bipartite entanglement for an interval on the line.


Introduction
The bipartite entanglement associated to a spatial bipartition has been largely studied during the past three decades in quantum field theories, quantum many-body systems and quantum gravity (see e.g. the reviews [1][2][3][4][5][6][7]). Consider a quantum system in a state characterised by the density matrix ρ and the bipartition A ∪ B of the space given by a region A and its complement B. Assuming that the Hilbert space of the system can be factorised as H = H A ⊗ H B , the reduced density matrix of the subsystem A is ρ A ≡ Tr H B ρ, where the normalisation condition Tr H A ρ A = 1 is imposed (hereafter we enlighten the notation by using Tr(· · · ) = Tr H A (· · · )).
The entanglement entropies provide an important set of quantities to study in order to understand the bipartite entanglement of the system in the state ρ for the bipartition A ∪ B. They are the entanglement entropy, the Rényi entropies and the single copy entanglement. The entanglement entropy S A is the von Neumann entropy of the reduced density matrix which can be obtained also through the analytic continuation α → 1 of the Rényi entropies (replica limit), defined in terms of the moments of ρ A as follows where α = 1 is a real and positive parameter. The replica limit in (1.1) naturally leads to identify S A ≡ S A . The single copy entanglement S (∞) A can be defined as the limit α → +∞ of the Rényi entropies (1.2) [8][9][10]. Since Tr(ρ α A ) = j λ α j , with λ j ∈ (0, 1) being the eigenvalues of ρ A , it is straightforward to realise that S (∞) A = − log(λ max ), where λ max is the largest eigenvalue of ρ A . Among these entanglement entropies, S A is the most important quantity because it measures the bipartite entanglement when ρ is a pure state.
Many fundamental results have been obtained for the entanglement entropies in relativistic quantum field theories in d + 1 spacetime dimensions. In this class of quantum field theories and when ρ is the ground state, the entanglement entropies are divergent as → 0, where is the ultraviolet (UV) cut off. The leading divergence exhibits the celebrated area law behaviour S (α) A ∝ Area(∂A)/ d−1 + · · · as → 0, where the dots correspond to the subleading terms [11][12][13][14]. Our analyses are restricted to d = 1 translation invariant quantum field theories on the line and to the bipartition given by an interval A; hence we are allowed to set A = [−R, R] without loss of generality. A notable exception to the area law behaviour are the conformal field theories in d = 1, where S (α) A = c 6 (1 + 1 α ) log(2R/ ) + · · · [15][16][17], with c being the central charge of the model, which occurs in the Virasoro algebra [18]. Another important result for the relativistic field theories in d = 1 involves the function C ≡ R ∂ R S A , constructed from the entanglement entropy of the interval A on the line when the entire system is in its ground state. Similarly to the Zamolodchikov's C function [19], the quantity C is UV finite and decreases monotonically along a renormalization group (RG) flow connecting a UV fixed point to an infrared (IR) fixed point [20]; hence it is often called entropic C function. The proof of this monotonic behaviour of C is based both on the strong subadditivity property of the entanglement entropy and on the relativistic invariance [20]. In relativistic quantum field theories, the effect of the finite density on the entanglement entropies has been also explored [21][22][23][24][25][26].
In order to gain some new insights about the relation between the spacetime symmetry and the characteristic features of entanglement, it is worth investigating the entanglement entropies in non-relativistic quantum field theories, which have also close connections with the quantum many-body systems. The properties of the bipartite entanglement quantifiers depend on whether the quantum field theory model displays a relativistic or a non-relativistic invariance. For instance, free Fermi systems at finite density exhibit a well known logarithmic violation of the area law of the entanglement entropy due to the occurrence of a Fermi surface [27,28] (for numerical results in lattice models, see e.g. [29,30]). Furthermore, an entropic C functions for non-relativistic quantum field theories in d = 1 is not known [26,31].
In the Wilsonian approach to quantum field theory, a fixed point corresponds to a scale invariant model and the scaling symmetry may not act on space and time in the same way. Under the assumption of spatial isotropy, the Lifshitz scale transformation is defined by t → χ z t and x → χ x for any spatial position vector x ∈ R d , where the parameter z > 0 is the Lifshitz exponent or dynamical critical exponent [32]. The Poincaré algebra, which characterises the relativistic field theories, has z = 1. The Schrödinger algebra has z = 2 [33][34][35][36][37][38].
We focus on d = 1 and consider the fermionic free Schrödinger field theory at finite density µ and on the infinite line, which is a free non-relativistic quantum field theory with z = 2. This model describes the dilute spinless Fermi gas in d = 1 [39]. The fermionic spinfull model in d = 1 in the presence of the quartic interaction has been studied through renormalization group methods [40][41][42]. Free fermionic spinless models in d = 1 with positive integer values of z have been also considered [43]. Let us remark that, in our analysis of the entanglement entropies for this free model, we do not approximate the dispersion relation with a linear dispersion relation at the Fermi points (Tomonaga's approximation) [44,45].
In d = 2, an interesting model with z = 2 called quantum Lifshitz model has been introduced in [46]: it is a free bosonic quantum field theory with the symmetries of a Lifshitz critical point with z = 2 and its bipartite entanglement has been studied in [47][48][49], finding an area law behaviour. This model belongs to a class of Lifshitz theories in d 2 having z = d [50], whose bipartite entanglement entropy for the ground state has been investigated in [51]. The entanglement entropy in non-relativistic Schrödinger models in d 2 has been computed also through heat kernel methods [52]. In d = 1, numerical studies of the entan-glement entropy in bosonic lattice models with various integer z have been reported e.g. in [53,54]. In the context of the gauge/gravity correspondence, gravitational backgrounds dual to Lifshitz spacetimes have been introduced [55][56][57][58] (see the review [59]) and the holographic entanglement entropy for z = 1 has been computed in various settings [23,[60][61][62][63][64][65][66][67].
In this manuscript we study the entanglement entropies of an interval A = [−R, R] on the line for the free fermionic spinless Schrödinger field theory at zero temperature and finite density µ. When µ = 0, we have S A = 0, as observed in [43,68,69].
At finite density µ > 0, we find that the entanglement entropies are finite functions of one variable given by the dimensionless parameter η ≡ R k F , where p F ≡ k F is the Fermi momentum. This parameter is proportional to the area of a limited rectangular region in the phase space that can be naturally identified from the interval A and the Fermi momentum. The finiteness of the entanglement entropies is proved by exploiting the properties of the solution of the spectral problem associated to the sine kernel in the interval A, first reported in a series of seminal papers by Slepian, Pollak and Landau [70][71][72][73] in terms of the prolate spheroidal wave functions (PSWFs) of order zero (see also the overview [74] and the recent book [75]). Efficient algorithms for the numerical evaluation of these functions have been developed (see [75] and references therein). The role of PSWFs in the context of the entanglement for free fermions has been studied in [76], which has inspired our work.
We find analytic results for the expansions of entanglement entropies in the regimes of small and large values of η. These results are obtained by employing the method proposed in [77,78] for the entanglement entropies in some spin chains (further developed in [79,80] to include the subleading terms) and the asymptotic expansions of the sine kernel tau function [81][82][83]. The latter expansions extend previous results [84][85][86][87][88] and have been found by adapting to the Painlevé V equation the method introduced in the seminal paper [89] to write the general solution of the Painlevé VI equation.
By applying the results of Fisher-Hartwig [90][91][92] and Widom [93,94] to the corresponding matrix spectral problem, the asymptotic behaviour of the entanglement entropies for various models and in different space dimensions have been extensively studied [77-80, 95, 96]. For the Schrödinger field theory on the line considered in this manuscript, the spectral problem in the continuum has an explicit solution in terms of the PSWFs and this allows to analyse the entropies in the whole range of parameters and not only in the asymptotic regime.
We study also the entanglement entropies of the interval A in the hierarchy of free Lifshitz fermions on the line considered in [43], in the massless case at zero temperature and finite density. These models are labelled by their integer Lifshitz exponent z 1. The special cases z = 1 and z = 2 correspond respectively to the relativistic chiral fermion and to the Schrödinger field theory introduced above.
The outline of the manuscript is as follows. In Sec. 2 the free fermionic spinless Schrödinger field theory on the line at finite density µ is briefly described. In Sec. 3 we present the solution of the spectral problem associated to the sine kernel in the interval and some properties of its spectrum. The spectrum of this kernel is employed in Sec. 4 to evaluate the entanglement entropies of the interval on the line. In Sec. 5 we discuss the behaviour of the entanglement entropies and of the quantity analogue to the relativistic entropic C function along the flow generated by the dimensionless parameter η. In Sec. 6 we explore the entanglement entropies of an interval for a class of d = 1 free fermionic models with integer Lifshitz exponents z 1. The expansions of the entanglement entropies in the asymptotic regimes of small η and large η are investigated in Sec. 7 and Sec. 8 respectively. In Sec. 9 we explore the Schatten norms of the sine kernel. In Sec. 10 we draw some conclusions. The derivations of some formulas and further auxiliary results are reported in the Appendices A, B, C, D, E, F and G.
2 Free Schrödinger field theory at finite density 2

.1 The model
We consider the non-relativistic spinless complex fermion field ψ(t, x) of mass m on the line, which evolves with the Schrödinger Hamiltonian and satisfies the equal-time canonical anticommutation relations ψ(t, x 1 ) , ψ(t, x 2 ) = ψ * (t, x 1 ) , ψ * (t, x 2 ) = 0 (2.2) where * stands for Hermitian conjugation. The solution of the equation of motion ∂ 2 x ψ(t, x) = 0 (2.4) which satisfies also (2.2) and (2.3) is [39,97,98] ψ(t, x) = ∞ −∞ a(k) e −iω(k)t e ikx dk 2π ω(k) = 2m k 2 (2.5) where {a(k), a * (k) : k ∈ R} generate the canonical anticommutation relation (CAR) algebra A given by a(k) , a(p) = a * (k) , a * (p) = 0 a(k) , a * (p) = 2π δ(k − p) . (2.6) In order to determine the state space of the system, a Hilbert space representation of the algebra A must be fixed. Before doing that, let us recall the symmetry content of the model. The system has an internal U (1) symmetry that implies the current conservation where (t, x) and j(t, x) are the particle density and current respectively, which are written through the field as follows (t, x) = ψ * (t, x) ψ(t, x) j(t, x) = i 2m (∂ x ψ * )(t, x) ψ(t, x)−ψ * (t, x) (∂ x ψ)(t, x) . (2.8) The space-time symmetries of (2.4) form the Schrödinger group [33,34,36], whose Lie algebra is generated by the momentum P , the Galilean boost G, the dilatation D and the special conformal transformation K. In terms of (2.8), these generators read (2.10) Using the canonical anti-commutation relations (2.2) and (2.3), one finds [36] [ where N in (2.12) is the particle number operator This operator commutes with all Schrödinger generators and defines a central extension of the Schrödinger algebra associated to a non-trivial cocycle [33,34]. The local form of the central term in the r.h.s. of the commutator [G, P ] in (2.12) reads c(t, y) ≡ m (t, y) . (2.14) For any representation of the CAR algebra A in (2.6), which is generated by a time and space invariant cyclic vector Ω, the expectation value is a t and x-independent dimensionless parameter (in our convention, t and x are measured in the same units) which characterises not only the central extension of the Schrödinger algebra, but also the representation of A.
When µ > 0 and at equal time t 1 = t 2 ≡ t, the correlators (2.28) and (2.29) become respectively [76] ψ * (t, where the Fermi wave number k F is defined in terms of the Fermi momentum p F as follows The two-point function (2.30) is a projector on the line; indeed, it satisfies +∞ −∞ sin(k F x 12 ) πx 12 sin(k F x 23 ) πx 23 dx 2 = sin(k F x 13 ) πx 13 (2.33) which tells us that the state providing the correlators (2.30) and (2.31) is pure [100]. The expectation value (2.30) provides the kernel of the integral operator whose spectrum is employed in the evaluation of the entanglement entropies.
In the zero temperature limit, the Gibbs representation with µ > 0 is characterised by In particular, c µ is proportional to the Fermi velocity v F = p F /m, i.e. the velocity at the Fermi surface, which is dimensionless in our conventions.

The parameter η
In this manuscript we consider the bipartition of the infinite line given by an interval A and its complement. The translation invariance allows to choose A = [−R , R] without loss of generality.
The mean value of the particle number operator N A ≡ A (t, x) dx in A is [102] where we have introduced the dimensionless parameter η ≡ R k F (2.38) and also its rescalings Notice thatã is the area of the following region in the phase space obtained through a space-momentum limiting process, which is the phase space analogue of the time-frequency limiting process considered for signals in the seminal papers [70][71][72][73][74]. In this analogy, the Heisenberg principle corresponds to the impossibility of the simultaneous confinement of a signal and of the amplitude of its frequency spectrum. Thus, non-zero bandlimited signals (i.e. signals whose frequencies belong to the finite band (−W, W), where W is called bandwidth) cannot be also time-limited (i.e. non-vanishing only in a finite interval of time). However, band-limited signal can be observed for a finite amount of time T . Within this parallelism, p F and R are the analogues of W and T respectively. Following this terminology, we call (2.40) limited phase space and its area isã = 4 η.
The dimensionless parameter η (or a, equivalently) plays a crucial role in the analysis of the entanglement entropies of A described in the subsequent sections. We always keep η finite and non-vanishing. The limits η → 0 and η → ∞ are taken only in the final results.

The spectral problem for the sine kernel
In this section we describe the spectral problem corresponding to the sine kernel in the interval A = [−R, R], whose spectrum determines the entanglement entropy of the bipartition of the line given by A and its complement, as discussed in Sec. 4.
The spectral problem defined by the correlator (2.30) restricted to the interval A = [−R, R] is equivalent to (3.1); indeed, it reads where x ∈ [−R, R] and the dimensionless parameter η has been introduced in (2.38). In particular, the same eigenvalues occur in (3.1) and (3.3) and the factor 1/ √ R in the eigenfunctions of (3.3) has been introduced to satisfy the normalisation condition imposed by the standard inner product of L 2 [−R, R].
In the two limits η → 0 and η → ∞, the kernel (3.2) degenerates to the vanishing kernel and to the Dirac delta respectively.
We remark that, in the relativistic model given by the massless Dirac field in one spatial dimension, the kernel defined by the two-point function of the massless Dirac field restricted to an interval (see (6.3) for z = 1) is not a compact operator. As a consequence, its spectrum is not discrete and spans the set [0, 1] (see e.g. [106][107][108][109]).
The eigenvalues and eigenfunctions in (3.1) can be expressed in terms of the prolate spheroidal wave functions (PSWFs), which have been introduced as solutions of the Helmholtz wave equation in spheroidal coordinates [110][111][112]. In particular, the eigenvalues can be written through the radial PSWFs of zero order R 0n as follows [70,75,105] γ n (η) = 2η π R 0n (η, 1) 2 n ∈ N 0 (3.8) and the eigenfunctions in terms of the angular PSWFs of zero order S 0n as The functions R 0n (η, x) and S 0n (η, x) in (3.8) and (3.9) can be studied e.g. through Wolfram Mathematica, where for the radial and angular PSWFs correspond to the built-in symbols R mn (x 1 , x 2 ) = SpheroidalS1[n, m, x 1 , x 2 ] and S mn (x 1 , x 2 ) = SpheroidalPS[n, m, x 1 , x 2 ] respectively [113], where m is the order of the PSWF and a particular normalisation is adopted. In (3.8) and (3.9), the normalisation chosen by Wolfram Mathematica for R 0n (η, x) and S 0n (η, x) is compatible with the normalisation induced by the standard inner product of L 2 [−1, 1] for the eigenfunctions f n (η; x), which has been imposed above. Notice that different normalisations for the PSWFs have been introduced in the literature [70,105].
The eigenvalues (3.8) are arranged in decreasing order. The eigenfunction f n (η; x) corresponding to γ n has a definite parity under x → −x, which is equal to the parity of n. Furthermore, f n (η; x) has n simple roots in (−1, 1) (see Theorem 2.3 in [75]).
In our numerical analysis we have employed an optimised Fortran code provided to us by Vladimir Rokhlin 1 , which is based on the results discussed in [75,114,115] regarding the numerical evaluation of the PSWFs. This code is faster than Mathematica and it also provides reliable results for the spectrum in the large η regime (η 350).
The arrangement of the eigenvalues in the decreasing order leads to identify a natural partition of the spectrum in the following three subsequent sets (see e.g. [116,117]): I a finite slow evolution region corresponding to 2η/π − n log(η) where γ n ≈ 1; II a finite plunge region for |n − 2η/π| log(η) where the eigenvalues rapidly change; III an infinite fast decay region corresponding to n − 2η/π log(η) where γ n ≈ 0 and γ n → 0 as n → +∞ at a super-exponential rate.
In Fig. 1 we show the first 31 eigenvalues of the spectra for four values of η such that the three parts I, II and III of the spectrum can be neatly identified. The partition of the spectrum in its three subsets I, II and III naturally leads to introduce the critical index n 0 as follows (see Eq. (3.9) of [75]) where x denotes the integer part of x. We remark that n 0 is a function of η. In a numerical inspection based on 50000 values of η ∈ [0.02, 1000], we have observed that γ n > 1/2 for all n < n 0 and γ n < 1/2 for all n > n 0 , while whether γ n 0 is greater or smaller than 1/2 depends on η. In this analysis we have also found that the eigenvalue closest to 1/2 is γ n 0 for the vast majority of the values of η explored, while for few cases it is given by From these observations and the discreteness of the spectrum, one realises that the behaviour of γ n 0 in terms of η is oscillatory. This is shown in Fig. 2, which displays the sawtoothed profile of γ n 0 in a typical finite range of η.

Entanglement entropies
In this section we study the entanglement entropies of the interval A = [−R, R] on the line for the Schrödinger field theory at zero temperature and finite density described in Sec. 2.
The Schrödinger field theory is a free fermionic model; hence its entanglement entropies (1.2) and (1.1) can be evaluated from the eigenvalues of the spectral problem (3.3) (or (3.1)) as follows [1,2,118] with finite α > 0 and α = 1. The limit S as α → +∞ defines the single-copy entanglement S (∞) A = − log(λ max ), which gives the largest eigenvalue λ max ∈ (0, 1) of the reduced density matrix [8][9][10]. For the model that we are considering, from (4.1) and (4.2), one obtains Since the spectrum of the sine kernel depends only on the dimensionless parameter η in (2.38), also the entanglement entropies in (4.1) and (4.3) depend only on η > 0. Let us remark that, denoting by {γ n ; n ∈ N 0 } the spectrum corresponding to the correlator (2.30) restricted to A (see (3.3)), the spectrum obtained from the correlator (2.31) restricted to A is {1 − γ n ; n ∈ N 0 }. Both these spectra provide the same entanglement entropies in (4.1) and (4.3) because the functions of x defined in (4.2) and (4.3) are invariant under x ↔ 1 − x.
Following [118], the single-particle entanglement energies ε n can be introduced as γ n = 1 e εn + 1 (4.4) i.e. ε n ≡ log(1/γ n − 1), which are negative for a finite number of n's. This allows to write the entanglement entropy in (4.1) in the following suggestive form [119] These expressions have the same form of the internal energy and of the thermodynamic potential of a quantum ideal Fermi gas respectively. Notice that Ω A < 0; hence S A > E A .
In the Schrödinger field theory at finite density, the entanglement entropies in (4.1) and (4.3) are finite. This crucial feature is due to the fast decay of γ n → 0 + as n → +∞. Indeed, for the eigenvalues γ n of the sine kernel, the following upper bound has been obtained (see Theorem 3.20 in [75]) 0 γ n γ n (4.7) whereγ n ≡g n η 2n+1g n ≡ 1 2 (n!) 2 (2n)! Γ(n + 3/2) (4.8) By employing the Stirling's approximation formula in (4.8), for the asymptotic behaviour of γ n as n → +∞ one finds (see [117,120] and Eq. (3.76) in [75]) which tells us that γ n vanishes as n → +∞ with a super-exponential decay rate. In the following we show that this important feature of the spectrum leads to finite entanglement entropies (4.1) and (4.3).
From (4.2) and (4.3), one first realises that 0 This sequence of inequalities tells us that the finiteness of S

(α)
A with α > 1 and the single-copy entanglement entropy S (∞) A . Considering a Rényi entropy with index α ∈ (0, 1), the corresponding function s α (x) is increasing for x ∈ [0, 1/2]; therefore we have where n 0 is the critical index (3.10) andγ n is the upper bound introduced in (4.9). In order to employ the ratio test for the convergence of a series for the one occurring in the last expression of (4.11), we consider the following limit where the first equality comes from (4.9), while the second equality has been obtained by using thatγ ∞,n → 0 as n → ∞ and that s α (x) = x α 1−α 1 + o(1) as x → 0 + . According to the ratio test, the vanishing of the limit in (4.12) implies the finiteness of the series occurring in the last expression of (4.11); hence (4.11) tells us that also S (α) A is finite. The finiteness of the Rényi entropies S
In Sec. 3, we mentioned the natural partition of the spectrum in the regions I, II and III. From a numerical inspection, we have observed that the contribution to the entanglement entropy coming from each part of this partition grows logarithmically for large values of η. In the Appendix A we discuss bounds for these three terms.
The entanglement entropy S A is a concave function of η. This important feature is a consequence of the strong subadditivity property of the entanglement entropy [121] which holds for any choice of two spatial regions A 1 and A 2 . Since the model that we are considering is defined in one spatial dimension and invariant under spatial translations, the entanglement entropy of an interval is a function of the length of the interval, for any given value of k F . Choosing A 1 = (−R, R) and A 2 = (−R + , R + ) for finite > 0, we have S A 1 = S A 2 because of the translation invariance; hence (4.13) becomes [122,123] (4.14) Assuming that the limit → 0 + of this inequality exists, one finds which implies that S A is a concave function of η. Figure 3: Entanglement entropies: The data points have been obtained from (4.1), while the dashed lines correspond to the small interval expansion (see (7.3)) and to the large interval expansion (see (8.14), (8.16) and (8.20)).
The data points displayed in Fig. 3 are numerical results for the entanglement entropies of the interval A = (−R, R) for some values of α and the data for the entanglement entropy correspond to α = 1. They have been obtained by using (4.1) and (4.3), with γ n computed through the Fortran code developed by Rokhlin and by evaluating the finite sums made by the contributions corresponding to the first n max = 5(n 0 + 2) eigenvalues, where n 0 is the critical value (3.10). We checked numerically that the entanglement entropies do not change significantly by increasing n max . This is due to the fact that n max is linear in n 0 , while the width of the region II depends logarithmically on n 0 (see (A.7)). This criterion for the truncation of the infinite sums have been adopted also for in the numerical determination of other quantities considered in this manuscript.
An analytic formula for the entanglement entropies as functions of η is not known in the entire domain η > 0; hence it is worth studying analytic expressions for their expansions in the regimes of small and large η. These analyses are discussed in Sec. 7 and Sec. 8 respectively. We study these expansions through two different approaches: one based on the expansions of the eigenvalues γ n (employed only in the regime of small η) [75,124] and another one based on the tau function of the sine kernel [81,82,87,125], which allows to obtain results both for η → 0 and for η → ∞.
The approach involving the sine kernel tau function is based on the method first employed in [77,78] for the entanglement entropies in some spin chains, where the spectrum of a Toeplitz matrix is involved. This method exploits the possibility to compute the entanglement entropies in (4.1) and (4.3) as specific contour integrals in the complex plane. Since 0 < γ n < 1 in our case, the closed path to consider must encircle the interval [0, 1] on the real axis. A natural choice is C = C 0 ∪ C − ∪ C 1 ∪ C + , where C 0 and C 1 are two arcs of radius /2 centered in 0 and 1 respectively, while C ± are the segments belonging to the horizontal lines x ± iδ with x ∈ R and intersecting C 0 and C 1 (see e.g. in Fig. 1 of [77], where a similar path is shown). Thus, the closed path C is parametrised by the infinitesimal parameters and δ. The entanglement entropies in (4.1) and (4.3) (we assume that α = 1 corresponds to the entanglement entropy) can be written as where the last expression has been obtained by using that s α (0) = 0 for the functions s α (x) defined in (4.2) and (4.3). By exchanging the summation with the integration along C in (4.18), one finds that the entanglement entropies can be computed as follows where s α (x) are the functions in (4.2) and (4.3) and τ is the sine kernel tau function (I denotes the identity operator) which corresponds to the Fredholm determinant of the sine kernel K in (3.2).
The expression (4.3) is obtained for holomorphic functions s α (z). This is not the case for the function s ∞ (z) in (4.3) providing the single copy entanglement. However, we apply (4.3) also for S (∞) A and then we check the outcomes with the corresponding numerical results.
The sine kernel tau function (4.20) has a long history and its expansions have been widely studied [84][85][86][87][88]. Its auxiliary function known as σ-form is the solution of a particular Painlevé V differential equation [84], as discussed in Sec. 5.2. Recent important advances in the analysis of the solutions of the Painlevé equations started with [89] have lead to find expansions of the sine kernel tau function to all orders [81][82][83]. In Sec. 7 and Sec. 8, we employ these results to obtain analytic expressions for the expansions in the small and large η regimes which correspond to the dashed curves in Fig. 3 (see (7.3) for small η and the combination of (8.14), (8.16) and (8.20) for large η).
The curves in Fig. 3 indicate that S

(α)
A vanish when η → 0. This is expected from (2.35), as emphasised in [43,68,69]; hence the limit k F → 0 and the evaluation of S (α) A commute. We find it worth anticipating that the leading term of the large η expansion for the entanglement entropy is S (α) A = 1 6 1+ 1 α log(η)+. . . as η → ∞ (see (8.14) and (8.16)), in agreement with the one dimensional case of the general result found in [27], obtained for fixed k F and R → ∞.
In Fig. 3 one observes that S

(α)
A with α = 1 display oscillations, while the entanglement entropy does not oscillate. The origin of the oscillatory behaviour can be identified with the sawtoothed behaviour of γ n 0 as function of η (see Fig. 2), which provides the largest contribution to the entanglement entropies. This argument leads to expect oscillations in the entanglement entropy as well, which, instead, are not observed. By numerical inspection, we noticed that this lack of oscillations in the entanglement entropy occurs only when all the infinite sum in (4.1) is taken into account. Indeed, by considering s(γ n 0 ) or finite sums p k=−p s(γ n 0 +k ) for some finite integer p 1 we observe oscillating curves and only for large values of p the oscillations disappear. It would be insightful to explore this feature further through the properties of the PSWFs.
Oscillations in the Rényi entropies have been investigated earlier also for a relativistic free massive fermion at finite density [26,126]. They could arise from localised terms on the defect that defines the Rényi entropy as a partition function on a sheeted Riemann manifold [126]. A further interpretation of this phenomenon has been provided in [26] through the defect operator product expansion in a relativistic setting. It would be interesting to apply the same method also in our non-relativistic scenario. However, we expect that the oscillations in the large η regime discussed in this manuscript (see Sec. 8) are recovered in the non-relativistic limit considered in [26].
In the Appendix B we discuss also the cumulants of the entanglement spectrum, which constitute an alternative to the moments of the reduced density matrix. We find that these quantities are finite and display an oscillatory behaviour (see Fig. 15), similarly to the Rényi entropies.
It is worth investigating the limits k F → 0 and k F → ∞. In these limits, the correlators in (2.30) and (2.31) become (2.35) and (2.36) respectively, and both these results give S A → +∞ respectively. The discrepancy of the limits in the large k F regime tells us that the limiting procedure k F → ∞ and the evaluation of the entanglement entropies do not commute. In other words, by writing S (α) we have that the two limits N → ∞ and η → ∞ do not commute. Indeed, by taking η → ∞ first and then N → ∞, one finds the expected result for k F → ∞ mentioned at the beginning of this paragraph because S A,N → 0 as η → ∞ for any finite N . This tells us that entanglement entropies in (4.1) and (4.3) are not uniformly convergent. It is worth considering the limit → 0 of our results (see also the recent analysis in [127]), which is determined by the behaviour of the dimensionless parameter η in this limit. From (2.38) and (2.32), one realises that different results can be obtained for entanglement entropies when → 0, depending on the quantities that are kept constant in this limit. For instance, for fixed m and R, the limit → 0 depends on the behaviour of µ/ 2 . In particular, if µ = O( 2 ) then η remains finite and non-vanishing [128].

Entanglement along the flow generated by η
Quantifying the loss/gain of information along a flow in the space of parameters (masses, coupling constants, etc.) characterising a given model is a challenging task. For instance, along a RG flow a loss of information is expected, in some heuristic sense [129]. In this section we discuss this issue in the specific non-relativistic free field theory that we are exploring and for the flow parameterised by η, which is not a RG flow. The information quantifiers that we consider are the entanglement entropy (Sec. 5.1) and the quantity analogue to the entropic C function introduced in [20] for the relativistic field theories in d = 1 (Sec. 5.2).

Entanglement entropy loss for decreasing η
In the free fermionic Schrödinger field theory on the line, at zero temperature and finite density, we found that the entanglement entropy S A of an interval A is finite (see Sec. 4). This crucial property makes S A a natural candidate to quantify the amount of information shared by the two parts of this bipartition. Let us recall that, for any quantum system in a pure state and for any bipartition A ∪ B of the space, S A = S B and it measures the bipartite entanglement associated to the state and to the bipartition.
Consider two finite and non-vanishing values of η, assuming 0 < η 1 < η 2 without loss of generality. The area of the limited phase space (2.40), given in (2.39), is proportional to the dimensionless parameter η and, in particular, it decreases as η decreases.
As for the entanglement entropy, the curves corresponding to α = 1 in Fig. 3 and Fig. 4 show that S A increases monotonically as η increases. This important feature is not observed for entanglement entropies with α = 1.
For a given interval, changing µ modifies the area of the limited phase space. In this case, S A increases as µ increases. For the states at zero temperature and different finite density that we are exploring, this behaviour is observed also for the mean energy E(t, x) ∞,µ , the mean density (t, x) ∞,µ and the central term c µ of the Schödinger algebra, given in (2.34).
For the translational invariant system we are considering, one can prove that S A is a strictly increasing function of η. The argument [130] is based on the concavity of S A (see (4.15), which is a consequence of the strong subadditivity (4.13) and of the translation invariance) and on the logarithmic growth of S A for η → ∞ (see e.g. Fig. 4 and Sec. 8). The proof is by contradiction. Assuming that S A is not strictly increasing, there exist two points 0 η 1 < η 0 such that S A (η 1 ) S A (η 0 ). Because of the logarithmic divergence at infinity, there exists also η 2 > η 0 such that S A (η 2 ) > S A (η 0 ). Let us fix the parameter t by requiring that η 0 = t η 1 + (1 − t)η 2 . Then one has 0 < t < 1 and concavity implies therefore that This contradiction concludes the argument and proves that S A is strictly increasing with η.
This feature of S A , combined with its differentiability (that we assume here), provides some  where γ n ≡ ∂ η γ n and the single-particle entanglement energies (4.4) have been used. The last two expressions in (5.2) have been found by employing the following remarkable relation occurring between any eigenvalue γ n and the corresponding eigenfunction f n (η; x) in the sine kernel spectral problem (3.1) (see Eq. (3.51) in [75]) The negative terms in the series of S A in (5.2) correspond to the eigenvalues γ n > 1/2, which have ε n < 0 and are a finite number.
An interesting inequality involves the relative entropy of the two reduced density matrices ρ A,2 ≡ e −K A,2 and ρ A,1 ≡ e −K A,1 (normalised to Tr(ρ A,2 ) = Tr(ρ A,1 ) = 1) associated to the same interval A and to two different values µ 1 < µ 2 , where the operator K A ≡ − log(ρ A ) is the entanglement (or modular) Hamiltonian [1,107,118,[131][132][133][134]. The relative entropy is [135,136] We remind that S(ρ A,2 |ρ A,1 ) 0, where the equality holds if and only if ρ A,1 = ρ A,2 . This inequality and the assumption ∆S A > 0 for η 2 > η 1 imply respectively that where the first inequality is the Bekenstein bound in the form discussed in [137]. In (5.8) both ∆S A with S(ρ A,2 |ρ A,1 ) are strictly positive and it would be interesting to compare them.
In the model we are exploring, it is natural to explore the majorization condition between two reduced density matrices ρ A,2 and ρ A,1 whose spectra are associated to two values η 1 and η 2 , which can be chosen η 2 > η 1 without loss of generality. This corresponds to either two different intervals with R 2 > R 1 at fixed k F or to the reduced density matrices of the same interval for two different states of the entire system characterised by Fermi momenta k F η 2 > k F η 1 . The single copy entanglement S (∞) A = − log(λ max ) shown in Fig. 4 (see also Fig. 7, Fig. 9 and Fig. 11) allows to explore the possible occurrence of a majorization relation between these two reduced density matrices. In particular, the oscillatory behaviour of S (∞) A rules out a majorization relation between ρ A,2 and ρ A,1 for two generic η 1 and η 2 because the validity of the inequality (5.9) for p = 0 depends on the specific choice of η 1 and η 2 .

The analogue of the relativistic entropic C function
In the class of the d = 1 relativistic quantum field theories, Zamolodchikov [19] constructed a finite C function that monotonically decreases along the RG flow and takes finite values at its fixed points equal to the central charges of the corresponding conformal field theories.
Considering the entanglement entropy associated to an interval of length in the infinite line when the system is in its ground state, Casini and Huerta [20] introduced a different function C ≡ ∂ S A , proving that it is UV finite, it takes finite values proportional to the central charge at the fixed points and it monotonically decreases along the RG flow; hence this function is usually called entropic C function. The proof of the monotonicity of this entropic C function is based on the relativistic invariance and on the strong subadditivity property of the entanglement entropy. We refer the interested reader to the review [145] for further discussions and references about this entropic C function and its generalisations to higher dimensions [24,[146][147][148]. Explicit examples of entropic C functions have been studied for relativistic free massive boson and Dirac fermion where M R is the dimensionless parameter generating the flow, with M being the relativistic mass of the field [2,149,150]. Both the Zamolodchikov's C function and the entropic C function are constructed through the ground state of the model along the RG flow.
In the non-relativistic model that we are exploring, it is worth considering the analogue of the entropic C function introduced for the relativistic models. This has been done in other settings e.g. in [26,151].
The entanglement entropies discussed in Sec. 4 allow us to introduce where C is the analogue of the relativistic entropic C function in our model. From (4.1) and the fact that γ n depend only on η, it is straightforward to find that the functions (5.10) can be written respectively as By employing (5.3), these expressions become respectively In the Appendix C we show that these series are well defined functions of the dimensionless parameter η. In our numerical analyses of C and C α , we have employed (5.12) by truncating the infinite sums as discussed in Sec. 4.
The expression (4.19) allows us to write the functions C and C α in (5.10) in the form in terms of the auxiliary function σ associated to a tau function σ ≡ η ∂ η log(τ ) (5.14) and in our case τ is the sine kernel tau function (4.20).
It is well known that the auxiliary function σ associated to the sine kernel tau function is the solution of the following Painlevé V differential equation (written in the σ-form) [84][85][86][87] equipped with the following boundary condition In Fig. 5 we show the quantity C defined in (5.10) in terms of η. The empty circles denote the data points obtained from (5.12), while the coloured dashed lines correspond to the asymptotic behaviours of C when η → 0 (green line) and when η → ∞ (red and magenta lines), which are analytic expressions derived respectively from (7.4) (see also (7.8)) and from (8.14), (8.16) and (8.20). In particular, we find that C = − 2 π η log(η) + . . . as η → 0 and that These two asymptotic behaviours and the assumption that C is a continuous function imply that C must possess at least one local maximum; hence it cannot be monotonous. This analysis cannot determine the number of local maxima. From the numerical data points in Fig. 5 we observe that C has only one local maximum. It would be insightful to find a proof for this numerical result.
In Fig. 6 we show the quantities C α introduced in (5.10) (see also Fig. 12). The numerical data points (empty circles) are obtained through (5.12) and the coloured dashed lines correspond to the asymptotic results derived in Sec 7 (see (7.4)) and Sec 8 (see (8.14), (8.16) and (8.20)). We remark that, while C does not oscillate, C α with α > 1 display an oscillatory behaviour, which tends to a sawtoothed curve as α → +∞, as discussed also at the end of Sec. 8.

Integer Lifshitz exponents
In this section we consider the family of Lifshitz fermion fields ψ(t, x) whose time evolution is given by In addition to (6.1), the equal time canonical anticommutation relations (2.2) and (2.3) are imposed on ψ(t, x). For z = 1, Eq. (6.1) gives the familiar relativistic equation of motion of a chiral fermion ψ(x − t). Instead, the case z = 2 corresponds to the Schrödinger equation (2.4). We find it worth considering the hierarchy of models corresponding to z ∈ N in a unified way by introducing the dispersion relation The solution of (6.1) satisfying the anticommutation relations (2.2) and (2.3) is still given by (2.5) with the substitution ω(k) −→ ω z (k). We observe in this respect that for even z one usually employs the alternative basis b(k) ≡ a(k) , c(k) ≡ a * (−k) : k 0 , where b(k) and c(k) are interpreted as annihilation operators of particles and antiparticles respectively. Performing the substitution ω(k) −→ ω z (k) in (2.24) and (2.25), one obtains the two point Lifshitz correlation functions in the Gibbs representation at temperature β and chemical potential µ.
Hereafter, it is convenient to distinguish between even and odd values of z. In the zero temperature limit β → ∞ one gets where n ∈ N and the Fermi momentum is Thus, p F,z=1 = µ, while p F,z → 2m as z → +∞. The special case µ = 0 has been already discussed in [43] and in the following we consider the case of non-vanishing µ.
The entanglement entropies of the interval A = [−R, R] for the Lifshitz fermions (6.1) can be obtained from (4.1), (4.3) and the spectrum of the kernel (6.3) restricted to A.
When z = 2n, it is straightforward to observe that the solution of this spectral problem is simply obtained by replacing η −→ R k F,2n in the solution of (3.3). Hence all the results obtained in this manuscript for the entanglement entropies in the model with z = 2 can be easily extended to the Lifshitz models with z = 2n through this simple replacement.
When z = 2n + 1, from (6.3) one must analyse the spectral problem the spectral problem (6.5) simplifies to which has the same eigenvalues as (6.6). The solution of the spectral problem (6.7) is well known [107,152,153] and, in particular, its eigenvalues read Contrary to case of even z, for odd z the spectrum (6.8) is continuous. This feature reflects the behavior of the dispersion relation (6.2), which is bounded from below for even z and unbounded for odd z. Another substantial difference between the eigenvalues (3.8) and (6.8) is that the latter ones do not depend on the Fermi momentum (6.4) and therefore they are independent of the Lifshitz exponent. As a consequence, all Lifshitz fermions with odd z have the same entanglement entropies, which coincide with the ones of the relativistic massless chiral fermion, i.e. [16,17] S (α) We remark that the independence of the spectrum on the Fermi momentum leads to a well known logarithmic ultraviolet divergency, which induces the presence of the UV cut off in (6.9).
Summarising, the entanglement entropies for the Lifshitz hierarchy of models given by (6.1) are ultraviolet divergent and µ-independent for odd z, while they are finite and µ-dependent for even z. This tells us that the entanglement entropies are heavily influenced by the global form of the dispersion relation and not only by its behaviour close to the Fermi momentum.
In the models characterised by odd values of z, we can find entanglement quantifiers that explicitly depend on the Lifshitz exponent. An important example is the entanglement Hamiltonian K A (also known as modular Hamiltonian) [131], which provides the reduced density matrix ρ A ∝ e −K A . For fermionic free models, this operator can be studied through the Peschel's formula [1,118], which has been largely explored [2,108,[153][154][155][156][157][158], also in its bosonic version [2,108,109,[159][160][161]. For Lifshitz fermion (6.1) with odd z and in the Gibbs state at zero temperature and finite density, we find that K A for the interval A ⊂ R is (the derivation is reported in the Appendix D) where k F,z is defined in (6.4), (t, x) is the particle density introduced in (2.8) and Notice that the two operators T tt (0, x) and (0, x) in (6.10) are normal ordered in the basis of oscillators given by a(k) , a * (k) : k ∈ R . In the special case of z = 1, the expression (6.10) becomes the entanglement Hamiltonian of the relativistic massless chiral fermion ψ(x − t) at zero temperature and finite density [21].
The reduced density matrix generates the one-parameter family of unitary operators {ρ iτ A : τ ∈ R}, which defines an automorphism on the operator algebra known as modular flow [131]. The modular flow of the field is defined as where ψ(x) is the initial configuration at t = 0. In Appendix D, by adapting the analysis described in [107,157], we find where β loc (x) is defined in (6.10) and (6.14) Notice that the solution (6.13) satisfies the initial condition ψ(0, x) = ψ(x), as expected.
Thus, both K A and the modular flow (6.13) explicitly depend on the Lifshitz exponent.
The expression for the field along the modular flow in (6.13) allows us to construct the corresponding correlation functions. For instance, we have where ξ j ≡ ξ(τ j , x j ) with j ∈ {1, 2}, ξ 12 ≡ ξ 1 − ξ 2 and τ 12 ≡ τ 1 − τ 2 . The correlator in the r.h.s. of (6.15) must be evaluated by employing (6.3) for odd values of z. Notice that (6.17) depends on the difference τ 12 of the modular parameters, as expected.
We remark that (6.17) has the structure identified in [154,157] for the modular correlators in other translation invariant cases. As a consequence, the correlator (6.17) satisfies the Kubo-Martin-Schwinger (KMS) condition [131] which is also a non-trivial consistency check of the entanglement Hamiltonian (6.10). In the Appendix D the partial differential equation satisfied by the correlator (6.15) is also reported.

Small η expansion
In this section we discuss the expansion of the entanglement entropies for small values of η. In Sec. 7.1 we employ the approach based on the tau function of the sine kernel (see (4.19) and (4.20)), while in Sec. 7.2 the expansion is obtained by exploiting the properties of the eigenvalues of the sine kernel.

Tau function approach
The expression (4.19), which is employed here also for the single copy entanglement, tells us that the expansion of S

(α)
A as η → 0 can be studied through the expansion of the sine kernel tau function (4.20) in this regime.
The small distance expansion of the tau function of a Painlevé V has been found in [81,82]. In the Appendix E.1 we specialise this result to the case of the sine kernel tau function (3.2), finding [81] where G(x) is the Barnes G function. This expansion corresponds to Eq. (5.13) of [81] written in the notation given in (E.5). The Taylor expansions of the functions B n (η) as η → 0 have been reported in the Appendix B of [81] for some values of n (they extend the earlier result reported in Eq. (8.114) of [88], which improves the previous expansions given in [85,86]). In (E.6)-(E.9) we have reported only the terms of these expansions employed in our analyses. We find it worth highlighting that the expansion (7.1) can be written also in terms of the area a of the limited phase space (see (2.39) and (2.40)).
Approximate analytic expressions for the entanglement entropies (4.19) are obtained by truncating the series (7.1) to a finite sum. In (7.1) a double series occurs because each B n (η) can be written through its Taylor expansion as η → 0. Given a positive integer N 1, the truncation condition of keeping all the term up to O(η N ) included leads to truncate also the series in n to a sum of N terms, where N satisfies N (N + 1) 2 − 1. Let us denote byτ N ,N the resulting finite sum, where N = N (N ). Since o(1/z N ) terms in (7.1) have been neglected, we have thatτ N ,N = P N,N (z)/z N , where P N,N (z) is a polynomial of order N whose coefficients are polynomials in η of different orders, up to order N included. Furthermore, since B n (0) = 1 for all the values of n that we consider (see [81] and the Appendix E.2), the polynomial P N,N (z) is monic. These observations lead to where z i ∈ P N,N are the zeros of P N,N (z), which are non-trivial functions of η whose explicit expressions depend on N .
The Abel-Ruffini theorem states that the roots of a polynomial of degree five or higher cannot be written through radicals. In our analysis, this implies that approximations corresponding to N 25 can be studied only numerically because N 5 is required in those cases. This leads us to consider N 4.
Plugging the finite sum (7.2) into (4.19) and exploiting the fact that s α (0) = 0, one obtains the following approximate result for the entanglement entropies where only the zeros of P N,N (z) belonging to [0, 1] contribute to this finite sum. In the Appendix E.2 we report the analytic expressions for the zeros of P N,N in terms of η, for various values of N 24. By introducing the N dependent parameter η * through the condition that at least one zero does not lie in (0, 1) or it has a non-vanishing imaginary part, for η ∈ (0, η * ) all the zeros belong to (0, 1) and therefore contribute to S

In this case (7.3) simplifies to S (α)
A;3 = s α (a 1 η). By expanding this result for η → 0, we obtain a leading term O η log(η) . This is confirmed by the numerical results reported in Fig. 7 and also by the expansion obtained through the PSWF in Sec. 7.2. Notice that, instead, expanding log(τ 3,1 ) as η → 0 first and then employing the resulting expansion in (4.19) leads to a wrong leading term O(η).
Improved approximations corresponding to N > 3 require N > 1; hence two or more terms can occur in (7.3). All the improved approximations characterised by N 24 are discussed in the Appendix E.2.
In Fig. 7 the numerical results for some entanglement entropies (black data points) are compared with the corresponding approximate analytic expressions S A; N reproduce the numerical data points becomes wider as N increases. The coloured dashed curves in Fig. 7 ends at some finite value η = η * where at least a zero of P N,N (z) lies outside the interval [0, 1]. The value η * depends on N and on the coefficients of the polynomial P N,N (z). In the simplest case we have P 1,3 (z) = z − 2η/π, hencez 1 = 2η/π, which leads to η * = π/2. For higher values of N we have determined η * numerically. The best approximation considered in this manuscript corresponds to N = 24 (red dashed curves). In this case, S

(α)
A; 24 perfectly agree with the corresponding numerical data for η 2.5, capturing also the second local maximum of S (α) A . In Fig. 8 we compare the numerical data points for C

(α)
A obtained numerically from (5.11) with the approximate analytic expressions C

(α)
A; N defined in (7.4), which hold for small η. The     in the small η regime. Different values of α are considered in the different panels. The numerical results correspond to the black empty markers and the coloured curves to the analytic approximate expressions as follows: the dashed curves are obtained from (7.3) for increasing values of N (see Sec. 7.1) and the solid blue lines from (7.8), (7.12) and (7.14) (see Sec. 7.2).    A; 24 , which is the best approximation considered in our analysis. The numerical data agree with the analytic results for the expansions as η → 0 + . At leading order we have C = O(η log(η)) when α = 1 and C α = O(η) when α > 1 (see also (7.8) and (7.13) respectively). Notice that C α with α > 1 displays an oscillatory behaviour whose amplitude grows with η. When α → ∞, this curve becomes sawtoothed (see also Fig. 6).

PSWF approach
The expansions of the entanglement entropies for small η can be studied also from (4.1) and (4.3) by employing the expansion of the eigenvalues as η → 0.
For the generic eigenvalue γ n of the spectral problem (3.1), it has been found that [124] whereγ n has been defined in (4.8). This tells us that γ n → 0 as η → 0, for any n 0.
Higher order terms in the expansion (7.7) can be written by including more terms in the expansion of s(x) = − x log(x) + x − x 2 /2 − x 3 /6 − x 4 /12 + O(x 5 ) and also taking into account the exponential correction occurring in the r.h.s. of (7.5). This leads to which has been obtained by considering only the terms coming from γ 0 and γ 1 . The expansion (7.8) is the best approximation allowed by (7.5). Indeed, s(γ n ) = O η 5 log(η) for n 2 and s(γ 0 ) contain a term of order O η 5 log(η) that cannot be evaluated because the O(η 4 ) term in the exponent of (7.5) has been neglected.
In order to study the expansion of the Rényi entropies with finite index α = 1 as η → 0, first we rewrite the function s α (x) in (4.2) as and then expand these two terms separately, obtaining By employing (7.5), (7.10) and (7.11), one finds that s α (γ n ) = O η min{5α,5} for n 2. Thus, by considering only the contributions coming from γ 0 and γ 1 in (4.1), for the Rényi entropies with finite α = 1 we obtain where the ordering of the terms based on their relevance depends on α. For instance, we have Some terms reported in (7.12) could be of order O η min{5α,5} , depending on α.
For the single-copy entanglement (4.3), by using s (7.14) In Fig. 7, the solid blue lines correspond to (7.8), (7.12), and (7.14). The width of the range of η where these curves agree with the numerical data (black data points) increases with α. We remark that better approximate expressions are obtained through the approach based on the sine kernel tau function. Indeed, the solid blue lines in Fig. 7 do not capture the first local maximum of the numerical data corresponding to α > 1, while the dashed lines nicely reproduce it when N is large enough. The curves corresponding to the analytic expressions obtained by applying the differential operator η ∂ η to (7.8), (7.12), and (7.14) have not been included in Fig. 8 to make this figure readable; but the ranges of η where they reproduce the numerical data are the same ones of Fig. 7, for any given value of α.
The results presented above can be compared with the ones discussed in Sec. 7.1. Considering e.g. the entanglement entropy, by expanding (7.3) with α = 1 as η → 0, one obtains a series that coincides with (7.8) up to a certain order which depends on the value of N chosen in (7.3). In order to obtain all the terms reported in (7.8), the expansion of (7.3) with N 4 must be considered. A similar analysis can be performed for the other entanglement entropies.

Large η expansion
The form (4.19) of the entanglement entropies (4.1) and (4.3) allows to study their expansions in the regime of large η by employing the expansion of the sine kernel tau function (4.20) in this regime.

Tau function
The asymptotic expansion of the sine kernel tau function for large η is a special case of the large distance expansion of the tau function for the general Painlevé V found in [82,83]. The details of this analysis are discussed in the Appendix F.1 and the result reads The functions D k (ν ) available in the literature are [82,83] D 0 (ν ) = 1 D 1 (ν ) = 4 ν 3 D 2 (ν ) = 8 ν 6 + 10 ν 4 . In our analyses we employ D k (ν ) with k ∈ {0, 1, 2, 3}.
We remark that the area a of the limited phase space (see (2.39) and (2.40)) is a natural variable for the expansion (8.1). We also stress that, differently from the small η expansion of the tau function given by the convergent series (7.1), the large η expansion (8.1) is asymptotic.
An approximate expression for (8.9) is obtained by considering only the summands corresponding to n = 1 and n = −1 in the second and third terms respectively. This gives Since the complex parameter ν takes the values (F.17) in the computation of the entanglement entropies, the approximation (8.10) allows to obtain their expansion up to O(1/η 4 ), as discussed in Sec. 8.2 and in the Appendix F.3 (see (8.14) and (8.20)).
By using the identities (2πz) 2 i 4ν = −[π/ sin(πν )] 2 and π sin(πw) = Γ(1−w) Γ(w), the factors multiplying the last two series in (8.10) can be written as follows where the term 1 in the last expression corresponds to k = 0 term of the first series.

Entanglement entropies
The expansion of the entanglement entropies for large η can be studied by employing (8.6) and (8.13) into (4.19). This leads to the decomposition A,∞ provide the oscillatory terms observed e.g. in Fig. 4. The term S (α) A,∞ in (8.14) can be computed by adapting to the case that we are considering in the continuum the analysis of [77,78] for the entanglement entropies of a block made by consecutive sites in the infinite one-dimensional spin- 1 2 Heisenberg XX chain in a magnetic field, which is based on the Fisher-Hartwig conjecture [90][91][92]. The result of this calculation, whose details are reported in the Appendix F.2, is where the argument of the logarithm is the area of the limited phase space (2.40) when = 1 and the constant term E α is defined as follows [77,80] In the limits α → 1 and α → +∞, this constant becomes respectively The result (8.16) can be obtained also by employing a result of Slepian [124] in this asymptotic regime, as shown in [76,163]. A rigorous derivation of the leading logarithm term in (8.16) has been provided in [165].
In order to study the expansion of the subleading term S

(α)
A,∞ in (8.14) as η → ∞, we adapt to our case in the continuum the analysis performed [80] in the spin- 1 2 Heisenberg XX chain in a magnetic field, which provides the terms subleading to the ones found in [77,78] by employing the generalised Fisher-Hartwig conjecture [91,92,166]. This analysis, described in the Appendix F.3.1, gives the following expansion for S A,∞,N for N ∈ {0, 1, 2, 3} are functions of η which can be conveniently written as The constant term S 12y 2 − 1 3 − 10 i j y 48y 4 − 16y 2 + 1 + 220y 4 − 114y 2 + 37 12 (8.27) which must be evaluated on the following points belonging to the imaginary axis  A at a given order as η → ∞. Notice that it is natural to write also (8.23) in terms of the area a of the limited phase space (see (2.39) and (2.40)).
The subleading contribution to the entanglement entropy corresponds to the limit α → 1 of S (α) A,∞ given by (8.20), (8.21), (8.22) and (8.23). It is remarkable that the oscillating quantity S A,∞ for large η has been studied, finding which can be written also in terms of the area a of the limited phase space (see (2.39) and (2.40)) and improves the expansion obtained in [80] in the lattice model up to term. We remark that in (F.88) all terms up to order O(1/η) are reported.
For the sake of completeness, in the appendix G we briefly discuss the lattice model where the Fisher-Hartwig formula and its generalisation have been applied and the double scaling limit providing the results in the continuum.
In Fig. 4 the entanglement entropies S

(α)
A are shown in a large domain η ∈ (0 , 100) which includes the one considered in Fig. 3. The filled circles are the data points obtained numerically, as discussed in Sec. 4, while the dashed lines correspond to the approximate analytic expressions derived in the regime of large η and given by (8.14), (8.16) and (8.20), which perfectly agree with the numerical results for η > 10. In this figure it is evident the logarithmic growth of S (α) A described by (8.16), the oscillatory behaviour for S

(α)
A when α > 1 (see (8.23)) and also the lack of such oscillations for the entanglement entropy. The approximate analytic expressions for large η mentioned above have been employed also in Fig. 3 and they nicely agree with the numerical data points for η > 1, which is quite remarkable.

The oscillations of the entanglement entropies S (α)
A with index α = 1 occurring in Fig. 4 have the same period equal to π/2, as one can observe from (8.23). In Fig. 9 we focus on the range η ∈ (990 , 1000) and show the change of these oscillations for different values of α in the regime of large η. The numerical data points in Fig. 9 (filled circles) are nicely reproduced by the dashed curves, obtained from (8.14), (8.16) and (8.20) with o(1/η 2 ) terms neglected. For small values of α > 1 harmonic oscillations are observed. As α increases, this behaviour changes and the local maxima of these oscillations become singular in the extreme case of the  single-copy entanglement.
Including more terms in the analytic expressions for the large η expansion leads to improved approximations for the entanglement entropies, as expected. This is shown in Fig. 10 for the case of α = 20, where we compare the numerical data (filled circles) with the coloured dashed curves obtained from (8.14), (8.16), (8.21) and with the sum in (8.23) truncated at different orders by neglecting the o(1/η r ) terms (in Fig. 10 and Fig. 11, "up to o(1/η r )" means that all the terms proportional to 1/η b with b r have been included, while "up to O(1/η r )" indicates that all the terms proportional to 1/η b with b < r have been considered). For instance, the almost horizontal grey curve is found by neglecting the o(1) terms; hence it corresponds to (8.16), which is the contribution given by the logarithmic and the constant term. Similarly the dashed cyan curve has been obtained by neglecting in (8.23) all the o(1/η 2 ) terms, finding a result which is almost indistinguishable from the numerical data in the range η ∈ (100 , 103). In our numerical analyses, we have also observed that, for a given range of large values of η, the agreement between the numerical data points and the curves obtained from the analytic results corresponding to a certain approximation improves as α decreases.
In Fig. 11 we show the single copy entanglement S approximates the numerical data points better than the expansion of the single copy entanglement. This may happen because (8.21) contains higher order terms as η → ∞ with respect to (8.29). In the inset of Fig. 11, we consider the domain η ∈ [996, 1000] and show that adding more terms in the expansion (8.29) does not necessarily improves the approximation of the numerical data point. Indeed, the best approximation corresponds to the truncation of (8.29) up to the term proportional to 1/(log η) 2 included (dashed red curve).
By applying the differential operator η ∂ η to the analytic expressions of the expansion obtained from (8.14), (8.16) and (8.21) and from (8.29), we obtain the large η asymptotics for the quantities C α introduced in (5.10). Although we do not report the explicit expressions here, we have employed the resulting analytic expressions (where O(1/η 2 ) terms have been discarded) to draw the dashed coloured curves in Fig. 12, which are compared with the numerical data points obtained from (5.12) (empty markers).
For the single copy entanglement, only the numerical data have been reported because in this range of η the analytic expression coming from the expansion at large η does not agree with the data points, as shown by the dashed red curve in Fig. 11. The oscillatory behaviour of C α with α = 1 about the constant value 1 6 1+ 1 α is due to the term η ∂ η S

(α)
A,∞,N,b coming from (8.23). The growing amplitudes of the oscillations of C α with α > 2 can be easily explained. Indeed, S

Schatten norms
In this section we study the p-th power of the Schatten p-norm K p of the sine kernel operator for integer p 1, which is defined as follows These series are convergent; hence they are well defined functions of η. For instance, for p = 1 we have P (1) A = 2η/π (see (3.7)). In the other cases, the convergence of (9.1) can be proved as follows. When p 1, we have that 0 < P A because the inequalities in (3.6) tell us that 0 < γ p n γ n . For 0 < p < 1, we observe that (4.7) provides an upper bound which is a convergent series, as one can show by using (4.9) and the ratio test, as done in (4.11) and (4.12) for the entanglement entropies.
We find it useful to write the p-th power of the Schatten p-norm of the sine kernel in terms of the sine kernel tau function (4.20). This can be done by adapting to (9.1) the procedure to obtain the contour integral (4.19) for the entanglement entropies, discussed in Sec. 4, and the result reads P (p) where we focus on integer values p 1, for simplicity.
(9.8) Figure 13: The p-th power of the Schatten p-norm of the sine kernel operator (see (9.1)) for some values of p > 0, in the regime of small η. When p = 1, the exact result is given by (3.7) for every η. The dashed red line corresponds to (3.7) (see also (9.9)) and the dashed black lines to (9.10-9.13).
The expansions of the functions B n (η) as η → 0 have been obtained in the Appendix B of [81] for n 5 and up to a certain order in η. We have reported them in the Appendix E.2 (see (E.6)-(E.9)), truncated to a certain order as discussed in Sec. 7.1.
By employing (9.7) into (9.6), the resulting integral can be evaluated by first expanding the integrand z p ∂ z log(τ ) as z → ∞ and then applying the residue theorem, which tells us to select the coefficient of the term corresponding to 1/z in the expansion of the integrand. For the first values of p, we obtain P (1) (9.13) in terms of the functions introduced in (9.8). Notice that (9.9) is equal to the mean value The p-th power of the Schatten p-norm of the sine kernel operator (see (9.1)) for some values of p > 0, where the linear divergence (9.16) for large η has been subtracted to highlight the subleading correction (see (9.17)). The black dashed straight line is −2η/π.
of the particle number operator given in (2.37) and to the trace of the sine kernel given in (3.7), as expected. We observe that the leading behaviour of P

(p)
A as η → 0 is O(η p ), which is determined by the term B 1 (η) p .
In Fig. 13 we show (9.1) in the regime of small η. The red dashed line is the exact result (3.7), while the black dashed lines correspond to the analytic expressions in (9.10)-(9.13), which hold as η → 0 and are polynomials up to O(η 31 ). Notice that the expansion of B 5 (η) in Appendix B of [81] has been employed for P (5) A and that B p (η) for p > 5 are not available in the literature.
In the large η regime, we can employ the factorised form (8.6), where the leading terms are contained inτ ∞ , defined in (8.7). Plugging this factorisation into (9.6), we find that the leading terms of (9.1) as η → ∞ are given by By using (F.6) with ν (z) defined in (8.2), the integrand of (9.14) reads in terms of the digamma function ψ(z) (see (F.7)).
The leading term of (9.14) comes from the linear term in η occurring in (9.15) and it can be easily evaluated by using the residues theorem, finding P (p) This linear divergence can be observed already in the range of η considering in Fig. 13. When p = 1, the exact result in (2.37) and (3.7) tells us that further corrections do not occur.
The linear divergence (9.16) can be understood also by observing that the largest contribution as η → ∞ comes from the eigenvalues γ n having n n 0 , where n 0 is the critical index (3.10). Furthermore, in the regime of large η, an upper bound for the contribution coming from the regions I and II of the spectrum (see Sec. 3) can be obtained by using γ p n < 1 and the Landau-Widom counting formula (A.12), where the coefficient of log(η) is positive.
The evaluation of the subleading term in (9.14), which corresponds to the term containing log(4η) in (9.15), is less straightforward. Indeed, splitting the contour C as discussed in Sec. 4, one finds that, while the integral over C 0 vanishes as → 0, the integral over C 1 diverges like log( ) + o(1). This logarithmic divergence cancels with the same divergence coming from the integral over C + ∪ C − , which gives a result proportional to − log( ) − γ E + ψ(p) . Combining these contributions, we obtain where γ E 0.577 is the Euler-Mascheroni constant. Notice that γ E + ψ(1) = 0, as expected. Furthermore, γ E + ψ(p) is a rational number when p ∈ N.
In order to highlight the subleading corrections for large η, in Fig

Conclusions
We have investigated the entanglement entropies S

(α)
A of an interval A in the free fermionic spinless Schrödinger field theory on the line, at finite density µ and zero temperature.
This problem can be studied without introducing approximations because the spectral problem (3.3) associated to the sine kernel in the interval A has been solved long ago by Slepian, Pollak and Landau in the seminal papers [70][71][72][73][74], which have generated a vast literature afterwards, in various directions [75]. The eigenvalues γ n of this spectral problem (see (3.8)) depend only on the dimensionless parameter η defined in (2.38), which is proportional to the area of the limited phase space (2.40). By employing the fact that γ n vanish superexponentially as n → ∞ (see (4.7) and (4.9)), we proved that the entanglement entropies

S (α)
A are finite functions of η. These functions are displayed e.g. in Fig. 3 and Fig. 4 for some values of α, where the data points have been generated through an efficient code optimised to evaluate numerically the PSWFs (kindly given to us by Vladimir Rokhlin).
In Sec. 5.1 it is shown that S A is a function that monotonically increases with the area of the limited phase space (see also the numerical results for S (α) A in Fig. 3 and Fig. 4). This property does not hold for the entanglement entropies with index α = 1. We proved that the analogue of the entropic C function for the d = 1 relativistic models [20], defined in (5.10), is not a monotonous function of η (see Fig. 5). Notice that, in the context of the gauge/gravity correspondence, non-monotonic holographic entropic C functions [170,171] have been also found when Lorentz invariance is broken [172].
We have studied also the entanglement entropies of an interval for a class of free fermionic massless Lifshitz models on the line (see (6.1)) at zero temperature and finite density, which are labelled by their integer Lifshitz exponent z 1. This class includes the relativistic massless chiral fermion (z = 1) and the above mentioned free fermionic spinless Schrödinger field theory (z = 2). Important qualitative differences in the entanglement entropies are observed, depending on the parity of z. For instance, the models with even z have finite and µ-dependent S

(α)
A , while the models with odd z have UV divergent and µ-independent S (α) A . For the subclass of models having odd Lifshitz exponents, we have computed the entanglement Hamiltonian (6.10), the modular flow of the field (6.13) and the corresponding correlation function (6.15), finding that, instead, these entanglement quantifiers explicitly depend both on µ and on the Lifshitz exponent through the Fermi momentum (6.4).
Finally, we have employed the method of [77,78] and the asymptotic expansions for small and large η of the sine kernel tau function (see (7.1) and (8.1) respectively) found in [81][82][83] to write the expansions of the entanglement entropies in these limiting regimes (see Sec. 7 and Sec. 8 respectively). The analytic expressions approximating the entanglement entropies obtained from these expansions give the dashed curves in all the figures from Fig. 3 to Fig. 12, which display a remarkable agreement with the numerical data points. These analytic results have also allowed to prove the non-monotonicity of the function C (see Sec. 5). Thus, our analysis provides a new application of the results for the general solution of the Painlevé V equation obtained in [81][82][83], specialised to the simple case of the sine kernel tau function.
Our results can be extended in many interesting directions. Since we have described an explicit example where the entanglement entropy is finite and monotonically increasing along the η-flow, it would be insightful to find whether this interesting feature occurs also in other models and whether it provides some new insights about the RG flows in non-relativistic field theories. For instance, it would be instructive to find a RG flow involving non-relativistic models where the entanglement entropy plays the role of the entropic C function.

B Cumulants of the entanglement spectrum
The reduced density matrix ρ A introduced in Sec. 1 (normalised by the condition Trρ A = 1) naturally leads to define the entanglement Hamiltonian K A as ρ A = e −K A . The spectrum of K A is called entanglement spectrum and its relevance has been discussed e.g. in [1,196]. Important information about the entanglement spectrum can be obtained by considering the moments of K A , i.e. Tr(ρ n A ), or, equivalently, the cumulants of the entanglement spectrum. In this appendix we evaluate the first cumulants of the entanglement spectrum corresponding to the interval A = [−R, R] ⊂ R for the free fermionic Schrödinger field theory at zero temperature and finite density considered in this manuscript.
The moments of the reduced density matrix are The moments of K A can be obtained through the following analytic continuation The logarithm of the generating function M A , provides the generating function of the cumulants of the entanglement spectrum C (p)   which are the connected correlators of K A . The entanglement entropy is the expectation value of K A ; hence it corresponds to p = 1, both in (B.2) and in (B.3). The second cumulant, which is (B.3) for p = 2 and gives the variance of K A , is called capacity of entanglement [197][198][199].
In the free fermionic Schrödinger field theory that we are exploring, S

(α)
A are (4.1) and (4.2); hence the cumulants of the entanglement spectrum (B.3) read For the first integer values of p, we obtain where τ is the sine kernel tau function (4.20) and C is the closed path in the complex plane introduced in Sec. 4. Thus, by adapting the analyses performed in Sec. 7 and Sec. 8, analytic expressions for the expansions of the cumulants either as η → 0 or as η → ∞ can be found.
In Fig. 15 we show the cumulants C A with p ∈ {2, 3, 4} display an oscillatory behaviour. From Fig. 15 we also observe that C A → 0 as η → 0 can be obtained by combining the fact that γ n (η) → 0 as η → 0 uniformly in n ∈ N and that c p (x) → 0 as x → 0.

C Finiteness of C α
In this appendix we show that the functions (5.10) are well defined functions of η, i.e. that the series (5.12) are convergent.

D Modular Hamiltonian and modular flow for odd Lifshitz exponents
In Sec. 6 we have observed that the entanglement entropies of an interval for the Lifshitz fermion fields with odd Lifshitz exponents are UV divergent and independent of k F,z , which contains both the density and the Fermi momentum. In this appendix we show that, instead, the modular Hamiltonian (6.10) and the corresponding modular flow of the field, defined in (6.12), depend on k F,z .
The free fermionic Lifshitz models that we are considering are quadratic field theories; hence the modular Hamiltonian K A of the interval A = [−R, R] ⊂ R can be written in the following quadratic form where : · · · : denotes the normal product in the oscillator algebra A introduced in Sec. 2. The kernel H A (x, y) in (D.1) can be found through the Peschel's formula [1,2,118] where C A is the kernel defined by the two point function (6.3) restricted to A.
When the Lifshitz exponent z is odd, the spectral problem is (6.5) and from (6.6), its eigenfunctions can be easily obtained from the eigenfunctions of the spectral problem (6.7), discussed in [107,152]. They are 3) The Peschel's formula (D.2) tells us that the kernels C A and H A have the same eigenfunctions and that the eigenvalues of H A are 2πs = log(1/γ s − 1) with s ∈ R, where γ s are the eigenvalues of C A , given in (6.8). Thus, the spectral representation of the kernel H A reads where (6.6) has been employed. The integral in the last expression of (D.4) does not contain k F,z ; hence H A (x, y) explicitly depends on k F,z . Plugging the explicit expression ofφ s (x) (given in (D.3)) into (D.4), we find Finally, the modular Hamiltonian (6.10) can be obtained by inserting this kernel into (D.1).
The modular flow of the field has been defined in (6.12). It is the solution the following partial differential equation where K A is (6.10), for a given initial field configuration ψ(x) when τ = 0. The explicit form of (D.6) reads This partial differential equation has the following structure where ψ(x) corresponds to the initial configuration of the field. The solution of (D.8) has been discussed e.g. in the Appendix B of [157] and it can be written as where ξ(t, x) is defined as follows which satisfies ξ(0, x) = x.
The partial differential equation (D.7) belongs to the class of partial differential equations defined by (D.8); indeed, it corresponds to In this special case, the functions and ξ(τ, x) and Φ(τ, x), defined respectively in (D.10) and (D.9), become respectively (6.14) and Finally, the expression (6.13) for the modular flow of the field is (D.9) specialised to the case given by (6.14), (D.11) and (D.12).

E On the small η expansion
In this appendix we derive the analytic expressions for the expansion of the entanglement entropies as η → 0. They are reported in Sec. 7.1 and employed in various figures.

E.1 Sine kernel tau function from the Painlevé V tau function
We find it instructive to report the derivation of (7.1) as a special case of the expansion of the solution of the general Painlevé V found in [81,82]. In the following, a given quantity indicated through a certain notation in [81,82] is denoted by the same symbol with the subindex .
We remark that the limit σ → 0 of (E.1) is not straightforward because G(m) = 0 for m ∈ Z and m 0. From (E.3) and the asymptotic behaviour of the Barnes G function (see e.g. Eq. (A.3) of [81]), we find In order to obtain the tau function for the sine kernel, we have to evaluate (E.1) in it (see Eq. (5.10) of [81]) and in the limit σ → 0. The function B 0 (n; it ) provides both the factor that cancels i n 2 in the denominator of (E.4) and B 0 (n; it ). In the Appendix B of [81] the function B 0 (n; it ) is called B sine (n; t ) and its expansions as t → 0 are reported, for various values of n (see also (E.6)-(E.9), where only the terms employed in our analysis are shown). By using (E.4), we can take the limit σ → 0 of (E.1) with N 0 = 1. Finally, we obtain (7.1) by introducing the following change of notation and by setting B 0 (n; it ) ≡ B n (η).

E.2 Approximate entanglement entropies
In Sec. 7.1 the approximate expression (7.3) for the entanglement entropies as η → 0 has been obtained from the small distance expansion of the sine kernel tau function given in (7.1). In the following we derive these approximate expressions up to O(η N ) included, for N 24. In order to consider terms up to O(η N ) included in the expansion of S

(α)
A as η → 0, we first write the finite sumτ N ,N = P N,N (z)/z N obtained by truncating (7.1), where P N,N (z) is a monic polynomial of order N . Then, the zeros of P N,N (z) are needed in (7.3) and they can be computed analytically through radicals only for N 4; hence we consider only N 4.
To writeτ N ,N , the expansions of the functions B n (η) in (7.1) as η → 0 must be taken into account up to the proper order. In the Appendix B of [81], the expansions of B n (η) up to a certain order in η for any n 5 have been obtained. In the following, by using also the first expression of (E.5), we report these expansions truncated to the order employed in this manuscript, where N 4. We use that The expression ∆ 3 is the discriminant of the cubic equation, whose three roots are real and distinct whenever ∆ 3 > 0. In our case ∆ 3 is a polynomial in η of high order, hence its positivity can be checked numerically. Specialising (7.3) to this case, one obtains the best approximation for the entanglement entropies when N = 3, namely which provides the curves corresponding to N = 15 in Fig. 7 as discussed in Sec. 7.1.

E.2.3 N = 4
The approximations corresponding to N = 4 have 16 N 24. Also in these cases we can obtain analytic expressions for the approximate entanglement entropies (7.3). Instead, for higher order approximations with N 25, where N 5 is required, the zeros of P N,N (z) cannot be found through analytic expressions involving radicals (Abel-Ruffini theorem). Thus, the best approximation accessible through analytic expressions corresponds to N = 24, which requires N = 4. In this case, by truncating (7.1) to these orders we find c ≡ c 9 η 9 − c 11 η 11 + c 13 η 13 − c 15 η 15 The optimal approximation for the entanglement entropies when N = 4, which corresponds to N = 24, is obtained by employing (E.43) into (7.3) specialised to N = 4. This gives the best approximation that we consider when η → 0 (as discussed in Sec. 7.1), namely which provides the curves corresponding to N = 24 in Fig. 7.

F On the large η expansion
In this appendix we discuss the derivation of the analytic results reported in Sec. 8 for the expansion of the entanglement entropies in the regime of large η, which have been employed in various figures.

F.1 Large η expansion of the sine kernel tau function
The asymptotic expansion (8.1) for the sine kernel tau function is a special case of the large η expansion for the tau function of the general Painlevé V, found in [82]. As done in the Appendix E.1, also in the following a certain quantity having the subindex corresponds to the same quantity (without this subindex) in [82] and in their notation, if not otherwise specified.
When θ 0 = θ t = θ * = 0 and in the limit σ → 0 (see (E.3)), the large distance expansion of the tau function given in Eq. (1.12a) of [82] simplifies to in terms of the parameters ρ and ν .

F.2 Contribution fromτ ∞
Consider the first integral in (8.15). From the definition ofτ ∞ in (8.7), we have where ν is defined in (8.2), which implies ∂ν ∂z = 1 2πi z(z−1) , and ψ(z) is the digamma function By using (F.6), the first integral in (8.15) can be written as the sum of a term proportional to η, a term proportional to log(η) and a constant term.
The leading term in the first integral in (8.15), proportional to η, vanishes because s α (0) = s α (1) = 0. Instead, the contribution proportional to log(η) is non-vanishing. It reads where in (F.9) we used that the integrals along C 0 and C 1 vanish and the following functions have been introduced in the remaining two integrals The result (F.10) provides both the logarithmic divergence and a contribution to the constant term. From (F.6), we have that the remaining contribution to the constant term in the first integral in (8.15) reads where the integral over C has been decomposed as done in (F.8) and we used the fact that, in the limits δ → 0 + and → 0 + , the contributions corresponding to C 0 and C 1 vanish and (F.11) holds. In the second step of (F.12), the integration variable y ≡ 1 2π log(1/x − 1) has been employed; hence (F.11) and (4.2) give respectively ν ± = −iy ± 1 2 and s(x(y)) = log(1 + e 2πy ) 1 + e 2πy + log(1 + e −2πy ) 1 + e −2πy (F.14) s α (x(y)) = 1 1 − α log 1 + e α2πy − α log 1 + e 2πy (F.15) which leads to s ∞ (x(y)) = log 1 + e −2π|y| . (F.16) The final expression (F.13) can be found by exploiting the identity ψ(1 + z) = ψ(z) + 1/z. Then, by employing in (F.13) the integral representation of the digamma function (see (F.7)) and exchanging the order of the two integrations, the integral in y can be performed, finding the integrand occurring in (8.17). Combining this result with (F.10), the expression (8.16) is obtained.
When α = 1, the function ∂ yŝ (y) has poles of order two at y k ≡ i(k + 1/2) with k ∈ Z. Instead, for finite α > 1, the term containing tanh(πy) in ∂ yŝα (y) has poles of order one at y k = i(k + 1/2) with k ∈ Z whose residue is α/(α − 1). The term containing tanh(απy) provides poles of order one with residue equal to 1/(1 − α) at y k ≡ i α k + 1 2 k ∈ Z (F.21) By introducing the following notation Ω(z) ≡ Γ(1 + z) Γ(−z) 2 (F. 22) for T ± ∞ in (F.17), which are obtained by evaluating (8.13) in ν = ν ± , we find (F.24) The expansions of these expressions for large η take the following form where R ± k are functions of η. By using (F.25) and introducing 26) the expressions occurring in the integrand of (F.19) can be expanded as We remark that our analysis is based on the approximate expression (8.10), obtained by considering only the summands corresponding to n ∈ {−1, 0, 1} in (8.8). Including also terms corresponding to |n| 2 leads to O(1/η 4 ) terms in the entanglement entropies. Hence, in order to be consistent with the approximation made in (8.10), we truncate (F.23) and (F.24) by considering only the terms up to O(1/η 3 ) included.
The above analysis tells us that (F.46) vanishes when α = 1 and that for finite α > 1 it can be written as follows , we find that D + 1 − D − 1 = 1 − 12 y 2 , which is an even function. Combining this observation with the fact that the functions in (F.20) are odd, we have that the integral in the first line of (F.51) vanishes because its integrand is an odd function. The integral in the second line of (F.51) can be evaluated by adapting the analysis made above to obtain (F.47) and (F.48). In particular, this integral gives a vanishing contribution when α = 1. Instead, for positive α = 1 we have where we used that P − 1 = − P + 1 = 1 − 12 y 2 and P − 1 |ỹ k = 1 + 3(2k + 1) 2 /α 2 , which have been obtained from (F.33), (F.17) and (F.21) and do not depend on the index j. (F.60) By using the fact thatỹ k is purely imaginary in (F.58), it is straightforward to observe that D − 2,j−1 |ỹ k is real; hence (F.60) is real. This is a real function because (F.65) evaluated along the imaginary axes is real.

F.3.2 Single copy entanglement
The large η asymptotic expansion of the single copy entanglement S (∞) A requires a separate discussion because the function (F.16) occurring in the integrand of (F.18) is not entire. This function has cusp singularities along the whole line Re(y) = 0; hence this case cannot be considered a special case of the above analysis, which employs the residue theorem.

G A double scaling limit of the lattice results
In this appendix we briefly mention some lattice results related to the quantities that have been studied in the main text through a particular double scaling limit.
The Hamiltonian of the free fermionic chain on the infinite line is [1] whereĉ i describe spinless fermionic degrees of freedom satisfying the anticommutation relation {ĉ i ,ĉ j } = δ i,j and h is the chemical potential. The ground state of (G.1) is a Fermi sea with Fermi momentum κ F = arccos(|h|) ∈ [0, π]. A Jordan Wigner transformation maps the Hamiltonian (G.1) into the Hamiltonian of the spin-1 2 Heisenberg XX chain in a magnetic field h. The two-point correlator of this lattice model reads Considering the bipartition of the infinite chain given by a block A made by L consecutive sites and its complement, many numerical analyses can be performed by considering the L×L reduced correlation matrix C A , whose generic element is (G.2) with i, j ∈ A, which turns out to be a Toeplitz matrix.
It is insightful to study the lattice results in the following double scaling limit [76,102,200] The full counting statistic generating function χ(ζ) allows us to study the cumulants of the particle number operator N A in the block A and its Toeplitz determinant representation reads [102,201] χ(ζ) ≡ e iζN A = det 1 + (e iζ − 1) where 1 is the L × L identity matrix and ζ = 2πν , with ν = ν (z) being defined in (8.2). The logarithm of (G.4) gives log χ(ζ) = log e iζN A = Tr log 1 + (e iζ − 1) C A (G.5) which generates the cumulants. The first cumulants (mean value, variance and skewness) read respectively N A = ∂ iζ log(χ) ζ=0 = Tr(C A ) (G.6) Since C A is a Toeplitz matrix, the Fisher-Hartwig conjecture [90] (proved in [91]) and its generalisation [80,91,166] can be employed to study the leading and the subleading terms respectively of the expansion of the determinant (G.4) as L → ∞. In particular, the double scaling limit (G.3) of the Fisher-Hartwig conjecture gives (8.7). Its generalised version include also some subleading corrections: for instance, the double scaling limit (G.3) of the expansion reported in Eq. (84) of [80] gives (8.13) with the infinite sums in k truncated to the finite sums including only the terms corresponding to k ∈ {0, 1, 2}. We emphasise that all order corrections occur in (8.1), which holds in the double scaling limit.
The cumulants generated by (G.5) have been studied by employing the Fisher-Hartwig conjecture and its generalisation e.g. in [102,168].
The entanglement entropies of the block A are obtained through the eigenvalues of the reduced correlation matrix C A [1,118]. Analytic results for these entanglement entropies for large L have been obtained by employing the Fisher-Hartwig conjecture for the leading terms [77,78] and the generalised Fisher-Hartwig conjecture combined with further computations for the subleading terms [80]. Taking the double scaling limit (G.3) in the result for the leading terms of S (α) A reported in [77] for large η, one obtains (8.16). Similarly, for the subleading terms of S (α) A , we have that the double scaling limit (G.3) of the lattice result in Eq. (10) of [80] gives (8.21) for N ∈ {0, 1, 2}. In order to obtain the N = 3 term in (8.21) through the double scaling limit (G.3), higher order terms must be computed in the lattice analysis, along the lines discussed in [80].
Beside the analysis involving the Fisher-Hartwig conjecture, in [77] lattice results have been reported that can be compared with the small η expansions discussed in this manuscript. As for the entanglement entropy, considering the double scaling limit (G.3) of the results of [77] for small η, we have that the first term in the last expression of (7.7) agrees with Eq. (5) of [77] for α = 1 up to missing factors of 2, whose absence there seems just a typo. As for the Rényi entropies, the expansion (7.13) agrees with Eq. (5) of [77] for α = 1. We find it worth mentioning also that the sine kernel tau function (4.20) evaluated at z = 1 provides the emptiness formation probability [202][203][204] of the XX chain and in the double scaling limit [205][206][207][208]. More recently, the approach of [77] has been extended to compute the EE of two disjoint intervals on the XX chain separated by a single site [181]. This problem has not exactly a well defined continuum limit, since in such a limit, the gap between the two intervals goes to zero giving a trivial subsystem. However, in the same lines, it would be very interesting to study the EE for two (or more) intervals separated by a finite distance in the continuum QFT. This problem has been addressed in [165] for the leading term. Subleading corrections to that result are still unknown.