<figure><img src="image_1.png" /></figure>

HAL
open science

# Inverse participation ratios in the XXZ spin chain

Grégoire Misguich, Vincent Pasquier, Jean-Marc Luck

► To cite this version:

Grégoire Misguich, Vincent Pasquier, Jean-Marc Luck. Inverse participation ratios in the XXZ spin chain. Physical Review B, 2016, 94, pp.155110. 10.1103/PhysRevB.94.155110. cea-01377658

HAL Id: cea-01377658

https://cea.hal.science/cea-01377658

Submitted on 18 Apr 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
Inverse participation ratios in the XXZ spin chain

Grégoire Misguich, Vincent Pasquier, and Jean-Marc Luck

Institut de Physique Théorique, Université Paris Saclay, CEA, CNRS, F-91191 Gif-sur-Yvette, France

(Received 5 July 2016; published 6 October 2016)

We investigate numerically the inverse participation ratios in a spin-1/2 XXZ chain, computed in the “Ising” basis (i.e., eigenstates of $\sigma_i^z$). We consider in particular a quantity $T$, defined by summing the inverse participation ratios of all the eigenstates in the zero-magnetization sector of a finite chain of length $N$, with open boundary conditions. From a dynamical point of view, $T$ is proportional to the stationary return probability to an initial basis state, averaged over all the basis states (initial conditions). We find that $T$ exhibits an exponential growth, $T \sim \exp(aN)$, in the gapped phase of the model and a linear scaling, $T \sim N$, in the gapless phase. These two different behaviors are analyzed in terms of the distribution of the participation ratios of individual eigenstates. We also investigate the effect of next-nearest-neighbor interactions, which break the integrability of the model. Although the massive phase of the nonintegrable model also has $T \sim \exp(aN)$, in the gapless phase $T$ appears to saturate to a constant value.

DOI: 10.1103/PhysRevB.94.155110

# I. INTRODUCTION

In recent years there has been a renewed interest in the dynamical properties of out-of-equilibrium isolated quantum systems [1–3]. For instance, in a “quantum quench” setup, the system is initially prepared in the ground state of some Hamiltonian $H_0$ and then evolved (unitarily) with a Hamiltonian $H \neq H_0$. A fruitful point of view is to consider a small spatial region of a large system and to investigate the long-time limit of observables defined in this subsystem. Some typical questions one may then address are whether these expectation values have well-defined long-time limits and if they are described by some statistical ensemble (either thermal or not) [4]. Such problems still represent an active subject of research. Clearly, the steady states that may be reached crucially depend on the structure of the eigenstates of $H$. In particular, it has been shown that integrable systems generally fail to reach thermal states [5]. Somewhat analogously, disordered many-body systems in the “many-body localized phase” also fail to “thermalize,” a phenomenon which has attracted a lot of attention recently [6].

Following Ref. [7], we adopt here a different point of view and consider how the system evolves as a whole, starting from a set of simple states $\{|a\rangle\}$ forming a “preferential” basis for the model. In the case of a particle system on a lattice, for instance, a natural choice for the preferential basis is the set of product states where the particles have fixed positions in real space. In a spin system (see below), one may choose the spin configurations that are eigenstates of all the $S_i^z$ operators for some choice of the quantization axis $z$. A quantity of interest is the typical time-average probability to return to the initial basis state. As will be explained below, these probabilities are related to the inverse participation ratios (IPRs) of the eigenstates, computed in the preferential basis. In the present study we explore numerically these quantities in a particular many-body problem, the spin-1/2 XXZ spin chain. The preferential basis is chosen to be the set of “Ising configurations” which are eigenstates of the $z$ component of the on-site magnetization. This spin-chain Hamiltonian depends on an anisotropy parameter $\Delta$. The main results of the paper concern the scalings of these IPRs in the gapped

phase ($|\Delta| > 1$) and in the gapless phase ($|\Delta| \le 1$) of the model.

The plan of this paper is the following. Section II is a review of isolated quantum systems, motivating the study of the IPR $t_n$ of individual energy eigenstates and of their sum $T$. In Sec. III we recall the definition of the XXZ spin chain and comment on the question of degeneracies. We then present and discuss our numerical results on various observables in the gapped phase and in the gapless phase in the main section (Sec. IV). Section V contains numerical results on a nonintegrable spin chain with second-neighbor interactions, and we briefly discuss our findings in Sec. VI.

# II. GENERALITIES

We consider an isolated quantum system with a finite-dimensional Hilbert space of dimension $D$. The eigenstates of the Hamiltonian $H$ are denoted by $|n\rangle$, and we assume for simplicity that they are nondegenerate. In the preferential basis $\{|a\rangle\}$, the IPR of an eigenstate is by definition

$$t_n = \sum_{a=1}^{D} |\langle a | n \rangle|^4. \qquad (1)$$

The maximum value of this quantity is reached when the eigenstate coincides with a single basis state. This is the completely localized case and gives $t_{max} = 1$. On the other hand, the minimum value of $t_n$ is reached for eigenstates which are uniform superpositions of all the basis states, with the same modulus $|\langle a | n \rangle| = 1/\sqrt{D}$. This maximally delocalized limit gives $t_{min} = 1/D$. These IPRs have been extensively used to measure the localization properties of a single-particle wave function [8,9] (for a review in the context of the Anderson localization see, for instance, Ref. [10]). They have also recently proved to be useful in the context of many-body localization [11,12]. They can be used as well to extract some universal long-distance properties from (clean) many-body ground-state wave functions [13].
Inverse participation ratios in the XXZ spin chain

Grégoire Misguich, Vincent Pasquier, and Jean-Marc Luck

Institut de Physique Théorique, Université Paris Saclay, CEA, CNRS, F-91191 Gif-sur-Yvette, France

(Received 5 July 2016; published 6 October 2016)

We investigate numerically the inverse participation ratios in a spin-1/2 XXZ chain, computed in the “Ising” basis (i.e., eigenstates of $\sigma_i^z$). We consider in particular a quantity $T$, defined by summing the inverse participation ratios of all the eigenstates in the zero-magnetization sector of a finite chain of length $N$, with open boundary conditions. From a dynamical point of view, $T$ is proportional to the stationary return probability to an initial basis state, averaged over all the basis states (initial conditions). We find that $T$ exhibits an exponential growth, $T \sim \exp(aN)$, in the gapped phase of the model and a linear scaling, $T \sim N$, in the gapless phase. These two different behaviors are analyzed in terms of the distribution of the participation ratios of individual eigenstates. We also investigate the effect of next-nearest-neighbor interactions, which break the integrability of the model. Although the massive phase of the nonintegrable model also has $T \sim \exp(aN)$, in the gapless phase $T$ appears to saturate to a constant value.

DOI: 10.1103/PhysRevB.94.155110

# I. INTRODUCTION

In recent years there has been a renewed interest in the dynamical properties of out-of-equilibrium isolated quantum systems [1–3]. For instance, in a “quantum quench” setup, the system is initially prepared in the ground state of some Hamiltonian $H_0$ and then evolved (unitarily) with a Hamiltonian $H \neq H_0$. A fruitful point of view is to consider a small spatial region of a large system and to investigate the long-time limit of observables defined in this subsystem. Some typical questions one may then address are whether these expectation values have well-defined long-time limits and if they are described by some statistical ensemble (either thermal or not) [4]. Such problems still represent an active subject of research. Clearly, the steady states that may be reached crucially depend on the structure of the eigenstates of $H$. In particular, it has been shown that integrable systems generally fail to reach thermal states [5]. Somewhat analogously, disordered many-body systems in the “many-body localized phase” also fail to “thermalize,” a phenomenon which has attracted a lot of attention recently [6].

Following Ref. [7], we adopt here a different point of view and consider how the system evolves as a whole, starting from a set of simple states $\{|a\rangle\}$ forming a “preferential” basis for the model. In the case of a particle system on a lattice, for instance, a natural choice for the preferential basis is the set of product states where the particles have fixed positions in real space. In a spin system (see below), one may choose the spin configurations that are eigenstates of all the $S_i^z$ operators for some choice of the quantization axis $z$. A quantity of interest is the typical time-average probability to return to the initial basis state. As will be explained below, these probabilities are related to the inverse participation ratios (IPRs) of the eigenstates, computed in the preferential basis. In the present study we explore numerically these quantities in a particular many-body problem, the spin-1/2 XXZ spin chain. The preferential basis is chosen to be the set of “Ising configurations” which are eigenstates of the $z$ component of the on-site magnetization. This spin-chain Hamiltonian depends on an anisotropy parameter $\Delta$. The main results of the paper concern the scalings of these IPRs in the gapped

phase ($|\Delta| > 1$) and in the gapless phase ($|\Delta| \le 1$) of the model.

The plan of this paper is the following. Section II is a review of isolated quantum systems, motivating the study of the IPR $t_n$ of individual energy eigenstates and of their sum $T$. In Sec. III we recall the definition of the XXZ spin chain and comment on the question of degeneracies. We then present and discuss our numerical results on various observables in the gapped phase and in the gapless phase in the main section (Sec. IV). Section V contains numerical results on a nonintegrable spin chain with second-neighbor interactions, and we briefly discuss our findings in Sec. VI.

# II. GENERALITIES

We consider an isolated quantum system with a finite-dimensional Hilbert space of dimension $D$. The eigenstates of the Hamiltonian $H$ are denoted by $|n\rangle$, and we assume for simplicity that they are nondegenerate. In the preferential basis $\{|a\rangle\}$, the IPR of an eigenstate is by definition

$$t_n = \sum_{a=1}^{D} |\langle a | n \rangle|^4. \qquad (1)$$

The maximum value of this quantity is reached when the eigenstate coincides with a single basis state. This is the completely localized case and gives $t_{max} = 1$. On the other hand, the minimum value of $t_n$ is reached for eigenstates which are uniform superpositions of all the basis states, with the same modulus $|\langle a | n \rangle| = 1/\sqrt{D}$. This maximally delocalized limit gives $t_{min} = 1/D$. These IPRs have been extensively used to measure the localization properties of a single-particle wave function [8,9] (for a review in the context of the Anderson localization see, for instance, Ref. [10]). They have also recently proved to be useful in the context of many-body localization [11,12]. They can be used as well to extract some universal long-distance properties from (clean) many-body ground-state wave functions [13].
Following Ref. [7] (see also [14]), we sum $t_n$ over all the eigenstates to get

$$T = \sum_{n} t_{n} = \sum_{a,n=1}^{D} |\langle a|n \rangle|^{4}. \qquad (2)$$

This quantity measures how much the eigenstates are localized in the preferential basis. It can range from $T_{min} = 1$ (all the eigenstates spread maximally over the whole basis) to $T_{max} = D$ (each eigenstate matches exactly a single basis state). This quantity was originally introduced [7] from a dynamical point of view, as the trace of the matrix $Q$ whose entries $Q_{ab}$ are the time-average probabilities to go from state $|a\rangle$ at the initial time to state $|b\rangle$ at time $t$. If the system is prepared in state $|a\rangle$ at $t = 0$, the probability to observe it in state $|b\rangle$ at time $t$ is $P_{ab}(t) = |\langle b|\psi(t)\rangle|^2$, i.e.,

$$P_{ab}(t) = \sum_{m,n} e^{i(E_n - E_m)t} \langle b|m\rangle \langle m|a\rangle \langle a|n\rangle \langle n|b\rangle. \qquad (3)$$

In the absence of degeneracies in the spectrum, the time-average probability reads

$$Q_{ab} = \lim_{t \to \infty} \frac{1}{t} \int_{0}^{t} P_{ab}(t') dt' = \sum_{n} |\langle a|n \rangle|^2 |\langle b|n \rangle|^2. \qquad (4)$$

In this dynamical picture, $T/D$ measures the stationary return probability to an initial basis state, averaged over all the basis states (initial conditions). The minimum value $T_{min}/D = 1/D$ is reached if the dynamics connects any initial basis state to all the other basis states, a limit of perfect “equilibration.” On the other hand, the maximum value $T_{max}/D = 1$ is reached when $\mathcal{H}$ is diagonal in the preferential basis, so that the system does not evolve at all if prepared at $t = 0$ in a basis state.

This quantity $T$ was studied in detail in the case of a single particle (tight-binding model) in a one-dimensional random potential [7]. In the present work, we analyze the quantity $T$ for a simple many-body problem without disorder, the spin-$\frac{1}{2}$ XXZ spin chain.

### III. XXZ SPIN CHAIN

#### A. Hamiltonian

We consider the spin-$\frac{1}{2}$ XXZ chain with open boundary conditions, with a Hamiltonian given by

$$\mathcal{H} = \sum_{i=1}^{N-1} (S_i^x S_{i+1}^x + S_i^y S_{i+1}^y + \Delta S_i^z S_{i+1}^z). \qquad (5)$$

As is well known, the system has a gapped spectrum for $|\Delta| > 1$, with long-range order and spontaneous symmetry breaking in the thermodynamic limit (“Ising” phase). On the other hand, the system has a gapless spectrum and displays an algebraic decay of spin-spin correlations at zero temperature for $|\Delta| \le 1$ [the so-called Tomonaga-Luttinger liquid (TLL) phase [15]]. Hereafter, we consider the “Ising configurations” $|a\rangle = |\uparrow\uparrow\downarrow\cdots\rangle, |\uparrow\downarrow\uparrow\cdots\rangle, \cdots$ (eigenstates of all the $S_i^z$) as a natural basis for this problem. We use this basis to define the IPR $t_n$ of the eigenstates.

#### B. Remarks on degeneracies

In the case of a degenerate multiplet of eigenstates $|n_1\rangle, \dots, |n_d\rangle$, the corresponding contribution to $T$ (as derived from the dynamical point of view [7]) should be $\sum_a (\sum_{\alpha=1}^d |\langle a|n_\alpha\rangle|^2)^2$. For open boundary conditions and a generic value of $\Delta$, it turns out that the only degeneracies in the spectrum of Eq. (5) are the twofold degeneracy ($S_{tot}^z \leftrightarrow -S_{tot}^z$) of the eigenstates with a nonzero magnetization. But since we focus here on the $S_{tot}^z = 0$ sector, we generically get nondegenerate eigenstates which are either even or odd under a global spin flip, as well as even or odd under the spatial (left-right) symmetry.

We finally note that some additional degeneracies occur at the free-fermion point ($\Delta = 0$). For simplicity, we avoid this special point when numerically computing the IPR.

### IV. NUMERICAL RESULTS

#### A. The quantity $T$

The eigenstates of Eq. (5) were obtained using a full numerical diagonalization for even sizes up to $N = 20$ in the $S_{tot}^z = 0$ sector. The total Hilbert space dimension is therefore $D = (N/2) \sim 2^N \sqrt{2/(\pi N)}$. Using the spin-flip parity and the spatial parity with respect to the center of the chain, the total space is decomposed into four blocks, each one being labeled by two parities. The results for $T$ are summarized in Fig. 1. Note that $T$ is independent of the sign of $\Delta$. This follows from the fact that changing the sign of $\Delta$ is equivalent to changing the sign of the $xy$ terms, and the latter can be undone by some unitary transformation ($\pi$ rotation about the $z$ axis on every second site) which only affects the sign of the wave functions in the Ising basis, not their modulus. We can therefore restrict ourselves to $\Delta \ge 0$.

#### 1. Gapped phase

We observe (top panel of Fig. 1) that, for sufficiently large $\Delta$, $\ln(T)$ is approximately proportional to the number of sites, i.e.,

$$\ln(T) \approx a(\Delta) N, \qquad (6)$$

with $0 < a(\Delta) < \ln(2)$. We can deduce that the average IPR $\bar{t} = T/D$ behaves as $\ln(\bar{t}) \approx -[\ln(2) - a(\Delta)]N$. We conjecture that these scalings hold in the whole gapped (massive) region ($\Delta > 1$).

At $\Delta = \infty$, any Ising configuration (basis state) becomes an eigenstate of $\mathcal{H}$, so that one may expect $T = D$ and $\ln(T)/N \approx \ln(2)$ in that limit. This simple reasoning is, however, not correct since many energy levels become degenerate in this limit. Thus, as soon as $\Delta$ is not strictly infinite, the actual eigenstates are nontrivial superpositions of several basis states due to the effect of the $S_i^+ S_{i+1}^- + \text{H.c.}$ terms (according to degenerate perturbation theory). We indeed have $\ln(T)/N \approx a(\infty)$ for very large $\Delta$, but the numerically obtained value in this regime, $a(\infty) \approx 0.25$ (see Fig. 1), is significantly smaller than $\ln(2)$. A very similar phenomenon has already been observed in the much simpler situation of a tight-binding particle in a random potential with binary disorder [7], where the hybridization of degenerate molecular states
Following Ref. [7] (see also [14]), we sum $t_n$ over all the eigenstates to get

$$T = \sum_n t_n = \sum_{a,n=1}^{D} |\langle a|n \rangle|^4. \qquad (2)$$

This quantity measures how much the eigenstates are localized in the preferential basis. It can range from $T_{min} = 1$ (all the eigenstates spread maximally over the whole basis) to $T_{max} = D$ (each eigenstate matches exactly a single basis state). This quantity was originally introduced [7] from a dynamical point of view, as the trace of the matrix $Q$ whose entries $Q_{ab}$ are the time-average probabilities to go from state $|a\rangle$ at the initial time to state $|b\rangle$ at time $t$. If the system is prepared in state $|a\rangle$ at $t=0$, the probability to observe it in state $|b\rangle$ at time $t$ is $P_{ab}(t) = |\langle b|\psi(t)\rangle|^2$, i.e.,

$$P_{ab}(t) = \sum_{m,n} e^{i(E_n - E_m)t} \langle b|m\rangle \langle m|a\rangle \langle a|n\rangle \langle n|b\rangle. \qquad (3)$$

In the absence of degeneracies in the spectrum, the time-average probability reads

$$Q_{ab} = \lim_{t \to \infty} \frac{1}{t} \int_0^t P_{ab}(t') dt' = \sum_n |\langle a|n \rangle|^2 |\langle b|n \rangle|^2. \qquad (4)$$

In this dynamical picture, $T/D$ measures the stationary return probability to an initial basis state, averaged over all the basis states (initial conditions). The minimum value $T_{min}/D = 1/D$ is reached if the dynamics connects any initial basis state to all the other basis states, a limit of perfect “equilibration.” On the other hand, the maximum value $T_{max}/D = 1$ is reached when $\mathcal{H}$ is diagonal in the preferential basis, so that the system does not evolve at all if prepared at $t=0$ in a basis state.

This quantity $T$ was studied in detail in the case of a single particle (tight-binding model) in a one-dimensional random potential [7]. In the present work, we analyze the quantity $T$ for a simple many-body problem without disorder, the spin-$\frac{1}{2}$ XXZ spin chain.

### III. XXZ SPIN CHAIN

#### A. Hamiltonian

We consider the spin-$\frac{1}{2}$ XXZ chain with open boundary conditions, with a Hamiltonian given by

$$\mathcal{H} = \sum_{i=1}^{N-1} (S_i^x S_{i+1}^x + S_i^y S_{i+1}^y + \Delta S_i^z S_{i+1}^z). \qquad (5)$$

As is well known, the system has a gapped spectrum for $|\Delta| > 1$, with long-range order and spontaneous symmetry breaking in the thermodynamic limit (“Ising” phase). On the other hand, the system has a gapless spectrum and displays an algebraic decay of spin-spin correlations at zero temperature for $|\Delta| \le 1$ [the so-called Tomonaga-Luttinger liquid (TLL) phase [15]]. Hereafter, we consider the “Ising configurations” $|a\rangle = |\uparrow\uparrow\downarrow\cdots\rangle, |\uparrow\downarrow\uparrow\cdots\rangle, \cdots$ (eigenstates of all the $S_i^z$) as a natural basis for this problem. We use this basis to define the IPR $t_n$ of the eigenstates.

#### B. Remarks on degeneracies

In the case of a degenerate multiplet of eigenstates $|n_1\rangle, \dots, |n_d\rangle$, the corresponding contribution to $T$ (as derived from the dynamical point of view [7]) should be $\sum_a (\sum_{\alpha=1}^d |\langle a|n_\alpha \rangle|^2)^2$. For open boundary conditions and a generic value of $\Delta$, it turns out that the only degeneracies in the spectrum of Eq. (5) are the twofold degeneracy ($S_{tot}^z \leftrightarrow -S_{tot}^z$) of the eigenstates with a nonzero magnetization. But since we focus here on the $S_{tot}^z = 0$ sector, we generically get nondegenerate eigenstates which are either even or odd under a global spin flip, as well as even or odd under the spatial (left-right) symmetry.

We finally note that some additional degeneracies occur at the free-fermion point ($\Delta = 0$). For simplicity, we avoid this special point when numerically computing the IPR.

### IV. NUMERICAL RESULTS

#### A. The quantity $T$

The eigenstates of Eq. (5) were obtained using a full numerical diagonalization for even sizes up to $N=20$ in the $S_{tot}^z = 0$ sector. The total Hilbert space dimension is therefore $D = (N/2) \sim 2^N \sqrt{2/(\pi N)}$. Using the spin-flip parity and the spatial parity with respect to the center of the chain, the total space is decomposed into four blocks, each one being labeled by two parities. The results for $T$ are summarized in Fig. 1. Note that $T$ is independent of the sign of $\Delta$. This follows from the fact that changing the sign of $\Delta$ is equivalent to changing the sign of the $xy$ terms, and the latter can be undone by some unitary transformation ($\pi$ rotation about the $z$ axis on every second site) which only affects the sign of the wave functions in the Ising basis, not their modulus. We can therefore restrict ourselves to $\Delta \ge 0$.

#### 1. Gapped phase

We observe (top panel of Fig. 1) that, for sufficiently large $\Delta$, $\ln(T)$ is approximately proportional to the number of sites, i.e.,

$$\ln(T) \approx a(\Delta) N, \qquad (6)$$

with $0 < a(\Delta) < \ln(2)$. We can deduce that the average IPR $\bar{t} = T/D$ behaves as $\ln(\bar{t}) \approx -[\ln(2) - a(\Delta)]N$. We conjecture that these scalings hold in the whole gapped (massive) region ($\Delta > 1$).

At $\Delta = \infty$, any Ising configuration (basis state) becomes an eigenstate of $\mathcal{H}$, so that one may expect $T = D$ and $\ln(T)/N \approx \ln(2)$ in that limit. This simple reasoning is, however, not correct since many energy levels become degenerate in this limit. Thus, as soon as $\Delta$ is not strictly infinite, the actual eigenstates are nontrivial superpositions of several basis states due to the effect of the $S_i^+ S_{i+1}^- + \text{H.c.}$ terms (according to degenerate perturbation theory). We indeed have $\ln(T)/N \approx a(\infty)$ for very large $\Delta$, but the numerically obtained value in this regime, $a(\infty) \approx 0.25$ (see Fig. 1), is significantly smaller than $\ln(2)$. A very similar phenomenon has already been observed in the much simpler situation of a tight-binding particle in a random potential with binary disorder [7], where the hybridization of degenerate molecular states
<figure><img src="image_1.png" /><figcaption>FIG. 1. Quantity T [defined in Eq. (2)], plotted as a function of the system size N (up to N = 20 spins) and of the anisotropy parameter Δ. Different scalings are used: ln(T)/N in the top panel and T/N in the middle. In the top panel, the crosses are the results of a simple extrapolation to N = ∞ using ln(T) ≈ a(Δ)N + b(Δ) + c(Δ)/N and three different system sizes (N = 16, 18 and 20). The bottom panel shows T as a function of N for three values of Δ in the gapless phase. The lines are guides to the eye and show that T is compatible with T ~ N in this regime. The data are computed using all the eigenstates in the $S_{tot}^z = 0$ sector.</figcaption></figure>

has been shown to result in a nontrivial asymptotic return probability $Q \approx 0.373$ in the limit of an infinitely strong disorder.

## 2. Gapless phase

In the gapless phase, T appears to scale approximately as N (bottom panel in Fig. 1). With the present finite-size data, it is, however, not possible to decide whether the exponent is exactly 1. In any case, this implies that $a(\Delta) = 0$ and that the average IPR $\bar{t} = T/D$ scales as $\ln(\bar{t}) \approx -N \ln(2)$ in the whole gapless phase.

The gapped and gapless phases are well known for having different behaviors at zero temperature (correlation functions, etc.). This sharp distinction between $|\Delta| \le 1$ and $|\Delta| > 1$, however, becomes a smooth crossover at finite temperature. Since Eq. (2) involves a sum over all the eigenstates, which is reminiscent of an infinite-temperature quantity, the

observation that T shows qualitatively distinct behavior in both phases is unexpected and remarkable. As discussed in Sec. IV C, it implies a qualitative change in the distribution of the IPR of highly excited eigenstates.

## 3. Remarks concerning the free-fermion point

Finding the scaling of T is a nontrivial question even at the free-fermion point ($\Delta = 0$).

Let us, however, notice that the IPR $t_0$ of the ground state of a periodic chain of N sites (with N even) at $\Delta = 0$ can be computed exactly. It is, indeed, related to the partition function of the Dyson-Gaudin gas at inverse temperature $\beta = 4$ [16], as well as to the Shannon-Rényi entropy for the special value $n = 2$ of the Rényi parameter [17]. The result reads

$$t_0 = \frac{N!}{(N/2)!(2N)^{N/2}} \quad (7)$$

and gives for large N

$$\ln(t_0) \approx -\frac{N}{2} + \frac{\ln(2)}{2}. \quad (8)$$

We are thus facing an explicit example of an eigenstate with $\alpha = 1/2$ (see below).

From another perspective, an approximate determination of the quantity T at the free-fermion point can be obtained by replacing the two matrices defining the Bogoliubov transformation which diagonalizes the quadratic fermionic Hamiltonian with two independent random orthogonal matrices. Skipping every detail, let us mention that this approach yields $T_{SO(N)} = \frac{1}{2}(\det(1 + \Omega)^2)$, where the brackets denote an average over the orthogonal matrix $\Omega$ with the uniform (Haar) measure on SO(N). The known value of the above average [18] leads to the simple expression

$$T_{SO(N)} = N + 1. \quad (9)$$

In spite of its approximate character, this result corroborates the observed scaling $T \sim N$ in the gapless phase.

## B. Distribution of the IPR of individual eigenstates

In order to understand the different scalings observed in the previous section, it is instructive to analyze how the IPRs $t_n$ are distributed over the energy spectrum. This information is represented in Figs. 2 and 3 for different values of $\Delta$. Figure 2 corresponds to a 12-site chain, and each individual eigenstate is represented by a cross. Figure 3 corresponds to 18 spins, and there the density of states is represented by a color scale. For $\Delta$ close to zero it appears that $t_n$ is weakly correlated with energy and that the IPRs of all the eigenstates are of the same order of magnitude. For $N = 16$ and $\Delta = 0.1$, $t_n$ range from 0.00022 to 0.0016, which corresponds to only a factor 7.

For a generic delocalized state, the IPR $t_n$ is expected to scale as $\exp(-\alpha_n N)$. For instance, Eq. (8) gives an example of a state with $\alpha = 1/2$. It therefore seems natural to introduce the quantity

$$\alpha_n = -\ln(t_n)/N \quad (10)$$

in order to compare the distributions of $t_n$ for systems with different sizes. In the example above ($\Delta = 0.1$),
<figure><img src="image_3.png" /><figcaption>FIG. 1. Quantity T [defined in Eq. (2)], plotted as a function of the system size N (up to N = 20 spins) and of the anisotropy parameter Δ. Different scalings are used: ln(T)/N in the top panel and T/N in the middle. In the top panel, the crosses are the results of a simple extrapolation to N = ∞ using ln(T) ≈ a(Δ)N + b(Δ) + c(Δ)/N and three different system sizes (N = 16, 18 and 20). The bottom panel shows T as a function of N for three values of Δ in the gapless phase. The lines are guides to the eye and show that T is compatible with T ~ N in this regime. The data are computed using all the eigenstates in the $S_{tot}^z = 0$ sector.</figcaption></figure>

has been shown to result in a nontrivial asymptotic return probability $Q \approx 0.373$ in the limit of an infinitely strong disorder.

## 2. Gapless phase

In the gapless phase, T appears to scale approximately as N (bottom panel in Fig. 1). With the present finite-size data, it is, however, not possible to decide whether the exponent is exactly 1. In any case, this implies that $a(\Delta) = 0$ and that the average IPR $\bar{t} = T/D$ scales as $\ln(\bar{t}) \approx -N \ln(2)$ in the whole gapless phase.

The gapped and gapless phases are well known for having different behaviors at zero temperature (correlation functions, etc.). This sharp distinction between $|\Delta| \le 1$ and $|\Delta| > 1$, however, becomes a smooth crossover at finite temperature. Since Eq. (2) involves a sum over all the eigenstates, which is reminiscent of an infinite-temperature quantity, the

observation that T shows qualitatively distinct behavior in both phases is unexpected and remarkable. As discussed in Sec. IV C, it implies a qualitative change in the distribution of the IPR of highly excited eigenstates.

## 3. Remarks concerning the free-fermion point

Finding the scaling of T is a nontrivial question even at the free-fermion point ($\Delta = 0$).

Let us, however, notice that the IPR $t_0$ of the ground state of a periodic chain of N sites (with N even) at $\Delta = 0$ can be computed exactly. It is, indeed, related to the partition function of the Dyson-Gaudin gas at inverse temperature $\beta = 4$ [16], as well as to the Shannon-Rényi entropy for the special value $n = 2$ of the Rényi parameter [17]. The result reads

$$t_0 = \frac{N!}{(N/2)!(2N)^{N/2}} \quad (7)$$

and gives for large N

$$\ln(t_0) \approx -\frac{N}{2} + \frac{\ln(2)}{2}. \quad (8)$$

We are thus facing an explicit example of an eigenstate with $\alpha = 1/2$ (see below).

From another perspective, an approximate determination of the quantity T at the free-fermion point can be obtained by replacing the two matrices defining the Bogoliubov transformation which diagonalizes the quadratic fermionic Hamiltonian with two independent random orthogonal matrices. Skipping every detail, let us mention that this approach yields $T_{SO(N)} = \frac{1}{2}(\det(1 + \Omega)^2)$, where the brackets denote an average over the orthogonal matrix $\Omega$ with the uniform (Haar) measure on SO(N). The known value of the above average [18] leads to the simple expression

$$T_{SO(N)} = N + 1. \quad (9)$$

In spite of its approximate character, this result corroborates the observed scaling $T \sim N$ in the gapless phase.

## B. Distribution of the IPR of individual eigenstates

In order to understand the different scalings observed in the previous section, it is instructive to analyze how the IPRs $t_n$ are distributed over the energy spectrum. This information is represented in Figs. 2 and 3 for different values of $\Delta$. Figure 2 corresponds to a 12-site chain, and each individual eigenstate is represented by a cross. Figure 3 corresponds to 18 spins, and there the density of states is represented by a color scale. For $\Delta$ close to zero it appears that $t_n$ is weakly correlated with energy and that the IPRs of all the eigenstates are of the same order of magnitude. For $N = 16$ and $\Delta = 0.1$, $t_n$ range from 0.00022 to 0.0016, which corresponds to only a factor 7.

For a generic delocalized state, the IPR $t_n$ is expected to scale as $\exp(-\alpha_n N)$. For instance, Eq. (8) gives an example of a state with $\alpha = 1/2$. It therefore seems natural to introduce the quantity

$$\alpha_n = -\ln(t_n)/N \quad (10)$$

in order to compare the distributions of $t_n$ for systems with different sizes. In the example above ($\Delta = 0.1$),
<figure><img src="image_2.png" /><figcaption>FIG. 2. IPR $t_n$ versus the eigenstate energy $E_n$ for all the eigenstates with $S_{tot}^z = 0$ in a 12-site chain. Top to bottom: $\Delta = 0.1$, 1, and 3. The energy is weakly correlated with $t_n$ in the gapless phase. On the other hand, for large $\Delta$, the states at the edges of the spectrum appear to be more localized (larger $t_n$) than those in the middle of the spectrum.</figcaption></figure>

<figure><img src="image_3.png" /><figcaption>FIG. 3. Color plot of the density of states as a function of the IPR $t_n$ and the eigenstate energy $E_n$ for all the eigenstates with $S_{tot}^z = 0$ in an 18-site chain. Top to bottom: $\Delta = 0.1$, 1, and 3.</figcaption></figure>

$\alpha$ goes from $\alpha_{min} = -\ln(0.0016)/16 \approx 0.402$ to $\alpha_{max} = -\ln(0.00022)/16 \approx 0.526$.

In order to get an idea of the distribution of $t_n$, we propose sorting them in decreasing order and look for the number of states $D_{1/2}$ one has to include in the sum (starting from the largest $t_n$) in order to get one half of $T$:

$$ \sum_{n=1}^{D_{1/2}} t_n \approx \frac{T}{2}. \qquad (11) $$

We find $D_{1/2}/D \approx 0.42$ for $\Delta = 0.1$ and 0.36 for $\Delta = 0.5$ (see Fig. 4). Importantly, these ratios are relatively stable when varying the system size from $N = 12$ to $N = 20$. For this reason we conjecture that the ratio $D_{1/2}/D$ remains finite in the thermodynamic limit in the whole gapless regime.

The situation turns out to be qualitatively different for large $\Delta$ (see the bottom panels of Figs. 2 and 3, corresponding to $\Delta = 3$). There, the eigenstates with extremal energies (close to the ground state or to the highest excited state) are significantly more “localized” (larger $t_n$) than those living in the middle of

<figure><img src="image_4.png" /><figcaption>FIG. 4. Proportion $D_{1/2}/D$ of states needed to get $T/2$ [see Eq. (11)] as a function of $\Delta$ and for different system sizes $N = 12, 14, 16, 18$, and 20. The calculation is restricted to $S_{tot}^z = 0$. This proportion weakly depends on $N$ for $\Delta < 1$, whereas it decreases significantly with $N$ in the massive phase.</figcaption></figure>
<figure><img src="image_4.png" /><figcaption>FIG. 2. IPR $t_n$ versus the eigenstate energy $E_n$ for all the eigenstates with $S_{tot}^z = 0$ in a 12-site chain. Top to bottom: $\Delta = 0.1$, 1, and 3. The energy is weakly correlated with $t_n$ in the gapless phase. On the other hand, for large $\Delta$, the states at the edges of the spectrum appear to be more localized (larger $t_n$) than those in the middle of the spectrum.</figcaption></figure>

<figure><img src="image_5.png" /><figcaption>FIG. 3. Color plot of the density of states as a function of the IPR $t_n$ and the eigenstate energy $E_n$ for all the eigenstates with $S_{tot}^z = 0$ in an 18-site chain. Top to bottom: $\Delta = 0.1$, 1, and 3.</figcaption></figure>

$\alpha$ goes from $\alpha_{min} = -\ln(0.0016)/16 \approx 0.402$ to $\alpha_{max} = -\ln(0.00022)/16 \approx 0.526$.

In order to get an idea of the distribution of $t_n$, we propose sorting them in decreasing order and look for the number of states $D_{1/2}$ one has to include in the sum (starting from the largest $t_n$) in order to get one half of $T$:

$$ \sum_{n=1}^{D_{1/2}} t_n \approx \frac{T}{2}. \qquad (11) $$

We find $D_{1/2}/D \approx 0.42$ for $\Delta = 0.1$ and 0.36 for $\Delta = 0.5$ (see Fig. 4). Importantly, these ratios are relatively stable when varying the system size from $N = 12$ to $N = 20$. For this reason we conjecture that the ratio $D_{1/2}/D$ remains finite in the thermodynamic limit in the whole gapless regime.

The situation turns out to be qualitatively different for large $\Delta$ (see the bottom panels of Figs. 2 and 3, corresponding to $\Delta = 3$). There, the eigenstates with extremal energies (close to the ground state or to the highest excited state) are significantly more “localized” (larger $t_n$) than those living in the middle of

<figure><img src="image_6.png" /><figcaption>FIG. 4. Proportion $D_{1/2}/D$ of states needed to get $T/2$ [see Eq. (11)] as a function of $\Delta$ and for different system sizes $N = 12, 14, 16, 18$, and 20. The calculation is restricted to $S_{tot}^z = 0$. This proportion weakly depends on $N$ for $\Delta < 1$, whereas it decreases significantly with $N$ in the massive phase.</figcaption></figure>
the spectrum (where the density of states is maximal). For $N = 16$ and $\Delta = 3$, $t_n$ range from 0.0016 to 0.47, which now corresponds to a factor 300. For this value of $\Delta$, the proportion of states required to get one half of $T$ is $D_{1/2}/D = 0.19$ for $N = 12$, 0.14 for $N = 16$, and 0.10 for $N = 20$. (see Fig. 4). This steady decrease suggests that the ratio $D_{1/2}/D$ vanishes in the thermodynamic limit, so that an infinitesimal proportion of all the eigenstates accounts for the major contribution to $T$. We conjecture that this is the case in the whole massive phase.

### C. Entropy function of the IPR

#### 1. Generalities on entropy function

The number $W(\alpha, \delta)$ of states with a value of $\alpha$ in the interval $[\alpha, \alpha + \delta]$ is expected to scale exponentially with the system size, in analogy with a conventional microcanonical density of states. This leads us to define an entropy $s(\alpha)$, such that

$$W(\alpha, \delta) \sim e^{Ns(\alpha)}. \qquad (12)$$

For sufficiently large $N$, the entropy defined in this way should no longer depend on the “bin size” $\delta$. As we will see, this function is the proper way to describe the distribution of $t_n$ in the thermodynamic limit. This entropy function can also be viewed as a large deviation function and is similar to the quantities used in multifractal analysis.

Some numerical results for this quantity are plotted in Fig. 5. Although we observe significant finite-size effects, the picture that emerges is the following:

(i) The entropy converges to a well-defined intensive function $s(\alpha)$ on an interval $[\alpha_{min}, \alpha_{max}]$ in the thermodynamic limit. A precise determination of this support is, however, not simple from the available data due to the small number of eigenstates close to the extremal values of $\alpha$.

<figure><img src="image_5.png" /><figcaption>FIG. 5. Entropy $s(\alpha)$ [see Eq. (12)] for two system sizes and two values of $\Delta$. Here the bin size is 1/200. The dotted lines have slope 1 and are guides to the eye to locate the saddle point $s'(\alpha_0) = 1$. In the thermodynamic limit we expect the maximum of the entropy to converge to $s = \ln(2)$. In addition, we conjecture that, in the thermodynamic limit and in the gapless phase, this maximum is attained at $\alpha_0 = \ln(2)$, where a cusp should develop (see text and Fig. 6).</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>FIG. 6. Schematic representation of the entropy function, as conjectured from the numerical analysis. Top: $|\Delta| > 1$ (gapped phase). Bottom: $|\Delta| < 1$ (gapless phase).</figcaption></figure>

(ii) The entropy $s(\alpha)$ has a maximum for some $\alpha = \alpha_{typ}$, with $\alpha_{min} < \alpha_{typ} < \alpha_{max}$, where these three quantities depend on $\Delta$. From a saddle-point evaluation of the total number of states

$$D \sim \int e^{Ns(\alpha)} d\alpha \sim 2^N, \qquad (13)$$

we predict that the maximum $s(\alpha_{typ})$ of the entropy converges to $\ln(2)$ for sufficiently large systems. The numerical values for $N = 18$ and $N = 20$ (see Fig. 5) are smaller than $\ln(2)$, but they significantly increase with $N$. Moreover, $\alpha_{typ}$ corresponds to the typical IPR: $t_{typ} \sim \exp(-\alpha_{typ}N)$.

The observations made previously on $T$ can be related to properties of $s(\alpha)$. The quantity $T$ can be estimated by evaluating the integral

$$T \sim \int e^{N(s(\alpha) - \alpha)} d\alpha \qquad (14)$$
the spectrum (where the density of states is maximal). For $N = 16$ and $\Delta = 3$, $t_n$ range from 0.0016 to 0.47, which now corresponds to a factor 300. For this value of $\Delta$, the proportion of states required to get one half of $T$ is $D_{1/2}/D = 0.19$ for $N = 12$, 0.14 for $N = 16$, and 0.10 for $N = 20$. (see Fig. 4). This steady decrease suggests that the ratio $D_{1/2}/D$ vanishes in the thermodynamic limit, so that an infinitesimal proportion of all the eigenstates accounts for the major contribution to $T$. We conjecture that this is the case in the whole massive phase.

### C. Entropy function of the IPR

#### 1. Generalities on entropy function

The number $W(\alpha, \delta)$ of states with a value of $\alpha$ in the interval $[\alpha, \alpha + \delta]$ is expected to scale exponentially with the system size, in analogy with a conventional microcanonical density of states. This leads us to define an entropy $s(\alpha)$, such that

$$W(\alpha, \delta) \sim e^{Ns(\alpha)}. \qquad (12)$$

For sufficiently large $N$, the entropy defined in this way should no longer depend on the “bin size” $\delta$. As we will see, this function is the proper way to describe the distribution of $t_n$ in the thermodynamic limit. This entropy function can also be viewed as a large deviation function and is similar to the quantities used in multifractal analysis.

Some numerical results for this quantity are plotted in Fig. 5. Although we observe significant finite-size effects, the picture that emerges is the following:

(i) The entropy converges to a well-defined intensive function $s(\alpha)$ on an interval $[\alpha_{min}, \alpha_{max}]$ in the thermodynamic limit. A precise determination of this support is, however, not simple from the available data due to the small number of eigenstates close to the extremal values of $\alpha$.

<figure><img src="image_7.png" /><figcaption>FIG. 5. Entropy $s(\alpha)$ [see Eq. (12)] for two system sizes and two values of $\Delta$. Here the bin size is 1/200. The dotted lines have slope 1 and are guides to the eye to locate the saddle point $s'(\alpha_0) = 1$. In the thermodynamic limit we expect the maximum of the entropy to converge to $s = \ln(2)$. In addition, we conjecture that, in the thermodynamic limit and in the gapless phase, this maximum is attained at $\alpha_0 = \ln(2)$, where a cusp should develop (see text and Fig. 6).</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>FIG. 6. Schematic representation of the entropy function, as conjectured from the numerical analysis. Top: $|\Delta| > 1$ (gapped phase). Bottom: $|\Delta| < 1$ (gapless phase).</figcaption></figure>

(ii) The entropy $s(\alpha)$ has a maximum for some $\alpha = \alpha_{typ}$, with $\alpha_{min} < \alpha_{typ} < \alpha_{max}$, where these three quantities depend on $\Delta$. From a saddle-point evaluation of the total number of states

$$D \sim \int e^{Ns(\alpha)} d\alpha \sim 2^N, \qquad (13)$$

we predict that the maximum $s(\alpha_{typ})$ of the entropy converges to $\ln(2)$ for sufficiently large systems. The numerical values for $N = 18$ and $N = 20$ (see Fig. 5) are smaller than $\ln(2)$, but they significantly increase with $N$. Moreover, $\alpha_{typ}$ corresponds to the typical IPR: $t_{typ} \sim \exp(-\alpha_{typ}N)$.

The observations made previously on $T$ can be related to properties of $s(\alpha)$. The quantity $T$ can be estimated by evaluating the integral

$$T \sim \int e^{N(s(\alpha) - \alpha)} d\alpha \qquad (14)$$
by means of a saddle-point approximation. The saddle point $\alpha_0$ is defined as the solution of $s'(\alpha_0) = 1$ and gives $T \sim \exp\{N[s(\alpha_0) - \alpha_0]\}$. In other words, $a(\Delta)$ [defined in Eq. (6)] is given by $a(\Delta) = s(\alpha_0) - \alpha_0$.

## 2. Gapless phase

The numerical data in Figs. 1 and 4 led us to conjecture that $a(\Delta) = 0$ and that $D_{1/2}/D$ reaches a finite limit in the gapless phase. These observations can be translated as follows in terms of $s(\alpha)$. By definition, to get one half of $T$, one needs to sum over the $D_{1/2}$ states with the highest $t_n$ (lowest $\alpha_n$). This amounts to integrating over $\alpha$ only up to some value $\tilde{\alpha}$. Within the saddle-point approximation, it is easy to see that we have, in fact, $\tilde{\alpha} = \alpha_0$. This implies $D_{1/2} \sim \exp[Ns(\alpha_0)]$. The ratio $D_{1/2}/D$ thus scales as $\exp\{N[s(\alpha_0) - \ln(2)]\}$. Thus, from the fact that $D_{1/2}/D$ appears to be finite in the gapless phase, we conclude that $s(\alpha_0) = \ln(2)$. Since $\ln(2)$ is the maximal value for the entropy function, we also have $\alpha_0 = \alpha_{typ}$. However, by definition, $s'(\alpha_0) = 1$. Thus, the entropy function has to have a cusp at its maximum. The fact that $\alpha_{typ} = \ln(2)$ also implies that the IPR of typical eigenstates is close to that of the maximally delocalized states. Finally, this also implies that $a(\Delta) = 0$, which is consistent with the scaling of $T$ as a power of $N$ in the gapless phase, as already discussed. The scenario

we propose for the gapless phases is summarized in the bottom panel of Fig. 6.

## 3. Gapped phase

Repeating the above arguments in the gapped phase, the observed scaling for $D_{1/2}/D$ leads to $s(\alpha_0) < \ln(2)$. The corresponding scenario for the entropy function is summarized in the top panel of Fig. 6, with a strictly positive value for $a(\Delta) = s(\alpha_0) - \alpha_0$ accounting for the exponential growth of $T$. It is interesting to note that, although $T$ grows exponentially with the system size, virtually all eigenstates have an exponentially small IPR, that is, a delocalized character.

## V. BEYOND THE XXZ CHAIN

In order to test the robustness of the results we have obtained for the XXZ chain, we consider some perturbation of the model by including second-neighbor interactions. The Hamiltonian is now defined as follows:

$$ \mathcal{H} = \sum_{i=1}^{N-1} (S_i^x S_{i+1}^x + S_i^y S_{i+1}^y + \Delta S_i^z S_{i+1}^z) + J_2 \sum_{i=1}^{N-2} (S_i^x S_{i+2}^x + S_i^y S_{i+2}^y + \Delta S_i^z S_{i+2}^z). \quad (15) $$

<figure><img src="image_7.png" /><figcaption>FIG. 7. Same as Fig. 1 for a spin chain with second-neighbor interactions [see Eq. (15)]. From left to right: $J_2 = 0.05, 0.15$, and $0.2$. The crosses in the top panels represent some extrapolations based on $N = 16, 18$, and $20$, as in Fig. 1.</figcaption></figure>
by means of a saddle-point approximation. The saddle point $\alpha_0$ is defined as the solution of $s'(\alpha_0) = 1$ and gives $T \sim \exp\{N[s(\alpha_0) - \alpha_0]\}$. In other words, $a(\Delta)$ [defined in Eq. (6)] is given by $a(\Delta) = s(\alpha_0) - \alpha_0$.

## 2. Gapless phase

The numerical data in Figs. 1 and 4 led us to conjecture that $a(\Delta) = 0$ and that $D_{1/2}/D$ reaches a finite limit in the gapless phase. These observations can be translated as follows in terms of $s(\alpha)$. By definition, to get one half of $T$, one needs to sum over the $D_{1/2}$ states with the highest $t_n$ (lowest $\alpha_n$). This amounts to integrating over $\alpha$ only up to some value $\tilde{\alpha}$. Within the saddle-point approximation, it is easy to see that we have, in fact, $\tilde{\alpha} = \alpha_0$. This implies $D_{1/2} \sim \exp[Ns(\alpha_0)]$. The ratio $D_{1/2}/D$ thus scales as $\exp\{N[s(\alpha_0) - \ln(2)]\}$. Thus, from the fact that $D_{1/2}/D$ appears to be finite in the gapless phase, we conclude that $s(\alpha_0) = \ln(2)$. Since $\ln(2)$ is the maximal value for the entropy function, we also have $\alpha_0 = \alpha_{typ}$. However, by definition, $s'(\alpha_0) = 1$. Thus, the entropy function has to have a cusp at its maximum. The fact that $\alpha_{typ} = \ln(2)$ also implies that the IPR of typical eigenstates is close to that of the maximally delocalized states. Finally, this also implies that $a(\Delta) = 0$, which is consistent with the scaling of $T$ as a power of $N$ in the gapless phase, as already discussed. The scenario

we propose for the gapless phases is summarized in the bottom panel of Fig. 6.

## 3. Gapped phase

Repeating the above arguments in the gapped phase, the observed scaling for $D_{1/2}/D$ leads to $s(\alpha_0) < \ln(2)$. The corresponding scenario for the entropy function is summarized in the top panel of Fig. 6, with a strictly positive value for $a(\Delta) = s(\alpha_0) - \alpha_0$ accounting for the exponential growth of $T$. It is interesting to note that, although $T$ grows exponentially with the system size, virtually all eigenstates have an exponentially small IPR, that is, a delocalized character.

## V. BEYOND THE XXZ CHAIN

In order to test the robustness of the results we have obtained for the XXZ chain, we consider some perturbation of the model by including second-neighbor interactions. The Hamiltonian is now defined as follows:

$$ \mathcal{H} = \sum_{i=1}^{N-1} (S_i^x S_{i+1}^x + S_i^y S_{i+1}^y + \Delta S_i^z S_{i+1}^z) + J_2 \sum_{i=1}^{N-2} (S_i^x S_{i+2}^x + S_i^y S_{i+2}^y + \Delta S_i^z S_{i+2}^z). \quad (15) $$

<figure><img src="image_9.png" /><figcaption>FIG. 7. Same as Fig. 1 for a spin chain with second-neighbor interactions [see Eq. (15)]. From left to right: $J_2 = 0.05, 0.15$, and $0.2$. The crosses in the top panels represent some extrapolations based on $N = 16, 18$, and $20$, as in Fig. 1.</figcaption></figure>
The parameter $J_2$ allows one to tune the strength of the perturbation.

For some detailed discussion about the rich phase diagram of this model, we refer the reader to Refs. [19,20] and references therein, but we simply focus here on the regime $0 \le J_2 \lesssim 0.24$. There, the ground state turns out to be a gapless TLL for $\Delta \le 1$, while it is a gapped state (Ising phase) for $\Delta > 1$. In other words, for sufficiently small $J_2$, the zero-temperature phase diagram is similar to the $J_2 = 0$ case. Note, however, that the model is not integrable for $J_2 \ne 0$.

In the TLL phase it is conventional to parametrize the long-distance properties of the correlation functions by the so-called Luttinger parameter $K$. Using a bosonization approach, one can show that the staggered part of $\langle S_0^+ S_r^- \rangle$ then behaves as $\sim (-1)^r r^{-1/2}$, and that of $\langle S_0^z S_r^z \rangle$ behaves as $\sim (-1)^r r^{-2K}$ [15]. $K$ is a priori some nontrivial function of the microscopic parameters of the model, here $\Delta$ and $J_2$. However, the SU(2) symmetry present at $\Delta = 1$ forces the correlations to be isotropic, so that, in fact, $2K = \frac{1}{2K}$ and $K = \frac{1}{2}$. Next, as discussed by Haldane [21], the perturbation (umklapp terms) that drives the transition from the TLL to the Ising phase is marginal (in the renormalization-group sense) when $K = \frac{1}{2}$. For this reason, the transition line between the TLL and the Ising phase lies exactly at $\Delta = 1$.

The numerical results concerning the scaling of $T$ for $J_2 = 0.05, 0.15$, and $0.2$ are displayed in Fig. 7. We first note that the “perturbed” model displays an exponential growth of $T[T \sim \exp(aN)$ with $a > 0]$ in the gapped phase, as for the $J_2 = 0$ case. This is rather clear from the plots in the top panels and from the associated extrapolation points (crosses).

The situation, however, seems different for $\Delta < 1$. The bottom panels of Fig. 7 indeed show that $T$ reaches a maximum and then slightly decreases for larger system sizes. Since $T$ cannot vanish (it is larger than 1 by construction), these data strongly suggest that $T$ is finite in the thermodynamic limit when $J_2 > 0$. How can this be reconciled with the observation that $T$ grows linearly with $N$ when $J_2 = 0$ and $|\Delta| < 1$ (Fig. 1)? The first possibility is that the $N \to \infty$ value of $T$ is finite for the nonintegrable models but diverges as $J_2$ approaches zero. This interesting scenario would imply a strong effect of the integrability of the model on the IPR distribution. We note that the integrability has been shown in Ref. [22] to have an impact on participation ratios, even though that study focused on different basis choices.

An alternative, albeit less probable, possibility is that $T$ also remains finite in the gapless phase of the integrable model ($J_2 = 0$). Although this is not what the data of Fig. 1 (bottom panel) suggest, one cannot exclude that $T$ converges to a finite limit for larger systems. If that were the case, studying larger chains would be needed to observe some saturation of $T$ for the unperturbed XXZ chain. Still, performing a full diagonalization on open chains beyond $N = 20$ sites would be numerically quite challenging. Indeed, for 22 spins at $S_{tot}^z = 0$, the Hilbert space dimension is larger than $1.7 \times 10^5$ in each of the four symmetry sectors (using space inversion and global spin flip).

## VI. CONCLUSIONS

We have investigated the IPRs (inverse participation ratios) of individual energy eigenstates in a preferential Ising basis for a spin-1/2 XXZ chain without disorder by means of an exact diagonalization of the Hamiltonian on finite chains of length $N$ (up to 20), with open boundary conditions. We have considered in particular the sum $T$ of all $t_n$, which has a dynamical interpretation: it yields the stationary return probability to a typical initial state of the preferential basis. Our main finding is the observation of a qualitatively different behavior of the latter quantities in the gapped and gapless phases. In the gapped phase ($|\Delta| > 1$), $T$ grows exponentially with $N$, and the entropy function $s(\alpha)$ describing the distribution of $t_n$ has a smooth maximum. In the gapless phase ($|\Delta| < 1$), $T$ seems to scale linearly with $N$, whereas the entropy function is singular at its maximum. We have also investigated the effect of next-nearest-neighbor interactions, which break the integrability of the model. Although $T$ still grows exponentially in the gapped phase, it now appears to saturate to a constant value in the gapless phase.

In some future work it would be very useful to make some progress concerning the free-fermion point ($\Delta = 0$) as well as in the limit $\Delta \to \infty$. As we explained, as far as the ground-state wave function(s) are concerned, these two limits are rather simple. However, computing $T$ for these models appears quite challenging. More generally, we need some deeper understanding of the qualitatively different scalings observed in the gapped and gapless phases. From this point of view, a Bethe ansatz formulation of the eigenstates seems a promising route to explore, as well as some continuum-limit approach to the problem in the gapless phase. The latter may also shed light on the possible existence of universal terms in the quantity $T$.

[1] A. Polkovnikov, K. Sengupta, A. Silva, and M. Vengalattore, *Rev. Mod. Phys.* **83**, 863 (2011).

[2] S. Goldstein, J. L. Lebowitz, R. Tumulka, and N. Zanghi, *Eur. Phys. J. H* **35**, 173 (2010).

[3] J. Eisert, M. Friesdorf, and C. Gogolin, *Nat. Phys.* **11**, 124 (2015).

[4] M. Rigol, V. Dunjko, and M. Olshanii, *Nature (London)* **452**, 854 (2008).

[5] M. Rigol, V. Dunjko, V. Yurovsky, and M. Olshanii, *Phys. Rev. Lett.* **98**, 050405 (2007).

[6] R. Nandkishore and D. A. Huse, *Annu. Rev. Condens. Matter Phys.* **6**, 15 (2015).

[7] J. M. Luck, *J. Phys. A* **49**, 115303 (2016).

[8] R. J. Bell and P. Dean, *Discuss. Faraday Soc.* **50**, 55 (1970).

[9] J. T. Edwards and D. J. Thouless, *J. Phys. C* **5**, 807 (1972).

[10] F. Evers and A. D. Mirlin, *Rev. Mod. Phys.* **80**, 1355 (2008).

[11] A. De Luca and A. Scardicchio, *Europhys. Lett.* **101**, 37003 (2013).

[12] D. J. Luitz, N. Laflorencie, and F. Alet, *Phys. Rev. B* **91**, 081103 (2015).
The parameter $J_2$ allows one to tune the strength of the perturbation.

For some detailed discussion about the rich phase diagram of this model, we refer the reader to Refs. [19,20] and references therein, but we simply focus here on the regime $0 \le J_2 \lesssim 0.24$. There, the ground state turns out to be a gapless TLL for $\Delta \le 1$, while it is a gapped state (Ising phase) for $\Delta > 1$. In other words, for sufficiently small $J_2$, the zero-temperature phase diagram is similar to the $J_2 = 0$ case. Note, however, that the model is not integrable for $J_2 \ne 0$.

In the TLL phase it is conventional to parametrize the long-distance properties of the correlation functions by the so-called Luttinger parameter $K$. Using a bosonization approach, one can show that the staggered part of $\langle S_0^+ S_r^- \rangle$ then behaves as $\sim (-1)^r r^{-1/2}$, and that of $\langle S_0^z S_r^z \rangle$ behaves as $\sim (-1)^r r^{-2K}$ [15]. $K$ is a priori some nontrivial function of the microscopic parameters of the model, here $\Delta$ and $J_2$. However, the SU(2) symmetry present at $\Delta = 1$ forces the correlations to be isotropic, so that, in fact, $2K = \frac{1}{2K}$ and $K = \frac{1}{2}$. Next, as discussed by Haldane [21], the perturbation (umklapp terms) that drives the transition from the TLL to the Ising phase is marginal (in the renormalization-group sense) when $K = \frac{1}{2}$. For this reason, the transition line between the TLL and the Ising phase lies exactly at $\Delta = 1$.

The numerical results concerning the scaling of $T$ for $J_2 = 0.05, 0.15$, and $0.2$ are displayed in Fig. 7. We first note that the “perturbed” model displays an exponential growth of $T$ [$T \sim \exp(aN)$ with $a > 0$] in the gapped phase, as for the $J_2 = 0$ case. This is rather clear from the plots in the top panels and from the associated extrapolation points (crosses).

The situation, however, seems different for $\Delta < 1$. The bottom panels of Fig. 7 indeed show that $T$ reaches a maximum and then slightly decreases for larger system sizes. Since $T$ cannot vanish (it is larger than 1 by construction), these data strongly suggest that $T$ is finite in the thermodynamic limit when $J_2 > 0$. How can this be reconciled with the observation that $T$ grows linearly with $N$ when $J_2 = 0$ and $|\Delta| < 1$ (Fig. 1)? The first possibility is that the $N \to \infty$ value of $T$ is finite for the nonintegrable models but diverges as $J_2$ approaches zero. This interesting scenario would imply a strong effect of the integrability of the model on the IPR distribution. We note that the integrability has been shown in Ref. [22] to have an impact on participation ratios, even though that study focused on different basis choices.

An alternative, albeit less probable, possibility is that $T$ also remains finite in the gapless phase of the integrable model ($J_2 = 0$). Although this is not what the data of Fig. 1 (bottom panel) suggest, one cannot exclude that $T$ converges to a finite limit for larger systems. If that were the case, studying larger chains would be needed to observe some saturation of $T$ for the unperturbed XXZ chain. Still, performing a full diagonalization on open chains beyond $N = 20$ sites would be numerically quite challenging. Indeed, for 22 spins at $S_{tot}^z = 0$, the Hilbert space dimension is larger than $1.7 \times 10^5$ in each of the four symmetry sectors (using space inversion and global spin flip).

## VI. CONCLUSIONS

We have investigated the IPRs (inverse participation ratios) of individual energy eigenstates in a preferential Ising basis for a spin-1/2 XXZ chain without disorder by means of an exact diagonalization of the Hamiltonian on finite chains of length $N$ (up to 20), with open boundary conditions. We have considered in particular the sum $T$ of all $t_n$, which has a dynamical interpretation: it yields the stationary return probability to a typical initial state of the preferential basis. Our main finding is the observation of a qualitatively different behavior of the latter quantities in the gapped and gapless phases. In the gapped phase ($|\Delta| > 1$), $T$ grows exponentially with $N$, and the entropy function $s(\alpha)$ describing the distribution of $t_n$ has a smooth maximum. In the gapless phase ($|\Delta| < 1$), $T$ seems to scale linearly with $N$, whereas the entropy function is singular at its maximum. We have also investigated the effect of next-nearest-neighbor interactions, which break the integrability of the model. Although $T$ still grows exponentially in the gapped phase, it now appears to saturate to a constant value in the gapless phase.

In some future work it would be very useful to make some progress concerning the free-fermion point ($\Delta = 0$) as well as in the limit $\Delta \to \infty$. As we explained, as far as the ground-state wave function(s) are concerned, these two limits are rather simple. However, computing $T$ for these models appears quite challenging. More generally, we need some deeper understanding of the qualitatively different scalings observed in the gapped and gapless phases. From this point of view, a Bethe ansatz formulation of the eigenstates seems a promising route to explore, as well as some continuum-limit approach to the problem in the gapless phase. The latter may also shed light on the possible existence of universal terms in the quantity $T$.

[1] A. Polkovnikov, K. Sengupta, A. Silva, and M. Vengalattore, *Rev. Mod. Phys.* **83**, 863 (2011).

[2] S. Goldstein, J. L. Lebowitz, R. Tumulka, and N. Zanghi, *Eur. Phys. J. H* **35**, 173 (2010).

[3] J. Eisert, M. Friesdorf, and C. Gogolin, *Nat. Phys.* **11**, 124 (2015).

[4] M. Rigol, V. Dunjko, and M. Olshanii, *Nature (London)* **452**, 854 (2008).

[5] M. Rigol, V. Dunjko, V. Yurovsky, and M. Olshanii, *Phys. Rev. Lett.* **98**, 050405 (2007).

[6] R. Nandkishore and D. A. Huse, *Annu. Rev. Condens. Matter Phys.* **6**, 15 (2015).

[7] J. M. Luck, *J. Phys. A* **49**, 115303 (2016).

[8] R. J. Bell and P. Dean, *Discuss. Faraday Soc.* **50**, 55 (1970).

[9] J. T. Edwards and D. J. Thouless, *J. Phys. C* **5**, 807 (1972).

[10] F. Evers and A. D. Mirlin, *Rev. Mod. Phys.* **80**, 1355 (2008).

[11] A. De Luca and A. Scardicchio, *Europhys. Lett.* **101**, 37003 (2013).

[12] D. J. Luitz, N. Laflorencie, and F. Alet, *Phys. Rev. B* **91**, 081103 (2015).
[13] J.-M. Stéphan, G. Misguich, and V. Pasquier, Phys. Rev. B **84**, 195128 (2011).

[14] F. Dukesz, M. Zilbergerts, and L. F. Santos, New J. Phys. **11**, 043026 (2009).

[15] T. Giamarchi, *Quantum Physics in One Dimension* (Clarendon, Oxford, 2004).

[16] M. Gaudin, J. Phys. (Paris) **34**, 511 (1973).

[17] J.-M. Stéphan, S. Furukawa, G. Misguich, and V. Pasquier, Phys. Rev. B **80**, 184421 (2009).

[18] J. P. Keating and N. C. Snaith, Commun. Math. Phys. **214**, 91 (2000).

[19] T. Hikihara, M. Kaburagi, and H. Kawamura, Phys. Rev. B **63**, 174430 (2001).

[20] S. Furukawa, M. Sato, S. Onoda, and A. Furusaki, Phys. Rev. B **86**, 094417 (2012).

[21] F. D. M. Haldane, Phys. Rev. B **25**, 4925 (1982).

[22] L. F. Santos and M. Rigol, Phys. Rev. E **81**, 036206 (2010).
[13] J.-M. Stéphan, G. Misguich, and V. Pasquier, Phys. Rev. B **84**, 195128 (2011).

[14] F. Dukesz, M. Zilbergerts, and L. F. Santos, New J. Phys. **11**, 043026 (2009).

[15] T. Giamarchi, *Quantum Physics in One Dimension* (Clarendon, Oxford, 2004).

[16] M. Gaudin, J. Phys. (Paris) **34**, 511 (1973).

[17] J.-M. Stéphan, S. Furukawa, G. Misguich, and V. Pasquier, Phys. Rev. B **80**, 184421 (2009).

[18] J. P. Keating and N. C. Snaith, Commun. Math. Phys. **214**, 91 (2000).

[19] T. Hikihara, M. Kaburagi, and H. Kawamura, Phys. Rev. B **63**, 174430 (2001).

[20] S. Furukawa, M. Sato, S. Onoda, and A. Furusaki, Phys. Rev. B **86**, 094417 (2012).

[21] F. D. M. Haldane, Phys. Rev. B **25**, 4925 (1982).

[22] L. F. Santos and M. Rigol, Phys. Rev. E **81**, 036206 (2010).