Entanglement entropy of the long-range Dyson hierarchical model

TL;DR

This paper investigates the entanglement entropy of a long-range quantum spin chain with hierarchical interactions, revealing an unexpected area law at criticality, contrasting with typical long-range models.

Contribution

It introduces a real-space renormalisation group approach to analyze entanglement in a hierarchical long-range model, uncovering novel area law behavior at criticality.

Findings

Entanglement entropy follows an area law at criticality.

Analytical results for large and small decay exponents $\sigma$.

Numerical validation of recursive equations.

Abstract

We study the ground state entanglement entropy of the quantum Dyson hierarchical spin chain in which the interaction decays algebraically with the distance as $r^{-1-\sigma}$. We exploit the real-space renormalisation group solution which gives the ground-state wave function in the form of a tree tensor network and provides a manageable recursive expression for the reduced density matrix of the renormalised ground state. Surprisingly, we find that at criticality the entanglement entropy obeys an area law, as opposite to the logarithmic scaling of short-range critical systems and of other non-hierarchical long-range models. We provide also some analytical results in the limit of large and small $\sigma$ that are tested against the numerical solution of the recursive equations.

Equations

$$S_A=-\text{Tr}[{\hat{\rho }}_A\ln {\hat{\rho }}_A].$$

$$\hat{H}$$

$$\hat{H}=\sum _{k=0}^{n-1}{\hat{H}}^{(k)},$$

$${\hat{H}}^{(k=0)}=-h\sum _{i=1}^L{\hat{\sigma }}_i^z-\sum _{i=1}^{L/2}{J}_0[\mu {\hat{\sigma }}_{2i-1}^x][\mu {\hat{\sigma }}_{2i}^x]=\sum _i^{L/2}{\hat{H}}_{2i}^{(0)}.$$

$$|{\lambda }_{\text{s/a}}^{+}{⟩}_{2i}=a_1|+{⟩}_{2i-1}|+{⟩}_{2i}\pm b_1|-{⟩}_{2i-1}|-{⟩}_{2i},{\lambda }_{\text{s/a}}^{+}=\mp \sqrt{(J_0{\mu }^2{)}^2+4h^2},$$

$$|{\lambda }_{\text{s/a}}^{-}{⟩}_{2i}=\frac{1}{\sqrt{2}}|-{⟩}_{2i-1}|+{⟩}_{2i}\pm \frac{1}{\sqrt{2}}|+{⟩}_{2i-1}|-{⟩}_{2i},{\lambda }_{\text{s/a}}^{-}=\mp J_0{\mu }^2,$$

$$|\pm {⟩}_{2i}^{[1]}\equiv |{\lambda }_{\text{s}}^{\pm }{⟩}_{2i}.$$

$${\hat{P}}_1\equiv \prod _i{\hat{P}}_{1,2i}=\prod _i|+⟩⟨+{|}_{2i}^{[1]}+|-⟩⟨-{|}_{2i}^{[1]},$$

$${\hat{H}}_R^{(0)}+{\hat{H}}_R^{(1)}=-\frac{L}{2}{e}_1{\hat{P}}_1-h_1\sum _i{\hat{\sigma }}_{2i}^{z[1]}-\sum _{i=1}^{L/4}{J}^{(1)}[{\mu }_1{\hat{\sigma }}_{4i-2}^{x[1]}][{\mu }_1{\hat{\sigma }}_{4i}^{x[1]}],$$

$$K_k\equiv \frac{J_k{\mu }_k^2}{h_k},\text{with}{K}_0=\frac{J_0{\mu }^2}{h}=\frac{J}{h}.$$

$$h_1$$

$$e_1$$

$${\mu }_1$$

$$|+{⟩}_i^{[k]}$$

$$|-{⟩}_i^{[k]}$$

$$a_k=\frac{1}{\sqrt{2}}{\left(1+\frac{2}{\sqrt{K_k^2+4}}\right)}^{1/2},b_k=\frac{1}{\sqrt{2}}{\left(1-\frac{2}{\sqrt{K_k^2+4}}\right)}^{1/2},$$

$$e_k=-h_k\frac{K_k+\sqrt{K_k^2+4}}{2}.$$

$$h_{k+1}=\frac{2h_k}{K_k+\sqrt{K_k^2+4}},{\mu }_{k+1}={\mu }_k\sqrt{2}\sqrt{1+\frac{K_k}{\sqrt{K_k^2+4}}},$$

$$K_k=\frac{K_{k-1}}{2^{1+\sigma }}\frac{{\left(K_{k-1}+\sqrt{K_{k-1}^2+4}\right)}^2}{\sqrt{K_{k-1}^2+4}}\equiv R[K_{k-1}].$$

$${\hat{H}}_R=-\sum _{k=1}^n\frac{L}{2^k}{e}_k{\hat{P}}_k-h_n{\hat{\sigma }}^{z[n]},$$

$$|{\psi }_{\text{GS}}⟩=|+{⟩}^{[n]},$$

$$|{\psi }_{\text{GS}}⟩=|+{⟩}^{[n]}=a_{n-1}|++{⟩}_1^{[n-1]}+b_{n-1}|--{⟩}_1^{[n-1]},$$

$$⟨{\hat{\sigma }}_i^x{\hat{\sigma }}_{i+r}^x⟩\sim \frac{1}{r^{2x_{\sigma }}},$$

$$x_{\sigma }=\frac{1-\sigma }{4}+\frac{\ln (K_c^2+4)}{8\ln 2}.$$

$${\hat{\tau }}_0^{[k]}\equiv {\text{tr}}_B|+⟩⟨+{|}^{[n]},{\hat{\tau }}_1^{[k]}\equiv {\text{tr}}_B|-⟩⟨-{|}^{[n]},$$

$${\hat{\tau }}_2^{[k]}\equiv {\text{tr}}_B|+⟩⟨-{|}^{[n]},{\hat{\tau }}_3^{[k]}\equiv {\text{tr}}_B|-⟩⟨+{|}^{[n]},$$

$${\hat{\tau }}_0^{[k]}=\left(\begin{array}{cc}{c}_k^{+}& 0\\ 0& d_k^{+}\end{array}\right),{\hat{\tau }}_1^{[k]}=\left(\begin{array}{cc}{c}_k^{-}& 0\\ 0& d_k^{-}\end{array}\right),{\hat{\tau }}_2^{[k]}=\left(\begin{array}{cc}0& e_k\\ f_k& 0\end{array}\right),{\hat{\tau }}_3^{[k]}={\left[{\hat{\tau }}_2^{[k]}\right]}^{\text{T}},$$

$$\{\begin{array}{l}{c}_k^{\pm }=a^2{c}_{k-1}^{\pm }+\frac{1}{2}{d}_{k-1}^{\pm }\\ d_k^{\pm }=b^2{c}_{k-1}^{\pm }+\frac{1}{2}{d}_{k-1}^{\pm },\end{array}\{\begin{array}{l}{e}_k=\frac{a}{\sqrt{2}}{e}_{k-1}+\frac{b}{\sqrt{2}}{f}_{k-1}\\ f_k=\frac{b}{\sqrt{2}}{e}_{k-1}+\frac{a}{\sqrt{2}}{f}_{k-1},\end{array}$$

$$\{\begin{array}{l}{c}_k^{+}=\frac{1+2b^2{e}^{-k\alpha }}{1+2b^2}\\ d_k^{+}=2b^2\frac{1-e^{-k\alpha }}{1+2b^2}\end{array},\{\begin{array}{l}{c}_k^{-}=\frac{1-e^{-k\alpha }}{1+2b^2}\\ d_k^{-}=\frac{2b^2+e^{-k\alpha }}{1+2b^2}\end{array},\{\begin{array}{l}{e}_k=\frac{1}{2}\left(e^{-\beta k}+e^{-\gamma k}\right)\\ f_k=\frac{1}{2}\left(e^{-\beta k}-e^{-\gamma k}\right)\end{array},$$

$$\alpha \equiv \ln \sqrt{K^2+4},\beta \equiv \ln \frac{\sqrt{2}}{a-b},\gamma \equiv \ln \frac{\sqrt{2}}{a+b}.$$

$$S_A=-c_k^{+}\ln c_k^{+}-d_k^{+}\ln d_k^{+}.$$

$$S_A=S_{\ast }+c_{\ast }{e}^{-k\alpha }+O(e^{-2k\alpha }),$$

$$S^{\ast }=\ln (1+2b^2)-\frac{2b^2}{(1+2b^2)}\ln 2b^2,c^{\ast }=\frac{2b^2\ln 2b^2}{1+2b^2}.$$

Full text

1 Introduction

Long-range interacting systems, because of their unconventional static and dynamical properties, have been extensively studied in classical statistical mechanics since many decades [1]. In the quantum domain, atomic, molecular and optical systems (AMO), as well as synthetic ones, are often described by long-range interactions. Yet, they have been a focus of great attention only relatively recently, as a consequence of their experimental implementation in different physical setups of many-body AMO systems, such as Rydberg atoms [2, 3], dipolar molecules [4, 5] or trapped ions [6, 7, 8, 9]. In particular trapped ions, which interact via ferromagnetic Ising long-range interactions [10, 11], represent one of the most promising platforms for quantum computation. These experimental advances motivated an increasing interest towards the understanding of quantum systems with long-range interactions both in thermodynamic equilibrium [12, 16, 13, 14, 15, 18, 17, 19, 30, 25, 29, 21, 22, 23, 24, 20, 26, 27, 28, 32, 31] and in several non-equilibrium situations [33, 34, 36, 37, 39, 38, 40, 41, 35, 42, 43, 44, 45, 46, 48, 50, 49, 52, 53, 54, 51, 47, 55]. Anyhow, there are still only few results concerning the equilibrium properties of quantum spin chains with ferromagnetic power-law decaying interactions, see e.g. Refs. [13, 14, 15] for some mean-field and perturbative solutions and Refs. [49, 56, 57, 58] for some numerical simulations. Even for the long-range ferromagnetic Ising model, the precise knowledge of the phase-diagram and universal features at the quantum phase transition are still lacking.

At the same time, it is by now well established that a large amount of information about many-body systems can be inferred by their entanglement properties [59, 60, 61]. The most useful quantity is surely the celebrated entanglement entropy, defined as follows. Let us consider a quantum system in a pure state with density matrix $\hat{\rho}=|\psi\rangle\langle\psi|$ and a bipartition of the Hilbert space $\mathcal{H}=\mathcal{H}_{A}\otimes\mathcal{H}_{B}$, the entanglement between subsystems $A$ and $B$ may be measured by the Von Neumann entropy of the reduced density matrix $\hat{\rho}_{A}=\text{Tr}_{B}{\hat{\rho}}$, i.e.

[TABLE]

The entanglement entropy turned out to be a very useful quantity, especially when $A$ is a block of $\ell$ contiguous lattice sites in a one dimensional chain. In fact, the scaling entanglement entropy is tightly connected to the presence of quantum phase transitions for short-range systems. Indeed, while for gapped systems with short-range interaction, $S_{\ell}$ obeys an area law [62, 63, 64] (i.e. it scales with surface of the boundary between the two partitions), it grows like $\ln\ell$ at critical point, with a proportionality constant given by the central charge of the underlying conformal field theory [65, 66, 67, 68].

However, for long-range quantum interactions conformal invariance is broken as a consequence of the lack of Lorentz symmetry (rotational one in Euclidean spacetime) as explicitly manifested by a dynamical critical exponent $z\neq 1$ [69]. Hence unusual ground state entanglement properties are expected [17, 18, 19, 21, 22, 23, 24, 20], as e.g. shown numerically for antiferromagnetic long-range one-dimensional systems, which violate area-law scaling also in the gapped phase [17, 26, 27, 28]. On the other hand, the ferromagnetic mean-field Ising model has a critical entanglement entropy that scales with the logarithmic of the number of spins in $A$ [70, 71]. Beyond mean field, the scaling of entanglement entropy for ferromagnetic long-range interaction is still an open issue, also in one spatial dimension.

In this work, we address such a problem, by considering the quantum version of the ferromagnetic Dyson hierarchical model [72, 73]. Its classical version was introduced by Dyson in 1969 [74, 75] and it provided, via an exact renormalisation group (RG) approach, analytical insights about the critical behaviour of long-range interacting one-dimensional spin chains. Also the quantum version has been recently solved by Monthus [72] using real-space renormalisation and the critical properties of the ground state have been obtained by a recursive projection onto the low lying energy states.

Within this frame, we study the entanglement entropy of the renormalised ground state of the Dyson hierarchical model at criticality. The result depends on the choice of the subsystem because the model is not translationally invariant. Hence, we consider three types of partitions and, for each of them, we determine recursively the reduced density matrix. We finally analyse the scaling of the entanglement entropy at the critical point and in some limits we derive analytical expansions. Our main and surprising result is that, although the correlation functions decay algebraically, the entanglement entropy obeys an area law. This is due to the hierarchical structure of the renormalised ground state that makes the rank of the reduced density matrix finite (the ground state is a tree tensor network [76, 77, 78, 79] with finite bond dimension). In spite of this saturation at the critical point, the renormalised ground-state reduced density matrix reproduces the power-law decay of the two-point correlation function with the exact critical exponent.

The paper is organised as follows. In Sec. 2, we describe the model and the main steps of the real-space RG procedure. Sec. 3 contains all our results, including a description of the three considered partitions and a first discussion of our findings. This section is divided as: in Subsecs. 3.1-3.3 we derive recursively the reduced density matrix and the entanglement entropy for the three considered partitions; we then show by an elementary approach that the rank of the reduced density matrix is finite (Subsec. 3.4) and that the correct power-law decay of the correlation functions follows from them (Subsec. 3.5). In the closing section we draw our conclusions and discuss some open questions. In the appendix we report some technical details of the calculations.

2 The model

The Dyson hierarchical model was introduced by Dyson in 1969 [74, 75] with the aim of studying phase transitions in one dimensional ferromagnetic Ising models with algebraically decaying interactions with exponent $1+\sigma$. A hierarchical Hamiltonian is defined with a tree structure on a spin chain of length $L=2^{n}$, where blocks of spins are subject to an interaction that decays as a power-law of their distance. In the classical case, the hierarchical model reproduces well the critical properties of the long-range Ising Hamiltonian [80]. In fact, for $\sigma<1/2$ the critical exponents of the former coincide with the mean-field of the latter. Furthermore, for $\sigma\geq 1/2$, the exponents of the two models coincide at first order in the $\epsilon=\sigma-1/2$ expansion [81].

The quantum version of this model has been introduced and studied with real-space RG by Monthus in [72, 73]. Let us describe the Hamiltonian which is pictorially represented in Fig. 1. At step $k=0$ each pair of spins at sites $2j-1$ and $2j$ interacts with a coupling $J_{0}=J$ as $J\hat{\sigma}_{2j-1i}^{x}\hat{\sigma}_{2j}^{x}$, where hereafter $\hat{\sigma}_{i}^{\alpha}$ denotes the Pauli matrices. The absence of interactions between the sites $2j$ and $2j+1$ is the main difference from standard Ising models. At the following step $k=1$, two blocks of two spins interact (always along $x$) with coupling $J_{1}=J/2^{1+\sigma}$, then at $k=2$ blocks of four spins interact with $J_{2}=J/2^{2(1+\sigma)}$ and so on. Hence, the total hamiltonian is

[TABLE]

where $J_{k}={J}{2^{-k(1+\sigma)}}$ and we introduced the transverse field $h>0$ and the magnetic moment $\mu$ (the latter will be set to 1, and it has been introduced for convenience in the renormalisation procedure). The real-space RG of Ref. [72] is constructed by projecting (at each step $k$) the Hamiltonian $\hat{H}^{(k)}$ onto the two lower energy states. Hence at each RG step, we have half of the spins. The renormalised transverse field $h_{k}$ and magnetic moment $\mu_{k}$ are obtained within this procedure.

2.1 Real-space RG for the Dyson hierarchical model

We report here the main ideas and the conceptual steps of the real-space RG approach. For all the details we refer to [72]. We rewrite the Hamiltonian (2) as a sum over the generations $k=0,1,..,n-1$

[TABLE]

where $\hat{H}^{(k)}$ for $k>1$ is read off Eq. (2), and in particular $k=0$ is

[TABLE]

The Hamiltonians $\hat{H}^{(0)}_{2i}$ are unrelated and each pair of spins may be diagonalised independently. It is convenient to work in the basis of eigenstates of $\hat{\sigma}^{z}_{i}$, i.e. $\hat{\sigma}^{z}_{i}|\pm\rangle_{2i}=\pm|\pm\rangle_{2i}$. Then, the four eigenstates and eigenvalues of $\hat{H}^{(0)}_{2i}$ are

[TABLE]

where “s/a” stands for the symmetric and anti-symmetric combination respectively and the coefficients are ${{a_{1}=\sqrt{1/2+{h}/\sqrt{(J_{0}\mu)^{2}+4h^{2}}}}}$, $b_{1}=\sqrt{1-a_{1}^{2}}$. In the RG approach, we keep only the two lowest energy states $|\lambda^{\pm}_{\text{s}}\rangle_{2i}$. Thus the RG transformation is the projection onto the symmetric subspace. The renormalised spins at level $k=1$ are then identified by the basis

[TABLE]

They are related to the physical spins by Eq. (5) and are pictorially represented in Fig. 1 by the green boxes. Hereafter, the superscript $[k]$ always refers to a renormalised spin at level $k$.

Eq. (6) allows us to define the projection operator

[TABLE]

which implements the RG transformation and selects the two lowest energy states of $\hat{H}^{(0)}$. When applied to $\hat{H}^{(0)}+\hat{H}^{(1)}$ in Eqs. (3-4), the projection gives the renormalised Hamiltonians

[TABLE]

where $h_{1}$, $\mu_{1}$ are the renormalised magnetic field and moment while $e_{1}$ is the contribution of the previous generation to the ground state energy. It is then clear that the RG flow should be written in terms of the control parameter

[TABLE]

Hence, at level $k=1$ the renormalised parameter are

[TABLE]

In Eq. (8) the operators $\sigma^{z/x\,[1]}_{2i}$ are Pauli matrices acting on the renormalised spins at level $k=1$ as $\sigma^{z\,[1]}_{2i}|\pm\rangle^{[1]}_{2i}=\pm|\pm\rangle^{[1]}_{2i}$ and $\sigma^{x\,[1]}_{2i}|\pm\rangle^{[1]}_{2i}=|\mp\rangle^{[1]}_{2i}$. The renormalised Hamiltonian (8), besides the projector $\hat{P}_{1}$, has exactly the same structure of $\hat{H}^{(0)}$ in Eq. (4), being a sum of independent two-spin Hamiltonians. Eqs. (10) are the renormalisation rules for the magnetic field, coupling, and energy.

The diagonalisation and projection transformations (5-8) are then applied again $n-1$ times on all the terms of the Hamiltonian (3). The $k$-th renormalised spin is the block of $2^{k}$ spins at level $k$. Its basis $|\pm\rangle_{i}^{[k]}$ is defined recursively as

[TABLE]

Here $|a\,b\,\rangle^{[k-1]}_{i}$ stands for the tensor product of two adjacent $(k-1)$-th renormalised spins as ${|a\,b\,\rangle^{[k-1]}_{i}=|a\rangle_{2i-1}^{[k-1]}\otimes|b\rangle_{2i}^{[k-1]}}$. The parameters $a_{k}$, $b_{k}$ are obtained from the diagonalisation and read

[TABLE]

while the contribution of each projection to the ground state energy is

[TABLE]

At every step, the renormalisation rules

[TABLE]

give explicitly the recursion relation for the control parameter (9) as

[TABLE]

2.2 The RG ground state

After renormalisation, the hierarchical Hamiltonian (3) becomes

[TABLE]

where $e_{k}$ is the energy contribution (13) of each renormalisation step, $\hat{P}_{k}$ is the projection operator onto the lowest energy eigenstates at level $k-1$ (i.e. the generalisation of Eq. (7)), $h_{n}$ is the renormalised transverse field (14) and $\hat{\sigma}^{z\,[n]}$ is the Pauli matrix at level $n$. Hereafter, we will not use the subscript for the Pauli matrices and eigenvectors at level $n$ because there is a single spin at that level. Since $h_{n}>0$, the ground state is simply

[TABLE]

and it can be written in the basis of the spins at the previous level using Eq. (11) as

[TABLE]

and the recursively up to writing it in the physical spins at level [math].

Note that not only the Hamiltonian has a binary tree form, but also the RG ground state has the structure of a tree tensor network where the isometries have all rank equal to two, see [76, 77, 78, 79]. In fact, the transformations of Eq. (11) automatically define the isometric tensors in the sense of Ref. [76].

2.3 The critical point

The critical point and all the universal exponents are determined from the study of the flow induced by the RG equation (15). The critical point $K_{c}=R[K_{c}]$ depends on $\sigma$. For $\sigma\ll 1$, $K_{c}\to 0$ (aka $h_{c}$ diverges) and all the spins align along the $z$ direction. Conversely, for $\sigma\gg 1$, $K_{c}\to\infty$ (aka $h_{c}\to 0$). Hence, as in the usual quantum Ising model [82], the critical point separates a ferromagnetic phase, for $h<h_{c}$ or $K>K_{c}$, from a paramagnetic one, for $h>h_{c}$ or $K<K_{c}$. At the critical point, the ground-state wave function simplifies considerably. Since $K_{k}=K_{c}$ for all $k$, the coefficients of Eq. (12) do not depend on $k$, i.e. $a_{k}=a\,,b_{k}=b$ for all $k$ and the recurrence relation (18) further simplifies, as we will exploit to calculate the entanglement entropy.

By linearising the flow equation (15) close to the fixed point, the critical exponents for the correlation length ($\nu_{\sigma}$), correlation function ($x_{\sigma}$), and the dynamical one ($z_{\sigma}$) have been determined [72]. We will be interested in the longitudinal correlation function between two spins at distance $r$, that at the critical point, behaves as

[TABLE]

where the exponent $x_{\sigma}$ is [72]

[TABLE]

We mention that in the limit $\sigma\ll 1$, all the critical exponents in [72] reproduce the mean-field results of [13, 14]. For larger $\sigma$, they are in good qualitative agreement with the values obtained numerically in a dissipative short-range model in the same universality class of the long-range Ising one[56, 57, 58].

3 The block entanglement entropy of the hierarchical model

In this section we present the main results of this paper about the entanglement entropy of the RG ground state of the Dyson hierarchical model at the critical point. We mention that a very large literature exists about the characterisation of the entanglement in many-body systems by real space RG, but its main focus is disordered systems and strong disorder RG [83, 84, 86, 85, 87, 88, 89, 92, 91, 90, 94, 93], see [95] as a recent review. Although the employed techniques are different, there are many qualitative similarities that also helped our understanding of the hierarchical model.

We consider bipartitions of the chain into a block of spins $A$ of length $\ell$ and the remainder $B$. The first step will be to reconstruct the reduced density matrix $\hat{\rho}_{A}$ of the ground state (17) as achieved by means of the recursive projection technique in the previous section. We then calculate the eigenvalues of $\hat{\rho}_{A}$ and consequently the entanglement entropy (1). Since local transformations within the subsystems do not alter the entanglement between $A$ and $B$, we will write $\hat{\rho}_{A}$ without tracing back all those renormalised spins that are not split by the bipartition: in tree tensor network language [76, 79], this means that all the isometries contained within $A$ do not contribute to the entanglement and cancel in the construction of the reduced density matrix. The size of the reduced density matrix is then determined by the number of the renormalised spins (11) that are split by the bipartition. As a simple example, let $A$ be the block composed of the first $\ell=L/2$ sites. This partition spits only the $n$-th effective spin at level $n$, as in the top panel of Fig. 2 for $k=1$. Hence, the reduced density matrix $\hat{\rho}_{A}$ is a $2\times 2$ matrix in the basis $|\pm\rangle^{[n-1]}_{1}$ and its eigenvalues will never change when rewriting in terms of the spins at lower levels up to the physical spins.

The same reasoning applies to all other kinds of intervals. However, since the model is not translationally invariant, the reduced density matrix depends not only on the length of the subsystems $A$, but also on its position in the chain. We consider the three following bipartitions as shown in Fig. 2:

the interval of length $\ell_{k}=L/2^{k}={2^{n-k}}$ (with $k=1,\dots n$) starting from left end of the chain as in the top panel of Fig. 2; 2. 2.

the interval of length $\ell_{pm}=(2^{p}-1)2^{m}+1$ (with $m,p\geq 1$) again starting from left end of the chain as in the central panel of Fig. 2; 3. 3.

the shifted interval of length $\ell=L/2$ starting from the second site of the chain as in the bottom panel of Fig. 2.

For these three cases, we will construct the reduced density matrix $\hat{\rho}_{\ell}$ recursively.

The partition (1), as previously discussed for $\ell=L/2$, splits only one renormalised spin at level $n-k$. Hence, the reduced density matrix is conveniently written in the basis at this level when it is again a $2\times 2$ matrix, whose elements will be obtained recursively.

Partition (2) instead cuts $p+1$ effective spins, hence it is a matrix of size $2^{p+1}\times 2^{p+1}$ in the basis given by the tensor product of the effective spins at corresponding levels (see below for details). Consequently the largest possible rank is obtained for $\ell=L/2-1$, i.e. $m=1$ and $p=n-2$.

The partition (3) is chosen in such a way to have maximum possible rank (as standard in tree tensor networks [76, 77, 78, 79]). Indeed, this partition cuts by construction $n$ effective spins and so $\hat{\rho}_{A}$ is, in principle, a $L\times L$ density matrix.

Numerical and analytical results for the entanglement entropy at criticality will be explicitly obtained in the following subsections for asymptotically large systems and subsystems. The most relevant and surprising result is that the reduced density matrix has alway a finite rank (at most equal to 16), even for the case (3) when the rank naively might have been maximal, i.e. $L$. The natural question at this point is how well the fixed point ground-state entanglement correctly describes the one of the real model, without integrating out the high-energy physics in the renormalisation procedure. To this aim, we compare the renormalised entanglement with the exact one. In Fig. 3 we report the comparison with the data from numerical exact diagonalisation for $L=16$ and partition (3) (because of the tree structure, the next possible size would be $L=32$ that is hard to exactly diagonalise; alternatively one can use tree tensor networks to get very accurate estimates for larger chain, but this is beyond the goal of this paper). We show the comparison for many values of $\sigma$. The two curves are qualitatively very similar, although $S_{A}$ is non-monotonic in $\sigma$. While for large and small $\sigma$ they compare very well, there are quantitative differences for $\sigma\sim 1$. Most probably this is just a finite size effect and indeed it is remarkable that for such small value of $L$, the two curves match so well. In the figure we report also the results from renormalisation group for larger $L$, showing that the saturation value is much larger for intermediate $\sigma$ as we shall discuss.

Since in the following we will only consider critical entanglement entropies, to lighten the notations we set $K_{c}=K$ to indicate the control parameter at the critical point and we set $J=1$.

3.1 Entanglement entropy for the partition (1).

We first consider the reduced density matrix for an interval $A$ of length $\ell_{k}=2^{n-k}$ starting from the chain boundary. The partial trace over $B$ (the complement of $A$, i.e. the last $L-\ell_{k}$ spins) may be obtained recursively by tracing away half of the spins at each step $k$ times, as shown in the top panel of Fig. 2. For later convenience, it is useful to define the auxiliary matrices $\hat{\tau}_{i}^{[k]}$, represented pictorially by the diagram $(a)$ in the inset of Fig. 2, as

[TABLE]

where the states $|\pm\rangle^{[n]}$ are defined in Eq. (11). Notice that the $k$ dependence is encoded in $B$. These matrix can be written explicitly in the basis of the first spin at level $n-k$ (i.e. in the basis $|\pm\rangle^{[n-k]}_{1}$) by recursion (see A.1) as

[TABLE]

with

[TABLE]

with initial conditions $c_{0}^{+}=1,d_{0}^{+}=0,c_{0}^{-}=0,d_{0}^{-}=1,e_{0}=1,f_{0}=0$. Here $a,b$ are the coefficients in Eq. (12) evaluated at the critical point $K_{k}=K$. These recurrence relations admit the exact solution

[TABLE]

with

[TABLE]

Being the ground state $|\psi\rangle_{\text{GS}}=|+\rangle^{[n]}$ as in Eq. (17), we just have $\hat{\rho}_{A}=\hat{\tau}^{[k]}_{0}$, which is already written in diagonal form. Hence, the exact fixed-point entanglement entropy for this bipartition is

[TABLE]

Note that, because of the tree structure of the ground state, $S_{A}$ is a function only of $k$ and not of $k$ and $n$ separately. For large $k$ we have the exact asymptotic expansion

[TABLE]

with

[TABLE]

Thus, the entanglement entropy saturates to a constant $S^{*}$, exponentially fast in $k$. Since $e^{-k\ln 2}=(\ell_{k}/L)$, Eq. (27) is an expansion for $\ell_{k}\ll L$.

In Fig. 4, we plot the asymptotic fixed-point value $S_{*}$ (28) as a function of $\sigma$. We also compare it with the exact diagonalisation results for $L=8$ and $L=16$ with all possible values of $\ell_{k}=1,2,4$ and $\ell=8$ for $L=16$. Although the size of the system is very small, the results for small $\ell$ lie very close to the saturation value $S^{*}$ (expected to be reached for $\ell_{k}\ll L$, but large $L$).

3.2 Entanglement entropy for the partition (2)

In this subsection, we obtain the renormalised reduced density matrix of the critical ground state (17) for a partition of length $\ell_{pm}=(2^{p}-1)2^{m}+1$. (We can also think to the bipartition as a binary number: the binary representation of $\ell_{pm}-1$ is just $p$ ones followed by $m$ zeros – we can place $\nu=n-p-m$ zeros in front of it – the result does not depend on $\nu$). The number of effective spins broken by the partition is $p+1$, as it should be clear also from Fig. 2. The rules to construct recursively the reduced density matrix are:

1. (a)

Trace away recursively half of the spins $\nu=n-m-p$ times to get the reduced density matrix $\hat{\rho}_{(a)}$ of the first $2^{p+m}$ spins. $\hat{\rho}_{(a)}$ is a $2\times 2$ matrix in the basis $|\pm\rangle^{[p+m]}_{1}$ as in Sec. 3.1. 2. (b)

Trace out the $(m-1)$-th effective spin on the right. $\hat{\rho}_{(b)}$ is a $2^{p+1}\times 2^{p+1}$ matrix in the basis $|\pm\rangle^{[m+p-1]}_{1}\otimes|\pm\rangle^{[m+p-2]}_{3}\otimes\dots\otimes|\pm\rangle^{[m]}_{2^{p}-1}\otimes|\pm\rangle^{[m-1]}_{2^{p+1}-1}$. This is the reduced density matrix of the first $2^{m+p}-2^{m-1}$ spins at level zero. 3. (c)

On the rightmost effective spin at level $m-1$, trace out all spins except one at level zero by means of the recurrence relations (21). This last step does not increase the size of the reduced density matrix because all involved effective spins have been already broken before. $\hat{\rho}_{(c)}$ is the desired $2^{p+1}\times 2^{p+1}$ density matrix in the basis $|\pm\rangle^{[m+p-1]}_{1}\otimes|\pm\rangle^{[m+p-2]}_{3}\otimes\dots\otimes|\pm\rangle^{[m]}_{2^{p}-1}\otimes|\pm\rangle^{[0]}_{\ell_{pm}}$.

These three rules are pictorially summarised in Fig. 5.

As an example, let us consider the case $p=1$ and $m=2$, i.e. an interval of length $\ell=5$ (in a chain of arbitrary size $L$): $\ell-1=4\to\dots 0100$, see also Fig. 5 (bottom). (a) Compute first the $2\times 2$ reduced density matrix of the first spin at level $p+m=3$, i.e. $\hat{\rho}_{(a)}=a^{2}\,|+\rangle\langle+|^{[3]}_{1}+b^{2}|-\rangle\langle-|^{[3]}_{1}$ . (b) Trace away the rightmost effective block-spin at level one. This yields the $4\times 4$ matrix $\hat{\rho}_{(b)}$, in the basis $|\pm\rangle^{[2]}_{1}\otimes|\pm\rangle^{[1]}_{3}$. (c) Trace out the sixth physical spin from the block spin at level one. Hence $\hat{\rho}_{\ell=5}$ is then a $4\times 4$ matrix in the basis $|\pm\rangle^{[2]}_{1}\otimes|\pm\rangle^{[0]}_{5}$.

Step (a) has been already discussed in Sec. 3.1. For the ground state (17), Eqs. (21-22) evaluated at $k=\nu$ read

[TABLE]

where $c_{\nu}^{\pm}$ are given in Eq. (24). The reduced density matrix is a $2\times 2$ matrix in the basis of the effective spin at level $n-\nu=p+m$.

Now we perform step (b) and construct $\hat{\rho}_{{(b)}}$. We trace away the rightmost $(m-1)$-th effective spins from the reduced density matrix in Eq. (29). Since the partial trace occurs at level $m$, it leaves untouched all the $p$ effective spins on the left. The trace of the $(m-1)$-th right effective spin from the $m$-th yields

[TABLE]

Let us now define for convenience the following matrices for a generic $p$

[TABLE]

They are obtained recursively as

[TABLE]

with the initial conditions given by Eq. (3.2). At the end of the recursive construction, the reduced density matrix is a $2^{p+1}\times 2^{p+1}$ block-diagonal matrix with symmetric and anti-symmetric components

[TABLE]

as it follows from the fact that the Hamiltonian (4) is block-diagonal. The $\hat{\rho}_{S/A}^{[p]}$ are $2^{p}\times 2^{p}$ matrices given directly in terms of Eqs. (24-3.2).

We are ready for the final step (c). The reduced density matrix $\hat{\rho}_{(c)}$ is obtained by tracing away the remaining $2^{m}-1$ physical spins from (33). This is achieved by means of the procedure (1) applied to the $m$-th effective spin $|\pm\rangle^{[m]}_{2^{p+1}-1}$. The final reduced density matrix then reads

[TABLE]

where ${\hat{\rho}}^{[p,m]}_{i}$ are given by the recurrence relations (3.2) with initial conditions given by ${\hat{\rho}}^{[0,m]}_{i}=\hat{\tau}^{[m]}_{i}$ in Eq. (22).

At this point, we have a general recursive form for the reduced density matrix whose elements can be easily constructed. We diagonalise the reduced density matrix (34) to yield the entanglement entropy. The results of this procedure are shown in Fig. 6 for various values of $\sigma$. We plot the entanglement entropy as function of $\ell$ for several values of $p$ and $m$ (as we stressed, the value of $L$ does not matter at the fixed point). For fixed $p$, $S_{A}$ initially grows with $m$ up to a saturation value depending on $p$. The same trend holds when $p$ is varied at fixed $m$. The scaling with $\ell$ approaching its asymptotic value depends on $\sigma$. In fact, for small $\sigma$ the entanglement entropy initially grows quickly and saturates after a given $\ell^{*}(\sigma)$. This saturation value $\ell^{*}(\sigma)$ increases as $\sigma$ gets larger. Conversely, at large $\sigma$, we observe a very small growth, roughly compatible with a logarithmic behaviour that persists for many decades (more than 15) before saturation. We will get an analytic understanding of this behaviour from the expansion at large $\sigma$. For fixed $\ell$ and $L$, the dependence of the entanglement entropy on $\sigma$ is very similar to the one in Fig. 3: starting from zero, it grows with $\sigma$ up to a maximum and then decreases.

3.2.1 Some analytical expansions.

We derive an analytical expression of the entanglement entropy in the case $p=1$, which corresponds to an interval of length $\ell_{m}=2^{m}+1$. In this case, the reduced density matrix is just a $4\times 4$ density matrix in the basis $|\pm\rangle^{[m]}_{1}\otimes|\pm\rangle^{[0]}_{\ell_{m}}$, whose elements are determined recursively via Eq. (34). These elements are easily worked out analytically and so are their expansion for $\sigma\ll 1$ and $\sigma\gg 1$. See A.2 for the derivation and the details.

For $\sigma\ll 1$, the critical coupling $K\to 0$, the expansion of the matrix elements (34) at order $\mathcal{O}(K^{3})$ provides the following asymptotic expression for large $\ell_{m}$

[TABLE]

Here we first take the limit of small $K$ and only after large $\ell$: the two limits do not commute. The entanglement entropy saturates for large $\ell$ to a constant value proportional to $K^{2}\ln K$. On the left panels of Fig. 7, we compare the entanglement entropy in Eq. (35) with the one from the recurrence relations (34), finding a perfect match of the results.

In the opposite limit $\sigma\gg 1$, i.e. for critical coupling $K\to\infty$, we find two different regimes for large $\ell_{m}$. At fixed $\ell_{m}$, for $K^{-1}\ll\ell_{m}\ll 2^{K^{2}}$, the entanglement entropy has the following asymptotic expansion

[TABLE]

In this regime, the entanglement entropy grows logarithmically with the length of the interval, but with a very small prefactor proportional to $K^{-2}$, as it was clear already from the data in Fig. 6. This prediction is checked against the exact data from the recurrence relations in Fig. 7 (right panels), finding perfect agreement. Furthermore, the results in Fig. 6 for $\sigma\gg 1$ show that the same logarithmic growth, numerically with the same prefactor, appears for all other lengths with $p\neq 1$, although we do not have an analytic handle on them. It is worth mentioning that this behaviour is reminiscent of the standard one for critical short-range systems [66]. However, the right panels of Fig. 7, show also that the exact data start deviating from this logarithmic scaling for $\ell_{m}\sim 2^{K^{2}}$. Indeed, at this value of $\ell_{m}$, saturation toward the asymptotic value starts taking place and Eq. (36) fails: in other words the limits $K\to\infty$ and $m\to\infty$ do not commute. In turn, when $\ell_{m}\gg 2^{K^{2}}$, the entanglement entropy saturates to

[TABLE]

3.3 Entanglement entropy for the partition (3)

We now consider the partition (3) as in Fig. 2 (bottom). This partition cuts $n$ effective spins, giving rise to a $L\times L$ reduced density matrix. The ground-state density matrix can be written as a function of the $(n-1)$-th effective spins as in Eq. (18)

[TABLE]

where the coefficients $a,b$ are given by Eq. (12) We trace away the sub-system $A$ and write the density matrix in the basis represented pictorially by the green circles in Fig. 2 (bottom-right). This trace over $A$ is relatively simple, because it can be performed separately on the first and on the second $(n-1)$-th effective spin, see Fig. 2. On the left, the partial trace is equivalent to the procedure (1) applied to $|\pm\rangle\langle\pm|_{1}^{[n-1]}$. It results in the $2\times 2$ matrices ${\hat{\tau}_{i}^{[n-1]}}$ defined in Eq. (21) in the basis $|\pm\rangle^{[0]}_{1}$ of the first physical spin. On the second block spin, the result obtained from the trace corresponds to Eq. (3.2) with $p=n-2$ in the basis of $|\pm\rangle^{[0]}_{2^{n-1}+2}\otimes|\pm\rangle^{[1]}_{2^{n-2}+2}\otimes\dots|\pm\rangle^{[n-2]}_{4}$. Accordingly, the reduced density matrix is

[TABLE]

The resulting entanglement entropy reported in Fig. 8 as a function of $\ell=L/2$ has the same qualitative behaviour of the partition (2): it initially grows with the system size and then saturates to a finite value.

For large $\sigma$ we can make also a more quantitive comparison. We focus on the partition (2) with $m=1$ and $p=n-2$ in (34), i.e. $\ell=L/2-1$ which is the choice that scales more similar to the partition of length $\ell=L/2$ considered in this subsection. These two entanglement entropy are compared in Fig. 9. We have analytic result for partition (2) only for $p=1$ and in the regime $K^{-1}\ll\ell\ll 2^{K^{2}}$, when

[TABLE]

with $f(y+1)=(1/4K^{2}\ln 2)\ln y\left[1-\ln(\ln y/4K^{2}\ln 2)\right]$). We have noticed in Fig. 6 that the same growth with $f(\ell)$ is approximately valid for all values of $p$, and indeed also in Fig. 9 the data for $m=1$ and $p=n-2$ are well described by Eq. (40). For partition (3), the same form describes very well the data, just multiplying $f(\ell)$ by 2. This factor is likely related to the number of boundaries in the partitions, in analogy to what happens for systems with short-range interaction, also at the critical point [66, 68].

3.4 Finite rank reduced density matrix.

We show now by a very elementary argument that the rank of the reduced density matrix is finite. We focus here on the partition (3), but the same reasoning applies to the partition (2) and also others that we did not consider. At criticality the reduced density matrix of partition (3) can always be written from Eq. (39) as

[TABLE]

where $\vec{c}=(a^{2},b^{2},ab,ab)$, $\hat{\tau}^{[n-1]}_{i}$ are $2\times 2$ matrices, whose coefficients have been determined exactly (cf. Eqs. (22-24)), and $\hat{\rho}^{[n]}_{i}$ are defined recursively in Eq. (3.2). The structure of $\hat{\rho}^{[n]}_{i}$ is such that for all $n$

[TABLE]

From the sub-additivity of the rank it follows that $\text{Rank}[\hat{\rho}_{A}]\leq 16$. For partition (3) this bound can be improved to $\text{Rank}[\hat{\rho}_{A}]=8$, but for other partitions is tight. In A.3, we provide a proof of Eq. (42) based on linear algebra. An analogous proof might be obtained using tensor networks.

The real-space RG procedure projects the ground state onto a finite entanglement state. In Fig. 10, we compare the spectrum of the reduced density matrix -aka the entanglement spectrum- of the fixed-point RG ground state (42) to the one computed from the ground state of from exact diagonalisation. In the RG procedure, we have only eight non-zero eigenvalues. These match extremely well the largest eigenvalues of the exact reduced density matrix. Consequently, the relative difference between the exact entanglement entropy and the one computed with eight eigenvalues is small for all values of $\sigma$ and it is at most of order $\Delta S/S\leq 10^{-5}$. This suggests that the portion of the entanglement spectrum not captured by the RG procedure is irrelevant for the exact entanglement entropy in the thermodynamic limit.

3.5 Power-law decaying correlation functions.

The most striking and surprising aspect of our result is that we have an area law state which supposedly capture power-law correlation functions. This fact is against many common beliefs in the literature. Such peculiar behaviour is due to the hierarchical structure of the ground state and, a fortiori, to the lack of translational invariance. Yet, one can be very suspicious whether this is really possible. For this reason, from the finite rank reduced density matrix in Eq. (41) we reconstruct the two-point correlation function

[TABLE]

that at the critical point should scale as $C(r)\sim{r^{-2\,x_{\sigma}}}$, with exponent given in Eq. (20). We test this behaviour numerically on the fixed-point reduced density matrix obtained recursively. The $i$-th physical spin is equally correlated with all the physical spins belonging to the same block and hence, we evaluate the correlation function between effective spins at distance $r=2^{k}$, being the reduced density matrix (41) written in the basis of the $k$-th effective spins. The resulting correlation function is reported in Fig. 11 and displays the expected power-law decay with exponent $x_{\sigma}$ predicted analytically by RG (20).

4 Conclusions and perspectives

In this paper, we presented a detailed RG analysis of the ground-state entanglement entropy of a spin-block in the quantum Dyson hierarchical model. Our main goal was to get an analytical insight on the entanglement of long-range interacting spin systems beyond mean-field approximation. We found, surprisingly, that entanglement entropy obeys the area law also at criticality when correlation functions decay algebraically. This peculiar behaviour is due to the particular simple structure of the RG ground state: it is a tree tensor network with finite bond dimension and hence it has a finite-rank reduced density matrix. We must mention that unusual scalings of the entanglement in critical (but non conformal) ground states have been already observed in other models [19, 97, 20, 96, 98], and so our results represent yet another example of anomalous scaling in the absence of conformal invariance. Although it is unlikely that the true long-range ferromagnetic Ising model obeys the area law at criticality, it is desirable to check this expectation from direct numerical simulations.

Because of its simplicity, the Dyson hierarchical Hamiltonian is an interesting playground to explore also other entanglement properties of long-range interacting spin chains. For example, it would be worthy to characterise better the simple structure of this RG state as a tree tensor network, since it undergoes a quantum phase transition with finite entanglement, as those introduced in [99]. It would be interesting also to study the entanglement entropy and negativity between two disjoint blocks (as in [100, 101] for short range critical systems) to investigate whether there is some form of long distance entanglement. Finally, the RG approach also gives access to the low-lying excited states and consequently to their entanglement. In the short-range models, the entanglement of low-lying states captures many interesting physical features (see, e.g., [102]) and it is natural to wonder whether the same is true for long-range systems.

Finally, questions about the real-time dynamics of long-range systems can in principle be addressed using this model. In particular, in mean-field long-range systems, the entanglement entropy after a quench grows logarithmically with time [54], rather than linearly (i.e. ballistically) as for short-range Hamiltonians [103, 104, 105, 106, 107]. A number of studies, e.g. [52, 53], have reported numerical evidence of a sub-linear growth also in one spatial dimension. Hence, this model appears as the suitable setup to explore the entanglement entropy dynamics with non-equilibrium techniques such as strong disorder RG [108, 109, 111, 110].

We acknowledge useful discussions with Amir Mohammad Aghaei, Vincenzo Alba, Rosario Fazio, Saverio Pascazio, Luca Tagliacozzo, Andrea Trombettoni, and Xhek Turkeshi. In particular we thank Luca Tagliacozzo for sharing with us his results [79] before publication. We acknowledge support from the European Research Council (ERC) under the EU Horizon 2020 research and innovation program grant agreements No 694925 (GP) and No 771536 (PC). Part of this work has been carried out during the workshop “Entanglement in quantum systems” at the Galileo Galilei Institute (GGI) in Florence. SP thanks Boston University’s Condensed Matter Theory Visitors program for support.

Appendix A Details of calculations

In this appendix we report the explicit calculations of some results presented in the main text. In A.1 we derive the recursive structure of the reduced density matrix for the interval $\ell_{k}=2^{n-k}$ (22) in partition (1). In A.2, we derive the asymptotic expansions of the entanglement entropy for $\ell=2^{m}+1$ in Eqs. (35) and (36). In A.3 we show the finiteness of the rank of the reduced density matrix (cf. Eq. (42)) by linear algebra methods.

A.1 Recurrence relations for the reduced density matrix in partition (1).

Here we derive explicitly the recurrence relations in Eqs. (22) and (23). We first compute the reduced density matrix of half-system in the ground state $|+\rangle\langle+|^{[n]}$ (cf. Eq. (11))

[TABLE]

where $\hat{\tau}_{0}^{[k]}$ is defined in Eq. (21) and $|\pm\rangle_{1}^{[n-1]}$ is the basis of the first effective spin at level $n-1$. The reduced density matrix for $\ell_{2}=L/4$ is obtained first by rewriting $\hat{\tau}_{0}^{[1]}$ in terms of the effective spins at the level $n-2$ as

[TABLE]

and then tracing out the second spin at level $n-2$, obtaining

[TABLE]

which is Eq. (22) with $k=2$. Iterating this procedure, one easily gets the recurrence relation for $\hat{\rho}_{L/2^{k}}$. At the critical point the coefficients in Eq. (12) are constant for all $k$. All the recursion relations simplify and admit as solution Eq. (24).

A.2 Asymptotic expansions of the entanglement entropy for partition (2)

Here we work out the asymptotic expressions of the entanglement entropy for $\ell_{m}=2^{m}+1$ in Eqs. (35-36) in the limits $\sigma\gg 1$ and $\sigma\ll 1$. In this case the reduced density matrix (34) for $p=1$, i.e.

[TABLE]

is a $4\times 4$ block-diagonal matrix in the basis $|\pm\rangle^{[m]}_{1}\otimes|\pm\rangle^{[0]}_{\ell_{m}}$, which makes all calculations straightforward. This equation can be solved using the recurrence relations (3.2) with the initial conditions (24), obtaining

[TABLE]

with $\hat{\rho}^{s/a}_{\nu,m}$ block matrices of the form

[TABLE]

where $A^{s/a}_{\nu},\cdots E_{\nu}^{s/a}$ are read from Eqs. (24) and (34). We take now the limit $n\to\infty$ and, since $\nu=n-p-m$, the elements $c^{+}_{\nu},\,d^{+}_{\nu}$ are constant and given by the limit for $k\to\infty$ of Eq. (24). Defining $\xi\equiv\sqrt{K^{2}+4}$, the coefficients in (48) are

[TABLE]

The reduced density matrix can be diagonalised and from its four eigenvalues $\{\lambda_{i}\}$, the entanglement entropy $S_{A}=-\sum_{i=1}^{4}\lambda_{i}\ln\lambda_{i}$ is obtained as a function of $K$. In the limits $\sigma\gg 1$ and $\sigma\ll 1$, we derive the asymptotic expression of $S_{\ell_{m}}$ in Eqs. (35-36).

Asymptotic result for $\sigma\ll 1$. We consider the limit $\sigma\ll 1$, equivalent to $K\to 0$, see Eq. (15). At order $\mathcal{O}(K^{3})$ the matrix elements (48) read

[TABLE]

With these elements, the eigenvalues of the symmetric and anti-symmetric matrices (48) are

[TABLE]

hence the asymptotic entanglement entropy as a function of $m$ reads

[TABLE]

that, as function of $2^{m}=\ell_{m}-1$, is the same as Eq. (35).

Asymptotic result for $\sigma\gg 1$. Here we consider the limit $\sigma\gg 1$, i.e. for $K\to\infty$, and we expand the matrix elements (48) at $\mathcal{O}(K^{-3})$. The saturation value is obtained for $\alpha\,m,\beta m,\gamma\,m\gg 1$, corresponding to have only diagonal terms in Eq. (48) given by

[TABLE]

leading to the following saturation entanglement entropy

[TABLE]

The preasymptotic regime, is obtained by taking the limit $K\to\infty$ at fixed $m$. In this limit

[TABLE]

and the exponent $\beta\to 0$, so that we cannot neglect the non-diagonal parts of Eq. (48). Thus, there is another regime in which $m\gg\alpha^{-1},\gamma^{-1}$ but with $\beta m\ll 1$, i.e., in terms of $K$, $(\ln K)^{-1}\ll m\ll 2K^{2}$. In this regime, the diagonal terms $A_{m}^{s/a},B_{m}^{s/a}$ are the same as in Eqs. (53), but the non-diagonal ones become

[TABLE]

The diagonalisation of the reduced density matrix (47) finally yields

[TABLE]

that in terms of $2^{m}=\ell_{m}-1$ is Eq. (36) in the main text.

A.3 Rank of reduced density matrices.

Here we show Eq. (42) of the main text. All the matrices $\hat{\rho}_{i}^{[n]}$ in Eq. (3.2) can be written in a compact form as

[TABLE]

where

[TABLE]

while the real coefficients of the tensor $\gamma^{i}_{jk}$ can be obtained from Eq. (3.2). For fixed $i,j$ there is only one $k_{ij}$ such that $\gamma^{i}_{jk_{ij}}\neq 0$ (e.g. $\gamma^{0}_{0k}=a^{2}\,\delta_{k0}$). The coefficients $\gamma^{i}_{jk_{ij}}$ satisfy the following property

[TABLE]

Eq. (58) is written in the basis of the spin blocks at level $n$, see Fig. 2 (e.g. $\hat{\rho}_{i}^{[n=0]}$ is written in terms of the physical spins). For generic $n$, $\hat{\rho}_{i}^{[n]}$ is a matrix of size $2^{n+1}\times 2^{n+1}$. By construction, up to a change of basis, the matrices $\rho_{i}^{[n]}$ are either diagonal or anti-diagonal

[TABLE]

In order to prove Eq. (42), we show that $\forall n$ it holds

[TABLE]

Proof.

We will prove by induction. The step $n=0$ is obvious. Let us first see what happens for $n=1$ in order to start the induction. We construct the matrices $\hat{\rho}_{i}^{[1]}$ with $\gamma^{i}_{jk}$ given by Eq. (3.2). Let us consider

[TABLE]

and analogously for the other $\rho_{i}^{[1]}$. After a change of basis, this defines the recurrence for the matrices $A,B,C,D$ as

[TABLE]

These matrices have rank one if the determinant is zero (since they are not zero, the rank is not null). For the matrix $A_{0}^{[1]}$ we have

[TABLE]

Because of Eqs. (3.2-59-60), the same is true for all the matrices.

Let us now assume Eq. (62) to be true at $n-1$ and show that the same holds at level $n$. Eq. (64) holds at every step after a permutation of the basis elements. We now want to prove that they have rank one. Let us re-write the block matrices $A_{0}^{[n]}$ setting $\gamma^{0}_{0,k_{00}}/\gamma^{0}_{2,k_{02}}=\gamma^{0}_{3,k_{03}}/\gamma^{0}_{1,k_{01}}=w_{0}$. Up to a multiplication of the first and the second rows by a coefficient, we have

[TABLE]

where we multiplied the first row by the block matrix $B_{1}^{[n-1]}$ on the right and of the second row by the block matrix $C_{2}^{[n-1]}$ on the left. Since we assumed $A_{0}^{[m]}\,B_{1}^{[m]}=C_{2}^{[m]}D_{3}^{[m]}$ to hold at $n-1$, this means that the second row is linearly dependent on the first, hence the rank is the one of the first rows $\text{Rank}[A_{0}^{[n]}]=\text{Rank}[\left(w_{0}A_{0}^{[n-1]},B_{1}^{[n-1]}\right)]=\min\left(\text{Rank}[A_{0}^{[n-1]}],\text{Rank}[B_{1}^{[n-1]}\right)=1$. It is immediate to verify that ${A_{0}^{[n]}\,B_{1}^{[n]}=C_{2}^{[n]}D_{3}^{[n]}}$.

∎