Emulating Many-Body Localization with a Superconducting Quantum Processor

The law of statistical physics dictates that generic closed quantum many-body systems initialized in nonequilibrium will thermalize under their own dynamics. However, the emergence of many-body localization (MBL) owing to the interplay between interaction and disorder, which is in stark contrast to Anderson localization, which only addresses noninteracting particles in the presence of disorder, greatly challenges this concept, because it prevents the systems from evolving to the ergodic thermalized state. One critical evidence of MBL is the long-time logarithmic growth of entanglement entropy, and a direct observation of it is still elusive due to the experimental challenges in multiqubit single-shot measurement and quantum state tomography. Here we present an experiment fully emulating the MBL dynamics with a 10-qubit superconducting quantum processor, which represents a spin-$1/2$ $XY$ model featuring programmable disorder and long-range spin-spin interactions. We provide essential signatures of MBL, such as the imbalance due to the initial nonequilibrium, the violation of eigenstate thermalization hypothesis, and, more importantly, the direct evidence of the long-time logarithmic growth of entanglement entropy. Our results lay solid foundations for precisely simulating the intriguing physics of quantum many-body systems on the platform of large-scale multiqubit superconducting quantum processors.

Figure 1

Experimental setup. (a) Diagram of the 10-qubit superconducting quantum processor. The qubits, shown as atoms with spins initialized in alternate orientations, are arranged along a circular chain with the nearest-neighbor couplings represented by red wavy lines, which are calibrated in a separate measurement as discussed in the Supplemental Material (the coupling for $Q_5$ and $Q_6$ is smaller than others) [35]. The long-range, i.e., beyond nearest-neighbor, spin-spin interactions are enabled by the central bus resonator $\mathcal{R}$ that couples to each individual qubit as illustrated by curved connecting lines. (b) Pulse sequences for emulating MBL. All ten qubits are initialized at their respective idle frequencies by applying $\pi$ pulses (dark-color sinusoids) to the even-site qubits to prepare the initial state $|0101010101⟩$, following which the rectangular pulses are applied to quickly bias all qubits to nearby $\mathrm{\Delta }/2\pi =-650\mathrm{MHz}$. Individual $Q_i$ is offset from $\mathrm{\Delta }$ by a small amount of $\delta h_{i,k}\in [-\delta h,\delta h]$, where $\delta h_{i,k}$ is randomly chosen but fixed for the $k$th pulse sequence (see the bottom right panel for the enlarged view of the pulse segment enclosed by dashed lines), and the ensemble of the $k=1–30$ pulse sequences effectively emulates the random disordered potential $\delta h_i$. After the 10-qubit system evolves for a specific time from 0 to 1000 ns under the square pulses, all qubits are biased back to their respective idle frequencies for the 10-qubit joint read-out, which returns binary outcomes of the ten qubits: We run all the $k=1–30$ sequences, each being executed for 3000 times, to count $2^{10}$ probabilities of $\{P_{00\dots 0},P_{00\dots 1},\dots ,P_{11\dots 1}{\}}_k$ for $k=1–30.$ Mean of the $k=1–30$ probability ensemble captures the effect of the random disorder. If the $N$-qubit QST is necessary, we insert tomographic rotation pulses (light-color sinusoids) to the involved qubits before the $N$-qubit joint read-out to obtain all tomographic probabilities, $3^N$ more than aforesaid, to calculate the $N$-qubit density matrix for the $k$th sequence. (Top right) The circular spin chain arranged in one dimension for illustrating the quantity $\delta h_{i,k}$. (c) The spin-spin coupling matrix $J_{ij}$ for $\mathrm{\Delta }/2\pi =-650\mathrm{MHz}$. Only nearest-neighbor couplings are positive, while all other interaction terms are negative, which do not decay over distance.

Figure 2

Dynamics of imbalance for the system initialized in $|0101010101⟩$. (a),(b) The time evolutions of the excited-state ($|1⟩$-state) probabilities $P_i$ (color bar on the far right) for $Q_i$ (indexed as qubit number $i$) when there is no disorder with $\delta h\approx 0$ and a strong disorder with $\delta h/2\pi \approx 12\mathrm{MHz}$, respectively. All probabilities are after read-out correction [40]. (c) Time evolutions of the system imbalance $\mathcal{I}$ at different disorder strengths $\delta h$ as listed. (d) The quasi-steady-state imbalance $\mathcal{I}$ (top) and standard deviation $\delta n$ of the $|1⟩$-state probability distributions taken at 1000 ns as functions of $\delta h$, showing that the system thermalizes for no disorder, but starts to enter the MBL phase with increasing disorder strength. Error bars are 1 SD calculated from all probability data of the $k=1$–30 pulse sequences. The transition looks a bit continuous for our 10-qubit system, but numerical simulations involving up to 20 qubits show that it becomes more pronounced and the associated critical disorder strength can be identified when the qubit number increases.

Figure 3

Subsystem density matrices for test of ETH. Shown are the experimental initial density matrices in absolute value for the 2-qubit subsystems in (a) and 5-qubit subsystems in (b), in comparison with those probed at 1000 ns when the system evolves to either the thermalized state ($\delta h\approx 0$) or the MBL phase ($\delta h/2\pi \approx 12\mathrm{MHz}$). Density matrix data of the 1-qubit subsystem $Q_3$ at a few selected evolution times can be found in the Supplemental Material [35].

Figure 4

Half-chain entanglement entropy for the 5-qubit subsystem $\{Q_3,Q_4,Q_5,Q_6,Q_7\}$, of which the density matrices are obtained by QST. (a),(b) The entanglement entropy $S$ as functions of the evolution time, in linear scale and in logarithmic scale, respectively, at different disorder strengths $\delta h$ as labeled. (c) Site-averaged $S$ at around 1000 ns as functions of the number of sites (qubits) $N$ at different disorder strengths $\delta h$ as indicated, which are calculated by taking the average of the entanglement entropies of all $N$-site choices ($N\le 5$) out of the directly measured 5-qubit subsystem. Dots connected by dashed lines are experimental data. The solid line shows $S=N\ln 2$. (d) Entanglement entropy $S$ as functions of the evolution time for comparison between MBL and Anderson localization (AL) with the disorder strength $\delta h/2\pi \approx 12\mathrm{MHz}$. Dots are experimental data, and lines are numerical simulation results as indicated (see Supplemental Material [35]). Error bars are 1 SD calculated from all tomography data of the $k=1–30$ pulse sequences. Since there is no straightforward method for calculating the entanglement entropy $S$ of a mixed state [44], here we estimate $S$ using the subsystem density matrix with the assumption that the 10-qubit system remains in a pure state. The apparent increase of $S$ under decoherence is due to an imperfection of such a pure-state assumption. We further verify with numerical simulations that the residual excited state population in each qubit, typically 0.01 or less due to unwanted excitations [26], has minimal impact on the signature of the long-time logarithmic growth of $S$ in the MBL phase.