Average eigenstate entanglement entropy of the XY chain in a transverse   field and its universality for translationally invariant quadratic fermionic   models

Lucas Hackl; Lev Vidmar; Marcos Rigol; Eugenio Bianchi

arXiv:1812.08757·cond-mat.stat-mech·February 14, 2019

Average eigenstate entanglement entropy of the XY chain in a transverse field and its universality for translationally invariant quadratic fermionic models

Lucas Hackl, Lev Vidmar, Marcos Rigol, Eugenio Bianchi

PDF

TL;DR

This paper demonstrates that the average eigenstate entanglement entropy in the XY chain exhibits universal volume-law scaling with a subleading correction depending on the subsystem ratio, extending previous results from the quantum Ising chain.

Contribution

It analytically and numerically extends the universality of eigenstate entanglement entropy scaling to the XY chain and refines bounds on the volume-law coefficient.

Findings

01

Universal volume-law entanglement entropy in XY chain.

02

Subleading correction depends on subsystem ratio.

03

Explicit bounds for the volume-law coefficient.

Abstract

We recently showed [Phys. Rev. Lett. 121, 220602 (2018)] that the average bipartite entanglement entropy of all energy eigenstates of the quantum Ising chain exhibits a universal (for translationally invariant quadratic fermionic models) leading term that scales linearly with the subsystem's volume, while in the thermodynamic limit the first subleading correction does not vanish at the critical field (it only depends on the ratio $f$ between the volume of the subsystem and volume of the system) and vanishes otherwise. Here we show, analytically for bounds and numerically for averages, that the same remains true for the spin-1/2 XY chain in a transverse magnetic field. We then tighten the bounds for the coefficient of the universal volume-law term, which is a concave function of $f$ . We develop a systematic approach to compute upper and lower bounds, and provide explicit analytic…

Tables5

Table 1. Table 1: Coefficients in t n subscript 𝑡 𝑛 t_{n} . List of the coefficients α n ( m ) superscript subscript 𝛼 𝑛 𝑚 \alpha_{n}^{(m)} in t n ( f ) = ∑ m = n 2 n − 1 α n ( m ) f m subscript 𝑡 𝑛 𝑓 subscript superscript 2 𝑛 1 𝑚 𝑛 superscript subscript 𝛼 𝑛 𝑚 superscript 𝑓 𝑚 t_{n}(f)=\sum^{2n-1}_{m=n}\alpha_{n}^{(m)}\,f^{m} , defined in Eq. ( 18 ), for n ≤ 4 𝑛 4 n\leq 4 .

$t_{n}$	$α_{n}^{(1)}$	$α_{n}^{(2)}$	$α_{n}^{(3)}$	$α_{n}^{(4)}$	$α_{n}^{(5)}$	$α_{n}^{(6)}$	$α_{n}^{(7)}$
$n = 1$	$2$
$n = 2$		$\frac{16}{3}$	$- 4$
$n = 3$			$22$	$- 48$	$32$
$n = 4$				$\frac{1816}{15}$	$- \frac{2592}{5}$	$\frac{54640}{63}$	$- 544$

Table 2. Table 2: ( 2 , 2 , 2 ) 2 2 2 (2,2,2) -contractions. List of all ( 2 , 2 , 2 ) 2 2 2 (2,2,2) -contractions ℭ ℭ \mathfrak{C} , their multiplicity # ℭ # ℭ \#\mathfrak{C} under translational symmetry, and their integrated prefactor ∫ ℭ ℭ \int\mathfrak{C} . Summing all terms, weighted with their multiplicity, gives ∫ 𝒞 ( 2 , 2 , 2 ) = 11 𝔠 1 3 f 3 = f 3 subscript 𝒞 2 2 2 11 superscript subscript 𝔠 1 3 superscript 𝑓 3 superscript 𝑓 3 \int\mathcal{C}_{(2,2,2)}=11\,\mathfrak{c}_{1}^{3}\,f^{3}=f^{3} .

$# ℭ$	$\int ℭ$	$# ℭ$	$\int ℭ$
$2$	$1$	$3$	$\frac{1}{2}$
$3$	$1$	$1$	$\frac{1}{2}$
$6$	$\frac{2}{3}$

Table 3. Table 3: ( 2 , 2 , 2 , 2 ) 2 2 2 2 (2,2,2,2) -contractions. List of all ( 2 , 2 , 2 , 2 ) 2 2 2 2 (2,2,2,2) -contractions ℭ ℭ \mathfrak{C} , their multiplicity # ℭ # ℭ \#\mathfrak{C} under translational symmetry, and their prefactor ∫ ℭ ℭ \int\mathfrak{C} under integration. Summing all terms, weighted with their multiplicity, gives ∫ 𝒞 ( 2 , 2 , 2 , 2 ) = 908 15 𝔠 1 4 f 4 = 908 15 f 4 subscript 𝒞 2 2 2 2 908 15 superscript subscript 𝔠 1 4 superscript 𝑓 4 908 15 superscript 𝑓 4 \int\mathcal{C}_{(2,2,2,2)}=\tfrac{908}{15}\,\mathfrak{c}_{1}^{4}\,f^{4}=\tfrac{908}{15}\,f^{4} .

$# ℭ$	$\int ℭ$	$# ℭ$	$\int ℭ$
$8$	$1$	$4$	$1$
$2$	$1$	$8$	$\frac{2}{3}$
$8$	$\frac{2}{3}$	$8$	$\frac{2}{3}$
$4$	$\frac{2}{3}$	$8$	$\frac{11}{30}$
$8$	$\frac{11}{30}$	$4$	$\frac{11}{30}$
$2$	$\frac{11}{30}$	$8$	$\frac{1}{2}$
$8$	$\frac{1}{2}$	$8$	$\frac{1}{2}$
$8$	$\frac{1}{2}$	$4$	$\frac{9}{20}$
$4$	$\frac{2}{5}$	$1$	$\frac{2}{5}$

Table 4. Table 4: ( 2 , 2 , 4 ) 2 2 4 (2,2,4) -contractions. List of all ( 2 , 2 , 4 ) 2 2 4 (2,2,4) -contractions ℭ ℭ \mathfrak{C} , their multiplicity # ℭ # ℭ \#\mathfrak{C} under translational symmetry, and their prefactor ∫ ℭ ℭ \int\mathfrak{C} under integration. Summing all terms, weighted with their multiplicity, gives 𝒞 ( 2 , 2 , 4 ) = 648 5 𝔠 1 2 𝔠 2 f 5 = − 1296 5 f 5 subscript 𝒞 2 2 4 648 5 superscript subscript 𝔠 1 2 subscript 𝔠 2 superscript 𝑓 5 1296 5 superscript 𝑓 5 \mathcal{C}_{(2,2,4)}=\frac{648}{5}\,\mathfrak{c}_{1}^{2}\mathfrak{c}_{2}\,f^{5}=-\frac{1296}{5}\,f^{5} .

$ℭ$	$# ℭ$	$\int ℭ$
$,,$	$8$	$1$
	$4$	$1$
$,,,,$	$8$	$\frac{2}{3}$
$,,,,$	$8$	$\frac{2}{3}$
$,,,,$	$8$	$\frac{1}{2}$
$,$	$4$	$\frac{1}{2}$
$,,,,,$	$8$	$\frac{9}{20}$
	$4$	$\frac{4}{9}$
	$2$	$\frac{4}{9}$

Table 5. Table 5: ( 4 , 4 ) 4 4 (4,4) -contractions. List of all ( 4 , 4 ) 4 4 (4,4) -contractions ℭ ℭ \mathfrak{C} , their multiplicity # ℭ # ℭ \#\mathfrak{C} under translational symmetry, and their prefactor ∫ ℭ ℭ \int\mathfrak{C} under integration. Summing all terms, weighted with their multiplicity, gives 𝒞 ( 4 , 4 ) = 1454 63 𝔠 2 2 f 6 = 5816 63 f 6 subscript 𝒞 4 4 1454 63 superscript subscript 𝔠 2 2 superscript 𝑓 6 5816 63 superscript 𝑓 6 \mathcal{C}_{(4,4)}=\tfrac{1454}{63}\,\mathfrak{c}_{2}^{2}\,f^{6}=\frac{5816}{63}\,f^{6} .

$# ℭ$	$\int ℭ$	$# ℭ$	$\int ℭ$
$4$	$1$	$8$	$\frac{2}{3}$
$8$	$\frac{2}{3}$	$2$	$\frac{2}{3}$
$8$	$\frac{11}{20}$	$4$	$\frac{11}{20}$
$1$	$\frac{151}{315}$

Equations237

S_{A} (∣ ψ ⟩) = - Tr (\overset{ρ}{^}_{A} ln \overset{ρ}{^}_{A}), with \overset{ρ}{^}_{A} = Tr_{H_{B}} ∣ ψ ⟩ ⟨ ψ ∣,

S_{A} (∣ ψ ⟩) = - Tr (\overset{ρ}{^}_{A} ln \overset{ρ}{^}_{A}), with \overset{ρ}{^}_{A} = Tr_{H_{B}} ∣ ψ ⟩ ⟨ ψ ∣,

\hat{H}_{XY} = - J j = 1 \sum L [(1 + γ) \hat{S}_{j}^{X} \hat{S}_{j + 1}^{X} + (1 - γ) \hat{S}_{j}^{Y} \hat{S}_{j + 1}^{Y}] - h \hat{S}_{j}^{Z},

\hat{H}_{XY} = - J j = 1 \sum L [(1 + γ) \hat{S}_{j}^{X} \hat{S}_{j + 1}^{X} + (1 - γ) \hat{S}_{j}^{Y} \hat{S}_{j + 1}^{Y}] - h \hat{S}_{j}^{Z},

\hat{H}_{XY}^{\pm} = κ \in K^{\pm} \sum ε_{κ} (\overset{η}{^}_{κ}^{†} \overset{η}{^}_{κ} - \frac{1}{2}) .

\hat{H}_{XY}^{\pm} = κ \in K^{\pm} \sum ε_{κ} (\overset{η}{^}_{κ}^{†} \overset{η}{^}_{κ} - \frac{1}{2}) .

κ \in K^{+}

κ \in K^{+}

K^{-}

\overset{η}{^}_{κ} = \frac{1}{L} j = 1 \sum L (u_{κ} e^{i κj} \hat{f}_{j} - v_{κ}^{*} e^{- i κj} \hat{f}_{j}^{†}),

\overset{η}{^}_{κ} = \frac{1}{L} j = 1 \sum L (u_{κ} e^{i κj} \hat{f}_{j} - v_{κ}^{*} e^{- i κj} \hat{f}_{j}^{†}),

u_{κ} ε_{κ} = \frac{ε _{κ} + a _{κ}}{2 ε _{κ} ( ε _{κ} + a _{κ} )}, v_{κ} = \frac{i b _{κ}}{2 ε _{κ} ( ε _{κ} + a _{κ} )}, = h^{2} + 2 h J cos (κ) + J^{2} + (γ^{2} - 1) J^{2} sin (κ)^{2} .

u_{κ} ε_{κ} = \frac{ε _{κ} + a _{κ}}{2 ε _{κ} ( ε _{κ} + a _{κ} )}, v_{κ} = \frac{i b _{κ}}{2 ε _{κ} ( ε _{κ} + a _{κ} )}, = h^{2} + 2 h J cos (κ) + J^{2} + (γ^{2} - 1) J^{2} sin (κ)^{2} .

\hat{H}_{NI} = - j = 1 \sum L (\hat{f}_{j}^{†} \hat{f}_{j + 1} + H.c.) .

\hat{H}_{NI} = - j = 1 \sum L (\hat{f}_{j}^{†} \hat{f}_{j + 1} + H.c.) .

\mathrm{i}\mathbb{J}(|\psi\rangle)\equiv\left(\begin{array}[]{c|c}\langle\psi|\hat{f}_{j}^{\dagger}\hat{f}_{l}-\hat{f}_{l}\hat{f}_{j}^{\dagger}|\psi\rangle&\langle\psi|\hat{f}_{j}^{\dagger}\hat{f}_{l}^{\dagger}-\hat{f}_{l}^{\dagger}\hat{f}_{j}^{\dagger}|\psi\rangle\\[2.84526pt] \hline\cr\\[-7.8245pt] \langle\psi|\hat{f}_{l}\hat{f}_{j}-\hat{f}_{j}\hat{f}_{l}|\psi\rangle&\langle\psi|\hat{f}_{j}\hat{f}_{l}^{\dagger}-\hat{f}_{l}^{\dagger}\hat{f}_{j}|\psi\rangle\end{array}\right)\,,

\mathrm{i}\mathbb{J}(|\psi\rangle)\equiv\left(\begin{array}[]{c|c}\langle\psi|\hat{f}_{j}^{\dagger}\hat{f}_{l}-\hat{f}_{l}\hat{f}_{j}^{\dagger}|\psi\rangle&\langle\psi|\hat{f}_{j}^{\dagger}\hat{f}_{l}^{\dagger}-\hat{f}_{l}^{\dagger}\hat{f}_{j}^{\dagger}|\psi\rangle\\[2.84526pt] \hline\cr\\[-7.8245pt] \langle\psi|\hat{f}_{l}\hat{f}_{j}-\hat{f}_{j}\hat{f}_{l}|\psi\rangle&\langle\psi|\hat{f}_{j}\hat{f}_{l}^{\dagger}-\hat{f}_{l}^{\dagger}\hat{f}_{j}|\psi\rangle\end{array}\right)\,,

S_{A} (∣ ψ ⟩) = - Tr [(\frac{\mathds 1 _{A} + [ i J ] _{A}}{2}) ln (\frac{\mathds 1 _{A} + [ i J ] _{A}}{2})],

S_{A} (∣ ψ ⟩) = - Tr [(\frac{\mathds 1 _{A} + [ i J ] _{A}}{2}) ln (\frac{\mathds 1 _{A} + [ i J ] _{A}}{2})],

S_{A} (∣ ψ ⟩) = L_{A} ln 2 - n = 1 \sum \infty \frac{Tr [ i J ] _{A}^{2 n}}{4 n ( 2 n - 1 )} .

S_{A} (∣ ψ ⟩) = L_{A} ln 2 - n = 1 \sum \infty \frac{Tr [ i J ] _{A}^{2 n}}{4 n ( 2 n - 1 )} .

S_{A}^{m +} (∣ ψ ⟩)

S_{A}^{m +} (∣ ψ ⟩)

S_{A}^{m -} (∣ ψ ⟩)

\langle S_{A}\rangle=\frac{1\,}{2^{L}}\sum_{\{n_{\kappa}\}}S_{A}\big{(}\ket{\{n_{\kappa}\}}\big{)}\,,

\langle S_{A}\rangle=\frac{1\,}{2^{L}}\sum_{\{n_{\kappa}\}}S_{A}\big{(}\ket{\{n_{\kappa}\}}\big{)}\,,

L \to \infty f \to 0 lim ⟨ S_{A} ⟩ = L_{A} ln 2 .

L \to \infty f \to 0 lim ⟨ S_{A} ⟩ = L_{A} ln 2 .

s (f) = L \to \infty lim \frac{⟨ S _{A} ⟩}{L _{A} ln 2} .

s (f) = L \to \infty lim \frac{⟨ S _{A} ⟩}{L _{A} ln 2} .

s_{m}^{\pm} (f) = L \to \infty lim \frac{⟨ S _{A}^{m \pm} ⟩}{L _{A} ln 2},

s_{m}^{\pm} (f) = L \to \infty lim \frac{⟨ S _{A}^{m \pm} ⟩}{L _{A} ln 2},

t_{m} (f) = L \to \infty lim \frac{⟨ Tr [ i J ] _{A}^{2 m} ⟩}{L _{A}} .

t_{m} (f) = L \to \infty lim \frac{⟨ Tr [ i J ] _{A}^{2 m} ⟩}{L _{A}} .

\mathrm{i}\mathbb{J}\equiv\frac{1}{L}\sum_{{\kappa}\in\mathcal{K}^{\pm}}N_{\kappa}\left(\begin{array}[]{c|c}|u_{\kappa}|^{2}e^{\mathrm{i}{\kappa}(j-l)}-|v_{\kappa}|^{2}e^{\mathrm{i}{\kappa}(l-j)}&u^{*}_{\kappa}v^{*}_{\kappa}\left(e^{\mathrm{i}{\kappa}(l-j)}-e^{\mathrm{i}{\kappa}(j-l)}\right)\\ \hline\cr u_{\kappa}v_{\kappa}\left(e^{\mathrm{i}{\kappa}(l-j)}-e^{\mathrm{i}{\kappa}(j-l)}\right)&-|u_{\kappa}|^{2}e^{\mathrm{i}{\kappa}(l-j)}+|v_{\kappa}|^{2}e^{\mathrm{i}{\kappa}(j-l)}\end{array}\right)\,.

\mathrm{i}\mathbb{J}\equiv\frac{1}{L}\sum_{{\kappa}\in\mathcal{K}^{\pm}}N_{\kappa}\left(\begin{array}[]{c|c}|u_{\kappa}|^{2}e^{\mathrm{i}{\kappa}(j-l)}-|v_{\kappa}|^{2}e^{\mathrm{i}{\kappa}(l-j)}&u^{*}_{\kappa}v^{*}_{\kappa}\left(e^{\mathrm{i}{\kappa}(l-j)}-e^{\mathrm{i}{\kappa}(j-l)}\right)\\ \hline\cr u_{\kappa}v_{\kappa}\left(e^{\mathrm{i}{\kappa}(l-j)}-e^{\mathrm{i}{\kappa}(j-l)}\right)&-|u_{\kappa}|^{2}e^{\mathrm{i}{\kappa}(l-j)}+|v_{\kappa}|^{2}e^{\mathrm{i}{\kappa}(j-l)}\end{array}\right)\,.

⟨ Tr [i J]_{A}^{2} ⟩ = \frac{1}{2} (⟨ Tr [i J]_{A}^{2} ⟩_{+} + ⟨ Tr [i J]_{A}^{2} ⟩_{-}) .

⟨ Tr [i J]_{A}^{2} ⟩ = \frac{1}{2} (⟨ Tr [i J]_{A}^{2} ⟩_{+} + ⟨ Tr [i J]_{A}^{2} ⟩_{-}) .

⟨ Tr [i J]_{A}^{2} ⟩_{\pm} = \frac{1}{2 ^{L - 1}} {N_{κ}}_{\pm} \sum Tr [i J]_{A}^{2},

⟨ Tr [i J]_{A}^{2} ⟩_{\pm} = \frac{1}{2 ^{L - 1}} {N_{κ}}_{\pm} \sum Tr [i J]_{A}^{2},

\displaystyle\begin{split}\langle N_{\kappa}\rangle_{\pm}=0\,,\qquad\langle N_{\kappa_{1}}N_{\kappa_{2}}\rangle_{\pm}=\delta_{\kappa_{1}\kappa_{2}}\,,\\ \langle N_{\kappa_{1}}\cdots N_{\kappa_{n}}\rangle_{\pm}=\left\{\begin{array}[]{ll}1&\text{$\kappa_{i}$ appear in pairs}\\ 0&\text{else}\end{array}\right.\,.\end{split}

\displaystyle\begin{split}\langle N_{\kappa}\rangle_{\pm}=0\,,\qquad\langle N_{\kappa_{1}}N_{\kappa_{2}}\rangle_{\pm}=\delta_{\kappa_{1}\kappa_{2}}\,,\\ \langle N_{\kappa_{1}}\cdots N_{\kappa_{n}}\rangle_{\pm}=\left\{\begin{array}[]{ll}1&\text{$\kappa_{i}$ appear in pairs}\\ 0&\text{else}\end{array}\right.\,.\end{split}

⟨ Tr [i J]_{A}^{2} ⟩_{\pm} = 2 \frac{L _{A}^{2}}{L} - \frac{2}{L ^{2}} j, l = 1 \sum L_{A} k \in K^{\pm} \sum 4∣ u_{κ} ∣^{2} ∣ v_{κ} ∣^{2} cos 2 κ (j - l) .

⟨ Tr [i J]_{A}^{2} ⟩_{\pm} = 2 \frac{L _{A}^{2}}{L} - \frac{2}{L ^{2}} j, l = 1 \sum L_{A} k \in K^{\pm} \sum 4∣ u_{κ} ∣^{2} ∣ v_{κ} ∣^{2} cos 2 κ (j - l) .

⟨ Tr [i J]_{A}^{2} ⟩ = 2 \frac{L _{A}^{2}}{L} - Z (L_{A}, L),

⟨ Tr [i J]_{A}^{2} ⟩ = 2 \frac{L _{A}^{2}}{L} - Z (L_{A}, L),

Z

Z

= \frac{1}{L ^{2}} j, l = 1 \sum L_{A} k = 1 - L \sum L - 1 4 ∣ u_{π k / L} ∣^{2} ∣ v_{π k / L} ∣^{2} cos [\frac{2 π k ( j - l )}{L}]

= \frac{1}{L ^{2}} k = 1 - L \sum L - 1 4 ∣ u_{π k / L} ∣^{2} ∣ v_{π k / L} ∣^{2} \frac{sin ^{2} ( \frac{π k L _{A}}{L} )}{sin ^{2} ( \frac{π k}{L} )},

k = 1 - L \sum L - 1 \frac{sin ^{2} ( \frac{π k L _{A}}{L} )}{sin ^{2} ( \frac{π k}{L} )} = 2 L_{A} (L - L_{A}),

k = 1 - L \sum L - 1 \frac{sin ^{2} ( \frac{π k L _{A}}{L} )}{sin ^{2} ( \frac{π k}{L} )} = 2 L_{A} (L - L_{A}),

0 \leq Z (L_{A}, L) \leq \frac{2 L _{A} ( L - L _{A} )}{L ^{2}} .

0 \leq Z (L_{A}, L) \leq \frac{2 L _{A} ( L - L _{A} )}{L ^{2}} .

t_{1} (f) = L \to \infty lim (2 \frac{L _{A}}{L} - \frac{Z ( L _{A} , L )}{L _{A}}) = 2 f .

t_{1} (f) = L \to \infty lim (2 \frac{L _{A}}{L} - \frac{Z ( L _{A} , L )}{L _{A}}) = 2 f .

z (κ)

z (κ)

= \frac{γ ^{2} J ^{2} sin ^{2} ( L _{A} κ )}{h ^{2} + 2 h J cos ( κ ) + ( γ ^{2} - 1 ) J ^{2} sin ^{2} ( κ ) + J ^{2}} .

Z (L_{A}, L) = \frac{1}{L ^{2}} k = 1 - L \sum L - 1 \frac{sin ^{2} ( \frac{π k L _{A}}{L} )}{2 + 2 cos ( \frac{π k}{L} )} = \frac{L _{A}}{L} - (\frac{L _{A}}{L})^{2},

Peer Reviews

No public reviews on file for this paper yet. If you reviewed it on a platform where reviews are public (OpenReview, ICLR, NeurIPS, ICML), you can paste yours below so the community can read it here.

Videos

No videos yet. Explain this paper in a talk, walkthrough, or lecture? Add one.

Full text

Average eigenstate entanglement entropy of the XY chain in a transverse field

and its universality for translationally invariant quadratic fermionic models

Lucas Hackl

Max Planck Institute of Quantum Optics, Hans-Kopfermann-Straße 1, D-85748 Garching bei München, Germany

Institute for Gravitation and the Cosmos, The Pennsylvania State University, University Park, PA 16802, USA

Department of Physics, The Pennsylvania State University, University Park, PA 16802, USA

Lev Vidmar

Department of Theoretical Physics, J. Stefan Institute, SI-1000 Ljubljana, Slovenia

Marcos Rigol

Department of Physics, The Pennsylvania State University, University Park, PA 16802, USA

Eugenio Bianchi

Institute for Gravitation and the Cosmos, The Pennsylvania State University, University Park, PA 16802, USA

Department of Physics, The Pennsylvania State University, University Park, PA 16802, USA

Abstract

We recently showed [Phys. Rev. Lett. 121, 220602 (2018)] that the average bipartite entanglement entropy of all energy eigenstates of the quantum Ising chain exhibits a universal (for translationally invariant quadratic fermionic models) leading term that scales linearly with the subsystem’s volume, while in the thermodynamic limit the first subleading correction does not vanish at the critical field (it only depends on the ratio $f$ between the volume of the subsystem and volume of the system) and vanishes otherwise. Here we show, analytically for bounds and numerically for averages, that the same remains true for the spin-1/2 XY chain in a transverse magnetic field. We then tighten the bounds for the coefficient of the universal volume-law term, which is a concave function of $f$ . We develop a systematic approach to compute upper and lower bounds, and provide explicit analytic expressions for up to the fourth order bounds.

I Introduction

Entanglement is a hallmark of quantum theory. The correlations between two entangled quantum systems that are in an overall pure state cannot be explained by a local classical theory epr_35 ; bell2001einstein . Different measures of entanglement have been extensively studied in the context of black hole physics Sorkin:2014kta ; Bombelli:1986rw ; srednicki_93 ; holzhey_larsen_94 , holography Ryu:2006bv ; Lewkowycz:2013nqa , ground states of condensed matter systems white_92 ; schollwoeck_05 ; schollwoeck_11 , and quantum information science nielsen_chuang_10 . Recently, the first measurements of an entanglement quantifier were done with ultracold atoms in an optical lattice islam_ma_15 ; kaufman_tai_16 .

Given a quantum system with Hilbert space $\mathcal{H}=\mathcal{H}_{A}\otimes\mathcal{H}_{B}$ , which can be decomposed as a tensor product of $\mathcal{H}_{A}$ and $\mathcal{H}_{B}$ , the bipartite von Neumann entanglement entropy (or, in short, the entanglement entropy) of a pure state $\ket{\psi}\in\mathcal{H}$ is defined as

[TABLE]

where a partial trace is taken over $\mathcal{H}_{B}$ to find the reduced (mixed) state $\hat{\rho}_{A}$ of $\ket{\psi}$ . Here, we are interested in universal features of the entanglement entropy of eigenstates of lattice Hamiltonians for a bipartition into two blocks, in which $L_{A}$ is the volume (number of lattice sites) of the smaller block. In particular, we investigate – both analytically and numerically – the average eigenstate entanglement entropy (over all eigenstates) of the spin-1/2 XY chain in a transverse magnetic field lieb61 , at and away the critical field for which the quantum phase transition occurs in the ground state.

The ground-state entanglement entropy of local quadratic Hamiltonians has been extensively studied in recent years amico_fazio_08 ; peschel_eisler_09 ; calabrese_cardy_09 ; eisert_cramer_10 ; audenaert_eisert_02 ; osterloh_amico_2002 ; osborne_nielsen_02 ; vidal_latorre_03 ; latorre_rico_04 ; calabrese_cardy_04 ; jin_korepin_04 ; calabrese_cardy_04 ; peschel_04 ; its_jin_05 ; wolf_06 ; gioev_klich_06 ; barthel_chung_06 ; li_ding_06 ; cramer_eisert_07 ; franchini_its_07 ; hastings_07 ; ding_brayali_08 , and several remarkable analytic results have been obtained. In contrast, despite growing interest alba09 ; moelter_barthel_14 ; storms14 ; beugeling15 ; lai15 ; nandy_sen_16 ; vidmar_hackl_17 ; riddell_mueller_18 ; zhang_vidmar_18 ; vidmar2018volume ; murciano_ruggiero_18 , an equivalent analytic understanding of the behavior of the entanglement entropy of excited eigenstates is still lacking. Recent analytic results include: (i) bounds for the average entanglement entropy over all eigenstates of translationally invariant quadratic fermionic models vidmar_hackl_17 , and (ii) the average entanglement entropy over eigenstates of the quadratic part of the Sachdev-Ye-Kitaev model liu_chen_18 . Those results have provided proofs that the leading term of the average entanglement entropy $\langle S_{A}\rangle$ in the corresponding quadratic models scales with the volume, $\langle S_{A}\rangle\propto L_{A}$ . In another recent study vidmar2018volume , it was shown that for the quantum Ising chain the leading term of $\langle S_{A}\rangle$ is independent of the Hamiltonian parameters, and it was conjectured to be universal for all translationally invariant quadratic fermionic Hamiltonians. Furthermore, it was shown that the subleading term of $\langle S_{A}\rangle$ in the quantum Ising chain depends on whether the system is at the critical field or not. Generalizing these results to the spin-1/2 XY chain in a transverse field is the first goal of this paper.

The bounds in Ref. vidmar_hackl_17 , as well as numerical results for noninteracting fermions storms14 ; vidmar_hackl_17 , have shown that the coefficient of the volume-law of the average eigenstate entanglement entropy of translationally invariant quadratic fermionic Hamiltonians is smaller than the theoretical maximum. The latter is attained by the (Haar measure) average over pure states in the Hilbert space page93 and by the average over eigenstates of quantum-chaotic Hamiltonians vidmar_rigol_17 . Motivated by questions about thermalization in generic isolated quantum systems dalessio_kafri_16 , the entanglement properties of excited eigenstates of quantum-chaotic Hamiltonians have been the focus of much recent work deutsch_10 ; santos_polkovnikov_12 ; hamma_santra_12 ; deutsch13 ; alba15 ; beugeling15 ; yang_chamon_15 ; vivo_pato_16 ; vidmar_rigol_17 ; garrison_grover_18 ; dymarsky_lashkari_18 ; nakagawa_watanabe_18 ; huang_19 ; luitz2016long ; lu2017renyi . In the context of thermalization, it is of particular interest to understand how entanglement properties differ between many-body energy eigenstates of quantum chaotic systems (systems that are expected to thermalize under unitary dynamics) and many-body energy eigenstates of integrable systems (systems that are not expected to thermalize under unitary dynamics).

An important property of the volume-law coefficient of the average eigenstate entanglement entropy in translationally invariant quadratic fermionic Hamiltonians (which are a class of particularly simple integrable models) is that it is a concave function of the ratio between the volume of the subsystem $L_{A}$ and the volume of the system $L$ (in short, the subsystem fraction $f=L_{A}/L$ ) for $f>0$ . To learn more about that concave function, the second goal of this paper is to sharpen the bounds reported in Ref. vidmar_hackl_17 for the volume-law coefficient.

The presentation is structured as follows. In Sec. II, we study average entanglement entropies of eigenstates of the spin-1/2 XY chain in a transverse field. Bounds are determined analytically, while exact averages are computed numerically. In Sec. III, we tighten the bounds reported in Ref. vidmar_hackl_17 for the conjectured universal leading term of the average entanglement entropy in translationally invariant quadratic fermionic Hamiltonians. In particular, we provide explicit analytical expressions for up to the fourth order bounds in our perturbative expansion for noninteracting fermions. A summary and discussion of our results is provided in Sec. IV. In Appendix B, we discuss details about the effect of eigenstates in the fermionic even/odd particle-number sectors, which are not eigenstates of the spin Hamiltonian, in calculations of averages. Appendix C provides a comprehensive review of our perturbative method, based on a generalized Wick’s theorem, and reports details of the calculations of up to the fourth order bounds.

Figure 1 reports the main results of Sec. III, which are among the main results of this work. In Fig. 1, we compare the average entanglement entropy, obtained numerically, with the fourth order upper and lower bounds, and with the exact asymptotic expression for the average for up to $f^{4}$ . That figure shows that the fourth order bounds, specially the lower one, are very tight, and that the asymptotic expression is accurate up to $f\approx 0.3$ .

II XY model in a transverse field

We consider the spin-1/2 XY model in a transverse field in one dimension, in systems with $L$ sites and periodic boundary conditions lieb61 . We denote sites as $j=1,\ldots,L$ , with site $L+1$ identified as site $1$ . For simplicity, we assume that $L$ is an even integer. The Hamiltonian can be written as

[TABLE]

where $\hat{S}^{\mathrm{X}}_{j}$ , $\hat{S}^{\mathrm{Y}}_{j}$ , and $\hat{S}^{\mathrm{Z}}_{j}$ are standard spin-1/2 operators, $J$ is the interaction strength, $h$ is the strength of the transverse field, and $0\leq\gamma\leq 1$ is the anisotropy parameter.

This model can be solved by applying first a Jordan-Wigner transformation and subsequently a Bogoliubov transformation. The resulting Hamiltonian can be written as sum $\hat{H}_{\mathrm{XY}}=\hat{H}_{\mathrm{XY}}^{+}\mathcal{P}^{+}+\hat{H}_{\mathrm{XY}}^{-}\mathcal{P}^{-}$ where $\mathcal{P}^{\pm}$ are projectors onto Hilbert sub spaces $\mathcal{H}^{\pm}$ with even ( $+$ ) and odd ( $-$ ) number of excitations. The individual terms

[TABLE]

are quadratic in fermionic creation and annihilation operators $\hat{\eta}_{k}^{\dagger}$ and $\hat{\eta}_{\kappa}^{\dagger}$ with

[TABLE]

We can express them in terms of localized creation and annihilation operators $\hat{f}_{j}^{\dagger}$ and $\hat{f}_{j}^{\dagger}$ at site $j$ via

[TABLE]

where the Bogoliubov coefficients and energies are

[TABLE]

A detailed derivation of this result can be found in Appendix A.

We will compare the results for the spin-1/2 XY model in a transverse field with those for noninteracting fermions, whose Hamiltonian has the form

[TABLE]

This is a special case of the more general class of quadratic fermionic Hamiltonians with translational invariance (see Eq. 41).

II.1 Eigenstate entanglement entropy

The eigenstates $\ket{\{n_{\kappa}\}}$ are fermionic Gaussian states peschel_03 ; peschel_eisler_09 ; dierckx_fannes_08 ; hackl_bianchi_17 . For simplicity, we denote them as $|\psi\rangle$ in what follows. We compute the entanglement entropy $S_{A}(\ket{\psi})$ , see Eq. (1), with respect to the Hilbert space decomposition $\mathcal{H}=\mathcal{H}_{A}\otimes\mathcal{H}_{B}$ by splitting $L$ sites into a connected piece of length $L_{A}$ and its complement with $L_{B}=L-L_{A}$ sites.

An efficient method of computing $S_{A}(\ket{\psi})$ is based on parametrizing Gaussian states in terms of the linear complex structure $\mathbb{J}(|\psi\rangle)$ , defined as

[TABLE]

where $j$ and $l$ run over all sites on the system, so that each matrix entry in Eq. (9) represents a block of size $L\times L$ , i.e., $\mathrm{i}\mathbb{J}(|\psi\rangle)$ is a $2L\times 2L$ matrix. For simplicity, in what follows we drop the state label $|\psi\rangle$ in $\mathrm{i}\mathbb{J}$ . The entanglement entropy of $\ket{\psi}$ , associated to the subsystem $\mathcal{H}_{A}$ , can be computed using the formula

[TABLE]

where $\mathds{1}_{A}$ and $[\mathrm{i}\mathbb{J}]_{A}$ are restrictions to the $2L_{A}\times 2L_{A}$ matrices with rows and columns associated to $j$ and $l$ in the subsystem only. The advantage of using this representation to evaluate Eq. (1) for Gaussian states relies on the reduction from working with $2^{L_{A}}\times 2^{L_{A}}$ many-body matrices to working with $2L_{A}\times 2L_{A}$ one-body matrices.

In Ref. vidmar_hackl_17 , we took a step forward in simplifying the evaluation of Eq. (10). It is based on the fact that the eigenvalues of $[\mathrm{i}\mathbb{J}]_{A}$ are limited to the interval $[-1,1]$ . As a result, $S_{A}(\ket{\psi})$ can be expressed as a converging series

[TABLE]

Since the traces of even powers of $[\mathrm{i}\mathbb{J}]_{A}$ are all positive, one can modify the series in Eq. (11) at any order $m$ to find upper and lower bounds $S^{m+}_{A}(\ket{\psi})$ and $S^{m-}_{A}(\ket{\psi})$ , respectively, given by

[TABLE]

We are interested in the average $\langle S_{A}\rangle$ over all eigenstates $\ket{\{n_{\kappa}\}}$ in the thermodynamic limit ( $L\to\infty$ ), with the average defined as

[TABLE]

namely, one first computes the entanglement entropy $S_{A}(\ket{\{n_{\kappa}\}})$ of each eigenstate $\ket{\{n_{\kappa}\}}$ , and then takes the average over all eigenstates. If the Hamiltonian has a non-degenerate spectrum, then the sum over the numbers $\{n_{\kappa}\}$ coincides with the sum over eigenstates of energy $E=\sum_{\kappa\in\mathcal{K}^{\pm}}n_{\kappa}\varepsilon_{\kappa}$ . If, on the other hand, the Hamiltonian has a degenerate spectrum, the sum over $\{n_{\kappa}\}$ resolves the degeneracy within each eigenspace by introducing the orthonormal basis of Gaussian states $\ket{\{n_{\kappa}\}}$ . This procedure is equivalent to removing the degeneracy by a vanishingly small perturbation of the $\varepsilon_{\kappa}$ ’s, such that no two states $\ket{\{n_{\kappa}\}}$ have the same energy.

In Ref. vidmar_hackl_17 , we proved that the leading term of $\langle S_{A}\rangle$ for translationally invariant quadratic fermionic models exhibits a volume-law scaling. We also proved that, in the thermodynamic limit when $f\rightarrow 0$ , where $f=L_{A}/L$ , the overwhelming majority of eigenstates exhibit the maximal entanglement entropy, namely,

[TABLE]

Here we are in interested in $\langle S_{A}\rangle$ when $L\to\infty$ for $f>0$ . We define the coefficient of its volume-law term as

[TABLE]

In Ref. vidmar2018volume we conjectured that, for quadratic fermionic Hamiltonians with translational invariance, $s(f)$ is a universal function of the subsystem fraction $f$ . $s(f)$ is the central object studied in Sec. III.2.

Accordingly, we also define the coefficients of the leading terms of the $m$ -th order bounds

[TABLE]

which, see Eqs. (12) and (13), depend on

[TABLE]

The average $\langle...\rangle$ in Eq. (18) is defined as in Eq. (14). As we show next, $t_{m}(f)$ can be computed from averages $\langle N_{\kappa_{1}}\cdots N_{\kappa_{2m}}\rangle$ over eigenstate occupation numbers.

II.2 First order bound: Leading term

Here we calculate the leading term in $t_{1}(f)$ [see Eq. (18)]. The first step is to calculate the matrix elements of $[\mathrm{i}\mathbb{J}]$ , see Eq. (9), for an eigenstate $\ket{\{n_{\kappa}\}}$ . The matrix $[\mathrm{i}\mathbb{J}]$ can be expressed in terms of $N_{\kappa}=(-1)^{n_{\kappa}}$ as

[TABLE]

Here the choice of the set $\mathcal{K}^{\pm}$ , introduced in Eqs. (4)–(5), depends on whether an eigenstate $\ket{\{n_{\kappa}\}}$ contains an even or an odd number of quasiparticles. We first compute the partial average $\langle\mathrm{Tr}[\mathrm{i}\mathbb{J}]_{A}^{2}\rangle_{\pm}$ , where we either take only eigenstates with an even ( $+$ ) or an odd ( $-$ ) number of quasiparticles, and then compute the average over all eigenstates

[TABLE]

It is easy to check that there is an equal number of eigenstates in the even and the odd sectors $\mathcal{H}^{\pm}$ of the Hilbert space.111Consider $0=(1-1)^{L}=\sum^{L}_{l=0}{L\choose l}(-1)^{l}(1)^{L-l}$ . Here, $l$ represents the number of excitations correctly weighted by the binomial coefficient. The positive terms represent the even sectors (even $l$ ) and the negative ones the odd sector (odd $l$ ). Clearly, both terms cancel each other, i.e., the even and odd sector both have Hilbert space dimension $2^{L-1}$ . The partial average is defined as

[TABLE]

where $\{N_{\kappa}\}_{\pm}$ corresponds to the sets of all $N_{\kappa}\in\{1,-1\}$ with $\kappa\in\mathcal{K}^{\pm}$ and $\prod_{\kappa\in\mathcal{K}^{\pm}}N_{\kappa}=\pm 1$ . We stress that we compute the average over the eigenstates of $\hat{H}^{+}_{\mathrm{XY}}$ with an even number of excitations and the average over the eigenstates of $\hat{H}^{-}_{\mathrm{XY}}$ with an odd number of excitations.

The partial average $\langle\mathrm{Tr}[\mathrm{i}\mathbb{J}]_{A}^{2}\rangle_{\pm}$ can be computed using the technique introduced in Ref. vidmar_hackl_17 , which is based on the partial average $\langle N_{\kappa_{1}}\cdots N_{\kappa_{n}}\rangle_{\pm}$ . The definition of the latter is analogous to the one in Eq. (21). As we are summing over all states, each number $N_{\kappa_{i}}$ takes values $\pm 1$ , subject to the constraint that the total number of $-1$ ’s is even or odd. However, this constraint does not affect the expectation values $\langle N_{\kappa_{1}}\cdots N_{\kappa_{n}}\rangle_{\pm}$ for $n<L$ , which are given by vidmar_hackl_17

[TABLE]

Only at order $n=L$ , these correlation functions are modified due to the above constraint, as explained in Appendix B. As we are taking the limit $L\to\infty$ , we can use the above expectation values for the bounds $s^{\pm}_{m}(f)$ at any order $m$ by assuming that $L>m$ .

Using Eq. (22), one can compute the partial average

[TABLE]

Therefore, taking the total average $\langle\mathrm{Tr}[\mathrm{i}\mathbb{J}]_{A}^{2}\rangle$ defined in Eq. (20) is equivalent to summing over all $\kappa\in\mathcal{K}^{+}\cup\mathcal{K}^{-}$ and dividing by $2$ . Thus, we find

[TABLE]

where

[TABLE]

in which the sum over $k$ runs only up to $k=L-1$ , because $u_{\pi}=0$ . We can use

[TABLE]

together with the inequality $|u_{\kappa}|^{2}|v_{\kappa}|^{2}\leq\frac{1}{4}$ , because $|u_{\kappa}|^{2}+|v_{\kappa}|^{2}=1$ , to find the bound

[TABLE]

This implies that $\lim_{L\to\infty}Z(L_{A},L)\leq{\rm const}$ . Using Eqs. (18) and (24), we arrive at

[TABLE]

This result agrees with the one obtained for translational invariant quadratic models vidmar_hackl_17 . Note, however, the XY chain in a transverse field has a nonquadratic boundary term. Our result for this model shows that the boundary term does not change the leading behavior of $t_{1}$ .

II.3 First order bound: Subleading term

Next, we calculate the first subleading term in $t_{1}(f)$ in the thermodynamic limit. While it is zero for noninteracting fermions vidmar_hackl_17 , it was found for the quantum Ising model ( $\gamma=1$ ) that it is a nonvanishing constant at the critical field $h=J$ and vanishes otherwise vidmar2018volume . Here we generalize that result to arbitrary values of $\gamma$ .

The subleading term of the first order bound is obtained by calculating the difference between the exact result $\langle{\rm Tr}[\mathrm{i}\mathbb{J}]_{A}^{2}\rangle$ and the leading term $2fL_{A}$ . This is exactly $Z(L_{A},L)$ , defined in Eq. (24). Following Eq. (25), one can write $Z(L_{A},L)=\frac{1}{L^{2}}\sum^{L-1}_{k=1-L}z(\pi k/{L})$ by introducing the function $z(\kappa)$ as

[TABLE]

We note immediately that, for $\gamma=0$ (noninteracting fermions), one has $Z(L_{A},L)=0$ . In order to bound the scaling of $Z(L_{A},L)$ as $L\to\infty$ for $\gamma>0$ , we need to distinguish two cases:

**(i) Critical line ( $h=J$ ):

**We first consider the simplest case, $\gamma=1$ , for which one can perform the sum analytically

[TABLE]

where the fact that $L_{A}$ is an integer allowed us to simplify the summation. This shows that the subleading term $Z(L_{A},L)$ does not vanish if $f>0$ , and only depends on $f$ . For $0<\gamma<1$ , we cannot compute the sum analytically, but we can estimate its scaling. Replacing $\sin^{2}(\kappa L_{A})$ by 1 leads to an upper bound

[TABLE]

which diverges at $\kappa=\pm\pi$ . However, the function $z(\kappa)$ does not diverge, it reaches its maximum

[TABLE]

where we used L’Hospital’s rule twice.

To show that, for a fixed $f$ , $\lim_{L\to\infty}Z(L_{A},L)$ reaches a nonzero value, we approximate $\frac{1}{L}\sum_{k=1-L}^{L-1}z(\pi k/L)\to\frac{1}{\pi}\int^{\pi}_{-\pi}z(\kappa)d\kappa$ . The function $z(\kappa)$ is shown in Fig. 2. The scaling of the integral can be predicted by bounding $z(\kappa)$ from below by the triangle function

[TABLE]

and from above either by the maximum value $L_{A}^{2}$ or by the asymptotic expansion $z^{+}(\kappa)=(\pi-|\kappa|)^{-2}$ as $|\kappa|\to\pi$ . The integral over both functions is proportional to $L_{A}$ . This implies that

[TABLE]

and proves that $\lim_{L\to\infty}Z(L_{A},L)\sim\mathcal{O}(L_{A}/L)$ for any $\gamma>0$ , i.e., that the first subleading correction does not vanish if $f>0$ , and only depends on $L_{A}$ through $f$ .

**(ii) Regime $h\neq J$ :

**We have $z(\kappa)\leq z^{+}(\kappa)$ using that $\sin^{2}(\kappa L_{A})\leq 1$ , where

[TABLE]

$z^{+}(\kappa)$ takes its maximal values at $\kappa=\pm\pi$ , given by

[TABLE]

Therefore

[TABLE]

This shows that, irrespective of the value of $f$ , the subleading term of the first order bound vanishes in the thermodynamic limit.

II.4 Average entanglement entropy

Here, we report numerical results for the average eigenstate entanglement entropy $\langle S_{A}\rangle$ in the spin-1/2 XY chain in a transverse field. They allow us to extract the leading and first subleading terms in the thermodynamic limit, thereby complementing the results for the quantum Ising model ( $\gamma=1$ ) reported in Ref. vidmar2018volume .

We define

[TABLE]

where $\langle S_{A}\rangle$ is the numerical average in a finite system, and $s(f)_{\rm NI}$ is an extrapolated value of $\langle S_{A}\rangle/(L_{A}\ln{2})$ for noninteracting fermions in the thermodynamic limit. For noninteracting fermions, we can accurately obtain the thermodynamic limit value with (at least) a four digit precision because finite-size effects decay exponentially with $L$ vidmar_hackl_17 . We find $s(f=1/2)_{\rm NI}=0.5378(1)$ and $s(f=1/4)_{\rm NI}=0.7939(8)$ .

Figure 3 shows scaling analyses of $\Delta s(f)$ at $f=1/2$ (left panels) and $1/4$ (right panels), for different values of $\gamma$ and $h/J$ . In all cases, we find $\lim_{L\to\infty}\Delta s(f)\to 0$ , as for the quantum Ising model vidmar2018volume . This provides further evidence, in addition to the one reported in Ref. vidmar2018volume , of the universality of the leading order term for quadratic fermionic Hamiltonians (and the spin models mappable to them) with translational invariance.

The finite-size scaling analyses of $\Delta s(f)$ in Fig. 3 also uncover the nature of the subleading term of $\langle S_{A}\rangle$ . They show that, for $\gamma>0$ and $f>0$ , $\Delta s(f)$ scales as $1/L$ at $h=J$ and as $1/L^{2}$ otherwise. This suggests that, in the spin-1/2 XY chain in a transverse field, $\langle S_{A}\rangle$ in the thermodynamic limit exhibits a subleading term that does not vanish in the critical line $h=J$ ( $\gamma>0$ ) for $f>0$ , and vanishes otherwise. Therefore, as found for the quantum Ising model in Ref. vidmar2018volume , the numerical results for the average and the analytical results for the first order bounds exhibit identical finite-size scalings of the leading and first subleading terms.

In Fig. 3 we also show results for the entanglement entropy averaged only over eigenstates with an even (or odd) number of quasiparticles. We define [see Eq. (21)]

[TABLE]

The results for $\Delta s(f)_{\pm}$ in Fig. 3 suggest that, generically, $\Delta s(f)_{+}$ and $\Delta s(f)_{-}$ exhibit a different subleading term. One may wonder what happens if one incorrectly takes the average over all eigenstates using a single set of $\kappa$ values (either $\kappa\in{\cal K}^{+}$ or $\kappa\in{\cal K}^{-}$ ) instead of the correct average involving the true eigenstates, one half of which come from the even sector ( $\kappa\in{\cal K}^{+}$ ) and the other half of which come from the odd sector ( $\kappa\in{\cal K}^{-}$ ). Figure 6 in Appendix B shows that the incorrect calculation produces an error that vanishes $\propto\exp(-\alpha L^{\beta})$ , where $\alpha$ depends, in general, on the Hamiltonian parameters, the sector used for the average, and on $f$ . Numerically, we find that $\beta=1$ for $f=1/4$ , and $\beta=1/2$ for $f=1/2$ .

III Universal volume-law coefficient

We focus next on tightening the bounds for the volume-law coefficient $s(f)$ defined in Eq. (16), which is a concave function of the subsystem fraction $f$ that was conjectured to be universal for quadratic fermionic Hamiltonians with translational invariance vidmar2018volume . In Sec. III.1, we provide details about the universality conjecture for $s(f)$ . In Sec. III.2, we report the fourth order bounds for $s(f)$ for noninteracting fermions.

III.1 Conjecture about universality

In Ref. vidmar2018volume , we reported numerical results supporting the conjecture that the volume-law coefficient $s(f)$ is universal for quadratic fermionic Hamiltonians with translational invariance. That conjecture is further supported by our earlier results here for the spin-1/2 XY chain in a transverse field, for which $s(f)$ was found to be identical to that of noninteracting fermions independently of the values of $\gamma$ , $h$ , and $f$ chosen.

We should stress that this universality cannot include all quadratic fermionic Hamiltonians. A simple counterexample is given by the ultra-local Hamiltonian

[TABLE]

with randomly selected $\epsilon_{j}$ ’s. Every eigenstate of this model has zero entanglement entropy. Another counterexample involving randomness is the quadratic part of the Sachdev-Ye-Kitaev model, for which the average eigenstate entanglement entropy was computed analytically liu_chen_18 and differs from our numerical results in Ref. vidmar_hackl_17 for translationally invariant quadratic systems.

The most general translationally invariant quadratic Hamiltonian in a one-dimensional lattice with periodic boundary conditions, which can be generalized to higher dimensions in a straightforward way, can be written as

[TABLE]

with $j,l\in\mathbb{Z}_{L}$ , namely, $j+L\equiv j$ . Hamiltonian (41) can be diagonalized by first applying a Fourier transformation and then performing a two mode squeezing. Explicitly, we define

[TABLE]

with $\mathcal{K}=\left\{\frac{2\pi k}{L}\,\big{|}\,k\in\mathbb{Z}\,,-\frac{L}{2}<k\leq\frac{L}{2}\right\}$ , giving rise to the Fourier transformed Hamiltonian

[TABLE]

with the parameters

[TABLE]

For the two mode squeezing, one determines the parameters $u_{\kappa}$ and $v_{\kappa}$ so that $\hat{\eta}_{\kappa}=u_{\kappa}\hat{c}_{\kappa}-v_{\kappa}^{*}\hat{c}_{-\kappa}^{\dagger}$ leads to the quadratic Hamiltonian

[TABLE]

where $E_{0}$ is an overall constant. Note that $u_{\kappa}$ and $v_{\kappa}$ satisfy the conditions

[TABLE]

which follow from the fact that one squeezes the mode pairs $(\kappa,-\kappa)$ . One always has $u_{0}=1$ and $v_{0}=0$ , and, for even $L$ , one also has $u_{\pi}=1$ and $v_{\pi}=0$ .

In order to test the universality of the volume-law coefficient $s(f)$ for translationally invariant quadratic Hamiltonians, in Ref. vidmar2018volume we studied a class of Hamiltonians defined by choosing the parameters $u_{\kappa}$ and $v_{\kappa}$ to be simple continuous and periodic functions on the interval $[0,\pi]$ .222A natural choice of regular limit functions are those continuous on the circle, i.e., periodic such that $u_{-\pi}=u_{\pi}$ and $v_{-\pi}=v_{\pi}=0$ . It is conceivable that one may also need to require additional regularity conditions on $u_{\kappa}$ and $v_{\kappa}$ . This is challenging to verify numerically and remains to be done. Finite-size scaling analyses showed that, in the thermodynamic limit, $s(f)$ in those Hamiltonians becomes identical to that of noninteracting fermions.

Based on those results, we conjectured that the volume-law coefficient $s(f)$ is universal for quadratic fermionic Hamiltonians with translational invariance vidmar2018volume . The conjecture can be formulated as follows.

Universality conjecture: For translationally invariant quadratic fermionic Hamiltonians in arbitrary dimensions, with $V$ sites and a regular subsystem region consisting of $V_{A}$ sites, the volume-law coefficient of the average (over all eigenstates) eigenstate entanglement entropy is a universal function in the thermodynamic limit:

[TABLE]

In the case of degenerate eigenspaces, the basis for the average should be that obtained in the limit of vanishingly small translationally invariant perturbations that break all degeneracies.

We also conjecture that $s(f)$ is typical for the entanglement entropy of eigenstates, which implies that:

[TABLE]

A few comments are in order about this conjecture.

(i) Translational invariance refers to a discrete symmetry of a regular lattice in $D$ dimensions. For definiteness, one can consider a hypercubic lattice with periodic boundary conditions and $L$ sites along each principal direction, so that the total number of sites is $V=L^{D}$ . This is the setup used to prove the universal first order bound with a Hamiltonian that is invariant under discrete lattice translations vidmar_hackl_17 .

(ii) Regular subsystem region refers to subsets consisting of $V_{A}$ sites that have a well-defined geometry in the thermodynamic limit. This intuition stems from our proof for the universality of the first order bound vidmar_hackl_17 . Even though we only considered subsystems defined as cuboids of length $L_{d}$ along the $d$ -th dimension, we expect that the results do not change if one deforms these regions as long as they have a well defined thermodynamic limit. This excludes fractal subsystems or other definitions that select subsystem sites based on some nongeometric criterion. The reason for this condition lies in the fact that already at the level of the first order bound, the generalization of Eq. (23) to arbitrary dimensions yields

[TABLE]

where one could select the points $\vec{x}_{i}$ in the subsystem as a function of the total system size to get specific contributions from the cosine function that would violate our bound. This is not possible if the subsystem approaches a regular geometric region in the thermodynamic limit.

(iii) Typicality is an additional part of the conjecture stating that the average volume-law coefficient $s(f)$ is typical for energy eigenstates. This means that, in the thermodynamic limit, all but a measure zero of the energy eigenstates have a leading term of the entanglement entropy described by $L_{A}\ln 2\,s(f)$ . In Ref. vidmar_hackl_17 we proved this for noninteracting fermions by showing that the variance $\Sigma_{s}(f)$ , defined in Eq. (48), vanishes at least as $1/\sqrt{L}$ with increasing system size.

III.2 Noninteracting fermions

Motivated by the conjecture discussed in Sec. III.1, and in order to learn more about $s(f)$ for quadratic fermionic Hamiltonians with translational invariance, we tighten the bounds reported in Ref. vidmar_hackl_17 for noninteracting fermions in one dimension [Eq. (8)]. The single-particle eigenstates of that model are plane waves, so this is the model for which the analytic calculation of the bounds is easiest.

III.2.1 Expansion up to fourth order

In the basis given by $\hat{\xi}^{a}\equiv(\hat{f}_{1},\cdots,\hat{f}_{L},\hat{f}_{1}^{\dagger},\cdots,\hat{f}_{L}^{\dagger})$ , the matrix $\mathrm{i}\mathbb{J}$ , see Eq. (9), is block diagonal and explicitly given by

[TABLE]

with $\mathcal{K}=\left\{\frac{2\pi k}{L}\,\big{|}\,k\in\mathbb{Z}\,,-\frac{L}{2}<k\leq\frac{L}{2}\right\}$ .

As we proved in Ref. vidmar_hackl_17 , the coefficients of the leading term of the first order bounds are universal for translationally invariant quadratic systems

[TABLE]

The same results follow for the spin-1/2 XY model in the presence of a transverse field, after substituting $t_{1}(f)$ from Eq. (28) in Eqs. (12) and (13).

Computing higher order bounds becomes increasingly complicated, and there are several subtleties to note. We develop a general procedure to compute averages of traces of $[\mathrm{i}\mathbb{J}]_{A}^{2n}$ based on generalized Wick contractions. This procedure is discussed in Appendix C. The main insight from our analysis is that the term $\langle{\rm Tr}[\mathrm{i}\mathbb{J}]_{A}^{2n}\rangle/L_{A}$ is a polynomial that, when $L\to\infty$ and $f\to 0$ , only contains powers from $f^{n}$ to $f^{2n-1}$ . Table 1 lists the coefficients of $f^{n}$ to $f^{2n-1}$ in $t_{n}(f)$ , defined in Eq. (18), for $n\leq 4$ (see Appendix C for the corresponding calculations).

Using $t_{n}(f)$ , one can write the volume-law coefficients of higher order bounds $s^{\pm}_{n}(f)$ [see Eqs. (12), (13), and (17)] as

[TABLE]

Substituting the coefficients from Table 1 in Eqs. (52) and (53), we find the following fourth order bounds

[TABLE]

where the function $t_{4}^{\mathrm{cor}}(f)$ is a piecewise defined polynomial, see Eq. (123) and its derivation in Appendix C.3. $t_{4}^{\mathrm{cor}}(f)$ arises because of a subtle nonanalytic contribution in the computation of the traces, which appears at fourth order for $f\geq\frac{1}{4}$ . In Fig. 4, we plot all the volume-law coefficients of the bounds $s_{n}^{\pm}(f)$ , from the first order to the fourth order, along with the numerically calculated volume-law coefficient of the average $s(f)_{\rm NI}$ . The bounds can be seen to be very tight, specially $s_{4}^{-}(f)$ .

III.2.2 Asymptotic analysis

As mentioned in Sec. III.2.1, and as discussed in Appendix C.3 and Ref. vidmar_hackl_17 ,

[TABLE]

where the requirement $f\in[0,\frac{1}{n}]$ ensures that nonanalytic contributions of the type $t_{n}^{\mathrm{cor}}(f)$ , seen already in Eqs (54) and (55), are absent. Using this, if one looks into Eqs. (52) and (53), one can see that the volume-law coefficient of the upper bound $s_{n}^{+}(f)$ describes the exact asymptotic expansion of $s(f)$ through order $f^{n}$ .

Therefore, the volume-law coefficient of the upper bound $s_{4}^{+}(f)$ , see Eq. (54), describes the exact asymptotics of $s(f)$ through order $f^{4}$

[TABLE]

We plot $s_{\rm asy}(f)$ in Fig. 1, and compare it to the numerical results for $s(f)_{\rm NI}$ and to the fourth order bounds, $s_{4}^{\pm}(f)$ , reported in Eqs. (54) and (55). The agreement between all the curves is remarkable up to $f\approx 0.3$ .

The fact that, at order $n$ , we generically find nonanalytic behavior for $f>\frac{1}{n}$ leads to the question of whether the limit function is nonanalytic. We are actually not aware of any simple analytic function whose asymptotic expansion around $f\to 0$ agrees with the terms computed so far. On the other hand, there is the possibility that the nonanalytic parts cancel each other in the total sum, so that question remains open at this point.

In light of the cumbersome expressions for the volume-law coefficients of the upper and lower bounds up to fourth order, we test their correctness by comparing our analytical results to those obtained using exact numerical calculations. We compute $\langle\mathrm{Tr}[\mathrm{i}\mathbb{J}]_{A}^{2n}\rangle/L_{A}$ numerically in finite systems, and then study the scaling of the difference between the numerical and the analytical (from Table 1) results with increasing the system size. Figure 5 shows that the difference

[TABLE]

vanishes as $1/L^{2}$ both for $t_{3}(f)$ and $t_{4}(f)$ . Of particular importance is the agreement at fourth order, which depends on the piecewise defined nonanalytic part $t_{4}^{\mathrm{cor}}(f)$ for $f>\frac{1}{4}$ .

IV Summary and Discussion

In Sec. II, we discussed analytical results for the first order bounds of the average eigenstate entanglement entropy $\langle S_{A}\rangle$ for the translationally invariant spin-1/2 XY chain in a transverse field. Those bounds were found to exhibit leading terms that scale linearly with the subsystem volume, which are independent of the model parameters ( $J$ , $\gamma$ , and $h$ ), and that are identical to the ones obtained for noninteracting fermions. The leading terms of the volume-law coefficients of the bounds only depend on the subsystem fraction $f$ . On the other hand, in the thermodynamic limit, the subleading terms of the bounds were found to be nonzero in the critical line $h=J$ ( $\gamma>0$ ) for $f>0$ , and vanish otherwise. The same behavior was obtained numerically for the average eigenstate entanglement entropy. These results generalize our findings in Ref. vidmar2018volume for the quantum Ising chain, which corresponds to $\gamma=1$ in Eq. (2).

In Sec. III, motivated by the conjectured universality of the leading order term of the average eigenstate entanglement entropy of translationally invariant quadratic fermionic Hamiltonians (see Sec. III.1), we tightened the bounds for the volume-law coefficient $s(f)$ for noninteracting fermions. $s(f)$ is a concave function of the subsystem fraction $f$ , and we computed up to its fourth order bounds. For this, we developed a series expansion method that can be applied to arbitrarily high orders in the expansion of $s(f)$ . About $f=0$ , our results determine the series expansion of the volume-law coefficient $s(f)$ through $f^{4}$ . Our analytical calculations reduce to the determination of combinatoric factors, which can be found algorithmically using the diagrammatic representation of the generalized Wick’s theorem presented in Appendix C.2.

Intriguingly, we found that nonanalytic contributions arise in volume-law coefficients that provide bounds for $s(f)$ . The nonanalytic contributions are piecewise defined polynomials of $f$ , and the series expansion allows us to compute them in each piecewise range of the subsystem fraction $f$ . The piecewise polynomial corrections to $s(f)$ for $f>\frac{1}{n}$ (where $n$ is the order of the series expansion) suggest that the function $s(f)$ may not be analytic. In principle, about $f=0$ , one can find the expansion in $f$ to arbitrarily high order, but for every interval $[0,\epsilon]$ , there exists an $n$ with $n>\frac{1}{\epsilon}$ , such that a nonanalytic correction appears. We cannot exclude the possibility that these piecewise defined terms ultimately cancel each other in the total sum, so the question of whether $s(f)$ is an analytic function of $f$ remains open. We should emphasize though that, even at fourth order, the piecewise correction represents an almost negligible change to the function value.

Acknowledgements.

We thank J. Eisert and I. Roth for discussions. We acknowledge support from a Mebus Fellowship (L.H.), the Max Planck Harvard Research Center for Quantum Optics (L.H.), the Slovenian Research Agency research core funding No. P1-0044 (L.V.), and the National Science Foundation Grant Nos. PHY-1707482 (M.R.) and PHY-1806428 (E.B.). The computations were done at the Institute for CyberScience at Penn State.

Appendix A Diagonalization of the XY model

For the XY model introduced in Eq. (2), we explain in detail the required transformations that lead to Eq. (3). We follow closely Ref. vidmar16 and draw attention to subtleties arising from the boundary terms.

First, expressing the spin operators $\hat{S}_{j}^{\mathrm{X}}$ and $\hat{S}_{j}^{\mathrm{Y}}$ in terms of ladder operators $\hat{S}^{\pm}_{j}=\hat{S}_{j}^{\mathrm{X}}\pm\mathrm{i}\hat{S}_{j}^{\mathrm{Y}}$ . Then, performing a Jordan-Wigner transformation by writing all spin operators in terms of fermionic creation and annihilation operators $\hat{f}_{j}^{\dagger}$ and $\hat{f}_{j}$ , namely,

[TABLE]

and $\hat{S}^{\mathrm{Z}}_{j}=\hat{n}_{j}-1/2$ , where $\hat{n}_{j}=\hat{f}_{j}^{\dagger}\hat{f}_{j}$ . Not only does this transformation provide an isomorphism between operators, but also preserves the bipartite entanglement, namely, the entanglement entropy associated to a region consisting of adjacent sites $(1,\,2,\ldots,j)$ is the same regardless of whether one uses the tensor product structure generated by the spin operators or by the fermionic creation and annihilation operators. This is a consequence of the remarkable fact that, despite the Jordan-Wigner transformation being nonlocal, the operator $\hat{S}^{\pm}_{j}$ only depends on the creation and annihilation operators $\hat{f}^{\dagger}_{l}$ and $\hat{f}_{l}$ in the range $1\leq l\leq j$ . The Hamiltonian expressed in terms of the fermionic operators is, up to a constant,

[TABLE]

where $\hat{P}=e^{\mathrm{i}\pi\hat{N}}$ is the parity operator, and $\hat{N}=\sum^{L}_{j=1}\hat{n}_{j}$ . The last term in Eq. (60) is a boundary term.

Because of the operator $\hat{P}$ , the Hamiltonian $\hat{H}_{\rm XY}$ in Eq. (60) is not quadratic in $\hat{f}_{j}^{\dagger}$ and $\hat{f}_{j}$ . The nonquadratic term containing $\hat{P}$ distinguishes the sectors with even and odd eigenvalues of the number operator $\hat{N}$ . To address this, the Hilbert space can be decomposed as direct sum $\mathcal{H}=\mathcal{H}^{+}\oplus\mathcal{H}^{-}$ , where $\mathcal{H}^{+}$ and $\mathcal{H}^{-}$ are the eigenspaces of the number parity operator $\hat{P}=e^{\mathrm{i}\pi\hat{N}}$ with eigenvalues $\pm 1$ . The projectors onto these eigenspaces are given by

[TABLE]

From here, we obtain exact eigenstates of $\hat{H}_{\mathrm{XY}}$ in both sectors. We diagonalize $\hat{H}_{\mathrm{XY}}$ over $\mathcal{H}^{\pm}$ individually by applying the Fourier transformations

[TABLE]

where $\kappa\in\mathcal{K}^{\pm}$ introduced in Eq. (5). The resulting Hamiltonian $\hat{H}_{\mathrm{XY}}=\hat{H}_{\mathrm{XY}}^{+}\hat{\mathcal{P}}^{+}+\hat{H}_{\mathrm{XY}}^{-}\hat{\mathcal{P}}^{-}$ is composed of the quadratic pieces

[TABLE]

with

[TABLE]

At this point, one only needs to perform individual fermionic two-mode squeezing transformations mixing the mode pairs $(\kappa,-\kappa)$ , with the form

[TABLE]

The transformation coefficients are

[TABLE]

As fermionic Bogoliubov coefficients, $u_{\kappa}$ and $v_{\kappa}$ satisfy $|u_{\kappa}|^{2}+|v_{\kappa}|^{2}=1$ . The cases $\kappa=0$ and $\pi$ , for which the Hamiltonian is already diagonal, are already included if one chooses the limits of $u_{\kappa}$ and $v_{\kappa}$ consistently. One has $u_{0}=0$ and $v_{0}=\mathrm{i}$ , while for $\kappa=\pi$ , one needs to distinguish

[TABLE]

(For $h=J$ , $\hat{c}_{\pi}$ and $\hat{c}_{\pi}^{\dagger}$ do not appear in the Hamiltonian.)

After performing the last transformation, the quadratic pieces take the diagonal form

[TABLE]

We use the eigenvalues of the quasiparticle number operators $\hat{n}_{\kappa}=\hat{\eta}_{\kappa}^{\dagger}\hat{\eta}_{\kappa}$ to write an orthonormal basis of eigenstates of $\hat{H}_{\mathrm{XY}}$ . The eigenstates $\ket{\{n_{\kappa}\}}$ satisfy

[TABLE]

with $n_{\kappa}\in\{0,1\}$ . Here, $\kappa\in\mathcal{K}^{\pm}$ depending on whether the sum $\sum_{\kappa\in\mathcal{K}^{\pm}}n_{\kappa}$ is even ( $+$ ) or odd ( $-$ ). For later use, we also define the variables

[TABLE]

which provide an equivalent characterization of eigenstates as the variables $n_{\kappa}$ . We emphasize that the quantity $N_{\kappa}$ is not related to the total number operator $\hat{N}=\sum^{L}_{j=1}\hat{f}_{j}^{\dagger}\hat{f}_{j}$ , which counts the total number of fermions.

The previous discussion shows that Hamiltonian (2), when written in terms of fermionic operators $\hat{f}_{j}$ and $\hat{f}^{\dagger}_{j}$ , becomes two quadratic fermionic Hamiltonians over direct addends of the fermionic Hilbert space in Eq. (3). Therefore, the overall Hamiltonian is almost quadratic, but not quite. One has two distinct notions of quasiparticles, described by $\hat{H}_{\mathrm{XY}}^{+}$ and $\hat{H}_{\mathrm{XY}}^{-}$ , depending on the number of quasiparticles. Namely, for an even number of quasiparticles one needs to use $\hat{H}_{\mathrm{XY}}^{+}$ while for an odd number of quasiparticles one needs to use $\hat{H}_{\mathrm{XY}}^{-}$ . In particular, the ground state, which has no quasiparticles, is given by the ground state of $\hat{H}_{\mathrm{XY}}^{+}$ , while the ground state of $\hat{H}_{\mathrm{XY}}^{-}$ is not an energy eigenstate of the system.

Appendix B Different averages over eigenstates

Here we study how the eigenstates of Hamiltonians $\hat{H}^{\pm}$ in Eq. (68), which are not eigenstates of Hamiltonian (60), change the average entanglement entropy. We define

[TABLE]

where $\langle S_{A}\rangle_{\pm}$ are the partial averages, as defined in Eq. (21), over eigenstates of the XY Hamiltonian (60) with only an even number of quasiparticles ( $+$ ) or only an odd number of quasiparticles ( $-$ ), respectively. On the other hand, $\langle S_{A}\rangle^{\mathrm{all}}_{\pm}$ are the averages over all eigenstates of Hamiltonians $\hat{H}^{\pm}$ in Eq. (68). The Hamiltonian $\hat{H}_{\mathrm{XY}}$ is not quadratic, but it can be written as the sum of two quadratic Hamiltonians $\hat{H}_{\mathrm{XY}}^{\pm}$ projected onto the even/odd sectors $\mathcal{H}^{\pm}$ of the Hilbert space. We are therefore exploring how the average eigenstate entanglement entropy of $\hat{H}_{\mathrm{XY}}$ compares to computing first the average entanglement entropy for $\hat{H}^{\pm}$ individually over all their eigenstates (rather than just over the even or odd sectors) and then averaging the two.

Results for $|\delta s_{\pm}(f)|$ are shown in Fig. 6 for $f=1/2$ (top panels) and $f=1/4$ bottom panels for different values of $h/J$ and $\gamma=1/2$ . Those results reveal that the differences between the correct and the incorrect calculations vanish in the thermodynamic limit. We can prove this observation analytically by studying which terms in our series expansion in Eq. (11) differ between the partial and the full average.

Let us consider a binomial distribution of $L$ independent random variables $(N_{1},\cdots,N_{L})$ that can take the values $\pm 1$ with equal probability, i.e., a Rademacher distribution. For simplicity, we still assume $L$ to be an even number. The $n$ -point correlation function $\langle N_{i_{1}}\cdots N_{i_{n}}\rangle$ is the one described in Eq. (22): it either vanishes, or it is $1$ iff each variable $N_{i}$ appears an even number of times. If we now restrict to the partial averages $\langle\dots\rangle_{\pm}$ over all possible values $N_{i}$ subject to the constraint $\prod^{L}_{i=1}N_{i}=\pm 1$ ( $L$ is still even), we find that the $n$ -point correlation function only changes in the very special case in which each random variable $N_{i}$ appears an odd number of times, in which case we find

[TABLE]

rather than $\langle N_{1}\dots N_{L}\rangle=0$ . If we consider the series expansion of the entanglement entropy $S_{A}$ from Eq. (11), the difference between the two averages is nonzero only for $n\geq L/2$ , i.e., we find

[TABLE]

This implies that $\lim_{L\to\infty}\delta s_{\pm}(f)=0$ , because the different averages agree to all order $n<L/2$ . With this, we have proved that the differences vanish in the thermodynamic limit, while for the exact exponential scaling we rely on the numerical results shown in Fig. 6.

Appendix C Computation of higher order bounds for noninteracting fermions

Here, we explain the method used to compute higher order bounds for noninteracting fermions [see Eq. (8)]. We also discuss the calculations of up to the fourth order bounds. They extend the results reported in Ref. vidmar_hackl_17 , where we computed the first and second order bounds.

C.1 Generalized Wick contractions and second order terms

In order to compute higher-order traces

[TABLE]

with $j(d_{i})=\frac{1}{L}\sum_{\kappa\in\mathcal{K}}e^{\mathrm{i}\kappa d_{i}}N_{\kappa}$ , we develop a systematic method of computing higher order correlation functions $\langle j(d_{1})\cdots j(d_{2n})\rangle$ . Technically, these are Fourier transformed correlation functions of the binomial distribution. They can be computed by adopting a strategy analogous to the one used for computing correlation functions of Gaussian distributions:

(i) The building blocks of the correlation functions are the so-called $2n$ -contractions given by

[TABLE]

where $\bar{\delta}(D)=1$ if $D=0\,(\mathrm{mod}\ L)$ and it is zero otherwise, i.e., $D=\sum^{2n}_{i=1}d_{i}$ is restricted to be an integer multiple of $L$ . The prefactors $\mathfrak{c}_{n}$ can be computed systematically as $\mathfrak{c}_{n}=L^{2n-1}\langle j(1)^{2n-1}j(1-2n)\rangle$ by using the explicit form $j(d_{i})=\frac{1}{L}\sum_{\kappa\in\mathcal{K}}e^{\mathrm{i}\kappa d_{i}}N_{\kappa}$ . Evaluating this formula for the first few cases, one realizes that the $\mathfrak{c}_{n}$ ’s are given by the tangent numbers with alternating signs, namely

[TABLE]

where $B_{2n}$ refers to the $2n$ -th Bernoulli number. Up to fourth order, we only need

[TABLE]

(ii) Once the $2n$ -contractions are known, one can compute a general correlation function as

[TABLE]

where each contraction consists of a product of different pairings, quadruplings, etc., of the $2n$ $j$ ’s. This is a generalization of Wick’s theorem because not only do we have $2$ -contractions, but also higher-order $2n$ -contractions.

To compute higher order traces up to order $n$ , we need to determine all prefactors up to $\mathfrak{c}_{n}$ . We can investigate the scaling in powers of $f=\frac{L_{A}}{L}$ , which appear in $\lim_{L\to\infty}\langle\mathrm{Tr}[\mathrm{i}\mathbb{J}]_{A}^{2n}\rangle/L_{A}$ . A general $2n$ correlation function is schematically given by

[TABLE]

To compute $\langle\mathrm{Tr}[\mathrm{i}\mathbb{J}]_{A}^{2n}\rangle/L_{A}$ , we sum over $x_{i}$ with $1\leq i\leq 2n$ , with $\sum^{2n}_{i=1}d_{i}$ automatically ensured to vanish, which implies that there is one redundant delta in each addend of Eq. (79).

The $2n$ -correlation functions consists of various products of contractions. We refer to such a type of product of contractions as $(2n_{1},\cdots,2n_{l})$ -contraction if it consists of $l$ contractions of types determined by $2n_{i}$ . Moreover, we refer to the full set or sum of $(2n_{1},\cdots,2n_{l})$ -contractions as $\mathcal{C}_{(2n_{1},\cdots,2n_{l})}$ . A specific contraction $\mathfrak{C}\in\mathcal{C}_{(2n_{1},\cdots,2n_{l})}$ gives $l$ delta functions $\bar{\delta}(D_{i})$ with $D_{i}$ being specific sums of the $d_{i}$ ’s according to the chosen contraction, namely,

[TABLE]

Provided that $D_{i}\leq L$ , we can ignore the $(\mathrm{mod}\ L)$ condition in $\bar{\delta}$ and replace it with a regular $\delta$ . This product of deltas appearing in the sum (74) then gives

[TABLE]

where we use the fact that each of the $l$ deltas reduces the dimension of the sum by one, except for the one that is redundant because $\sum^{l}_{i=1}D_{i}=\sum{2n}_{i=1}d_{i}=0$ . Therefore, the sum gives a polynomial of degree $2n-l+1$ . Applying this result to all contractions, we find

[TABLE]

as $f\to 0$ , where ${\cal O}(L_{A}^{n+1})$ refers to a polynomial of degree $n+1$ . Note that we included the requirement that all $D_{i}<L$ by considering this expression in the limit $f\to 0$ . Therefore, we conclude that

[TABLE]

i.e., that it is a polynomial function containing only powers from $f^{n}$ up to $f^{2n-1}$ . This result is instrumental in the analysis of the behavior of the entropy for vanishing subsystem fraction $f\to 0$ . It implies that it is sufficient to compute the expansion through order $n$ to get the dominant terms through order $f^{n}$ . However, note that in general we expect additional corrections due to the $(\mathrm{mod}\ L)$ condition in $\bar{\delta}$ for $f>\frac{1}{n}$ . Interestingly, these corrections are absent in the relevant interval $f\in[0,\frac{1}{2}]$ for $n\leq 3$ , but for $n=4$ we will see that additional terms do appear for $f>\frac{1}{4}$ .

Second order trace. To illustrate this method, let us apply it to the second order correction containing $\langle\mathrm{Tr}[\mathrm{i}\mathbb{J}]_{A}^{4}\rangle$ . We need to compute the $4$ -point correlation function $\langle j(d_{1})j(d_{2})j(d_{3})j(d_{4})\rangle$ . We find three $2$ -contractions and one $4$ -contraction given by

[TABLE]

Plugging this expression into the sum (74) for the trace gives

[TABLE]

as already found in Ref. vidmar_hackl_17 .333Note that there was a missing factor of $2$ in the supplemental material of Ref. vidmar_hackl_17 . It then follows that

[TABLE]

where the first term comes from the sum over the $2$ -contractions, while the second term comes from the $4$ -contraction for which $\bar{\delta}(\sum^{4}_{i=1}d_{i})=1$ due to $\sum^{4}_{i=1}d_{i}=(x_{1}-x_{2})+\cdots+(x_{4}-x_{1})=0$ . Note that, in Eq. (86), we assume $L_{A}\leq L/2$ to avoid the complication that $\bar{\delta}(d_{i_{1}}+d_{i_{2}})$ is also nonzero for $d_{i_{1}}+d_{i_{2}}=0\,(\mathrm{mod}\ L)$ .

C.2 Wick diagrams and third order terms

When computing higher order correlation functions, we need to sum a larger number of contractions and it is important to keep track of all cases. Contractions can be efficiently organized by taking symmetry into account. The evaluation of a contraction in which one shifts all $d_{i}$ ’s by an overall constant does not change the result. It is important to emphasize that the $d_{i}$ ’s should be thought of as arranged in a circle. Hence, representing the equivalence class of contractions under translations is best done on a circle. We refer to such a pictorial representation as a Wick diagram

[TABLE]

where we connect the $2n$ dots associated to the $d_{i}$ ’s accordingly to get a specific Wick diagram. Each diagram $\mathfrak{C}\in\mathcal{C}_{(n_{1},\cdots,n_{l})}$ comes with some multiplicity $\#\mathfrak{C}$ associated to translations that do not lead to the same diagram. We refer to the underlying product of Kronecker deltas in a given diagram $\mathfrak{C}$ as

[TABLE]

where we ignore any contributions from $(\mathrm{mod}\ L)$ in $\bar{\delta}(D_{i})$ . In the thermodynamic limit, we can evaluate the sum over $\delta_{\mathfrak{C}}$ as an integral

[TABLE]

by formally changing from Kronecker deltas to delta distributions. Note that we removed one of the Kronecker deltas due to the redundancy mentioned earlier, to avoid having an ill defined $\delta(0)$ . Moreover, we define the dimensionless quantity

[TABLE]

in which we remove the dependence on $L_{A}$ by making the replacement $x_{i}\to x_{i}\,L_{A}$ , which only affects the range of integration. As long as we can ignore the $(\mathrm{mod}\ L)$ condition in $\bar{\delta}$ , i.e., for $f\in[0,\frac{1}{n}]$ , we have

[TABLE]

Finally, the summation over all diagrams within a class $\mathcal{C}_{(n_{1},\cdots,n_{l})}$ is given by

[TABLE]

Representing contractions by Wick diagrams, and replacing sums by integrals, will prove to be very efficient to evaluate the third and fourth order traces in the thermodynamic limit. Let us first reconsider the simpler example of the second order trace.

Second order trace. When we computed the terms of order $2$ , we had four generalized Wick contractions. These correspond to only three Wick diagrams, namely

[TABLE]

because the second diagram has multiplicity $2$ . We can compute their prefactors as

[TABLE]

of which only the first diagram requires some work. An efficient way of evaluating the integrals involved is to introduce a dummy variable $y$ , leading to

[TABLE]

where we defined the piecewise polynomial function $I_{2}(y)$ over the intervals $[0,1]$ and $[1,2]$ . Distinguishing between the two cases where $y$ is in either of the two intervals gives

[TABLE]

Using this result, it is easy to compute $\frac{2}{3}$ as the prefactor. When adding up the contributions of the three diagrams, taking their multiplicity $\#\mathfrak{C}$ , $\mathfrak{c}_{1}^{2}=1$ , and the overall prefactor of $2$ into account, we find

[TABLE]

which reproduces the result in Eq. (86).

Third order trace. We can apply the representation of contractions as Wick diagrams to compute the third order terms. The main ingredients in the trace are

[TABLE]

•

$(2,2,2)$ **-contraction: $\int\mathcal{C}_{(2,2,2)}=11\,\mathfrak{c}_{1}^{3}\,f^{3}=11\,f^{3}$

**The $15$ distinct contractions fall into $5$ different symmetry classes with multiplicities ranging from one to six. We list a representative $\mathfrak{C}$ of each class in table 2, and note its multiplicity $\#\mathfrak{C}$ , and the prefactor $\int\mathfrak{C}$ resulting from the integration over the subsystem.

•

$(2,4)$ **-contraction: $\int\mathcal{C}_{(2,4)}=12\,\mathfrak{c}_{1}\mathfrak{c}_{2}\,f^{4}=-24\,f^{4}$

**Due to the redundancy in Kronecker deltas, we can ignore the delta from the $4$ -contraction and only focus on the $2$ -contraction. There are $15$ distinct $2$ -contractions, of which $6$ contract adjacent points, while the remaining $9$ do not. The prefactor is $1$ for the adjacent ones and $\frac{2}{3}$ for the remaining ones. We thus find $\int\mathcal{C}_{(2,4)}=(6\cdot 1+9\cdot\frac{2}{3})\,\mathfrak{c}_{1}\mathfrak{c}_{2}\,f^{4}=-24f^{4}$ .

•

$(6)$ **-contraction: $\int\mathcal{C}_{(6)}=\mathfrak{c}_{3}\,f^{5}=16\,f^{5}$

**This integration is trivial because the only contraction leads to a redundant delta leading to the prefactor $1$ .

Combining all the terms computed above, and including the prefactor of $2$ , we find

[TABLE]

C.3 Delta terms, fourth order terms, and

non-analytic contributions

An important subtlety arises starting with the fourth order bound, because it matters that our generalized Kronecker deltas are defined $(\mathrm{mod}\ L)$ . This means that we have $\delta(D)=1$ not only for $D=0$ , but also for any multiple of the full system size $L$ . We will see that this leads to surprising effects and turns the bounds into piecewise defined polynomials rather than a single analytic polynomial for $f\in[0,\tfrac{1}{2}]$ . Such a correction is expected to appear for $f>\frac{1}{n}$ for traces $t_{n}$ of order $n$ . This subtlety comes from the fact that our generalized Kronecker deltas are defined as

[TABLE]

where $D=\sum^{2j}_{i=1}d_{l_{i}}$ with $d_{l}=x_{l}-x_{l+1}$ and $x_{L+1}=x_{1}$ . Up to third order, the $(\mathrm{mod}\ L)$ condition played no role. The explanation for this consists of two parts:

(i) For a $2$ -contraction $\contraction{}{j}{(d_{l})}{j}j(d_{l})j(d_{m})=\frac{\mathfrak{c}_{1}}{L}\delta(x_{l}-x_{l+1}+x_{m}-x_{m+1})$ , the largest possible value of $D=d_{l}+d_{m}$ is given by $D_{\max}=2(L_{A}-1)$ , namely for $x_{l}=x_{m}=L_{A}$ and $x_{l+1}=x_{m+1}=1$ . However, we only compute the trace for $L_{A}\leq L/2$ leading to $D_{\max}\leq L-2<L$ . Thus for any $2$ -contraction, we can just ignore the $(\mathrm{mod}\ L)$ condition altogether.

(ii) Based on the first statement, we only need to take the $(\mathrm{mod}\ L)$ condition into account for $\delta(\sum^{2j}_{i=1}d_{l_{i}})$ when $j>2$ . Such generalized Kronecker deltas only appear for $4$ -contractions and higher orders. Moreover, we should recall that at bound order $n$ , we always have $\sum^{2n}_{i=1}d_{i}=0$ , which implies that we can always choose one generalized Kronecker delta to be redundant. In the second and third order bounds, we do have terms that contain $4$ - and $6$ -contractions, but they always appear alone or in a product with only $2$ -contractions. In those cases, we can always choose the $4$ - or $6$ -contraction to be the redundant one and we are left with $2$ -contractions, discussed in point (i) above.

To understand how to take the $(\mathrm{mod}\ L)$ condition into account, as well as its effects, let us consider the fourth order trace.

Fourth order trace. The fourth order bound is computed from

[TABLE]

Using again the generalized Wick’s theorem, the structure of the $8$ -point correlation function is given by

[TABLE]

where we grouped the $(4,4)$ -contractions together with the $(2,6)$ -contractions, because they contribute a power $L_{A}^{7}$ . Using straightforward combinatorics, one can count how many contractions there are per type,

[TABLE]

leading to a total of $414$ distinct contractions. As before, many contractions are related by translational symmetry, so it is sufficient to compute one representative per symmetry class.

•

$(2,2,2,2)$ **-contractions: $\int\mathcal{C}_{(2,2,2,2)}=\frac{908}{15}\,f^{4}$

**The $105$ distinct contractions fall into $18$ different symmetry classes, with multiplicities ranging from one to eight. We list a representative $\mathfrak{C}$ of each class in table 3 and note both, its multiplicity $\#\mathfrak{C}$ and its prefactor $\int\mathfrak{C}$ under integration over the subsystem.

•

$(2,2,4)$ **-contractions: $\int\mathcal{C}_{(2,2,4)}=-\frac{1296}{5}\,f^{5}$

**The $210$ distinct contractions fall into $29$ different symmetry classes, with multiplicities ranging from one to eight. We list a representative $\mathfrak{C}$ of each class in table 4 and note both, its multiplicity $\#\mathfrak{C}$ and its prefactor $\int\mathfrak{C}$ under integration over the subsystem.

•

$(4,4)$ **-contractions: $\int\mathcal{C}_{(4,4)}=\frac{5816}{63}\,f^{6}$

**The $70$ distinct contractions fall into $7$ different symmetry classes, with multiplicities ranging from one to eight. We list a representative $\mathfrak{C}$ of each class in table 5 and note both, its multiplicity $\#\mathfrak{C}$ and its prefactor $\int\mathfrak{C}$ under integration over the subsystem.

•

$(2,6)$ **-contractions: $\int\mathcal{C}_{(2,6)}=\frac{1024}{3}\,f^{6}$

**Similar to $\mathcal{C}_{(2,4)}$ , we can ignore the delta from the $6$ -contraction and only focus on the $2$ -contraction. There are $8\cdot 7/2=28$ distinct $2$ -contractions, of which $8$ contract adjacent points, while the remaining $20$ do not. As before, the prefactor for adjacent contractions is $1$ , and it is $\frac{2}{3}$ for the remaining ones. Overall, we thus find $\int\mathcal{C}_{(2,4)}=(8\cdot 1+20\cdot\frac{2}{3})\,\mathfrak{c}_{1}\mathfrak{c}_{3}\,f^{6}=\frac{1024}{3}\,f^{6}$ .

•

$(8)$ **-contractions: $\mathfrak{c}_{4}=\int\mathcal{C}_{(8)}=\mathfrak{c}_{4}=-272\,f^{7}$

**This integration is trivial because the only contraction gives a redundant delta leading to the prefactor $1$ .

Putting all terms together results in the fourth order trace

[TABLE]

However, Eq. (106) only holds for $r\leq\frac{1}{4}$ . The reason for this is exactly the subtlety mentioned before, namely, there are diagrams for which $\bar{\delta}(D_{i})$ actually gives contributions due to the $(\mathrm{mod}\ L)$ condition. Clearly, only $(4,4)$ -contractions can give such a contribution, because all other contractions only consist of $2$ -contractions plus $4$ - or $6$ -contractions that we can choose to be the redundant ones. Going through the list of $(4,4)$ -contractions in table 5, we find the following contractions that lead to $\bar{\delta}(D_{i})$ with $D_{i}=\pm L$ for $f<\frac{1}{2}$ :

[TABLE]

With our computation of $\int\mathfrak{C}$ , as listed in table 5, we already covered the case where $D_{i}=0$ . However, for the Wick diagrams in Eq. (107), we also need to compute the contributions from

[TABLE]

where $I_{4}(y)=\int_{0}^{1}d^{4}x\,\delta(x_{1}+x_{2}+x_{3}+x_{4}-y)$ . Here, we use again the trick of introducing a dummy variable $y$ and rescaling by $L_{A}$ , such that the $(\mathrm{mod}\ L)$ condition becomes $\mod{\frac{L}{L_{A}}}=\mod{f^{-1}}$ . Computing the function $I_{4}(f)$ requires some work because one needs to differentiate between several cases. Ultimately, we find a piecewise defined polynomial, with compact support on the interval $[0,4]$ , given by

[TABLE]

with $I^{(1)}_{4}(y)=\frac{y^{3}}{6}$ , and $I^{(2)}(y)=\frac{1}{6}(-3y^{3}+12y^{2}-12y+4)$ .

One can deal with the Wick diagrams in Eqs. (108) and (109) in an analogous way, we find

[TABLE]

where the integral of both diagrams can be identified by renaming the variables $x_{i}$ . We introduced $I_{3}(y)=\int_{0}^{1}d^{4}x\,\delta(x_{1}+x_{2}+x_{3}-y)$ , which is given by

[TABLE]

with $I^{(1)}_{3}(y)=\frac{y^{2}}{2}$ and $I^{(2)}_{3}(y)=-y^{2}+3y-\frac{3}{2}$ .

The integrands in Eqs. (110) and (112) are therefore given by the product of two piecewise defined polynomials, as illustrated in Fig. 7. The integrals can be computed analytically in each interval, and are given by the piecewise defined functions:

[TABLE]

For higher order $n$ , the non-analytic part of $t_{n}$ will contain integrals of the form

[TABLE]

The relevant functions $I_{n}$ will be compactly supported in the region $[0,n]$ , such that the integral in Eq. (121) will vanish for $mf^{-1}>n$ , i.e., $f<m/n$ . The smallest possible value is $m=1$ . Thus, we only expect non-analytical terms in $t_{n}$ for $f>\frac{1}{n}$ .

To compute the correction for $t_{4}(f)$ , let us recall that each of these terms comes with an additional multiplicity of $2$ due to the two different choices $\pm L$ in the $(\mathrm{mod}\ L)$ condition, which can be translated into two distinct choices of $\pm f^{-1}$ in our rescaled variables. Furthermore, we need to take the prefactors $\mathfrak{c}_{2}^{2}=4$ and the diagram multiplicities

[TABLE]

into account. With this in hand, we find the correction term $t^{\mathrm{cor}}_{4}(f)$ introduced in Eq. (55) of the main text:

[TABLE]

Note that the trace formula (74) contains an additional factor of $2$ , and we need to include again the power $f^{7}$ , which we separated from the integral by rescaling $x_{i}\to x_{i}\,L_{A}$ .

Bibliography68

The reference list from the paper itself. Each links out to its DOI / PubMed record.

1(1) A. Einstein, B. Podolsky, and N. Rosen, Can quantum-mechanical description of physical reality be considered complete?, Phys. Rev. 47 , 777 ( 1935 ). · doi ↗
2(2) J. S. Bell, Einstein-podolsky-rosen experiments, John S Bell on the Foundations of Quantum Mechanics , 74–83 (World Scientific, 2001).
3(3) R. D. Sorkin, On the Entropy of the Vacuum outside a Horizon, Proceedings, 10th International Conference on General Relativity and Gravitation: Padua, Italy, July 4-9, 1983 , volume 2, 734–736 (1984). eprint 1402.3589.
4(4) L. Bombelli, R. K. Koul, J. Lee, and R. D. Sorkin, A Quantum Source of Entropy for Black Holes, Phys. Rev. D 34 , 373 ( 1986 ). · doi ↗
5(5) M. Srednicki, Entropy and area, Phys. Rev. Lett. 71 , 666 ( 1993 ). · doi ↗
6(6) C. Holzhey, F. Larsen, and F. Wilczek, Geometric and renormalized entropy in conformal field theory, Nuc. Phys. B 424 , 443 ( 1994 ). · doi ↗
7(7) S. Ryu and T. Takayanagi, Holographic derivation of entanglement entropy from Ad S/CFT, Phys. Rev. Lett. 96 , 181602 ( 2006 ). · doi ↗
8(8) A. Lewkowycz and J. Maldacena, Generalized gravitational entropy, J. High Energ. Phys. 08 , 090 ( 2013 ). · doi ↗