From chaos to many-body localization: some introductory notes

TL;DR

This paper introduces the concept of many-body localization, exploring its quantum coherence origins, relation to classical chaos, and methods to identify it in interacting spin chains, highlighting its significance in quantum dynamics.

Contribution

It provides an introductory overview of many-body localization, connecting classical chaos, quantum coherence, and disorder effects, with practical identification methods for interacting systems.

Findings

Dynamical localization arises from quantum coherence in chaotic systems

Disorder-induced Anderson localization is discussed in the context of quantum systems

A practical method to identify many-body localization in spin chains is demonstrated

Abstract

Staring from the kicked rotator as a paradigm for a system exhibiting classical chaos, we discuss the role of quantum coherence resulting in dynamical localization in the kicked quantum rotator. In this context, the disorder-induced Anderson localization is also discussed. Localization in interacting, quantum many-body systems (many-body localization) may also occur in the absence of disorder, and a practical way to identify its occurrence is demonstrated for an interacting spin chain.

Equations

$$\displaystyle H\big{(}t\big{)}=\frac{p}{2}+k\cos\theta\sum\limits_{n}\delta\big{(}t-n\big{)}.$$

$${\theta }_{n+1}={\theta }_n+p_{n+1},$$

$$p_{n+1}=p_n+k\sin {\theta }_{n+1}.$$

$$\displaystyle\hat{H}\big{(}t\big{)}=-\frac{1}{2}\frac{\partial^{2}}{\partial\theta^{2}}+k\cos\theta\sum\limits_{n}\delta\big{(}t-n\big{)},$$

$$T_m{U}_m+\sum _{r\ne 0}{W}_r{U}_{m+r}=E{U}_m.$$

$$|{\psi }^{+}(t+1)⟩=\hat{F}|{\psi }^{+}(t)⟩,$$

$$|{\psi }^{-}(t)⟩=\exp (-i{\hat{H}}_0/\hslash )|{\psi }^{+}(t-1)⟩.$$

$$|{\psi }^{\pm }(t)⟩=\sum _n{\psi }_n^{\pm }(t)|n⟩,$$

$$|{\psi }^{\pm }(t)⟩=\underset{0}{\overset{2\pi }{\int }}d\theta {\psi }^{\pm }(\theta ,t)|\theta ⟩,$$

$${\psi }_n^{-}(t+1)=\exp (-i{n}^2/2){\psi }_n^{+}(t),$$

$${\psi }^{+}(\theta ,t)=\exp (-ik\hat{V}(\theta )){\psi }^{-}(\theta ,t).$$

$${\psi }_m^{+}(t+1)=⟨m|{\psi }^{+}(t+1)⟩=$$

$$⟨m|\exp (-ik\hat{V}(\theta )){\psi }^{-}(t+1)⟩=$$

$$\sum _n⟨m|\exp (-ikV(\theta ))|n⟩\exp (-i{n}^2/2){\psi }_n^{+}(t).$$

$$[1-i\hat{W}(\theta )]u(\theta )=\exp (i(\varphi -{\hat{H}}_0))[1+i\hat{W}(\theta )]u(\theta ),$$

$$i\frac{\exp (i(\varphi -{\hat{H}}_0))}{\exp (i(\varphi +{\hat{H}}_0))}u(\theta )-\hat{W}(\theta )u(\theta )=0.$$

$$T_m{u}_m+\sum _{r\ne 0}{W}_r{u}_{m+r}=E{u}_m$$

$$I=\frac{1}{2\pi }\oint pdq,\theta =\frac{\partial S(q,I)}{\partial I},$$

$$\displaystyle H\big{(}I,\theta\big{)}=H_{0}\big{(}I\big{)}+\varepsilon V\big{(}I,\theta\big{)}.$$

$$\displaystyle P(S)=\frac{\pi S}{2}\exp\bigg{(}-\frac{\pi}{4}S^{2}\bigg{)},$$

$$\displaystyle P(S)=\frac{32S^{2}}{\pi^{2}}\exp\bigg{(}-\frac{4S^{2}}{\pi}\bigg{)},$$

$$\displaystyle P(S)=\frac{2^{18}}{3^{6}\pi^{3}}S^{4}\exp\bigg{(}-\frac{64}{9\pi}S^{2}\bigg{)},$$

$$\displaystyle\hat{p}_{i}=\frac{1}{\bigg{(}N\sum\limits_{r}|\psi_{i}(r)|^{4}\bigg{)}}.$$

$${\hat{P}}_i=\sum _r|{\psi }_i(r){|}^4.$$

$$\displaystyle P^{(k)}(E)=\bigg{\langle}\sum\limits_{i,r}|\psi_{i}(r)|^{2k}\delta(E-e_{i})/\sum\limits_{i}\delta(E-e_{i})\bigg{\rangle},$$

$$P^{(k)}\approx L^{-{\tau }_k},{\tau }_k=d(k-1)+{\Delta }_k.$$

$$\displaystyle i\frac{\partial\psi}{\partial t}=-\frac{1}{2m}\bigg{(}\frac{\partial^{2}\psi}{\partial x^{2}}+\Lambda\delta(x-x_{0})\psi\bigg{)}$$

$$\psi (x)=a^{+}\exp (ikx)+a^{-}\exp (-ikx),x<x_0,$$

$$\psi (x)=b^{+}\exp (ikx)+b^{-}\exp (-ikx),x>x_0.$$

$$\psi (x_0^{+})=\psi (x_0^{-}),$$

$${\partial }_x\psi (x_0^{+})-{\partial }_x\psi (x_0^{-})=\Lambda \psi (x_0).$$

$$\displaystyle G=\frac{2e^{2}}{\hbar}\sum\limits_{n}T_{n}=\frac{2e^{2}}{\hbar}Tr\big{(}t^{{\dagger}}t\big{)}.$$

$$P(S)=c\exp [-\beta Trf(t^{†}t)],$$

$$\displaystyle P\big{(}\{T_{n}\}\big{)}=c\prod\limits_{i<j}|T_{i}-T_{j}|^{\beta}\times$$

$$\prod _k{T}_k^{-1+\beta /2}\exp [-\beta f(T_k)].$$

$$r_n=\min \left({\delta }_n,{\delta }_{n-1}\right)/\max \left({\delta }_n,{\delta }_{n-1}\right),$$

$$r=\frac{1}{N-2}\sum _{n=3}^N{r}_n$$

Full text

I Hamiltonian Chaos vs. random impurities

As put by Edward Lorenz, classical deterministic chaos is Lorenz : "when the present determines the future, but the approximate present does not approximately determine the future". The instability of phase trajectories is quantified by the Lyapunov exponent Zaslavsky : In phase space, the state of the $N$ dimensional dynamical system with $2N$ degrees of freedom $\mathbf{x}\equiv(\mathbf{p}_{n},\mathbf{q}_{n}),~{}n=1,..N$ is mapped to a single dot and the evolution of the system is described by a phase trajectory circumscribed by the system’s state. Let us track the time evolution of two different phase trajectories see Fig.(1) $\mathbf{x}(t)+\mathbf{\delta}(t),~{}~{}\mathbf{x}(t)$ emanated from slightly different initial conditions $\mathbf{x}_{0},~{}\mathbf{x}_{0}+\mathbf{\delta}_{0}$. During the time evolution the distance between phase trajectories \|\mathbf{\delta}(t)\|=\|\mathbf{\delta}(0)\|\exp\big{(}\lambda t\big{)} increases provided that the Lyapunov exponent is positive $\lambda>0$. Therefore, arbitrary small uncertainties in the initial conditions accumulate over a long period to the substantial error, calling thus for a statistical description when evaluating observable quantities.

I.1 Model of the minimal D=3/2 Hamiltonian chaos

The kicked rotator is a prototype of systems exhibiting chaos. The Hamiltonian of the one dimensional classical periodically kicked rotator reads

[TABLE]

Here $(p,\theta)$ is the canonical pair of momentum-angle variables and $n$ specifies the number of kicks. Solution of the classical kicked rotator problem is described by the Chirikov’s map

[TABLE]

The phase space of the Chirikov’s map consists of regions of regular and chaotic motions. If the strength of the kick exceeds the specific value $k>0.97$, regular islands disappear, and the chaotic sea covers the whole phase space see Fig.(2) leading to a diffusive growth of the momentum in time (in the number of kicks) \langle\big{(}p_{n}-p_{0}\big{)}^{2}\rangle=\frac{k^{2}}{2}n. The kicked rotator is realized in a number of different physical setting. For example, in Ref.Komnik the model was employed to investigate spin systems. Contrary to the classical case, the quantum counterpart of the model Eq.(1)

[TABLE]

shows a suppression of the diffusion due to the destructive role of quantum coherence Haake . Hence, instead of the diffusive growth localization occurs. We note, in both the classical and the quantum case the key issue is the discrete time evolution. The quantum kicked rotator can be mathematically mapped onto the one dimensional Anderson tight-binding model Anderson ;

[TABLE]

Here $W_{r}$ are the hopping amplitudes, the on-site energies $T_{m}$ are uniformly distributed in the interval $T_{m}\in[-T_{max}/2,T_{max}/2]$ and mimic the effect of randomly embedded impurities. Localization is indicated by the absence of mobility of electrons in the localized phase and the exponential decay of amplitudes of the wave functions. If the electron is initially localized around the site $n_{in}$, the probability to find the electron on the site $m$ exponentially decays U_{m}^{(n_{in})}\approx\exp\big{(}-|n_{in}-m|/l_{c}\big{)} with the distance $|n_{in}-m|$, and the localization length is a function of energy l_{c}\big{(}E\big{)}. The Floquet operator $\hat{F}=\exp(-ik\hat{V}(\theta))\exp(-i\hat{p}^{2}/2)$ describes the evolution of the kicked rotator during the one period of driving. The unperturbed Hamiltonian $\hat{H}_{0}=\hat{p}^{2}/2$ is diagonal $\hat{p}|n\rangle=n|n\rangle$ in the eigenbasis $|n\rangle=\frac{1}{\sqrt{2\pi}}\exp(in\theta)$ of the momentum operator. Therefore, for the wave functions before and after the n-th kick we use the following ansatz

[TABLE]

Here $|\psi^{+}(t+1)\rangle$ and $|\psi^{+}(t)\rangle$ are wave functions after the $n+1$-th and $n$-th kicks, while $|\psi^{-}(t)\rangle$ is the wave function before the $n$-th kick. For treating the free and kicked evolution parts, we utilize the following parametrization

[TABLE]

leading to the following recurrent relations for the expansion coefficients:

[TABLE]

Taking into account Eq.(I.1), Eq.(I.1) and Eq.(I.1) we deduce the quantum map for the expansion amplitudes in Eq.(I.1)

[TABLE]

In the case $\hat{V}=\cos\theta$, the matrix element $J_{nm}=\langle m|\exp(-ikV(\theta))|n\rangle$ becomes Bessel functions. To find the equivalence between the kicked rotator and the Anderson model Eq.(4) we use the eigenfunctions and eigenvalues of the Floquet operator $\hat{F}|u^{+}\rangle=\exp(-i\phi)|u^{+}\rangle$, and consider the following identity $u^{-}(\theta)=\exp(ikV(\theta))u^{+}(\theta)$ in the $\theta$ representation. Taking into account $\exp(-i\phi)\hat{F}u^{+}(\theta)=\exp(-ikV)\exp(i(\phi-\hat{H}_{0}))u^{+}(\theta)=u^{+}(\theta)$ we obtain $u^{-}(\theta)=\exp(i(\phi-\hat{H}_{0}))u^{+}(\theta)$ or in the eigenbasis of $\hat{H}_{0}$ $u^{-}_{n}=\exp(i(\phi-n^{2}/2))u^{+}_{n}$. Let us introduce the auxiliary operator $\hat{W}$ defined as follows \exp(-ik\hat{V})=\frac{1+i\hat{W}}{1-i\hat{W}}=\tan\big{(}\frac{k\hat{V}}{2}\big{)} and the vector |u\rangle=\frac{1}{2}\big{(}|u^{+}\rangle+|u^{-}\rangle\big{)}. Then, with the identity $u(\theta)=(1+i\hat{W})^{-1}u^{+}(\theta)=(1-i\hat{W})^{-1}u^{-}(\theta)$ we infer the operator equation:

[TABLE]

or

[TABLE]

Introducing the notations T_{m}=i\frac{1-\exp(i(\phi-m^{2}/2))}{1+\exp(i(\phi-m^{2}/2))}=\tan\big{(}(\phi-m^{2}/2)/2\big{)}, $u_{m}\langle m|u\rangle$, $E=-W_{0}$, $W_{r}=\langle m|W|m+r\rangle$ we rewrite Eq.(10) in the matrix form

[TABLE]

which demonstrates the mapping between the kicked rotator Eq.(3) and the Anderson tight-binding model Eq.(4). We note that the Anderson model, as well as the kicked rotator, describe single particle noninteracting systems. These models do not capture Many-body localization in interacting systems.

I.2 D>2 classical and quantum chaos, KAM theorem and random matrix theory

For nonlinear systems with a few degrees of freedom the coupling between different variables is decisive. At certain values of the coupling strength, the internal resonances overlap giving rise to Hamiltonian chaos and new phenomena such as Arnold’s diffusion in phase space Arnold . Let us introduce the canonical pair of action-angle variables:

[TABLE]

and present theHamiltonian of the system in the form:

[TABLE]

Here H_{0}\big{(}I,\big{)},~{}~{}\forall~{}I_{i}\equiv I,i\in(1,..N) is the integrable part of the Hamiltonian and V\big{(}I,\theta\big{)},~{}~{}\forall~{}\theta_{i}\equiv I,i\in(1,..N) is the perturbation. We presume that the Hamiltonian is not degenerate i.e. the Hessian $det|\frac{\partial\omega}{\partial I_{i}}|=det|\frac{\partial^{2}H_{0}}{\partial I_{i}\partial I_{j}}|\neq 0$ is not zero.

Following the Kolmogorov–Arnold–Moser (KAM) theorem one finds that the motion of an integrable system is confined to an invariant torus. If the system is subjected to a weak nonlinear perturbation, some of the invariant tori are deformed and survive.

The mechanism of Arnold’s diffusion is as follows: In the region of the destroyed invariant torus, in the vicinity of the internal resonances condition $\sum\limits_{j=1}^{N}n_{j}\omega_{j}=0,~{}~{}\forall~{}n_{j}\in\mathbb{N}$ a homocyclic structure emerges. The internal resonance condition determines a family of surfaces H_{0}\big{(}I_{R}\big{)}=E with corresponding resonant action variables $I_{R}$. Nodal crossing points of different resonances construct a stochastic net. The Arnold’s diffusion sets in the stochastic net and is a topological generalization of classical kicked rotator problem. The Arnold’s diffusion was explored for the Heisenberg model, e.g., in Ref.Chotorlishvili1 . The specific energy level statistics characterize the quantum counterpart of the classically chaotic nonlinear $D>2$ system Eq.(13). Depending on the symmetry of the Hamiltonian $\hat{H}$, in the region of destroyed phase space the invariant tori, according to the random matrix theory (RMT)Haake the mean inter-level energy distance $S=\langle E_{n+1}-E_{n}\rangle$ follows either a Gaussian orthogonal ensemble (GOE), a Gaussian unitary ensemble (GUE) or a Gaussian symplectic ensemble (GUE), respectively, i.e.

[TABLE]

while for the integrable case the statistics is Poissonian $P(S)=\exp(-S)$. We note that the termination of the quantum Arnold’s diffusion for the particle moving in a quasi-1D waveguide was studied in Demikhovskii . The case $D=3/2$ was addressed in Ref. Chotorlishvili2 . Classical chaotic Hamiltonian systems are insensitive to GOE, GUE, and GSE symmetries.

II Diagnostic tools of the localized state

Before analyzing MBL in more details, we briefly review diagnostic tools of a localized state.

II.1 Multifractality of the wave functions and scaling properties

a The inverse participation ratio averages the fourth power of the wave function. It is positive for localized states and vanishes for extended states in the thermodynamic limit. If the wave function $|\psi_{i}\rangle$ at site $i$ of a tight-binding model is normalized then partition ratio $\hat{p}_{i}$ (PR) is given by

[TABLE]

The inverse participation ratio (IPR) has been defined as

[TABLE]

If the wave function spreads over $l$ lattice sites with equal amplitude $|\psi_{i}(r)|^{2}=1/l$ and vanishes elsewhere one can deduce $\hat{p}_{i}=l/N,~{}~{}\hat{P}_{i}=1/l$. As stated by F. WegnerWegner , PR describes the proportion of the total number of atoms in a system which contribute effectively to an eigenstate, whereas IPR is the inverse number of orbitals contributing effectively to this state. For localized states, IPR is larger than zero whereas PR vanishes in the thermodynamic limit. The wavefunction of the system shows anomalous scaling properties in the vicinity of the transition point to the localized phase. In particular, the ensemble averaged quantity Evers

[TABLE]

shows fractal scaling properties with the characteristic size of the system $L$ and non-integer $\tau_{k}$. Here $\psi_{i}(r)$ is the amplitude of the eigenfunction $|i\rangle$ with the energy $e_{i}$ at site $r$ of the particle in a tight binding model and non-integer $\Delta_{k}\notin\mathbb{N}$.

II.2 RMT beyond Quantum chaos. The transfer-matrix methodBeenakker

Let us consider the simplest formulation of the scattering problem, a scattering of an electron on the $\delta$-like impurities located at $x=x_{0}$. The Schrödinger equation

[TABLE]

admits the solution

[TABLE]

and the continuity conditions

[TABLE]

define the elements of the scattering matrix: $c_{out}=\begin{bmatrix}b^{+}\\ b^{-}\end{bmatrix}$

$\begin{bmatrix}r&\acute{t}\\ t&\acute{r}\end{bmatrix}$ $\begin{bmatrix}a^{+}\\ a^{-}\end{bmatrix}$

$S\begin{bmatrix}a^{+}\\ b^{-}\end{bmatrix}=Sc_{in}.$

Generally, the channel can be arbitrarily large and the dimension of the $t$ matrix is equal to the number of scattering channels. The explicit expressions of the $t$ matrix elements is obtained after solving Eq.(II.2) and expressing the coefficients $b^{+},b^{-}$ in terms of $a^{+},a^{-}$. The Landauer-Büttiker conductanceBeenakker is given by:

[TABLE]

We note, a) When the system is in the localized phase, the eigenvalues of the transfer matrix are random numbers. However, there is a principle difference between RMT of quantum transport and chaotic nonlinear systems. The RMT describes quantum transport in terms of transmission eigenvalues of the open system (i.e., leads, impurities, etc.), while RMT is an internal property of nonlinear non-integrable quantum systems. The Thouless energy characterizes the energy scale of disordered conductors and is defined via the formula $E_{T}=\hbar D/L^{2}$, where $D$ is the diffusion constant, and $L$ is the characteristic size of the system. b) B. L. Altshuler and B. I. Shklovskii Altshuler ; eeanderson showed that for energy separations greater than the Thouless energy $|E^{{}^{\prime}}-E|>E_{T}$, the correlation function deviates from random-matrix theory. If the condition $|E^{{}^{\prime}}-E|<E_{T}$ holds, the distribution function of eigenvalues of the transfer matrix read

[TABLE]

Here P\big{(}\{T_{n}\}\big{)} describes the correlated distribution of eigenvalues, $S=\langle E_{n+1}-E_{n}\rangle$ is the mean distance between neighboring eigenvalues, $\beta=1$ for GOE, $\beta=2$ for GUE, and $\beta=4$ for GSE symmetries, respectively, and the function $f(T_{k})$ is determined through the average spectral density $f(T_{k})=\int dT^{{}^{\prime}}_{k}\sigma(T^{{}^{\prime}}_{k})\ln|T_{k}-T^{{}^{\prime}}_{k}|$.

III Energy level statistics for the MBL phase

MBL is a localization phenomenon in interacting system Basko ; Devakul ; Luitz ; Vasseur ; Bardarson ; Serbyn ; topical ; huse ; huse1 ; huse2 ; Bauer , distinct from the Anderson localization: We note, in absence of electron-electron interaction in low dimensional systems $d=1,2$ all states are localized even in the presence of arbitrary small disorder and the conductivity $\sigma_{d=1,2}(T)=0$. For $d=3$ the conductivity follows the Arrhenius law $\sigma_{d=3}(T)\sim\exp(-E_{c}/T)$. Here $E_{c}$ is the distance between the Fermi level and the mobility edge (mobility edge refers to the border between localized and extended bands). Therefore, for any finite $E_{c}$, the Anderson localization of electronic states leads to a metal-insulator transition at zero temperature. The conductivity becomes finite at any finite temperature. MBL can be viewed as Anderson like localization of many-body wave functions. The expectation value of an arbitrary quantum operator is given by formula $\langle\mathcal{\hat{O}}\rangle_{t}=\langle\psi(t)|\mathcal{\hat{O}}|\psi(t)\rangle$. Taking into account that $|\psi(t)\rangle=\exp(-i\hat{H}t/\hbar)|\psi(0)\rangle$ and utilizing the expansion over the basis of energy eigenstates $|\psi(0)\rangle=\sum\limits_{\alpha}c_{\alpha}|E_{\alpha}\rangle$ we deduce: $\langle\mathcal{\hat{O}}\rangle_{t}=\sum\limits_{\alpha,\beta}c_{\alpha}^{\ast}c_{\beta}\mathcal{O}_{\alpha\beta}\exp(-i(E_{\beta}-E_{\alpha})t/\hbar)$, where $\mathcal{O}_{\alpha\beta}=\langle E_{\alpha}|\mathcal{\hat{O}}|E_{\beta}\rangle$. The Eigenstate Thermalization Hypothesis(ETH) statesIlievski ; Pozsgay ; Mierzejewski ; ETH1 ; ETH2 ; ETH3 that for quantum ergodic systems, in the long time limit the expectation value is given by the following formula: $\langle\mathcal{\hat{O}}\rangle_{t}^{t\rightarrow\infty}=\sum\limits_{\alpha}|c_{\alpha}|^{2}\mathcal{O}_{\alpha\alpha}$. The MBL state violates ETH. Therefore, MBL can be referred to as a non-ergodic phase. In the ergodic phase, depending on the symmetry of the Hamiltonian $\hat{H}$, the energy level statistics is either GOE, GUE, or GSE, while in the MBL phase level statistics is Poissonian. The transition from the ergodic phase to the MBL phase occurs at specific critical strength of disorder. To infer precisely the MBL transition point, one needs to obtain and analyze a large number of statistics, i.e., consider different realizations of disorder and take ensemble average. A key quantity is the disorder average of the ratio

[TABLE]

where $\delta_{n}=E_{n}-E_{n-1}$ is the distance between two neighboring energy levels labelled by $n$ and $N$ is the number of eigenstates. In the ergodic phase for the Gaussian orthogonal ensemble GOE $r_{\texttt{{GOE}}}=0.5307$ is found, while in the MBL phase $r_{\texttt{{Poisson}}}=0.3863$. The disorder strength has a strong influence on the system’s spectral characteristics. To avoid finite size effects and other numerical artifacts calculations should be done for different sizes of the system and then the obtained data should be collapsed into a single universal size-independent result. The conventional finite-size scaling analysis works as followsLuitz : One selects the middle part of the spectrum of the system. The level statistics of $\left\langle r\right\rangle$ as a function of disorder $h$ are scaled with $f_{L}(h)=L^{1/\nu}(h-h_{c})$, where $\nu=1$ was assumed and $h_{c}$ defines the critical strength of disorder. For best fitting parameter $h_{c}$, the scaling procedure admits a minimizing function $w(h)$

[TABLE]

where $h_{1}$ and $h_{2}$ are defined by the common integration domain $h_{1}=\max\left(L_{i}^{\frac{1}{\nu}}(a-h_{c})\right),\;h_{2}=\min\left(L_{i}^{\frac{1}{\nu}}(b-h_{c})\right)$ where $L_{i}$ denote the system sizes to be analyzed and $a$, $b$ are boundaries limits for un-scaled data set.

IV Systems with complex symmetry properties

For systems with a mixed symmetry, $\hat{H}_{tot}=\hat{H}_{GOE}+\lambda\hat{H}_{GUE}$ the total Hamiltonian doesn’t belong to a particular symmetry class. In this case, finite-size scaling analysis may fail when GOE/GUE statistics show different scaling features. Depending on the value of the parameter $\lambda$ the spectrum of the Hamiltonian $\hat{H}_{tot}=\hat{H}_{GOE}+\lambda\hat{H}_{GUE}$ displays different features: In the limit of small $\lambda$ the level statistics obey GOE, for large $\lambda$ GUE, and what is interesting, in the transition area the system shows qualitatively different properties than $\hat{H}_{GOE}$ and $\hat{H}_{GUE}$. We need a new method to explore the MBL phase in systems with mixed complex symmetries.

In our recent work stefan we observed that in spite of the difference of the GOE and GUE statistics, the enhanced fluctuations at the MBL transition point bears universal physical character typical for the both GOE and GUE symmetry. Thus, our method can be implemented for systems with complex symmetry properties when standard finite-size scaling procedure fails. We explored the MBL problem Hamiltonian with a dynamical Dzyaloshinskii Moriya interaction (DMI) chiral ; Azimi

[TABLE]

Note, in the absence of the DMI term the Hamiltonian has GOE symmetry and the symmetry of the DMI term is GUE. The DMI constant $D$ combines the effect of an electric field and the magnetoelectric coupling. The nearest neighbor exchange interaction is ferromagnetic $J_{1}<0$ while the next nearest neighbor one is antiferromagnetic $J_{2}>0$ leading in general to a frustrated spin order.

Computationally we are able to deal with only small chains; experimentally relevant are for instance the Fe chains on the $(5\times 1)-$Ir(001) surface Wiesendanger . $\hat{H}$ is block-diagonal. Each block is identified via the conserved total spin component $\hat{S}^{z}=\sum^{L}_{i=1}\hat{S}_{i}^{z}$. Of a special interest is the largest subspace of states $|\Psi_{n}\rangle$ obeying $\hat{S}^{z}|\Psi_{n}\rangle=M|\Psi_{n}\rangle$ with $M=0$ for even $L$ or $M=1$ for odd $L$, respectively. A uniform magnetic field $B^{z}_{i}=B^{z}$ shifts the eigenvalues equally in each subspace and has no prominent effect on the interlevel distance $r_{n}$, while randomness incorporated in the magnetic field $B^{z}_{i}\in\left\langle-h,h\right\rangle$ can induce a qualitative change of the spectral properties from Wigner-Dyson to Poisson level spacing statistics. The strength of disorder is measured on a scale set by $J_{1}$. In what follows we work with dimensionless units such that $J_{1}=1$.

To formulate a possibly general criterion for MBL that is applicable in such cases, as well we analyze the full statistics for each realization $\alpha$ of the random magnetic fields $r^{(\alpha)}$. The histograms corresponding to a counting classification of $r^{(\alpha)}$ for a given disorder strength $h$ is presented on Fig. 3. As can be inferred, the histograms are narrow far away from the MBL transition, while the histograms become particularly broad close to the transition point mimicking the behavior of fluctuations near conventional phase transitions. The histograms become more and more pronounced for increasing $L$. Fig. 4 demonstrates the convergence of the histograms of for chains with different lengths. As we see already for $L=14$ counts histograms amalgamate underlying that an analysis of the histograms can serve as a further indicator in addition to finite-size scaling. The convergence of histograms even for relatively small systems endorses our method as less computationally demanding which is a major advantage for exact diagonalization approaches that are considered as well suited for MBL studies. As for the histograms of consecutive level spacing, Fig. 3 illustrates the broadening of the histograms when approaching the transition point between the ergodic and the MBL phases. As evident, the effect of broadening is even more prominent for systems with a larger size Fig. 3. We note that the observed phenomena is not related to a particular type of level statistics but it is rather akin to the transition regime. Away from the transition point on the ergodic side (GOE statistics), and on the MBL phase side (Poisson statistics) the width of histograms are narrower. The broadening is linked to the enhanced quantum fluctuations Fig. 5. This behavior is of a general character and is maintained even after adding the next nearest neighbor interaction and DMI terms.

Physically, the broadening of histograms is attributable to the enhanced fluctuations near phase transitions (cf. Figs. 5, 6). Hence, such broadening serves as a further indicator for approaching the MBL phase.

Disorder in the exchange coupling or in $D$ may also occur. The latter (cf. Eq. IV) can be viewed as random change in $\mathbf{E}$ or a random elastic energy change ($\mathbf{E}.{\bf P}=g_{ME}\mathbf{E}\sum\limits^{L}_{i=1}\langle\mathbf{e}_{x}\times(\mathbf{\hat{S}}_{i}\times\mathbf{\hat{S}}_{i+1})\rangle$), and thus, it is important for spin-phonon-coupled systems at finite temperatures. Calculations evidence the robustness of the MBL phase against randomizing $D$ within a physically reasonable range, an example is depicted in Fig. 7.

Acknowledgements.

We are indebted to Marcin Mierzejewski for numerous discussions and suggestions. We acknowledge financial support from the Deutsche Forschungsgemeinschaft through SFB 762, and BE 2161/5-1.