Entanglement critical length at the many-body localization transition

The study of disordered, isolated quantum systems has recently received a lot of attention due to both technological progress in the isolation and manipulation of mesoscopic quantum systems [1–3] and theoretical advances. Among the latter, the phenomenon of many-body localization (MBL) has taken center stage since its theorization [4] as the surprising breakdown of ergodicity due to quantum effects in a generic (i.e. non-integrable) disordered system [5, 6]. The MBL phase, which is the continuation of the celebrated Anderson localized phase [7] to interacting systems, has been characterized, among others, via the absence of transport of conserved quantities (particles or spin and energy) [4, 6, 8, 9], nature of the eigenstates [10–12], spreading or evolution of entanglement [13–17] and emergent local integrability [18–22]. MBL is emerging as a compelling and universal phenomenon and its implications for our understanding of quantum dynamics and applications to future technologies are becoming evident [23–25].

Entanglement is commonly used as a characterization of a quantum state and of its non-classical properties [26, 27]. As said before, it has been used to characterize the MBL phase as well [13–17, 28]. In disordered systems one should discuss the probability distribution over realizations of the quenched disorder of the various entanglement measures. In particular, the typical value of the entanglement entropy of a single eigenstate has been shown to be a good discriminator of the two phases: the system's entanglement follows a volume law in the ergodic phase (in line with the eigenstate thermalization hypothesis, see [29, 30] and the recent review [31]) and an area law in the MBL phase [10]. The sample-to-sample fluctuations of the entanglement entropy have been instead shown [32] to be good locators of the transition, and have been proposed as a signal of its first-order nature. The entanglement spectrum [33] has proved in various situations to contain a lot more information than the entanglement entropy alone, and, where it is the case, it shows signatures of topological order.

In the context of MBL, the entanglement spectrum has been studied in some recent papers [34–36] (but see also in a similar setting [37]); when focusing on the MBL phase, these works find power-law distributions for the entanglement spectrum (similar to those found in [11]), while in the ergodic phase the distribution is found similar to a Marchenko–Pastur distribution, which is the distribution of eigenvalues of the reduced density matrices of random states [38–40]; [35] instead focuses on the level statistics of the entanglement spectrum, revealing interesting details of the transition on the MBL side.

The aim of this paper is to show that deviations from Marchenko–Pastur of the probability distribution of the entanglement spectrum in the ergodic phase provide an important characterization of the ergodic phase of such a disordered system. Moreover, such deviations from Marchenko–Pastur can be used to define a correlation length and predict the location and finite-size scaling exponents of the MBL transition. We will show that the entanglement spectrum therefore appears to be a crucial quantity, being able to identify even the most subtle correlations that are present in the ergodic phase of a disordered quantum system. The tools developed in this work should be useful for future studies of disordered quantum systems in their delocalized phases.

We consider what has become the standard model of MBL in one dimension: the spin-1/2 Heisenberg chain with random fields

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle H=-J\sum_{i=1}^L \ \vec{s}_i\cdot \vec{s}_{i+1} - \sum_{i=1}^L h_i s^z_i, \nonumber \end{align} \tag{ 1 }$$

where we fix $J=1$ and the random fields h_i are drawn from a uniform distribution in $[-h, h]$, and $s^{xyz}=\frac{1}{2}\sigma^{xyz}$ ($\sigma^k$ being the Pauli matrices). Furthermore, we choose to work in the subspace with total magnetization $S^z=0$, whose dimension is ${\mathcal N}=\binom{L}{L/2}$. In the numerics we study $L=8$ to 18, and consider periodic boundary conditions. The localization transition happens, for these parameters, at $h=h_c\approx3.7$ according to various numerical investigations [10, 12, 41].

We consider a single eigenstate of H, $\renewcommand{\ket}[1]{\vert #1\rangle} \ket{\psi}$, computed through the exact diagonalization of H, so that the corresponding eigenvalue is in the center of the energy spectrum; this can be obtained by targeting an energy $E={\mathcal N}^{-1} {\rm Tr} H$ in a shift-invert diagonalization algorithm (this estimates the mean spectral energy and accounts for its finite size corrections). Not averaging over a number of eigenstates corresponding to the same disorder realization removes a source of correlation.

To study the entanglement properties we focus on half of the chain, of size $L/2$, indicated as region A, the complementary region being B. We compute the reduced density matrix of this (pure) state as

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\ketbra}[2]{\vert #1\rangle\langle#2\vert} \renewcommand{\ket}[1]{\vert #1\rangle} \displaystyle \rho_A={\rm Tr}_B\rho={\rm Tr}_B\ketbra{\psi}{\psi} \label{eq:rdm} \nonumber \end{align} \tag{ 2 }$$

and determine its $2^{L/2}$ eigenvalues $\lambda_i$; this set (or better, the set of the logarithms $\ln\lambda_i$) is called the entanglement spectrum. Because of the $S_z=0$ constraint, (2) is in a block-diagonal form. We obtain numerical data with 10⁶ disorder realizations up to $L=14$, $7.5\cdot 10^4$ realizations for $L=16$ and 2500 to 5000 realizations for $L=18$.

In the delocalized phase the eigenvalues are all of the same order. Because of the normalization condition $\sum_i\lambda_i=1,$ we have therefore $\lambda_i\sim 2^{-L/2}$. One can think of the eigenvalues of order $\newcommand{\e}{{\rm e}} 2^{-\ell}$ as originating from states that are (approximately) uniformly entangled over a region of $\newcommand{\e}{{\rm e}} \ell$ spins, so that $\log_2\lambda_i$ can be thought of as a length scale.

A few examples will clarify the main idea. Consider first a maximally entangled state, such that the reduced density matrix of half system A has all eigenvalues $\lambda_i = 2^{-L/2}$. The only scale in the system is clearly $L/2$, as entanglement is evenly spread all over subsystem A. As a second example, consider a state such that the maximum eigenvalue is $\lambda_1=\mu = O(1)$ and all the others are $\lambda_i= O(2^{-L/2})$ such that $\mu + \sum_{j\neq 1} \lambda_j = 1$; then there are two scales, one $O(\log_2 \mu)=O(1)$ and the other one "global", $O(L/2)$. This situation typically arises when one biases the purity of a random state, so that one eigenvalue 'evaporates' out of the continuum of the other ones $O(2^{-L/2})$ [40], carrying most entanglement. As a third example, consider a situation in which the eigenvalues split is two classes: $\lambda_i= O(1/M_1)=O(2^{-L/2})$ and $\mu_j= O(1/M_2)=O(2^{-\sqrt{L}})$, their multiplicity such that $\sum_j \mu_j + \sum_{i} \lambda_i = 1$. Then there are two scales, one more 'local', $\newcommand{\e}{{\rm e}} \ell_2 = O(\log_2 M_2)=O(\sqrt{L})$ and the other one 'global' $\newcommand{\e}{{\rm e}} \ell_2 = O(\log_2 M_1)=O(L)$; these two scales will be in general intertwined in the system, in the sense that there will presumably be correlated 'islands' of size $O(\sqrt{L})$ interspersed in the system size $O(L)$. The above examples clarify that the scales are not necessarily related to the physical qubits, but rather to given directions in the Hilbert spaces of system A (and B), given by the Schmidt decomposition.

The above intuition is centerpiece in the numerical methods which are based on the truncation of the reduced density matrix [42] (used in the context of MBL in [43]) and will be useful in the following. In the next section we will investigate the details of the distributions and, specifically, the effects that arise from increasing the system size L and approaching the critical value of the disorder h_c. As noticed in previous work [36], in the MBL phase ($h>h_c$) the existence of a many-body localization length $L_l(h)$ implies a power-law distribution for the $\lambda_i$. In this paper we show the signs of another correlation length $L_s(h)$ in the entanglement spectrum, that will appear in the sub-leading corrections for $h<h_c$, and which diverges at the critical point. The different lengths (to be) introduced in this Article are summarized in table 1.

Table 1. Lengths and scales.

System size	MBL length	Correlation length	Boundary length scale	Bulk length scale
L	$\newcommand{\e}{{\rm e}} L_\ell$	Ls	$\newcommand{\e}{{\rm e}} \ell_1$	$\newcommand{\e}{{\rm e}} \ell_2$

In order to identify the different length scales in the $P(\lambda)$ it is useful to consider the entanglement entropy, defined as

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\Tr}{{\rm Tr}} \displaystyle S_A=-\Tr{\rho_A\ln\rho_A}=-\sum_i\lambda_i\ln\lambda_i . \nonumber \end{align} \tag{ 3 }$$

In recent years a lot of progress has been done on understanding the properties of the entanglement entropy close to a quantum phase transition [44]. In our case, it will not be possible to get too close to the MBL phase transition due to prohibitive numerical difficulties, but we think it is useful to discuss the general setup in which these kind of questions are considered.

Using the replica trick [45] S_A can also be written as (we shall drop henceforth the subscript _A on entropies)

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\Tr}{{\rm Tr}} \displaystyle S=-\lim_{n\to 1}\frac{\partial \Tr(\rho_A^n)}{\partial n} \equiv -\lim_{n\to 1}\frac{\partial S^{(n)}}{\partial n}. \nonumber \end{align} \tag{ 4 }$$

Usually, one takes for $\renewcommand{\ket}[1]{\vert #1\rangle} \ket{\psi}=\ket{E_0}$ (assuming $E_0=0$ for simplicity) the ground state of some local Hamiltonian which can be obtained as

$$\begin{align} \newcommand{\e}{{\rm e}} \renewcommand{\ket}[1]{\vert #1\rangle} \renewcommand{\bra}[1]{\langle#1\vert} \displaystyle \label{eq:gsproj} \ket{E_0}\bra{E_0}=\lim_{\beta\to\infty} {\rm e}^{-\beta H}. \nonumber \end{align} \tag{ 5 }$$

Physically, one needs $\beta\gtrsim 1/\Delta_0$, where $\Delta_0$ is the gap on the ground state. One can say that $1/\Delta_0$ defines a characteristic (imaginary) time after which the properties of the ground state can be extracted from the (imaginary) time dynamics.

Since

$$\begin{align} \newcommand{\e}{{\rm e}} \renewcommand{\braket}[2]{\langle#1\vert #2\rangle} \renewcommand{\ket}[1]{\vert #1\rangle} \renewcommand{\bra}[1]{\langle#1\vert} \displaystyle \rho_A(s^{\prime}_A\vert s_A)=\sum_{s_B}\braket{{s^{\prime}}_A ; s_B}{E_0}\braket{E_0}{s_A ; s_B},\quad \ket{s_A ; s_B}\equiv \ket{s_1,\dots,s_{L_A};s_{L_A+1},...,s_L}, \nonumber \end{align} \tag{ 6 }$$

we can consider $\rho_A(s^{\prime}_A\vert s_A)$ as the transition amplitude from configuration $\{s_A\}$ to configuration $\{s^{\prime}_A\}$, for the imaginary time evolution ${\rm e}^{-\beta H}$. Following this construction, one can think of $\newcommand{\Tr}{{\rm Tr}} S^{(n)}=\Tr(\rho_A^n)$ as a sum over all possible configurations of spins $s_{i}$ arranged on a n-sheeted (lattice) surface with two conical singularities $\pm 2\pi n$ at the boundaries between A and B. At criticality when one of the parameters $g\to g_c$, a correlation length ξ diverges, and one sees scale invariance in the limit $1\ll L\lesssim \xi$. For a conformal field theory (CFT) $\xi\sim 1/\Delta\sim \vert g-g_c\vert ^{-\nu}$ where ν is some critical exponent characteristic of the theory.

Whether the continuum limit yields a CFT or not, we can define twist operators [44, 45], which create conical singularities at 0 and L_A in the continuum limit, and write the entanglement entropy as a correlation function of these operators:

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\Tr}{{\rm Tr}} \displaystyle S^{(n)}(L_A)= \Tr(\rho_A^n)=\langle O^n(0)O^{-n}(L_A)\rangle, \nonumber \end{align} \tag{ 7 }$$

and by taking derivatives wrt to n one obtains the entanglement entropy. Notice that by knowing $S^{(n)}$ for every n one can reconstruct the entanglement spectrum (and vice versa).

If the energy E of the eigenstate is extensively away from the ground state, we need to modify (5) as

$$\begin{align} \newcommand{\e}{{\rm e}} \renewcommand{\ket}[1]{\vert #1\rangle} \renewcommand{\bra}[1]{\langle#1\vert} \displaystyle \ket{E}\bra{E}=\lim_{T\to\infty}\int_{-T}^T \frac{{\rm d}t}{2T}{\rm e}^{-{\rm i}t(H-E)} + O\left(\frac{1}{T} \right), \nonumber \end{align} \tag{ 8 }$$

where now $T\gtrsim1/\Delta$, the gap around the energy E, as one can easily show by looking at the matrix elements $\renewcommand{\bra}[1]{\langle#1\vert} \renewcommand{\braket}[2]{\langle#1\vert #2\rangle} \braket{s^{\prime}}{E}\braket{E}{s}$ ⁹. Now however, in a generic many-body model $\Delta\sim {\rm e}^{-L \sigma(E)}$ is exponentially small in L ($\sigma(E)$ is the entropy per spin at energy E) so we see that the scaling region requires $\xi\sim \ln(1/\Delta)\sim \vert g-g_c\vert ^{-\nu}$, which is probably the origin of claims of a dynamical exponents $z=\infty$ for the MBL transition [10].

As one takes larger and larger subsystems L_A we have a change of behavior when $L_A\simeq \xi$, entering the scaling region. We now go back to the description in terms of the eigenvalues of $\rho_A$ and in the rest of the paper we will see how to identify the lengthscale ξ in the entanglement spectrum. As said before, in the delocalized phase the eigenvalue probability distribution will have a finite support and a mean value that scales as $O(2^{-L_A})$. We can then use the system-size rescaled eigenvalues $2^{L_A}\lambda_i\to \lambda_i$, so that the mean $\langle{\lambda}\rangle=1$ does not scale with the system size. Typically, one observes (scaled) distributions such as those displayed in figure 1 (obtained for $h=1$). If this scaling is observed, then it is easy to show that

$$\begin{align} \newcommand{\e}{{\rm e}} \newcommand{\Tr}{{\rm Tr}} \displaystyle \Tr(\rho_A^n)=S^{(n)}(L_A)=2^{-(n-1) L_A}\phi_n(L_A). \nonumber \end{align} \tag{ 9 }$$

Zoom In Zoom Out Reset image size

By taking derivatives, the entanglement entropy reads

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle \label{eq:SL} S=\ln 2\ L_A-s_{L_A}, \nonumber \end{align} \tag{ 10 }$$

where L_A is the length of the interval A. The second term is related to the function $\phi_n$ and contains non-trivial information on the correlations between the twist operators. In particular it may also contain a non-trivial length scale. Notice that $s_{L_A}\geqslant 0$ since the first term on the right-hand side is also the maximum possible value of S_A. As we have chosen $L_A=L/2$, we will denote all the quantities with their L (rather than L_A) dependence. We will also indicate with λ the rescaled eigenvalues of the reduced density matrix.

Previous works [34, 36] noticed that in the delocalized phase, a substantial part of entanglement spectrum is well approximated at low disorder by a Marchenko–Pastur (MP) distribution

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle P(\lambda)=\frac{1}{2 \pi} \sqrt{\frac{4 - \lambda}{\lambda}}. \label{eq:mp} \nonumber \end{align} \tag{ 11 }$$

For MP, the constant in equation (10), in the thermodynamic limit, is $s_\infty=1/2$. However, finite-L corrections are important, and, if we assume that the corrections follow the same pattern as those in the MP distribution we should find

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle \label{eq:sLMP} s_L\simeq s_\infty+c\ {\rm e}^{-m L/2}, \nonumber \end{align} \tag{ 12 }$$

where for MP $s_\infty=1/2$, $m=\ln 2$, $c\simeq 0.38(3)$ (m and c come from the exact diagonalization of random Wishart matrices: an explicit calculation should be possible using methods of random matrix theory [38, 40, 46] and is left for future work). Notice that equation (12) is of the form of the corrections expected from quantum field theory (see [47] and references therein) where m is the mass of the lightest excitations. Using this form for fitting the numerical data we find $s_\infty\simeq 0.66(1)$ and $m\simeq 0.20(1)$ (around $h=1$, with a small dependence on h), therefore making evident that MP is not the limit of the rescaled eigenvalues distribution. In order to collect more evidence about such difference and study the evolution of the entanglement spectrum at the MBL (critical) transition point, we need to look at the whole distribution.

Although equation (12) might suggest a fast approach to the thermodynamic limit, we remark that the asymptotic limit is achieved only when L is larger than the correlation length ξ discussed before, which diverges at the transition at $h=h_c=3.7$. We will see how to identify this length, which we will call $L_s(h)$, in deference to the appearance of a spinodal point in the entanglement spectrum, and show that $L_s\gtrsim 20$ already at $h=2.5,$ representing therefore the main roadblock to the observation of a true asymptotic region in current numerics. The width of the numerically accessible critical length demands that we understand more accurately the sub-leading corrections to the entanglement quantities in MBL, if we want some more information on the theory of the critical point.

The remainder of this paper is focused on the definition and study of this correlation length $L_s(h)$ in the delocalized region, by scrutinizing its divergence and the associated finite-size corrections.

To better characterize the distribution we consider its logarithmic derivative

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle \beta(\ln\lambda,h,L)=\frac{{\rm d}\ln P}{{\rm d}\ln \lambda}, \nonumber \end{align} \tag{ 13 }$$

as a function of $\ln\lambda.$ For the MP distribution one would have $\beta=-2/(4-\lambda)$ and, in the limits $\log\lambda\to-\infty, \, \ln 4$, one would get $\beta\to -1/2, \, -\infty$, respectively. However, one numerically observes the peculiar behavior in figure 2 (at $h=1$) that, we will argue, extends all the way to the localization transition point $h_c=3.7$. For fixed h and small L, β is monotonically decreasing from $-1/2$ (at small λ) to $-\infty$ (at large λ): see for example the $L=8$ curve in figure 2. As L increases, there appears a value $L=L_s(h)$ (approximately $L=10$ in figure 2) where a tricritical point is found and where a minimum and a maximum are born in β; for $L>L_s(h)$ these extremal points are located at $\newcommand{\e}{{\rm e}} \ln\lambda\equiv -\ell_{1, 2}$, and we denote their values with $\newcommand{\e}{{\rm e}} \beta(-\ell_{1, 2}, L, h)\equiv \beta_{1, 2}$, respectively (see inset in figure 2).

Zoom In Zoom Out Reset image size

The tricritical point L_s is the solution of the equations (the primes denote derivatives wrt $\ln\lambda$):

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle \beta^{\prime}(\ln\lambda,L,h)=0,\nonumber \\ \beta^{\prime\prime}(\ln\lambda,L,h)=0, \nonumber \end{align} \tag{ 14 }$$

yielding the solution $\newcommand{\e}{{\rm e}} (L_s, \ln\lambda=-\ell_{1}=-\ell_{2})$ as a function of h. These equations can be solved easily, once a polynomial interpolating function for β is obtained from the data. The main lesson to be learned from this observation is the existence of the length scales $\newcommand{\e}{{\rm e}} L_s, \ell_1, \ell_2$. Before going to the physical interpretation of these scales we look at the evolution at fixed h of the pairs $\newcommand{\e}{{\rm e}} (\ell_{1, 2}, \beta_{1, 2})$ when $L\to\infty$. As can be seen from figure 3, the maximum $\beta_2$ approaches the MP value, $\beta_2\to -1/2$, with its position $\newcommand{\e}{{\rm e}} \ell_2$ remaining of $O(1)$ (inset); at the same time, the minimum $\beta_1$ goes to a value which is dependent on h (for $h=1$, $\beta_1\to -0.8$), but its position $\newcommand{\e}{{\rm e}} \ell_1$ escapes to the arbitrarily small eigenvalue region, $\newcommand{\e}{{\rm e}} \ell_1\to\infty$ (inset). We offer no explanation for this value of $\beta_1$ but we notice that in our numerics, as $h\to h_c$, it seems that $\beta_1\to-1$ (similarly to [11]).

Zoom In Zoom Out Reset image size

The distribution P therefore gets closer to MP for $\newcommand{\e}{{\rm e}} -\ln\lambda\geqslant \ell_2$ (figure 2). Since, as explained before, $-\ln \lambda$ is to be considered a length scale, we propose that $\newcommand{\e}{{\rm e}} \ell_2$ represents a typical length scale of the system, that is related to local phenomena, and characterizes the bulk properties of the system, describing for example the dynamics of local excitations of the eigenvector, enabling us to distinguish it from the eigenvector of a random matrix. We also surmise that the second, larger length scale $\newcommand{\e}{{\rm e}} \ell_1$, which goes to $\infty$ as $L\to\infty$, defines the effects on the distribution due to the size of the cut $L_A=L/2$, hence boundary effects. The length scale $L_s(h)$ is the minimum system size at which these two properties can be separated, and therefore it can be identified as the distance after which disturbances decay, or the correlation length ξ of our system.

We plot L_s versus h in figure 4 and observe its divergence at the MBL phase transition. The fit (dashed red line) gives a result consistent with the data in [12]:

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle \label{eq:Lshfit} L_s(h)\simeq\frac{a}{(h_c-h)^{\nu}}, \nonumber \end{align} \tag{ 15 }$$

where $h_c=3.7\pm 0.4$ and $\nu=0.9\pm 0.2$. On the other hand, by fixing the transition point at $h_c=3.72$ in the fitting procedure we get $\nu=0.88$ with smaller errors (in both cases $a\approx 23$); see figure 4. Notice that L_s is a proper length, measured in lattice sites. Notice also that the critical exponents are close to 1. This value is again in agreement with what is found in the numerical literature, but inconsistent with a Harris-like bound in [48], which suggest $\nu\geqslant 2$, and renormalization group analyses [49, 50] which give $\nu\simeq 3$. We must remark however that, upon identification of the length $L_s(h)$ which defines the beginning of the asymptotic region, the true critical exponents should be extracted by analyzing data at $L\gg L_s(h)$ and $h\lesssim h_c$. Current technology allows $L\lesssim 24$, which means $h<2.8$, therefore subleading corrections to the critical scaling might be substantial (an analog situation, in a similar context is observed in [51]). In particular a fit of the form $L_s\left({h}/{h_c}\right)\simeq a\, {\left(1-\frac{h}{h_c}\right){\hspace{0pt}}^{-3}\left[1-b_1\, \left(1-\frac{h}{h_c}\right)+b_2\, \left(1-\frac{h}{h_c}\right){\hspace{0pt}}^2\right]}$ is a decent alternative to (15) even with $h_c=3.72$ fixed, yielding $b_1\sim 7$ (admittedly not a small correction). This situation is unfortunate but cannot be mended in absence of either a full theory of the critical point of MBL (providing with the nature of the subleading corrections due to irrelevant operators, in the language of renormalization group analysis) or some radically new numerical tools (e.g. a quantum computer), enabling the study of eigenstates/eigenvalues for significantly larger system sizes (say $L\sim 100$).

Zoom In Zoom Out Reset image size

We now return to the subleading terms s_L in the entanglement entropy, discussing the finite size effects in this quantity as well. For $h\lesssim h_c$, if $L\ll L_s(h)$ the system should show critical behavior. All indicators from previous works (like absence of level repulsion, statistics of eigenvalues, entanglement entropy etc [12, 50, 52, 53]) show that the critical behavior is very similar to the MBL phase and that a jump of the entanglement entropy occurs at the transition. An analogous result is obtained also in studies of the (single particle) Anderson model on the Bethe lattice, which shows a critical behavior indistinguishable from that of the localized region [54, 55]. So if $L\ll L_s(h)$ the system shows localization properties. Eigenstates in the MBL phase obey an area law [13–15]; thus, since in the localized region $S=o(L)$, in equation (10) $s_L=\frac{\ln 2}{2}L+o(L)$ and so its derivative is $\frac{\partial s_L}{\partial L}=\frac{\ln2}{2}>0$. On the other hand, upon increasing L, if the finite-size corrections in s_L follow the pattern of the corrections in equation (12) we should observe, for $L\gg L_s(h)$, a negative derivative, $\frac{\partial s_L}{\partial L}=-c\, \frac{m}{2}\, {\rm e}^{-m L/2}<0$; notice that this behavior and the sign of the derivative is independent of the precise value of the coefficients. Summarizing, we should see that, for fixed h (and hence fixed $L_s(h)$) and as we increase L, the derivative $\partial s_L/\partial L$ must change sign, going from positive to negative, and then approach zero from negative values. This is exactly what is observed in figure 5. The value of L at which $\partial s_L/\partial L=0$ is proportional to L_s (the constant of proportionality being very close to 1), as is the position of the minimum of $\partial s_L/\partial L$ (the constant here is close to $1.3$).

Zoom In Zoom Out Reset image size

We analyzed the distribution of the entanglement spectrum of a single eigenstate of a spin-1/2 Heisenberg chain with random magnetic fields, identifying the finite-size corrections and the phase transition to an MBL phase. Our central result is the identification of three length scales: the first one $L_s(h)$ determines the minimum system size at which one can separate bulk and boundary effects; the other two $\newcommand{\e}{{\rm e}} \ell_1$, $\newcommand{\e}{{\rm e}} \ell_2$ are born at $L=L_s(h)$, but while the first one ($\newcommand{\e}{{\rm e}} \ell_1$) diverges as $L\to\infty$ and defines the boundary effects, the second one ($\newcommand{\e}{{\rm e}} \ell_2$) remains finite in that limit and determines the entanglement properties of local dynamics. We have developed new precision tools for the study of the distribution of the entanglement spectrum and we foresee that they will be useful for the analysis of the numerics of many other models, both disordered and clean.

The study of the size dependence of the entanglement entropy of a finite region [44, 45, 47, 56] can unveil crucial details about the critical theory underlying the MBL transition, which is, currently, unknown. Although MBL transitions (in particular because of the infinite temperature) are not the natural setting where quantum field theory has been applied, an effective model which uses the emergent scale invariance at h_c probably has a natural description in terms of QFT. If this is the case, the study of entanglement in MBL can be the back-door to develop a QFT of the MBL transition.

We would like to thank B L Altshuler, P Calabrese, A De Luca, D A Huse, F Slanina and M Schiulaz for useful discussions. FP would also like to thank ICTP for hospitality during the completion of this work. SP was partially supported by INFN through the project 'QUANTUM'. This project has received funding from the European Research Council (ERC) under the European Union's Horizon 2020 research and innovation programme (grant agreement 694925).