Comparing many-body localization lengths via nonperturbative construction of local integrals of motion

Many-body localization (MBL), characterized by the absence of thermalization and the violation of conventional thermodynamics, has elicited much interest both as a fundamental physical phenomenon and for practical applications in quantum information. A phenomenological model which describes the system using a complete set of local integrals of motion (LIOMs) provides a powerful tool to understand MBL but can usually be computed only approximately. Here we explicitly compute a complete set of LIOMs with a nonperturbative approach by maximizing the overlap between LIOMs and physical spin operators in real space. The set of LIOMs satisfies the desired exponential decay of the weight of LIOMs in real space. This LIOM construction enables a direct mapping from the real-space Hamiltonian to the phenomenological model and thus enables studying the localized Hamiltonian and the system dynamics. We can thus study and compare the localization lengths extracted from the LIOM weights, their interactions, and dephasing dynamics, revealing interesting aspects of many-body localization. Our scheme is immune to accidental resonances and can be applied even at the phase transition point, providing a tool to study the microscopic features of the phenomenological model of MBL.

Figure 1

The flow diagram shows an example of the construction of a complete set of LIOMs in a system with $L=3$. Gray blocks represent undetermined matrix elements, and orange (green) blocks represent $+1$ ($-1$) matrix elements. The procedure works as follows: after diagonalizing the Hamiltonian, find the candidate LIOM ${\tau }_{j_M}$ that maximizes the overlap with the physical spin operator $\langle {\tilde{\tau }}_z^j{\sigma }_z^j\rangle$ ($j_M=2$ here). Divide the eight eigenstates into two sectors, each containing four states according to $\langle n\left|{\sigma }_z^2\right|n\rangle$, and assign ${\tau }_z^2={\tilde{\tau }}_z^2$. For each sector, find $j_M$ within the sector, divide into two sectors each containing two states, and assign ${\tau }_z^{j_M}={\tilde{\tau }}_z^{j_M}$. Repeat the step one more time, and then all LIOMs are determined.

Figure 2

(a) and (b) Median of the LIOM weights $|f_{n,k}^j|$ as a function of distance $n$ for two disorder strengths: (a) $W=20$, deep in the MBL phase, where the median decays exponentially to zero, and (b) $W=1$, in the ergodic phase, where the median saturates to a nonzero value. For each $j$, the median is taken over the index $k$ in $|f_{n,k}^j|$ as well as 20 different disorder realizations. The darker color represents LIOMs in the middle of the chain, $j=L/2$ [as shown in the bottom of (b)], and left (right) half of the LIOMs is represented by dashed blue curves. (c) and (d) Median of the interaction strength as a function of range $r$ for two disorder strengths. Dotted curves represent $l$-body interaction terms $|V_{ij}|,|V_{ijk}|,\cdots$ ($l=2,\cdots ,9$), where the median is taken over all indices $i,j,\cdots$, as well as 100 disorder realizations. The solid curve represents the median of all interaction terms for a given range $V\left(r\right)$, regardless of how many LIOMs are involved. $L=10$ in all subplots.

Figure 3

(a) Probability distribution of LIOM weights ${log}_{10}\left|f_{n,k}^j\right|$. For a given distance $n$, samples are taken from all possible $j$ and $k$ as well as 200 disorder realizations. The distribution shows one single Gaussian peak that shifts toward smaller weights with increasing distance $n$, signaling the localization of IOMs. (b) Probability distribution of the interaction strength ${log}_{10}\left(\right|V\left|\right)$. For a given range $r$, samples are taken from all terms in Eq. (3) as well as 10 000 disorder realizations. Two peaks can be observed: the left peak is due to the localized cases as it shifts to smaller interaction strengths for longer range; the right peak shows the delocalized cases (rare events) as it is independent of interaction range. $L=10$ and $W=20$ for both (a) and (b).

Figure 4

(a) Dephasing of the physical spin operator ${\sigma }_x^L$ (dark green dashed curve) and LIOM ${\tau }_x^L$ (green solid curve). The initial state is a product state with each spin pointing randomly in the $xy$ plane, i.e., $|\psi \left(0\right)\rangle ={\otimes }_{j=1}^L{\left(\right|+\rangle }_j+e^{i\varphi }{|-\rangle }_j)/\sqrt{2}$, with $\varphi$ randomly sampled in $[0,2\pi ), {\sigma }_z^j{|+\rangle }_j={|+\rangle }_j{\sigma }_z^j{|-\rangle }_j=-{|-\rangle }_j$ for the red curve and ${\tau }_z^j{|+\rangle }_j={|+\rangle }_j{\tau }_z^j{|-\rangle }_j=-{|-\rangle }_j$ for the blue curve. $L=10, W=20$. Averaging is performed over 20 different initial states and 20 disorder realizations. The error bar represents the standard deviation of all configurations. (b) and (c) Localization lengths as a function of disorder strength $W$ for $L=12$. The LIOM localization length $\xi$ is extracted by fitting $\mathrm{Tr}\left({\tau }_z^j{\sigma }_z^k\right)\sim exp(-|k-j|/\xi )$ as a function of $k$ (for $j=1$). The interaction localization length is obtained from the fit of the interaction as a function of range, $V\left(r\right)\sim exp(-r/\kappa )$. The dynamical localization length describes the LIOM dephasing shown in (a) and is obtained by the fit to $\langle {\langle {\tau }_x^L\rangle }^2\rangle \sim t^{-{\xi }^{\prime }ln2}$. Here the error bars derive from the fitting error. The dephasing curve used to fit ${\xi }^{\prime }$ is extracted from the median of 50 disorder realizations and 50 initial states. $\xi$ and $\kappa$ are extracted from the median of 5000 disorder realizations. (b) is a zoom of (c) near the transition point.

Figure 5

LIOM in the noninteracting model. (a) and (b) Median interaction strength in the LIOM basis $\left\{{\tau }_z^j\right\}$ vs range $r$. $J_z=0$ (blue dashed line) corresponds to a noninteracting model, and $J_z=0.1$ (red line) corresponds to a model with a small interaction. We note that for large disorder the noninteracting model has only nonzero $V$ for the range $r=1$, which denotes the single-particle Hamiltonian ${\xi }_j{\tau }_z^j$, and even for small disorder $V(r=1)$ is the dominant term. This is in contrast to the behavior of the interacting model. (c) and (d) Median overlap between LIOMs and physical spins $\mathrm{Tr}\left(\widehat{\mathcal{O}}{\sigma }_z^k\right)$ for the noninteracting model $J_z=1$, with $\widehat{\mathcal{O}}={\mathrm{\Sigma }}_z^j$ (solid lines) for single-particle LIOM and $\widehat{\mathcal{O}}={\tau }_z^j$ (dashed lines) for LIOMs obtained using the scheme proposed in this paper. Blue stands for smaller $j$, and red stands for larger $j$ except for the single-particle LIOMs in (d), where the color is randomly chosen because single-particle LIOMs are too delocalized to be ordered. In (a) and (c) $W=20$. $r>1$ interaction strength is below machine precision $\sim {10}^{-15}$. ${\mathrm{\Sigma }}_z$ and ${\tau }_z$ show little difference. In (b) and (d) $W=0.5$. ${\tau }_z^j$ is more localized at site $j$, but the interaction among LIOMs is not zero. $L=10$ and 500 disorder realizations are used in all plots.

Figure 6

Comparison of LIOMs and physical spin operators. (a) Blue dots show the Frobenius norm of the difference between same-site physical spin-1/2 Pauli matrices and local integral of motion $\left|\right|{\sigma }_z^j-{\tau }_z^j\left|\right|$ as a function of disorder strength in a size $L=10$ system. For disorder $W>W_c\sim 7$ (vertical black line), the norm scales as $1/W$ (green line). The red dashed line shows the norm of the difference between two LIOMs at different sites for comparison. (b) Weight of the first LIOM $f_{n,k}^1$ (blue curves) and overlap of the first LIOM with physical spins ${\sum }_j\mathrm{Tr}\left({\tau }_z^1{\sigma }_z^{1+n}\right)/L$ (red curves) as a function of distance $n$ for $W=20$ (solid curves) and $W=40$ (dashed curves). $L=10$, and the median is taken over $k$ and 100 disorder realizations.

Figure 7

(a) Median $\left|\mathrm{Tr}\right({\tau }_z^j{\sigma }_z^k\left)\right|$ and (b) normalized probability distribution of ${log}_{10}\left|V\left(r\right)\right|$ for a chain of size $L=14$ with disordered field only on sites 1 to 8 ($i_e=8$). In (a), green dashed (blue solid) curves represent the LIOMs in the disorder-free (disordered) region. In (b), blue curves show only the coupling terms within the disordered region, and red curves show the distribution of all coupling terms. $W=50$ and 100 random realizations are used in both (a) and (b).

Figure 8

(a) $\langle {\langle {\tau }_x^L\rangle }^2\rangle$ under real (green line) and two artificial Hamiltonians (blue and red lines). The dark green dashed curve shows $\langle {\langle {\sigma }_x^L\rangle }^2\rangle$ under the real Hamiltonian. $L=12, W=40$, so that the delocalized cases are negligible. All $\langle {\langle {\tau }_x^L\rangle }^2\rangle$ are averaged over 20 random initial states in the $xy$ plane and 20 disorder realizations. $\langle {\langle {\sigma }_x^L\rangle }^2\rangle$ is averaged over ten random initial states in the $xy$ plane and ten disorder realizations. (b) The red curve shows ${\kappa }^{-1}$ as in Fig. 4 of the main text. The green solid curve represents the dephasing localization length ${\xi }^{\prime \prime }$ extracted from the artificial Hamiltonian with $\left|V\left(r\right)\right|=e^{-(r/\kappa )}$. The green dashed curve shows ${\xi }^{\prime \prime }+ln2/2$, which overlaps with ${\kappa }^{-1}$ within the error bars.