Symmetry resolved entanglement in free fermionic systems

TL;DR

This paper analyzes the asymptotic behavior of symmetry-resolved entanglement entropies in one-dimensional free fermionic systems, providing exact constants and corrections, and extends the results to free Fermi gases.

Contribution

It applies the Fisher-Hartwig conjecture to derive detailed asymptotics of symmetry-resolved entanglement, including non-universal constants and oscillatory terms, in free fermionic models.

Findings

Exact asymptotic formulas for symmetry-resolved entanglement entropies.

Identification of oscillatory subleading contributions.

Extension of results to one-dimensional free Fermi gases.

Abstract

We consider the symmetry resolved R\'enyi entropies in the one dimensional tight binding model, equivalent to the spin-1/2 XX chain in a magnetic field. We exploit the generalised Fisher-Hartwig conjecture to obtain the asymptotic behaviour of the entanglement entropies with a flux charge insertion at leading and subleading orders. The o(1) contributions are found to exhibit a rich structure of oscillatory behaviour. We then use these results to extract the symmetry resolved entanglement, determining exactly all the non-universal constants and logarithmic corrections to the scaling that are not accessible to the field theory approach. We also discuss how our results are generalised to a one-dimensional free fermi gas.

Equations

$$S_{\mathrm{vN}}=-\mathrm{tr}{\rho }_A\ln {\rho }_A.$$

$$S_n=\frac{1}{1-n}\ln \mathrm{tr}{\rho }_A^n.$$

$$S_{\mathrm{vN}}=\frac{c}{3}\ln \ell ,S_n=\frac{c}{6}\frac{(n+1)}{n}\ln \ell ,$$

$${\rho }_A={\oplus }_q{\Pi }_q{\rho }_A={\oplus }_q\left[p(q){\rho }_A(q)\right],$$

$$S_{\mathrm{vN}}=\sum _{q}p(q)S_{\mathrm{vN}}(q)-\sum _{q}p(q)\ln p(q)\equiv S^c+S^f,$$

$$S_{\mathrm{vN}}(q)\equiv -\mathrm{tr}\left[{\rho }_A(q)\ln {\rho }_A(q)\right].$$

$$S_n(q)\equiv \frac{1}{1-n}\ln \mathrm{tr}{\left[{\rho }_A(q)\right]}^n.$$

$${\mathcal{Z}}_n(q)\equiv \mathrm{tr}\left({\Pi }_q{\rho }_A^n\right),$$

$$S_n(q)=\frac{1}{1-n}\ln \left[\frac{{\mathcal{Z}}_n(q)}{{\mathcal{Z}}_1(q{)}^n}\right]S_{\mathrm{vN}}(q)=-{\partial }_n{\left[\frac{{\mathcal{Z}}_n(q)}{{\mathcal{Z}}_1(q{)}^n}\right]}_{n=1}.$$

$$p(q)={\mathcal{Z}}_1(q).$$

$${\mathcal{Z}}_n(q)={\int }_{-\pi }^{\pi }\frac{d\alpha }{2\pi }{e}^{-iq\alpha }{Z}_n(\alpha ),Z_n(\alpha )\equiv \mathrm{tr}\left({\rho }_A^n{e}^{i{Q}_A\alpha }\right),$$

$${\mathcal{T}}_{n,\alpha }(x,\tau ){\varphi }_i(x^{\prime },\tau )=\{\begin{array}{l}{\varphi }_{i+1}(x^{\prime },\tau )e^{i\alpha {\delta }_{i,j}}{\mathcal{T}}_{n,\alpha }(x,\tau )(x<x^{\prime }),\\ {\varphi }_i(x^{\prime },\tau ){\mathcal{T}}_{n,\alpha }(x,\tau )\text{otherwise}.\end{array}$$

$$Z_n(\alpha )=⟨{\mathcal{T}}_{n,\alpha }(\ell ,0){\tilde{\mathcal{T}}}_{n,\alpha }(0,0)⟩.$$

$$h_{n,\alpha }=h_n+\frac{h_{\alpha }}{n},h_n=\frac{c}{24}\left(n-\frac{1}{n}\right),$$

$$Z_n(\alpha )=c_{n,\alpha }{\ell }^{-\frac{c}{6}\left(n-\frac{1}{n}\right)-2\frac{h_{\alpha }+{\bar{h}}_{\alpha }}{n}},$$

$$h_{\alpha }={\bar{h}}_{\alpha }=\frac{1}{2}{\left(\frac{\alpha }{2\pi }\right)}^{2}K.$$

$${\mathcal{Z}}_n(q)\simeq {\ell }^{-\frac{c}{6}\left(n-\frac{1}{n}\right)}\sqrt{\frac{n\pi }{2K\ln \ell }}{e}^{\frac{n{\pi }^2(q-⟨Q_A⟩{)}^2}{2K\ln \ell }}.$$

$$S_n(q)=S_n-\frac{1}{2}\ln \left(\frac{2K}{\pi }\ln \ell \right)+O({\ell }^0),S_{\mathrm{vN}}(q)=S_{\mathrm{vN}}-\frac{1}{2}\ln \left(\frac{2K}{\pi }\ln \ell \right)+O({\ell }^0).$$

$$H=-\sum _{i=-\infty }^{\infty }\left[c_i^{†}{c}_{i+1}+c_{i+1}^{†}{c}_i-2h\left(c_i^{†}{c}_i-\frac{1}{2}\right)\right],$$

$$\rho_{A}=\det C_{A}\exp\Big{(}\sum_{i,j}\left[\ln(C_{A}^{-1}-1)\right]_{ij}c_{i}^{\dagger}c_{j}\Big{)},$$

$$(C_A{)}_{i,j}=\frac{\sin k_F(i-j)}{\pi (i-j)},i,j\in A.$$

$$\mathrm{tr}{\rho }_A^n=\prod _{i=1}^{\ell }\left[{\left(\frac{1+{\nu }_i}{2}\right)}^n+{\left(\frac{1-{\nu }_i}{2}\right)}^n\right],$$

$$S_n=\sum _{i=1}^{\ell }{e}_n({\nu }_i),e_n(x)\equiv \frac{1}{1-n}\ln \left[{\left(\frac{1+x}{2}\right)}^n+{\left(\frac{1-x}{2}\right)}^n\right].$$

$$Z_n(\alpha )=\prod _{i=1}^{\ell }\left[{\left(\frac{1+{\nu }_i}{2}\right)}^n{e}^{i\alpha }+{\left(\frac{1-{\nu }_i}{2}\right)}^n\right].$$

$$Z_1(\alpha )\equiv \mathrm{tr}{\rho }_A{e}^{i{Q}_A\alpha }=\prod _{i=1}^{\ell }\left[\left(\frac{1+{\nu }_i}{2}\right)e^{i\alpha }+\left(\frac{1-{\nu }_i}{2}\right)\right],$$

$$Z_1(\alpha )=e^{i\alpha \ell /2}\prod _i\left[\left(\frac{1+{\nu }_i}{2}\right)e^{i\alpha /2}+\left(\frac{1-{\nu }_i}{2}\right)e^{-i\alpha /2}\right]=e^{i\alpha \ell /2}g(\alpha ),$$

$$Z_n(\alpha )=e^{i\alpha \ell /2}\prod _i\left[{\left(\frac{1+{\nu }_i}{2}\right)}^n{e}^{i\alpha /2}+{\left(\frac{1-{\nu }_i}{2}\right)}^n{e}^{-i\alpha /2}\right]=e^{i\alpha \ell /2}{g}_n(\alpha ),$$

$${\mathcal{Z}}_n(q)=\sum _{{\mathcal{S}}_q}\prod _{i\in {\mathcal{S}}_q}{\left(\frac{1+{\nu }_i}{2}\right)}^n\prod _{j\in \overline{{\mathcal{S}}_q}}{\left(\frac{1-{\nu }_j}{2}\right)}^n,$$

$$\ln Z_n(\alpha )=\sum _{i=1}^{\ell }{f}_n({\nu }_i,\alpha ),f_n(x,\alpha )=\ln \left[{\left(\frac{1+x}{2}\right)}^n{e}^{i\alpha }+{\left(\frac{1-x}{2}\right)}^n\right],$$

$$\ln Z_n(\alpha )=\frac{1}{2\pi i}\oint d\lambda f_n(\lambda ,\alpha )\frac{d\ln D_{\ell }(\lambda )}{d\lambda },$$

Full text

1* International School for Advanced Studies (SISSA) and INFN, Via Bonomea 265, 34136 Trieste, Italy

2International Centre for Theoretical Physics (ICTP), Strada Costiera 11, 34151 Trieste, Italy *

1. 1 Introduction

2. 2 Symmetry resolution and flux insertion

3. 2.1 Replica method and results from QFT

4. 3 Free fermions on a lattice: flux insertion and charged entropies

5. 3.1 Charged entropies via the generalised Fisher-Hartwig conjecture

6. 3.1.1 Leading term $(m=0)$

7. 3.1.2 Leading corrections $(m=\pm 1)$

8. 4 Free fermions on a lattice: symmetry resolved entropies

9. 4.1 $Q_{A}$-resolved moments via Fourier trasform

10. 4.2 Symmetry resolved Rényi and entanglement entropy

11. 5 Charged and symmetry resolved entanglement for the Fermi gas

12. 6 Conclusions

13. A Appendix A: Details of calculations

1 Introduction

Entanglement measures turned out to be fundamental tools for a finer description of the quantum world. Nowadays a lot has been understood in the context of low-dimensional many body quantum systems, both from the theoretical side [1, 2, 3, 4] and, more recently, also experimentally [5, 6, 7, 8]. For bipartite pure states described by a density matrix $\rho$, the most useful measure is surely the entanglement entropy, defined as the von Neumann entropy of the reduced density matrix (RDM) $\rho_{A}={\rm tr}_{B}\rho$ associated to one of the two subsystems (denoted by $A$ and $B$ respectively), i.e.,

[TABLE]

$S_{\rm vN}$ is the limit for $n\to 1$ of a larger family of entropies, known as Rényi entropies

[TABLE]

Within quantum field theory (QFT) both Eqs. (1.1) and (1.2) are usually obtained from the integer moments $Z_{n}={\rm tr}\rho_{A}^{n}$ of $\rho_{A}$, which are in turn easily expressed in the path integral formalism as partition function on suitable $n$-sheeted Riemann surfaces [9, 10]. This approach, when applied to critical systems, whose low energy physics is described by $(1+1)$ dimensional conformal field theory (CFT), leads to the famous scaling results [11, 9, 10, 12]

[TABLE]

for a subsystem $A$ made of a single interval of length $\ell$ embedded in an infinite one-dimensional system (and similarly for finite systems, systems at finite temperature, and other situations: it is sufficient to replace $\ell$ with the relevant length in the considered regime, see e.g. [10]).

However, a recent experiment, in the context of disordered systems [8], showed that it is also important to understand the “internal symmetry structure” of the entanglement as well. In particular, looking at systems possessing an internal global symmetry, entanglement turned out to have two different contributions, dubbed configurational and number or fluctuation entanglement [8]. These two contributions account for the entanglement within symmetry sectors and fluctuations thereof, respectively (see below for a precise definition).

At the same time, a new theoretical framework has been developed to address the problem of extracting the symmetry resolved contributions for different entanglement measures [13, 15, 16, 17]. Indeed, these contributions have been related to the moments of $\rho_{A}$ where twisted boundary conditions are imposed along the cuts of the Riemann surface: we will refer to them as charged moments. As we are going to see more in details, the twist can be implemented geometrically within field theory via threading an appropriate Aharonov-Bohm flux through the multisheeted Riemann surface [13]. The relation between twisted boundary conditions and flux insertion was actually previously explored, for example, in the context of free field theories [18, 20, 19]. Moreover, similar quantities have also been introduced in the holographic setting [21, 22] and in the study of entanglement in mixed states [15, 24, 25].

If on one hand the field theory approach is very powerful and versatile in order to provide the scaling limit of both the charged and symmetry resolved entanglement entropies, on the other, it does not give access to non-universal model-dependent pieces which are also very important to accurately characterise critical systems. As we are going to see, for the special case of free fermions, we can go beyond the field theory results. Indeed, we can rely on the (generalised) Fisher-Hartwig conjecture [26, 27, 28] to compute systematic expansions of these entropies which reproduce the field theory results and provide exact expressions for the non-universal terms. In fact, this method, which has already been explored for the standard entanglement and Rényi entropy [27, 28], can be simply generalised to the same quantities with a further flux insertion and therefore to their symmetry resolved analogue.

The paper is organised as follows. In Section 2 we carefully define all the quantities we will be dealing with and give an overview of the field theory results. Sections 3 and 4 are the core of this paper where we derive results for free fermions on a lattice for the charged and the symmetry resolved entanglement entropies, respectively. In Section 5 we show how all results derived for the lattice model can be directly adapted to a free Fermi gas. We conclude in Section 6 with some remarks and discussions. Some details of the calculations can be found in the appendix.

2 Symmetry resolution and flux insertion

Let us consider a many body quantum system with an internal $U(1)$ symmetry. Let $\rho$ be the density matrix in a given (pure) state and $Q$ the operator generating such symmetry. If the system is in a given representation of the charge $Q$, i.e., in an eigenstate of $Q$ corresponding to a definite eigenvalue, then $[\rho,Q]=0$.

We will be interested in a bipartition of the total system into two complementary spatial subsystems $A$ and $B$, with $\rho_{A}={\rm tr}_{B}\rho$ being the reduced density matrix of the subsystem $A$. Usually the operator $Q$ splits in the sum $Q=Q_{A}+Q_{B}$, meaning that $Q$ comes from local degrees of freedom within the two subsystems. Consequently, by taking the trace over $B$ of $[\rho,Q]=0$, we find that $[\rho_{A},Q_{A}]=0$. This implies that $\rho_{A}$ acquires a block diagonal form, in which each block corresponds to a different charge sector with a definite eigenvalue $q$ of $Q_{A}$, i.e.,

[TABLE]

where $\Pi_{q}$ is the projector on eigenspace of fixed value of $q$ in the spectrum of $Q_{A}$. In the last equality we factorised $p(q)={\rm tr}(\Pi_{q}\rho_{A})$, the probability of finding $q$ as the outcome of a measurement of $Q_{A}$. Note that in this way the density matrices $\rho_{A}(q)$ of different blocks are normalised as ${\rm tr}\rho_{A}(q)=1$.

Our goal is to understand how the entanglement is distributed in the different charge sectors. Focusing on the von Neumann entanglement entropy as a prototypical example, Eq. (2.1) implies the following decomposition

[TABLE]

where we defined the symmetry resolved entanglement entropy as the one associated to $\rho_{A}(q)$ in (2.1)

[TABLE]

The two different contributions in (2.2) are the configurational entanglement entropy, $S^{c}\equiv\sum_{q}p(q)S_{\rm vN}(q)$ [23, 8], measuring the total entropy due to each charge sector (weighted with their probability) and the fluctuation entanglement entropy $S^{f}=-\sum_{q}p(q)\ln p(q)$ [8], which instead takes into account the entropy due to the fluctuations of the value of the charge within the subsystem $A$.

Similarly, one defines also symmetry resolved Rényi entropies as

[TABLE]

In general evaluating such symmetry resolved quantities is a highly non-trivial problem, mainly due to the non local nature of the projector $\Pi_{q}$. As mentioned in the introduction, recently, this problem has been understood from a different perspective in [13, 16]. This new approach works as follows. Let us first define the (unnormalised) quantity

[TABLE]

which is related to the entanglement and Rényi entropies in (2.3) and (2.4) (respectively) through

[TABLE]

Also the probability $p(q)$ is read off ${\cal Z}_{n}$ as

[TABLE]

The key observation of Refs. [13, 16] is that (2.5) is given by the following Fourier transform

[TABLE]

where $Z_{n}(\alpha)$ are the charged moments mentioned in the introduction. Note that $Z_{n}(0)={\rm tr}\rho_{A}^{n}$. Therefore, we can access the symmetry resolved entanglement entropy by studying $Z_{n}(\alpha)$ (which, as explained below, are much easier to compute) and after Fourier transforming.

2.1 Replica method and results from QFT

In Ref. [13] a geometric approach in the framework of the replica trick has been introduced and it is applicable to generic (1+1)-dimensional QFT. The main idea is to insert an appropriate conjugate Aharonov-Bohm flux through a multi-sheeted Riemann surface $\mathcal{R}_{n}$, such that the total phase accumulated by the field upon going through the entire surface is $\alpha$. The result is that $Z_{n}(\alpha)$ is the partition function on such modified surface.

In QFT language, the insertion of the flux corresponds to a twisted boundary condition, which, as usually done in this context, can be implemented by the action of a local operator, acting at the boundary of the subsystem $A$. This operator is a modified twist field $\mathcal{T}_{n,\alpha}$ whose action, in operator formalism, is defined by [13]

[TABLE]

In this way one can further reformulate the problem in terms of a correlation function of twist fields [14]. In the simplest case of the subsystem consisting of a single interval $A=[0,\ell]$

[TABLE]

where $\tilde{\mathcal{T}}$ is the antitwist field. If we now specialise to (1+1) dimensional CFT, $\mathcal{T}_{n,\alpha}$ and $\tilde{\mathcal{T}}_{n,\alpha}$ behave as primary operators with conformal dimension given by [13]

[TABLE]

meaning that the phase shift is implemented by a composite twist field that can be written as $\mathcal{T}_{n,\alpha}=\mathcal{T}_{n}\cdot\mathcal{V}_{\alpha}$. This immediately implies

[TABLE]

where $c$ is the central charge of the CFT and $c_{n,\alpha}$ the unknown non-universal normalisation of the composite twist-field.

The focus of Ref. [13] was a free boson compactified on a circle of radius $R$, i.e., a Luttinger liquid with Luttinger parameter $K$. In this case the operator $\mathcal{V}_{\alpha}$ implementing the twisted boundary conditions is a vertex operator with (holomorphic and antiholomorphic) scaling dimensions

[TABLE]

From Eq. (2.12), the symmetry resolved moments are found by taking the Fourier transform as in Eq. (2.8). At leading order for large $\ell$, this reads [13]

[TABLE]

Notice that we set a posteriori the average number of the charge in the subsystem $\langle Q_{A}\rangle$, since it is a non-universal quantity, not encoded in the CFT. For a given microscopical model, its origin can be easily traced back, e.g. as a phase shift in the bosonisation rule.

Through Eq. (2.6), this leads to the following result at leading order for the Rényi and the von Neumann entropy

[TABLE]

This result has been dubbed equipartition of entanglement [16]: at leading order the entanglement is the same in the different charge sectors, just the probability $p(q)$ of being in a given sector varies.

3 Free fermions on a lattice: flux insertion and charged entropies

Eq. (2.15) provides the leading symmetry resolved entanglement entropies of all microscopical models with a $U(1)$ symmetry, that, at low energy, are described by a CFT. Indeed the results in Eq. (2.15) have been tested numerically both for free fermions [13, 16] and in interacting spin chains [16, 29]. In this Section we are going to provide an analytic derivation for the special case of a chain of free fermions, whose scaling limit is indeed described by a free compact boson with $K=1$. Our analysis will also provide the exact value of the non-universal constants, as well as the corrections to (2.15) for this specific model.

We consider the tight binding model in one dimension with hamiltonian

[TABLE]

where $c_{i}$ are free fermionic spinless degrees of freedom, satisfying the anticommutation relations $\{c_{i},c_{j}^{\dagger}\}=\delta_{ij}$ and $h$ is the chemical potential. $H$ is diagonal in momentum space and its ground state is a Fermi sea with Fermi momentum $k_{F}=\arccos|h|$. As it is clear from (3.1), the particle number $Q=\sum_{i}c_{i}^{\dagger}c_{i}$ is a conserved $U(1)$ charge of the model. It is also local and $Q=Q_{A}+Q_{B}$ for any spatial bipartition of the chain. By Jordan Wigner transformation, Eq. (3.1) is mapped to the XX spin chain in a magnetic field $h$ and the charge $Q$ becomes the spin in the $z$ direction.

In the following we will be interested in the bipartition where the subsystem $A$ is given by $\ell$ contiguous lattice sites. The corresponding RDM of the subsystem can then be written as [30, 31, 32]

[TABLE]

where the $\ell\times\ell$ matrix $(C_{A})_{ij}\equiv\langle c_{i}^{\dagger}c_{j}\rangle$ is the correlation matrix restricted to the subsystem $A$, that for the ground-state of an infinite chain has elements

[TABLE]

Note that $C_{A}$ is a Toeplitz matrix, meaning that its entries only depend on the difference $(i-j)$: this is a key point for what follows.

If one write the eigenvalues of the matrix $C_{A}$ as $(1+\nu_{k})/2$ (with $k\in[1,\ell]$), then simple algebra leads to the moments of $\rho_{A}$ as

[TABLE]

and, equivalently, the Rényi entropies read

[TABLE]

It has been first noticed in Ref. [13] that the $\alpha$-dependent moments $Z_{n}(\alpha)$, defined in (2.8), can be also easily written in terms of the eigenvalues of the correlation matrix with a simple modification of the above formulas, i.e.,

[TABLE]

The interpretation of this equation is straightforward: each particle carries a weight $e^{i\alpha}$ while the holes carry weight $1$. Eq. (3.6) provides a very simple method for the numerical computation of $Z_{n}(\alpha)$. Not only: as we are going to discuss next, it is also the right starting point to get the asymptotic analytic expressions of $Z_{n}(\alpha)$.

Before embarking into the study of $Z_{n}(\alpha)$ a quick recap of its properties and limits is necessary, also to provide useful consistency checks for our calculations. First, for $n=1$

[TABLE]

is the moment-generating function of $Q_{A}$. This quantity has been already studied in the literature [33, 34, 35, 36, 37] also because of its relation with the entanglement entropy itself. We will see that it is simply related to another quantity usually introduced in this context ($D_{\ell}(\lambda)$ of the next subsection). The first moment is just the average number of particle in $A$, i.e., $\langle Q_{A}\rangle=\ell k_{F}/\pi$. Hence $Z_{1}(\alpha)=1+i\alpha\langle Q_{A}\rangle+O(\alpha^{2})$. At half-filling $Z_{1}(\alpha)$ further simplifies as a consequence of the fact that, by particle-hole symmetry, for each $\nu_{i}$ there is a $\nu_{j}$ such that $(1-\nu_{i})=(1+\nu_{j})$. Thus we have

[TABLE]

with $g(\alpha)$ real and even. For general filling instead the odd cumulants are non-vanishing (the odd derivatives of $\ln Z_{1}(\alpha)$ are non zero) and $Z_{1}(\alpha)$ has no particular parity or reality properties. Indeed, Eq. (3.8) remains true for generic $n$

[TABLE]

with $g_{n}(\alpha)$ real and even. This symmetry of $Z_{n}(\alpha)$ at half-filling represents a cross check of our numerical and analytic calculations. Again, away from half-filling, $\ln Z_{n}(\alpha)$ has all non-zero derivatives.

We notice that by Fourier transforming (3.6) one easily gets

[TABLE]

where the sum is over all subset $\mathcal{S}_{q}$ of $\mathcal{S}={1,\cdots,\ell}$ of cardinality $q$ and $\bar{\mathcal{S}_{q}}$ denotes the complementary subset. Unfortunately Eq. (3.10) is not a very convenient way to get $\mathcal{Z}_{n}(q)$, since one has to sum over $\ell!/((\ell-q)!q!)$ terms and this is impossible already for moderate values of $\ell$. The most convenient way to extract $\mathcal{Z}_{n}(q)$ is by direct Fourier transform of $Z_{n}(\alpha)$.

3.1 Charged entropies via the generalised Fisher-Hartwig conjecture

The method that we employ takes advantage of the Toeplitz nature of the correlation matrix that can be handled with the (generalised) Fisher-Hartwig conjecture providing the asymptotics of determinant of Toeplitz matrices. This technique has been used already in the context of entanglement in free lattice models to derive the leading term and the corrections to entanglement entropies [27, 28, 38, 39, 40, 41, 42, 43]. We are going to show that the same technology applies also to the $\alpha$-dependent moments $Z_{n}(\alpha)$ and therefore, as a consequence of the discussion above, to their symmetry resolved equivalent $\mathcal{Z}_{n}(q)$.

The starting point of our analysis is to rewrite the logarithm of Eq. (3.6)

[TABLE]

as a contour integral

[TABLE]

where the contour of integration encircles the segment $[-1,1]$. Here we defined the characteristic polynomial of $C_{A}$ as the determinant

[TABLE]

where $\mathbb{I}_{A}$ is the identity matrix restricted to $A$. In the basis that diagonalises $C_{A}$, such determinant simply becomes $D_{\ell}=\prod_{i}(\lambda-\nu_{i})$ and therefore, by residue theorem, Eq. (3.12) is the same as (3.11). Notice that $D_{\ell}(\lambda)$ is related to the generating function $Z_{1}(\alpha)$ as

[TABLE]

In Refs. [27, 28] it has been exploited the fact that the matrix $(\lambda+1)\mathbb{I}_{A}-2C_{A}$ has a Toeplitz form. Therefore the asymptotics for large $\ell$ of the determinant $D_{\ell}$ in (3.13) is obtained by means of the generalised Fisher-Hartwig conjecture. The interested reader can find the derivation in Ref. [28], we just report here the final result which is [28]

[TABLE]

where $G(\cdot)$ is the Barnes $G$-function, $L_{k}=2\ell|\sin k_{F}|$ and

[TABLE]

For the moments ${\rm tr}\rho_{A}^{n}$, i.e., $Z_{n}(\alpha=0)$ in (3.11), the leading term in the sum for $D_{\ell}$ in Eq. (3.15) is the one with $m=0$, first evaluated in [27]. The next to leading contributions come from the terms with $m=\pm 1$ (at the same order) as shown in [28]. The situation for $\alpha\neq 0$ is slightly more complicated. For $-\pi<\alpha<\pi$ the leading term is always the one with $m=0$. Since $Z_{n}(\alpha)$ is periodic in $\alpha$ with period $2\pi$ we will restrict ourselves to $\alpha\in[-\pi,\pi]$, having in mind that, if required, the function can be extended outside of this interval by periodicity. Concerning the subleading contributions, the terms with $m=\pm 1$ have different power laws, but one of them is always dominating, as we shall see. Anyhow, for values of $\alpha$ close to $\pm\pi$, also next-to-next leading terms should be taken into account in order to get reasonable results for moderately large values of $\ell$. In the following we first compute the leading term and then we move to the calculations of the corrections.

3.1.1 Leading term $(m=0)$

For $\alpha\in[-\pi,\pi]$, the leading behaviour of Eq. (3.12) is given by term with $m=0$ in (3.15), i.e.,

[TABLE]

so that the integral (3.12) is $\ln Z_{n}(\alpha)=\ln Z_{n}^{(0)}(\alpha)+o(\ell^{0})$ with

[TABLE]

where

[TABLE]

are respectively the linear, the logarithmic and the constant term (in $\ell$) coming from $\ln D^{(0)}_{\ell}$ in Eq. (3.17). These three integrals are explicitly calculated in Appendix A with final result

[TABLE]

where

[TABLE]

in analogy with the definition $\Upsilon(n)$ in [27], which is recovered when $\alpha=0$. We stress that $\Upsilon{(n,\alpha)}$ is real for $\alpha$ real, even if not apparent from the formula. For future reference it is useful to write $\Upsilon{(n,\alpha)}$ as

[TABLE]

Eq. (3.18) contains several pieces of information. The linear term is just the mean number of particles in $A$, $\langle Q_{A}\rangle=k_{F}\ell/\pi$, as expected. Anyhow, this is the only term with an imaginary part up to order $O(1)$. We know that this is exactly true at half-filling ($k_{F}=\pi/2$), cf. Eq. (3.9). For generic filling, it is not true in general and we will observe in numerics tiny deviations at small $\ell$ in the imaginary part of $\ln Z_{n}(\alpha)$. The term $\propto\ln L_{k}$ provides the dimension of the modified twist field which comes out from the field theory calculation: the result agrees with the one found by CFT methods in (2.12) when specialised to a compact boson with $K=1$. It was important to test analytically this result that was already checked numerically in [13]. The constant term in Eq. (3.22) is probably the most interesting one, first because it is a result that was not known by other means (being non-universal cannot be fixed by field theory), and second because it provides few physical consequences. It is real and even in $\alpha$, a property that was guaranteed only at half filling. It is independent from $k_{F}$, as its limit for $\alpha=0$ [27]. Finally, it is very close to a parabola, but all the even terms in the series expansion are non zero, although $\epsilon(n,\alpha)$, cf. Eq. (3.24), is very small. In Figure 1 we report $\Upsilon{(n,\alpha)}$ as function of $\alpha$ for some $n$ and compare it with the quadratic approximation $\Upsilon{(n)}+\gamma_{2}(n)\alpha^{2}$. The closeness of the two curves shows that the quadratic approximation will be enough for most of the applications, as we shall explicitly show. The accuracy of the quadratic approximation is also evident from the plot of $\epsilon(n,\alpha)$ in the right panel of Figure 1. On passing, this precision of the quadratic approximation of $\Upsilon{(n,\alpha)}$ explains, a posteriori, the quality of the symmetry resolved spectrum obtained in Ref. [13] exploiting the method of Stieltjes transform [44] which implicitly assumes this approximation.

In Figure 2 we report the numerical data for Rényi entropies with the insertion of a flux $\alpha$ for several values of $n$ and $\alpha$ and with fillings $k_{F}=\pi/2$ (left) and $k_{F}=\pi/3$ (right). The theoretical prediction for the leading scaling in Eq. (3.22) is also reported for comparison. It is evident that the analytical result correctly describes the asymptotic data, but large and oscillating corrections to the scaling are present. The amplitude of these oscillations increase with $n$ and with $\alpha$. This peculiar $n$ dependence was already known at $\alpha=0$ [28, 45, 46, 48, 47]. In the following subsection we will explicitly consider these oscillations and work out their analytical description.

3.1.2 Leading corrections $(m=\pm 1)$

The leading correction to the determinant $D_{\ell}(\lambda)$ comes from the terms with $m=\pm 1$ in (3.15) and is given by [28]

[TABLE]

We define the difference

[TABLE]

that for large $L_{\kappa}$ is

[TABLE]

Changing variable to $\lambda=\tanh(\pi w)$, in the final integral we only need the discontinuities across the branch cut that for the two cases are

[TABLE]

where we have dropped terms of order $O(L_{k}^{-4})$ compared to the leading ones and we have defined

[TABLE]

Integrating by parts and using (A.4) we finally get

[TABLE]

This integral can be evaluated on the complex plane by residue theorem. For the first piece of the integral in square bracket, we should close the contour in the upper half plane, while for the second piece in the lower half plane. In principle we should sum over all residues inside the integration contour, but if we are interested in the limit of large $L_{k}$, we can limit ourself to consider the singularities closest to the real axis. For the first integral this is at $w=i/(2n)(1-\alpha/\pi)$ while for the second one it is at $w=-i/(2n)(1+\alpha/\pi)$. Summing up the two contributions we finally have

[TABLE]

Let us comment this result. First, it is obvious that the formula is valid only for $-\pi<\alpha<\pi$, else one of the power laws would blow up as a consequence of the fact that one of the terms for $m=\pm 1$ becomes the leading ones. Then we see that this correction is real only for $\alpha=0$ and at half-filling (when it should be real at all orders, cf. Eq. (3.9)). Away from half filling, there is generically a non-zero imaginary part. The two contributions have a different power-law decays (for $\alpha\neq 0$) and so only one of them is the leading correction depending on the sign of $\alpha$. However, for $\alpha$ close to zero, the two powers are too close in magnitude and they should be both taken into account in order to have an accurate description of the data for moderately large values of $\ell$. When $\alpha$ gets closer and closer to $\pm\pi$, Eq. (3.30) becomes accurate only for very large $\ell$ because the term with $m=0$ is about of the same order of magnitude as the one with $m=\pm 1$ (depending on the sign of $\alpha$). A better description of the asymptotic behaviour may be achieved using Eq. (3.15) without expanding as in Eq. (3.25). Finally, let us notice that, while in the absence of flux ($\alpha=0$) the oscillating corrections to the scaling vanish in the limit $n\to 1$ [28], for $\alpha\neq 0$ also the von Neumann entropy presents leading oscillating corrections described by Eq. (3.30).

In Figures 3 and 4 we report the difference $d_{2}(\alpha)$ as calculated numerically for $k_{F}=\pi/2,\pi/3,\pi/4$ and for $\alpha=1$ and $\alpha=2$ as function of $\ell$. The numerical data are compared with the leading prediction in Eq. (3.30) and the agreement is extremely good. We actually observe that this prediction works slightly worst for $\alpha=2$ (cf. Fig. 4) than for $\alpha=1$ (cf. Fig. 3). We indeed checked that the match becomes worst and worst when $\alpha$ moves close to $\pi$, when the leading term in the generalised Fisher-Hartwig changes. In principle it is possible to systematically analyse further corrections to $\ln Z_{n}(\alpha)$ by taking into account the known expansion of $D_{\ell}(\lambda)$ in powers of $\ell$ [28], but this is very cumbersome and far beyond the scope of this paper.

4 Free fermions on a lattice: symmetry resolved entropies

In this section we finally move to the symmetry resolved entropies and to their analysis.

4.1 $Q_{A}$-resolved moments via Fourier trasform

The first step toward the symmetry resolved entropies is to calculate $\mathcal{Z}_{n}(q)$, the Fourier transform of $Z_{n}(\alpha)$ as defined in Eq. (2.8). We will show that we may obtain a very accurate prediction by keeping only the $m=0$ term in (3.15), but with all non-universal pieces. Within this approximation the Fourier transform $\mathcal{Z}_{n}(q)$ is

[TABLE]

where we defined the “bare variance”

[TABLE]

As a first step we use Eq. (3.24) to rewrite $\mathcal{Z}_{n}(q)$ as

[TABLE]

where $Z_{n}^{(0)}(0)=Z_{n}^{(0)}(\alpha=0)=e^{\Upsilon(n)}L_{k}^{-\frac{1}{6}\left(n-\frac{1}{n}\right)}$, we defined the “renormalised variance”

[TABLE]

and $g_{n}(\alpha)\equiv e^{\epsilon(n,\alpha)}$. Up to this point, we only rewrote the starting expression (4.1). We now proceed by treating the integral, for large subsystem size $\ell$, by means of the saddle point approximation. When $\ell\gg 1$, the large parameter in (4.3) is $b_{n}$. Furthermore, we assume that $g_{n}(\alpha)=1$, because we have shown in the previous section that the function $\epsilon(n,\alpha)\ll 1$, cf. Fig. 1. Within this approximation we finally get

[TABLE]

where we defined $\bar{q}\equiv\langle Q_{A}\rangle=\ell k_{F}/\pi$.

Eq. (4.5) is one of the main results of this paper. Let us discuss its features. First, in the limit $\ell\to\infty$, we recover the CFT result (2.14) for $K=1$, but with the correct normalisation of $Z_{n}(0)$. Although this normalisation was not previously known rigorously (at least to the best of our knowledge), it could have been easily guessed from the results in the absence of flux (i.e., $\alpha=0$ of Ref. [27]). The mean of the gaussian term is $\bar{q}$ and it is not changed compared to the result (2.14). Consequently, the main new insight from Eq. (4.5) is the prediction for the constant term to add to $\ln L_{k}$ (or equivalently, the multiplicative scale for $L_{k}$ as in [16]). Although this non-universal constant is a correction to the leading behaviour (expanding for large $L_{k}$, it gives a term going like $(\ln L_{k})^{-1}$), it is very important: the $(\ln L_{k})^{-1}$ decay is so slow that it must be taken into account even for very large $L_{k}$ in order to quantitatively describe the data, as we shall see.

Given the importance that the quantity $\gamma_{2}(n)$ has in this analysis, we report its analytic expression

[TABLE]

as well as its explicit numerical value for some $n$: $\gamma_{2}(1)=-0.0799027$, $\gamma_{2}(2)=-0.0462208\dots$, and $\gamma_{2}(3)=-0.0319926\dots$ (in particular $\gamma_{2}(1)=-(1+\gamma_{E})/(2\pi^{2})$ with $\gamma_{E}$ the Euler constant, as anticipated in [16]). The importance of this constant in the description of the numerical data, was understood already in [16], where the authors define $g_{n}=2e^{-2\pi^{2}n\gamma_{2}(n)}$ and provide the analytic results for $n=1$, as well as a numerical estimate for $n=2$, i.e., $g_{2}\sim 12.39$ which is very close to the exact value that we have found $g_{2}=12.4022\dots$.

Let us briefly discuss the terms that have been neglected in the derivation of Eq. (4.5). The most relevant one comes from having approximated $g_{n}(\alpha)$ with $1$. By expanding this function in powers of $\alpha$, it is immediate to realise that the series coefficient $\alpha^{2k}$ (with $k\geq 2$) provides a correction to the leading term of the order $(\ln L_{k})^{-k}$. Anyhow, these factors influence little the final result because the amplitude of the various terms is very small. Another correction comes from the extremes of integration that we pushed to $\pm\infty$ instead of $\pm\pi$. Although their effect can be taken into account as done in Ref. [16], they provide corrections which decay as $e^{-\pi^{2}b_{n}^{R}}/b_{n}^{R}$, i.e., algebraically in $L_{k}$, and so negligible at this level. Also the corrections due to the terms with $m\neq 0$ in (3.15) decay as power laws in $L_{k}$ and can be safely neglected at this stage.

In Figure 5 we report the numerically calculated symmetry resolved partition sums ${\cal Z}_{n}(q)$. We compare the numerical data for $n=1,2,3$ with the CFT prediction without fixing the non universal constant as in Eq. (2.14). The qualitative agreement is reasonable, but quantitatively far. We also report the prediction for ${\cal Z}_{n}(q)$ at order $O(\ell^{0})$: the curves moves closer to the numerical data, but the match is still not perfect. Only when we use the complete Fisher-Hartwig prediction (4.5) (with the correct value of $\gamma_{2}(n)$), the data are perfectly reproduced. As anticipated, including the logarithmic corrections is fundamental to have an accurate description of the data. Also the $q$-dependence of ${\cal Z}_{n}(q)$ is perfectly captured by (4.5) as shown in the lower panels of Figure 5.

4.2 Symmetry resolved Rényi and entanglement entropy

We now use Eq. (4.5) to calculate the symmetry resolved Rènyi and the Von Neumann entropies. Let us start from the former. Eq. (2.6) implies

[TABLE]

The first ratio in (4.7) just gives the total Rényi entropy of order $n$, with the right additive constant (and indeed this is true at all orders). The other $q$-independent term is

[TABLE]

The constant $\delta_{n}$ has been introduced to resum partially the subleading corrections to the scaling and it is given by

[TABLE]

The last term is the ratio of the two Gaussian factors which is the only one depending on $q$. For this last contribution we have

[TABLE]

where the constant

[TABLE]

has been introduced, again, to resum partially the subleading corrections.

Putting together the three pieces we have

[TABLE]

This equation not only predicts the leading diverging behaviour for large $\ell$ which was already known from CFT [13, 16] (cf. Eq. (2.15)), but also the non-universal additive constant, as well as the some subleading corrections in $\ln L_{k}$. The latter are not only important to correctly describe the data, but are also the leading $q$-dependent contributions. So while the leading and finite terms satisfy the equipartition of entanglement [16], within our approach we are able to identify the leading term that breaks this equipartition.

Taking now the limit for $n\to 1$, we get the von Neumann entropy

[TABLE]

with $\ln\delta_{1}=2\pi^{2}\gamma^{\prime}_{2}(1)$ and $\gamma^{\prime}_{2}(1)=0.0545724$.

These Fisher-Hartwig calculations for the symmetry resolved entanglement are compared with the numerical data in Figure 6. It is evident in these figures that the results for different $q$ are not on top of each other although we reported $\ell$ as large as $3000$. Indeed their difference (that we know to go to zero as $(\ln\ell)^{-2}$) can be easily misinterpreted as a different additive constant if one would proceed with a fit of the numerical data. Only the exact knowledge of the asymptotic behaviour (4.12) and (4.13) allow us to correctly understand the data. In the figure we also report Eqs. (4.12) and (4.13) truncated at $o(1)$ (just for $q=\bar{q}$), showing that these leading curves are far from the data and that the distance between the two barely reduces. We stress that not only the prefactor of the logarithmic corrections are important, but also the precise values of the amplitudes (4.11) and (4.9), as it is easy to check. Finally we observe that increasing $(q-\bar{q})$ the corrections to the scaling that we neglected become more important.

We finish the section commenting about the double log contribution in Eqs. (4.12) and (4.13). It may seem rather awkward that all the symmetry resolved contributions have a double log correction, while the total entanglement entropy does not. Indeed, when calculating the total entanglement this double log cancels when summing to the fluctuation entanglement $S^{f}$ as in Eq. (2.2). Indeed, taking as a prototypical example the von Neumann entropy and using that the probability is $p(q)={\cal Z}_{1}(q)$, we have

[TABLE]

Note that both leading terms in $S^{f}$ in the above equation cancel exactly with the corresponding ones in the symmetry resolved entanglement in Eq. (4.13). The same is true for all Rényi entropies of arbitrary order.

5 Charged and symmetry resolved entanglement for the Fermi gas

In this section we derive the symmetry resolved entanglement entropy for a Fermi gas using the overlap matrix approach [49, 50]. This technique has been successfully applied to the calculation of entanglement in many different circumstances [58, 49, 50, 51, 52, 53, 54, 55, 56, 57].

The system we are going to study consists of a gas of $N$ free spinless non-relativistic fermions with some suitable boundary conditions in order to have a discrete energy spectrum. The many body wave functions $\Psi(x_{1},...,x_{N})$ is the Vandermonde determinant $\Psi(x_{1},...,x_{N})={\rm det}[\phi_{m}(x_{n})]/\sqrt{N!}$, built with the occupied single particle eigenstates with wave functions $\phi_{m}(x)$. The many body ground state is obtained by filling the $N$ levels with lowest energies. Given that there is no lattice, the particle number $N$ provides also the ultraviolet cutoff.

Our focus is again a subsystem, now taken to be a single interval $A$, and its entanglement with the rest of the system. The RDM $\rho_{A}$ in the subsystem $A$ is obtained as the continuous limit of Eq. (3.2). Therefore we can still exploit its gaussian nature. A crucial quantity is the overlap matrix associated to $A$ defined as

[TABLE]

In fact, it has been shown [50, 49] that the (continuous) correlation matrix and the (discrete) overlap matrix share the same eigenvalues $(1+\nu_{k})/2$ (with $k\in[1,N]$) and so, for example, the moments of the RDM may be written as

[TABLE]

In the case of a system with periodic boundary conditions in the interval $[0,L]$, the eigenfunctions are plane waves $\phi_{k}(x)=e^{2\pi ikx/L}/\sqrt{L}$ with integer wave-numbers $k$. When the subsystem is also an interval, say $A=[0,\ell]$, the overlap matrix is easily calculated and reads

[TABLE]

A crucial observation made in [50] is that such matrix is identical to the lattice correlation matrix, Eq. (3.3), upon identifying $k_{F}/\pi$ with $\ell/L$ and $N$ with $\ell$. As a consequence, this simple replacement allows to translate all the results from the lattice to the continuous model. In particular all formulas derived for the tight binding model are valid also for the gas, where now $L_{k}$ is not anymore $2\ell\sin k_{F}$, but

[TABLE]

This replacement allows to show rigorously the CFT scaling in finite size for this model.

For the symmetry resolved entanglement we can straightforwardly make our predictions for the gas. First of all, the entanglement Rényi entropies in the presence of a flux are

[TABLE]

where we only included the leading contributions at $m=0,\pm 1$ in the generalised Fisher-Hartwig conjecture. We tested this result against exact numerical computations and, as for the lattice model, we found that it provides a very accurate description as long as $\alpha$ is not close to $\pm\pi$.

Similarly, for the symmetry resolved entropies the predictions for the gas are obtained simply by plugging Eq. (5.4) into Eqs. (4.12) and (4.13) for Rényi and von Neumann entropy respectively (with $\bar{q}=\frac{\ell N}{L}$).

6 Conclusions

In this manuscript we derived exact formulas for the asymptotic behaviour of the symmetry resolved entanglement entropies in free fermion systems. First we obtained an exact expression for the charged entropies given by Eqs. (3.22) (asymptotic behaviour up to order $O(1)$) and (3.30) (leading corrections to the scaling). The leading logarithmic term in (3.22) perfectly match the CFT prediction, but we also determined the non-universal $O(1)$ contribution. The $o(1)$ corrections present interesting oscillatory behaviour and a power law decay with exponents that depend on the flux $\alpha$. We then moved to the truly symmetry resolved entropies given by the Fourier transform of the charged ones. The partition sums are given by Eq. (4.5), while Rényi anf von Neumann entropy by Eqs. (4.12) and (4.13) respectively. These equations agree in the limit of large $\ell$ with the CFT results, but we also determine a number of non-universal constants as well as logarithmic corrections to the scaling which are fundamental for an accurate description of the numerical data. Our analysis also provide the first term in the expansion for large $\ell$ which depend on the symmetry sector, hence breaking the equipartition of entanglement [16]. We also related the the double logarithmic correction to the fluctuation entanglement.

While we have considered the specific case of free fermions, many features we find are in fact universal. The CFT results [13, 16] (cf. Eq. (2.15)) shows how the leading term of the charged entropies get renormalised by the Luttinger liquid parameter $K$. A first natural question is how the exponent of the leading corrections to the scaling gets renormalised. It would be very interesting to adapt the field theoretical approach of Refs. [46, 47] to understand how this new universal exponent (equal to $2/n(1\pm\alpha/\pi)$ for free fermions) can be obtained in CFT. For the symmetry resolved entanglement we showed the presence of very large logarithmic corrections to the scaling. The natrual question here is whether they are universal and if they can be also understood within CFT. Furthermore, we find that many non-universal constants entering in these corrections are related to each other (e.g. Eqs. (4.9) and (4.11)): it is interesting to understand also the level of universality of these relations.

Finally there are few possible generalisations of the present calculations that can be done following the same logic as here; for example the case of an open system can be analysed exploiting the generalised Fisher-Hartwig results in [40], disjoint intervals using the approach in [42], and trapped gases can be studied by random matrix techniques [37, 57] to recover results from curved CFT [59].

Acknowledgments.

We thank Luca Capizzi, Sara Murciano and Hassan Shapourian for useful discussions and collaboration on closely related topics. All authors acknowledge support from ERC under Consolidator grant number 771536 (NEMO).

Appendix A Appendix A: Details of calculations

In this appendix we report the details of the calculation needed in Section 3.

The first integral is the linear term in $\ell$

[TABLE]

While the first piece is analytic in $\lambda$ (because $f_{n}(-1,\alpha)=0$), the second piece has a simple pole in $\lambda=1$, leading to

[TABLE]

The second integral is

[TABLE]

The only non zero contribution to the contour integral comes from the discontinuity at the cut of $\beta_{\lambda}$,

[TABLE]

and hence we finally have

[TABLE]

This integral may be considered as the final answer, but indeed it may be simply evaluated analytically. First one performs the change of variable

[TABLE]

to get

[TABLE]

At this point, integrate by part using

[TABLE]

to get

[TABLE]

Consider now the difference

[TABLE]

The derivative with respect to $\alpha$ of the last integral is

[TABLE]

Integrating back and using the boundary condition that the difference for $\alpha=0$ is zero by definition, we finally have

[TABLE]

that using the known value of $a_{1}|_{\alpha=0}$ [27] is equivalent to (3.22).

The last integral is

[TABLE]

Let us first integrate by parts

[TABLE]

Again the only discontinuity at the cut comes from the function $\beta_{\lambda}$, so that

[TABLE]

where we used $\Gamma(z)=G(z+1)/G(z)$. We then have

[TABLE]

Changing integration variable from $\lambda$ to $w$ as in (A.6) and using the derivative (A.8), we finally get

[TABLE]