# Observing emergent hydrodynamics in a long-range quantum magnet

M. K. Joshi  Institute for Quantum Optics and Quantum Information, Austrian
Academy of Sciences, Technikerstraße 21a, 6020 Innsbruck, Austria  F. Kranzl
Institute for Quantum Optics and Quantum Information, Austrian Academy of
Sciences, Technikerstraße 21a, 6020 Innsbruck, Austria  Institut für
Experimentalphysik, Universität Innsbruck, Technikerstraße 25, 6020 Innsbruck,
Austria  A. Schuckert  Department of Physics and Institute for Advanced Study,
Technical University of Munich, 85748 Garching, Germany  Munich Center for
Quantum Science and Technology (MCQST), Schellingstr. 4, 80799 München,
Germany  I. Lovas  Department of Physics and Institute for Advanced Study,
Technical University of Munich, 85748 Garching, Germany  Munich Center for
Quantum Science and Technology (MCQST), Schellingstr. 4, 80799 München,
Germany  C. Maier  Institute for Quantum Optics and Quantum Information,
Austrian Academy of Sciences, Technikerstraße 21a, 6020 Innsbruck, Austria
AQT, Technikerstraße 17, 6020 Innsbruck, Austria  R. Blatt  Institute for
Quantum Optics and Quantum Information, Austrian Academy of Sciences,
Technikerstraße 21a, 6020 Innsbruck, Austria  Institut für Experimentalphysik,
Universität Innsbruck, Technikerstraße 25, 6020 Innsbruck, Austria  M. Knap
Department of Physics and Institute for Advanced Study, Technical University
of Munich, 85748 Garching, Germany  Munich Center for Quantum Science and
Technology (MCQST), Schellingstr. 4, 80799 München, Germany  C. F. Roos
Institute for Quantum Optics and Quantum Information, Austrian Academy of
Sciences, Technikerstraße 21a, 6020 Innsbruck, Austria  Institut für
Experimentalphysik, Universität Innsbruck, Technikerstraße 25, 6020 Innsbruck,
Austria

###### Abstract

Identifying universal properties of non-equilibrium quantum states is a major
challenge in modern physics. A fascinating prediction is that classical
hydrodynamics emerges universally in the evolution of any interacting quantum
system. Here, we experimentally probe the quantum dynamics of 51 individually
controlled ions, realizing a long-range interacting spin chain. By measuring
space-time resolved correlation functions in an infinite temperature state, we
observe a whole family of hydrodynamic universality classes, ranging from
normal diffusion to anomalous superdiffusion, that are described by Lévy
flights. We extract the transport coefficients of the hydrodynamic theory,
reflecting the microscopic properties of the system. Our observations
demonstrate the potential for engineered quantum systems to provide key
insights into universal properties of non-equilibrium states of quantum
matter.

Universality in equilibrium asserts that microscopic details are irrelevant
for the emergent quantum phases of matter and their transitions. Rather,
symmetries and topology determine the essential macroscopic properties. By
contrast, all scales, from low to high energies, are relevant for quantum
systems which are driven far from their thermal equilibrium. Recent
experimental progress in engineering coherent and interacting quantum systems
made it possible to create and explore exotic non-equilibrium states, which
can exhibit unconventional relaxation dynamics [1, 2, 3, 4, 5, 6], dynamical
phases [7, 8, 9, 10, 11, 12], and transitions between them [13, 14].

Despite this wealth of observed quantum phenomena, a common anticipation is
that classical hydrodynamics of a few conserved quantities emerges universally
for any complex quantum system, as strong interactions entangle and
effectively mix local degrees of freedom [15, 16]. However, verifying this
assumption, and furthermore determining the non-universal transport
coefficients of the emergent hydrodynamic theory for specific systems is
challenging. Recently, enormous efforts have been devoted to detect
hydrodynamic transport in quantum gases [17, 18, 19, 20] and condensed matter
systems [21, 22, 23, 24]. While transport is generally expected to be
diffusive, a variety of largely unexplored classes of hydrodynamics have been
predicted theoretically including anomalous subdiffusive [25, 26] and
superdiffusive transport [27, 28, 29].

In this work, we experimentally probe the dynamics of a long-range quantum
magnet, realized in a quantum simulator of 51 individually controlled ions. We
develop a protocol to measure space- and time-resolved spin correlations in an
engineered infinite temperature state, enabling us to experimentally establish
that hydrodynamics emerges in the non-equilibrium quantum state. By tuning the
long-range character of the interactions, we observe a whole family of
hydrodynamical universality classes, ranging from normal diffusion to
anomalous superdiffusion.

The (pseudo-)spins of our quantum system are realized with two electronic
states of the 40Ca+ ion: $\ket{S_{1/2},m=+1/2}$ as $\ket{\downarrow}$ and
$\ket{D_{5/2},m=+5/2}$ as $\ket{\uparrow}$. The quantum state of individual
ions is controlled by a tightly focused, steerable laser beam capable of
addressing any ion in the string, in conjunction with a laser beam that
collectively interacts with the ions. A two-tone laser field realizes
approximately power-law decaying Ising interactions between the (pseudo-)spins
by off-resonantly coupling motional and electronic degrees of freedom of the
ion chain. Application of a strong transverse field energetically penalizes
spin non-conserving contributions [30]. The effective dynamics is then
described by the long-range XY model

$\hat{H}=\sum_{i<j}\frac{J}{|i-j|^{\alpha}}(\hat{\sigma}^{+}_{i}\hat{\sigma}^{-}_{j}+\hat{\sigma}^{-}_{i}\hat{\sigma}^{+}_{j}),$
(1)

where $J$ and $\alpha$ determine the strength and the range of the spin-spin
interaction, respectively [31]. They can be tuned by varying the amplitude and
frequency of the two-tone laser field. Further experimental details are given
in the supplementary materials [32], Sec. I-IV.

Figure 1: Emergent hydrodynamics in a long-range quantum magnet. A) In our
system, the dynamics of a locally created spin excitation is effectively
governed by a hydrodynamic theory of Lévy flights, which are random walks with
step sizes drawn from a distribution with algebraic tails. This occasionally
leads to long-distance jumps of the excitation. B) The family of hydrodynamic
universality classes can be tuned by the power-law exponent $\alpha$ of the
spin-spin interactions. C) We measure infinite-temperature correlations by
averaging over initial product states while preparing the central ion
deterministically in the same state (blue box); [32] Sec. II. Picture of three
exemplary initial states in a chain of 40Ca+ ions ($\ket{\uparrow}$ dark
spots, $\ket{\downarrow}$ bright spots). White squares (circles) indicate the
intended preparation of $\ket{\uparrow}$ ($\ket{\downarrow}$), achieved with
$99\%$ fidelity per ion. D) Measured auto-correlations (green squares) for 51
ions and $\alpha=1.1$ (gray line: guide to the eye). Error bars, denoting the
standard error of the mean are smaller than the symbols (Sec. V [32]). At
short times (red shading), the spin excitation quickly relaxes to a local
equilibrium state. At late times (blue shading) global conservation laws
constrain the relaxation of spin excitations, leading to a slow power-law
decay of the auto-correlations. Insets: Auto-correlations on a double
logarithmic scale for different values of $\alpha$ highlighting the tunable
transport. Gray dashed lines are power laws with the predicted exponent from
Lévy flights. Here, $J=248$ rad/s, $129$ rad/s, and $116$ rad/s for
$\alpha=0.9$ (51 ions), $\alpha=1.1$ (51 ions), and $\alpha=1.5$ (25 ions).
Figure 2: Spatial correlation profiles. Spatially resolved correlations for A)
$\alpha=1.5$ and B) $\alpha=1.1$. The spatial correlations are measured at
short times (left), intermediate times (middle), and in the hydrodynamic late-
time regime (right). STE: analytic short-time expansion, ED: exact
diagonalization accessible only for the shorter chain of $25$ ions. In the
hydrodynamic regime (right) the measured profiles are compared to predictions
from Lévy flights. The spatial profile in A) is compatible with a Gaussian
(dashed) and in B) with a Lorentzian (solid), as supported by the reduced
$\chi^{2}$ values of the fit: A) $\chi^{2}_{\mathrm{L}}=3.9$,
$\chi^{2}_{\mathrm{G}}=1.3$; B) $\chi^{2}_{\mathrm{L}}=1.1$,
$\chi^{2}_{\mathrm{G}}=3.6$ (obtained by fitting over the central 27 sites).

Conservation laws determine the macroscopic late-time dynamics of quantum
systems. In the long-range XY model the total magnetization
$\sum_{j}\hat{\sigma}_{j}^{z}$ is conserved, which constrains the relaxation
of spin excitations at late times. The long-range couplings enable transport
of spin excitations over many lattice sites. As a consequence, the effective
dynamics is expected to be governed by classical Lévy flights comprised of
long-distance jumps [29], illustrated in Fig. 1A. These long-distance jumps
become irrelevant when the variance of the step size of excitations is finite.
In our model this is the case for $\alpha>3/2$, implying conventional
diffusion. However, the variance diverges for $\alpha\leq 3/2$ [33] and
hydrodynamics is expected to be described by the anomalous diffusion equation

$\partial_{t}\langle\sigma^{z}_{j}\rangle=D_{\alpha}\nabla^{2\alpha-1}\langle{\sigma^{z}_{j}}\rangle,$
(2)

with $D_{\alpha}$ denoting the associated transport coefficient; [32], Sec.
VIII. Therefore, transport is modified by the long-range character of the
interactions, ranging from conventional diffusion over superdiffusion to a
regime in which the spins are so strongly connected that the lattice geometry
becomes irrelevant, and a mean-field regime is entered. The predicted
dynamical phase diagram is shown in Fig. 1B.

Figure 3: Contrasting the infinite temperature background with the spin
polarized background. A) The deterministically prepared excitation of the
central ion strongly interacts with all the other excitations of the infinite
temperature background and slowly spreads through the system following the
laws of classical hydrodynamics. B) By contrast, a single excitation on top of
the fully down-polarized state
$\ket{\downarrow\cdots\downarrow\uparrow\downarrow\cdots\downarrow}$ has no
other excitations to scatter off and therefore spreads freely, exhibiting
quantum interference patterns. Data is measured for $25$ ions with power-law
exponent $\alpha=1.5$.

Hydrodynamic transport is most directly probed by creating a spin excitation
at time $t=0$ in the center of the chain and tracking how it propagates in
space and time, as measured by the unequal-time correlation function

$C_{j}(t)=\langle{\hat{\sigma}^{z}_{j}(t)\hat{\sigma}^{z}_{0}}\rangle_{T=\infty}.$
(3)

Over the course of the quantum dynamics, the initial excitation has to scatter
off a highly excited background in order to quickly reach the hydrodynamic
regime. Therefore, we measure $C_{j}(t)$ at infinite temperatures $T$ with
magnetization $\sim 0$.

Both, measuring spin correlations at different times and preparing an infinite
temperature state is difficult in isolated, engineered quantum systems. We
overcome these challenges by expressing the infinite temperature expectation
value as an equally weighted trace over product states $\ket{\phi}$ in the
$\hat{\sigma}^{z}$ basis; Sec. X [32]. When preparing the central spin in the
same polarization, we have $\hat{\sigma}_{0}^{z}\ket{\phi}=+\ket{\phi}$; Fig.
1C. The correlation function can then be directly evaluated as
$C_{j}(t)\sim\braket{\phi_{1}}{\hat{\sigma}_{j}^{z}(t)}{\phi_{1}}+\braket{\phi_{2}}{\hat{\sigma}_{j}^{z}(t)}{\phi_{2}}+\ldots$.
It can be accurately obtained by sampling a finite number of initial product
states; Sec. X [32]. Yet, for a comparatively small and experimentally
accessible number of initial states, large statistical fluctuations are
expected in the measured correlations. We remove these fluctuations by
sampling pairs of conjugate product states, $\ket{\phi}$ and
$\ket{\overline{\phi}}$, where in the second configuration all spins are
flipped except for the central one. For each pair of product states, initial
correlations are unity in the center of the system while being zero elsewhere,
reproducing directly this property of the full trace; see first two initial
states in Fig. 1C. With this procedure, convergence is already achieved for a
small number of initial product states; Sec. X [32]. For $\alpha=0.9$ and
$1.1$ ($\alpha=1.5$), we create 60 (120) initial product state configurations
each of which is realized, evolved, and measured for 50 to 200 times.

The auto-correlation function $C_{0}(t)$ determines the residual excitation in
the center of the chain. We measure $C_{0}(t)$ for 51 ions; Fig. 1D. The
relaxation dynamics occurs in multiple stages. A local equilibrium is reached
already after a few collisions with the abundant excitations of the infinite
temperature state. As a consequence, at short times the auto-correlation
exhibits a rapidly damped oscillation. At later times the system enters the
hydrodynamic regime and eventually approaches a global equilibrium through the
slow rearrangement of spin excitations constrained by the conserved
magnetization. The power-law decay of the correlations becomes manifest on a
double logarithmic scale; see insets for various values of $\alpha$. Our
experimental data is consistent with the hydrodynamic theory of Eq. (2), which
predicts a power-law exponent of $1/(2\alpha-1)$. Exploiting the experimental
control over the long-range exponent $\alpha$, hence enables us to tune the
transport from normal diffusion to anomalous superdiffusion.

Figure 4: Extracting transport coefficients. At late times the correlations
exhibit a self-similar scaling, relating space and time with a microscopic
transport coefficient $D_{\alpha}$. Rescaled correlations $\tau_{x}C_{j}(t)$
as a function of rescaled space $j/\tau_{x}$, where
$\tau_{x}=(Jt)^{1/(2\alpha-1)}$. A) $\alpha=0.9$, 51 ions, B) $\alpha=1.1$, 51
ions, and C) $\alpha=1.5$, 25 ions. Correlations are shown for times $Jt>5$
(darker colors correspond to later times).

A more stringent test of the emergent hydrodynamics is obtained from the full
spatial correlation profile; Fig 2. We realize the long-range spin model for
$\alpha=1.5$ with 25 ions and for $\alpha=1.1$ with 51 ions. At early times
(3.1ms) the quantum dynamics is well-described by an analytic short-time
expansion of the equations of motion ([32] Sec. IX). At intermediate times
(14.1ms), the excitation starts spreading through the system, but some quantum
coherence remains, indicated by the spatial oscillations. For 25 ions, the
measured dynamics compares well with results obtained from exact
diagonalization, demonstrating the coherence of our quantum system; Sec VII
[32]. For the 51 ion chain a comparison to exact diagonalization is not
possible due to the exponential growth of the Hilbert space, reaching a
dimension of $2^{51}\approx 10^{15}$.

At the latest times shown, interactions have entangled local degrees of
freedom with the rest of the system. Quantum interference patterns are
averaged out and the hydrodynamic regime is entered. This becomes even more
apparent in the space-time correlations; Fig. 3A. By contrast, a single
excitation on top of a spin-polarized state cannot scatter, and the associated
correlations exhibit coherent space-time oscillation patterns instead [30];
Fig. 3B.

The theory of Lévy flights predicts the following scaling form for the spatio-
temporal profiles

$C_{j}(t)=(D_{\alpha}t)^{-\frac{1}{2\alpha-1}}F_{\alpha}\left(\frac{|j|}{(D_{\alpha}t)^{1/(2\alpha-1)}}\right),$
(4)

where $F_{\alpha}$ is given by the family of stable symmetric distributions
[33]. The scaling between space and time entering the distribution
$F_{\alpha}$ can be deduced directly from the generalized hydrodynamic
equation (2). The distribution $F_{\alpha}$ cannot be expressed in terms of
elementary functions, except for $\alpha=3/2$, where a Gaussian indicates
conventional diffusion, and $\alpha=1$, where a Lorentzian appears.

The measured hydrodynamic profiles are compatible with a Lorentzian for
$\alpha=1.1$ ($F_{\alpha=1.1}$ is still very close to a Lorentzian) and with a
Gaussian for $\alpha=1.5$, see $\chi^{2}$ analysis in Fig. 2 (right panels),
in agreement with Eq. 4.

One of the most striking predictions of hydrodynamics is the self-similar
scaling of the correlations, relating time and space by a universal dynamical
scaling exponent, given by $1/(2\alpha-1)$, and a non-universal constant, the
transport coefficient $D_{\alpha}$. In analogy with the universal scaling of
correlation functions near a second-order equilibrium phase transition, the
hydrodynamic scaling indicates the proximity of the dynamics to a thermal
fixed point. While the scaling form of the data is universal and can be
predicted from purely hydrodynamic reasoning, the transport coefficient
$D_{\alpha}$ depends on the full quantum many-body spectrum and is therefore
challenging to predict from analytical or numerical methods. From a scaling
collapse of the experimental data, shown in Fig. 4, we obtain a transport
coefficient of
$D_{\alpha}/J=0.5_{-0.1}^{+0.2},0.8_{-0.2}^{+0.3},2.6_{-0.7}^{+0.9}$ for
$\alpha=0.9,1.1,1.5$ (Sec. VI [32]).

Large, coherent quantum systems with local control provide key insights into
fundamental properties of non-equilibrium quantum states. In our experiment,
we measure the emergent macroscopic hydrodynamics in an infinite temperature
state of a strongly interacting spin chain; a regime that is notoriously
challenging to describe by controlled analytical or numerical calculations. By
measuring the full spatio-temporal profile of the hydrodynamic scaling
functions, we experimentally establish a tunable family of transport ranging
from conventional diffusion to anomalous superdiffusion.

An exciting prospect is to investigate transport in models with enhanced
symmetries, which can be realized in our setting by suitably designed Floquet
protocols [34]. Moreover, our tools for measuring high-temperature
correlations can be readily applied to other quantum devices with local
control, including quantum gas microscopes, Rydberg atom arrays, and
superconducting qubits. New classes of hydrodynamics can be realized and
probed in that way.

Note added: During the completion of this manuscript, we became aware of
related work demonstrating superdiffusive transport in an integrable
Heisenberg chain with nearest-neighbor superexchange interactions [35].

Acknowledgements. We acknowledge support from the Technical University of
Munich - Institute for Advanced Study, funded by the German Excellence
Initiative and the European Union FP7 under grant agreement 291763, the Max
Planck Gesellschaft (MPG) through the International Max Planck Research School
for Quantum Science and Technology (IMPRS-QST), the Deutsche
Forschungsgemeinschaft (DFG, German Research Foundation) under Germany’s
Excellence Strategy–EXC–2111–390814868. The project leading to this
application has received funding from the European Union’s Horizon 2020
research and innovation programme under grant agreement No 817482, from the
European Research Council (ERC) under the European Union’s Horizon 2020
research and innovation programme (grant agreement No. 851161, No. 741541 and
No. 771537), the Institut für Quanteninformation GmbH, and Austrian Science
Fund through the SFB BeyondC (F7110).

Author contributions. The research was devised by MJ, CM, MK, and CFR. AS, IL
and MK developed the theoretical protocols. MJ, FK, CM, RB, and CFR
contributed to the experimental setup. MJ and FK performed the experiments.
MJ, AS, IL, MK and CFR analyzed the data and MJ, AS, and CFR carried out
numerical simulations. MJ, AS, MK, and CFR wrote the manuscript. All authors
contributed to the discussion of the results and the manuscript. All
correspondence should be addressed to MK<EMAIL_ADDRESS>and CFR
<EMAIL_ADDRESS>
Competing interests. The authors declare no competing interests.

## Supplemental materials

##  I Experimental setup

The experiments are performed with a string of 40Ca+ ions confined in a linear
radiofrequency trap. For confining long ion strings, the trapping potential is
made strongly anisotropic, with oscillation frequencies of 2.93 MHz and 2.898
MHz in the radial plane, and an axial trapping frequency of 126.3 kHz (126.6
kHz) for trapping 25 ions (51 ions). At the above stated confining
frequencies, the 25-ion (51-ion) string has a length of 173 $\mu$m (247
$\mu$m). (Pseudo-)spins are encoded in superpositions of 4S${}_{1/2},m=1/2$
and 3D${}_{5/2},m=5/2$ electronic states and coherently manipulated by a
frequency-stable laser at 729 nm coupling to the quadrupole transition. A low-
noise magnetic field of 4.17 Gauss produced by an assembly of SmCo magnets is
used to lift the degeneracy of the Zeeman states.

At the beginning of each experiment, all transverse collective modes of the
ion string are cooled close to the ground state by sideband cooling on the
S${}_{1/2}\leftrightarrow\text{D}{}_{5/2}$ quadrupole transition. The axial
modes are polarization-gradient cooled to sub-Doppler temperatures and the
ions are prepared in the S${}_{1/2},m=+1/2$ state by frequency-resolved
optical pumping [36].

At the end of each experimental run, a quantum state measurement projecting
onto the (pseudo-)spin basis states is carried out by spatially resolved
measurements of the ions’ fluorescence emitted on the S1/2 to P1/2 transition.
The state assignment error probability is below $10^{-3}$ per ion.

Figure 5: Initial state preparation fidelities of the A) 25-ion and B) 51-ion
strings. Individual ions are prepared in spin up/down states with about $99\%$
fidelity. Blue bars indicate the magnetization of the center ion whereas red
bars indicate the total magnetization of all other ions.

##  II Initial state preparation

Preparation of the (pseudo-)spins in the desired initial product states is
achieved by using a spin-echo sequence that combines resonant laser pulses
coupling to all spins with nearly the same strength and off-resonant pulses
addressed to individual spins that are sandwiched between the global pulses.
For a 25-ion crystal, the state preparation starts with a set of composite
$\pi/2$ pulses that prepares all calcium ions in an equal superposition of the
4S${}_{1/2},m=1/2$ and 3D${}_{5/2},m=5/2$ electronic states. A tightly focused
laser beam, with beam radius of about 2 $\mu$m, is steered over the ion string
(ion spacing $>$ 4 $\mu$m) with an acousto-optic deflector. Pre-calibrated
addressing pulses are placed in a specific order within the two halves of the
spin-echo sequence such that errors resulting from cross-talk while addressing
neighbouring ions are minimized. The pulse sequence is completed by rotating
the unaddressed (addressed) ions to the S (D) state by the second global
$\pi/2$ pulse. For a 51-ion crystal, we perform the above described sequence
but on 4S${}_{1/2},m=1/2$ and 3D${}_{5/2},m=-3/2$ electronic states before
restoring the population to the desired qubit states. Due to the restriction
of the polarization angle and beam geometry, we only achieve optimum off-
resonant coupling, thus short interaction times, with the addressing beam
while working on 4S${}_{1/2},m=1/2$ and 3D${}_{5/2},m=-3/2$ states.

Composite laser pulses: Ions are collectively excited with an elliptically
shaped laser beam with a beam diameter of about 470 $\mu$m, which is shone
over the ion string of length 247 $\mu$m from the transverse direction of the
ion string. This leads to a non-uniform illumination, resulting in spatially
dependent rotation angles of the ions when performing the global single qubit
rotations. In our experiment, we reduce the impact of these unequal rotations
of the spins by employing Composite pulse sequences [37]. The Composite pulse
sequences for $\theta=\pi/2$ and $\theta=\pi$ angle rotations are given as

$\displaystyle
U^{x}_{\theta=\pi/2}=e^{-i\frac{\pi}{4}\sigma_{y}}e^{-i\frac{\pi}{4}\sigma_{x}},$
(5) $\displaystyle
U^{x}_{\theta=\pi}=e^{-i\frac{\pi}{4}\sigma_{x}}e^{-i\frac{\pi}{2}\sigma_{y}}e^{-i\frac{\pi}{4}\sigma_{x}}.$
(6)

Off-resonant rotations with a tightly focused beam: The tightly focused off-
resonant laser beam is used in conjunction with the global laser beam. Pulse
length and laser intensity of the addressing beam are chosen such that the
specific ions are prepared in the desired quantum states with a single ion
fidelity of about $99\%$. Addressing of 51(25) ions, required about 0.5(0.25)
ms of interaction time, which is much shorter than the coherence time ($\sim$
50 ms) of our system, hence the initial state is well protected against any
decoherence effects. Slow spatial drift of the ion string along the weak
confining axis, which impacts the quality of addressing operations, is
corrected on-the-fly by routinely measuring the ion positions and
consecutively applying compensating potentials to the trap electrodes. Laser
beam intensities of both global and addressing laser beams are stabilized via
a sample and hold intensity stabilization technique. Furthermore, the
addressing pulse lengths are routinely calibrated and corrected against any
unwanted slow variations.

Product state fidelities: Fig. 5 shows initial state preparation fidelities of
the 25-ion and 51-ion strings. Each spin is prepared in the desired up/down
basis states with about 98.8$\%(99.3\%)$ fidelity, which corresponds to a
fidelity of about $73$%$(70\%)$ of the desired product state for a 25(51) ion
string.

Figure 6: Spin-spin interaction matrix terms displayed for A) $\alpha=0.9$, B)
$\alpha=1.1$, and C) $\alpha=1.5$. Bars in red represent experimentally
realized interaction elements and bars in blue represent the interaction
elements that are approximated by a power law. Figure 7: A) Experimental
spread of excitation as a function of entangling laser-ion interaction time
for $\alpha=1.1$ (51 ions). Here we initialize the middle ion as spin-up and
remaining ions as spin-down. B) Numerical simulations performed for the same
experimental conditions.

##  III Engineered spin-spin interaction and many body dynamics

Spin-spin interactions are engineered with a two-tone laser beam inducing
Mølmer-Sørensen type interactions; i.e. the laser off-resonantly couples to
both the lower and upper vibrational sidebands of all $2N$ collective motional
modes describing the ion motion in the plane that is normal to the direction
of an $N$-ion string. The tunable long-range coupling is provided by a
specific choice of laser parameters that generates a long-range Ising-type
spin-spin interaction with a coupling strength approximately described by a
power-law [38, 39]. The laser beam is elliptically shaped and elongated along
the axis of the ion string to provide near-uniform illumination over the ion
string. Unwanted light shifts, which occur due to the presence of other
electronic transitions, are compensated with a third laser beam that is mixed
with the two-tone laser beam [30]. The laser beam propagates from the
transverse direction of the axis of the ion string and imparts a near flat
waveform over the ion string, which otherwise can have adverse effects in the
laser-ion interaction.

The engineered spin-spin interaction, for $B\gg J$ gives rise to a flip-flop
type interaction described by

$\hat{H}=\sum_{i<j}\frac{J}{|i-j|^{\alpha}}(\hat{\sigma}^{+}_{i}\hat{\sigma}^{-}_{j}+\hat{\sigma}^{-}_{i}\hat{\sigma}^{+}_{j}),$
(7)

where exponent $\alpha$ and strength of the spin-spin interaction $J$ can be
tuned by changing the detuning and intensity of the entangling laser beam. The
experimentally realized Hamiltonian slightly differs from the power law
expression given in equation 7. The exact spin-spin matrix connecting all
spins is given by

$J_{ij}=\frac{\Omega_{i}\Omega_{j}}{2}\sum_{m=1}^{2N}\frac{\eta_{im}\eta_{jm}}{\Delta_{m}},$
(8)

where $\eta_{jm}$ is the Lamb-Dicke parameter of the $j^{th}$ ion and the
$m^{th}$ motional mode. $\Delta_{m}$ is the detuning of the two-tone laser
field from the motional mode of interest and $\Omega_{i}$ is the Rabi
frequency of the $i^{th}$ ion. The detuning of the laser beams from the
vibrational sideband ranges from a few tens of kHz for the center-of-mass
(COM) mode to about 1.5 MHz for the lowest-frequency mode, resulting in
effective spin-spin interactions that decay with distance between the spins.
The structure of the coupling matrix can be experimentally controlled by
varying the laser detuning from the COM mode and the anisotropy of the
trapping potential. A comparison between the approximated power law and
experimentally realized interaction Hamiltonian is shown in Fig. 6 for all
three $\alpha$ values that are used in the main manuscript. The figure
displays a small mismatch between the power-law interaction of eq. 7 and the
effective interaction realized in the trapped ion systems [39], which is
inherently present in any the laser-ion interaction. Additionally, due to
technical reasons, a non-uniform illumination of the laser beam over the ion
chain alters the spin-spin interaction and slightly weakens the interactions
of the spins at the ends of the chain. The numerical simulations presented in
the main manuscript are performed for the experimentally realized spin-spin
interaction which includes the non-uniform laser intensity over the ion
string.

The engineered entangling interaction is experimentally tested by time
evolving an initial state where the middle ion (remaining ions) is initialized
in a spin-up (spin-down) state. In Fig. 7(a) we show dynamics of spin
excitation for $\alpha=1.1$. We compare the experimental results with
numerical simulations that we carry out in a subspace where the magnetization
is conserved, see Fig. 7(b) and find a good agreement at early times whereas
some discrepancy is visible at late times. This discrepancy is attributed to
decoherence effects that arise from non-perfect initial state preparation and
non-conservation of magnetisation, discussed in the following section.

##  IV Decoherence effects

While experimentally realizing the spin-spin Hamiltonian, the system undergoes
some unwanted decoherence effects that can alter the dynamics. In this
section, we describe the two major decoherence channels that are important for
the initial state preparation and dynamics with long ion strings. The first
one is dephasing noise, which arises from laser frequency/phase noise and
fluctuating magnetic fields. Jointly, they reduce the coherence time of our
spins. Most importantly, the superposition, which is an essential step for
individual control of the spins for preparing arbitrary product states, decays
due to the dephasing noise. The laser frequency noise is reduced by locking
the laser to an ultra-stable reference cavity. The trap setup is placed inside
a mu-metal shield to suppress the magnetic field noise. The overall coherence
time of the system is about 50 ms. Furthermore, we reduce the impact of the
dephasing noise by applying spin-echo sequences. While these effects play an
important role for the state preparation, the dynamics is not affected, as we
work in the decoherence-free-subspace in which the spins are protected against
the dephasing noise.

The second type of decoherence that spins experience can be described by
amplitude damping and depolarization channels. Due to the finite lifetime of
the metastable state D5/2 ($\tau\sim$1 s), in which our spin-up state is
encoded, decays to the S1/2 electronic state. In addition to this decay
channel, the spins experience unwanted excitation/decay due to the bichromatic
laser field coupling to the motional sidebands. These processes violate the
conservation of magnetization that would be expected from the flip-flop
interaction of eq. (7). We experimentally measure the spin-flip rates for the
experimental conditions that are used in the main manuscript. In Fig. 8 we
show magnetization variation as a function of interaction time for a 51-ion
chain. Here, we model two spin flip rates: spontaneous decay of the excited
state $\Gamma$ and the laser induced spin flips rate $\gamma_{\text{flip}}$,
similar to what has been described in reference [10]. The model assumes that
the spin-flip rate due to incoherent laser-ion coupling does not depend on the
spin orientations. The flip rates that we extract by fitting the experimental
data are $\gamma_{\text{flip}}=0.78$/s and $\Gamma=0.91$/s. These two
competing processes in the end give rise to a dominating decay of the spins
that are encoded as spin up over those are encoded as spin down. For the
measurements that are presented in this manuscript, we remove the biased decay
of the excitation to leading order by preparing a second set of initial
states, where the central spin is prepared in a $\ket{\downarrow}$-state and
average over the two distinct sets of initial states.

Figure 8: Magnetization variation as a function of interaction time. Filled
circles are experimental data points and solid lines are fitted to the
experimental data. Here we randomly sample the initial states while keeping
the total magnetization to be +1. The measurement protocol is identical to the
one we described in Fig. 1 of the main manuscript. Error bars are estimated to
be 1/$\sqrt{N_{m}}$, where $N_{m}$ is number of measurements.

The magnetization decay model, which we assumed for fitting of the
experimental data, considers that all the ions have the same spin flip rate.
However, we experimentally observed that this rate depends slightly on the
position of the ions in the string. The exact mechanism underpinning this
effect is still under investigation. We speculate that the spatially dependent
laser-ion interactions, coupling to the low-frequency motional modes of the
ion chain, could induce such adverse effects. Most importantly, cross-
frequency components that arise while producing the entangling laser beam
could excite motional modes which in an ideal scenario should be suppressed.
Fig. 9 displays one of these artifacts which is observed while time evolving a
spin state with all spins prepared in spin-down orientation. It may be noted
that the magnetization varies the most for the middle ions in comparison to
the outer ions. Ideally, the Hamiltonian does not couple this state to any
other states hence we expect the state to be unaltered. However, due to the
incoherent spin flips the total magnetization of the time evolved state
changes with the interaction time. We further investigate the impact of
spatially dependent incoherent processes in our current study of emergent
hydrodynamics. In Fig. 10, we show the correlation function that we
experimentally measure at $t=40$ ms. The measurement sequence is identical to
that has been described in Fig. 1 of the main manuscript, except that we
exchange the encoding of ion number 26${}^{\text{th}}$ with the first ion in
the 51-ion string. Note that the current measurements reveal no signature of
unwanted variation of the magnetization which signifies the robustness of the
measurement protocol against magnetization errors.

Figure 9: Position dependent variation of the magnetization of a 51-ion chain
after 40 ms of entangling operation. Here, all the ions are prepared in a
spin-down product state before switching on the entangling interaction.
Figure 10: Correlation measurements shown for 51 ions after 40 ms of
entangling operation, while storing the excitation at one of the edges of the
ion string. The experimental sequence is identical to that has been discussed
in the main text.

Figure 11: Spatio-temporal evolution of the correlations for A) $\alpha=0.9$
(51 ions), B) $\alpha=1.1$ (51 ions) and C) $\alpha=1.5$ (25 ions). Dashed
lines represent projection of each curve in the ‘XZ’ and ‘YZ’ planes.

##  V Error estimation

For the present measurements, the experimental error bars for an individual
ion are estimated from the quantum projection noise model [40]. In our
experiment, we repeat the measurements $N_{m}$ times for each initial state.
Each initial state is chosen from a uniform random distribution function that
samples the individual spin state with an equal probability of spin-ups and
spin-downs, while restricting the total magnetization of the whole system to
be $+1$. The measurements are performed $N_{u}$ times. The mean standard
deviation of measurements for $i^{\text{th}}$ ion, while preparing the initial
state $\ket{\Phi_{u}}=\bigotimes_{i=1}^{N}\ket{\phi_{u}^{i}}$, is given by

$\sigma_{u}^{i}=\sqrt{p_{u}^{i}(1-p_{u}^{i})/N_{m}},$ (9)

where $p_{u}^{i}$ is the probability of finding the $i^{th}$ ion, prepared in
$\ket{\phi_{u}^{i}}$ state, in the measurement basis. We further propagate
this formula for the auto-correlation function described in Eq. (3) of the
main text and estimate the mean standard deviation to be

$\sigma_{corr}=2\frac{\sqrt{\sum_{u=1}^{N_{u}}\sigma_{u}^{2}}}{N_{u}}.$ (10)

.

##  VI Transport coefficients

To extract transport coefficients from the experimental data, we perform a
least-square fitting with the spatio-temporal function

$C_{j}(t)=(D_{\alpha}t)^{-\beta}F_{\alpha}\left(\frac{|j|}{(D_{\alpha}t)^{\beta}}\right),$
(11)

where $D_{\alpha}$ represents the transport coefficient and the exponent
$\beta$ determines the rate at which auto-correlations decay with time.
$F_{\alpha}$ are the stable symmetric distributions defined in Eq. 19 below.

Fitting of scaling functions: In Fig. 2 of the main text we fitted a
normalized Lorentzian and Gaussian to the data points which are finite within
error bars. Hence, the full profile was fitted for $\alpha=1.5$, while the fit
was restricted to the central 27 sites for $\alpha=1.1$.

Extraction of diffusion coefficients:  In Fig. 4 of the main text, we used the
expected values for the exponent $\beta=1/(2\alpha-1)$ to rescale the data. We
then fitted the Ansatz in Eq. 11 to the rescaled data to extract the diffusion
coefficient $D_{\alpha}$. The error bars are deduced by taking into account
the uncertainty in extracting the $\alpha$ parameters.

Spatio-temporal evolution of correlations:  The unscaled spatio-temporal
evolution of correlations for all three $\alpha$ values are shown in Fig. 11.

##  VII Auto-correlation function

Figure 12: Auto-correlation function from experiment and exact diagonalization
(ED). In ED we used the experimental $J_{ij}$ matrix for $\alpha=1.5$ and
$L=25$ as well as the $240$ initial states used in the experiment.

In Fig. 12 we compare the experimentally measured auto-correlation function
with an exact simulation for the 25 ions chain, finding excellent agreement.

##  VIII Anomalous transport from Lévy flights

Transport in the long-range XY model (7) is described by Lévy flights [29].
The Hamiltonian connects states with an up spin at site $i$ and a down spin at
site $j$, $\ket{\uparrow_{i}\downarrow_{j}}$, with the opposite configuration
$\ket{\downarrow_{i}\uparrow_{j}}$. As these states have equal energy, we can
use Fermi’s golden rule to estimate the transition probability for a spin
exchange as

$W_{i\rightarrow
j}\equiv|\braket{\uparrow_{i}\downarrow_{j}}{\hat{H}}{\downarrow_{i}\uparrow_{j}}|^{2}=\frac{\lambda}{|i-j|^{2\alpha}},$
(12)

with $\lambda\propto J^{2}$. Moreover, $W_{i\rightarrow j}=W_{j\rightarrow
i}\equiv W_{ij}$. Hence, we can construct a classical Master equation for the
spin excitation probability density
$f_{i}\equiv(\braket{\hat{\sigma}^{z}_{i}}+1)/2$,

$\partial_{t}f_{i}=\sum_{j}W_{ij}(f_{j}-f_{i}),$ (13)

where the first and second term on the right hand side are the “gain” and
“loss” terms, respectively.

We can solve this equation by Fourier transforming to momentum $k$ and taking
both a continuum and infinite system size limit while retaining a short-
distance cutoff of unity, leading to

$\partial_{t}f_{k}=W_{k}f_{k},$ (14)

with the transition rates

$W_{k}=\left(-c_{\alpha}|k|^{2\alpha-1}+\frac{k^{2}}{3-2\alpha}\right)\lambda,$
(15)

where $c_{\alpha}=-2\Gamma(1-2\alpha)\sin(\alpha\pi)$ and $\Gamma$ is the
Gamma function. Solving Eq. (14) and Fourier transforming back for an
initially localized excitation $f_{j}=\delta_{j,0}$, we get

$f_{j}(t)=\int\frac{dk}{2\pi}e^{ikj+W_{k}t}.$ (16)

The exact form of the solution crucially depends on the value of $\alpha$. In
the following, we discuss the three cases shown in the dynamical phase diagram
in Fig. 1, main text.

Diffusion for $\alpha\geq 1.5$: In this case, the second term in Eq. (15)
dominates. Inserting $W_{k}$ into Eq. (16), we find

$f_{j}^{\alpha>1.5}(t)=(D_{\alpha}t)^{-1/2}G\left(\frac{|j|}{(D_{\alpha}t)^{1/2}}\right),$
(17)

with a diffusive Gaussian profile $G(y)=\exp(-y^{2}/4)/\sqrt{4\pi}$ and
diffusion coefficient $D_{\alpha}=\lambda/(2\alpha-3)$. The long-range
character of the transition rates only shows up in subleading corrections,
leading to algebraic tails of the distribution function [41, 29].

Superdiffusion for $0.5<\alpha<1.5$: In this case, we get

$f_{j}^{0.5<\alpha<1.5}(t)=(D_{\alpha}t)^{-1/(2\alpha-1)}F_{\alpha}\left(\frac{|j|}{(D_{\alpha}t)^{1/(2\alpha-1)}}\right),$
(18)

with $D_{\alpha}=\lambda c_{\alpha}$. The distribution function $F_{\alpha}$
is given by the stable symmetric distributions defined by

$F_{\alpha}(y)=\int\frac{dk}{2\pi}e^{iyk-|k|^{2\alpha-1}},$ (19)

which only has a known elementary solution for $\alpha=3/2$ (Gaussian) and
$\alpha=1$ (Lorentzian).

Mean-field regime for $\alpha<0.5$: For $\alpha<0.5$, the system enters the
mean-field (or unstable transport) regime. Here, a consistent hydrodynamic
description is no longer possible in the thermodynamic limit, as revealed by
the diverging transition rates in Eq. (15) for $k\rightarrow 0$. Instead, the
system forms a zero-dimensional “quantum dot” rather than a one-dimensional
chain, and relaxes to thermal equilibrium instantaneously [29, 42], due to the
non-local nature of the interactions. This behavior can be best understood
from the limit of all-to-all interactions $\alpha=0$, allowing to rewrite the
Hamiltonian in terms of a single global spin operator $\vec{\hat{\sigma}}_{\rm
global}=\sum_{i=1}^{L}\vec{\hat{\sigma}}_{i}$. In the thermodynamic limit, the
diverging spin eigenvalue associated with this operator gives rise to a
vanishing relaxation time.

##  IX Short time expansion

Here, we calculate the short time dynamics of

$\operatorname{Tr}\left[\hat{\sigma}^{z}_{i}(t)\hat{\sigma}^{z}_{j}(0)\right]$
(20)

for time evolution under the Hamiltonian

$\hat{H}=\sum_{i<j}J_{ij}(\hat{\sigma}^{+}_{i}\hat{\sigma}^{-}_{j}+h.c.).$
(21)

We expand $\hat{\sigma}^{z}_{i}(t)$ in the Heisenberg picture,

$\displaystyle\hat{\sigma}^{z}_{i}(t)=e^{i\hat{H}t}\hat{\sigma}^{z}_{i}e^{-i\hat{H}t}$
$\displaystyle=\hat{\sigma}_{i}^{z}+iJt\left[\frac{\hat{H}}{J},\hat{\sigma}_{i}^{z}\right]+\frac{(iJt)^{2}}{2!}\left[\frac{\hat{H}}{J},\left[\frac{\hat{H}}{J},\hat{\sigma}_{i}^{z}\right]\right]+\mathcal{O}((Jt)^{3}).$

When evaluating the commutators and taking the trace, the first order
contribution vanishes and we get

$\displaystyle\operatorname{Tr}[\hat{\sigma}^{z}_{i}(t)\hat{\sigma}^{z}_{j}(0)]=$
$\displaystyle t^{2}\sum_{k\neq
i}J_{ik}^{2}\left(\operatorname{Tr}[\hat{\sigma}_{k}^{z}\hat{\sigma}_{j}^{z}]-\operatorname{Tr}[\hat{\sigma}_{i}^{z}\hat{\sigma}_{j}^{z}]\right)$
$\displaystyle+$
$\displaystyle\operatorname{Tr}[\hat{\sigma}_{i}^{z}\hat{\sigma}_{j}^{z}].$

Furthermore,
$\operatorname{Tr}[\hat{\sigma}_{k}^{z}\hat{\sigma}_{m}^{z}]=\delta_{km}$ and
hence

$\operatorname{Tr}[\hat{\sigma}^{z}_{i}(t)\hat{\sigma}^{z}_{j}(0)]=\bigg{\\{}\begin{array}[]{lr}1-t^{2}\sum_{i\neq
k}J_{ik}^{2}&\text{for }i=j\\\ t^{2}J_{ij}^{2}&\text{for }i\neq j\\\
\end{array}.$ (22)

In Fig. 3 of the main text we used the explicit $J_{ij}$ matrix determined
from the laser parameters in the experiment to obtain the short time
expansion.

Figure 13: Comparing the evaluation of the full trace with sampling of product
initial states. Data obtained from exact diagonalization on an $L=19$ chain
with $\alpha=1$. For product state sampling, we used $M=240$ initial states,
for full trace we used $10$ Haar-random states.

##  X Product state sampling of the trace

The infinite temperature correlation function in Eq. (3) of the main text can
be written as

$\displaystyle C_{j}(t)$
$\displaystyle=\braket{\hat{\sigma}^{z}_{j}(t)\hat{\sigma}^{z}_{0}(0)}_{T=\infty}$
(23)
$\displaystyle=\frac{1}{Z}\mathrm{Tr}\left[e^{-\hat{H}/T}\hat{\sigma}^{z}_{j}(t)\hat{\sigma}^{z}_{0}(0)\right]\bigg{|}_{T=\infty}$
(24)
$\displaystyle=\frac{1}{Z}\mathrm{Tr}\left[\hat{\sigma}^{z}_{j}(t)\hat{\sigma}^{z}_{0}(0)\right],$
(25)

with the partition sum $Z=\mathrm{Tr}[e^{-\hat{H}/T}]|_{T=\infty}=2^{L}$ for
systems of size $L$. We expand the trace in the $\hat{\sigma}^{z}$ basis to
obtain

$\displaystyle C_{j}(t)$
$\displaystyle=\frac{1}{Z}\sum_{\sigma_{-\frac{L}{2}},\dots,\sigma_{\frac{L}{2}}}\braket{\sigma_{-\frac{L}{2}}\dots\sigma_{\frac{L}{2}}}{\hat{\sigma}^{z}_{j}(t)\hat{\sigma}^{z}_{0}(0)}{\sigma_{-\frac{L}{2}}\dots\sigma_{\frac{L}{2}}}$
(26)
$\displaystyle=\frac{1}{Z}\sum_{\sigma_{-\frac{L}{2}},\cdots,\sigma_{\frac{L}{2}}}\sigma_{0}\braket{\sigma_{-\frac{L}{2}}\cdots\sigma_{\frac{L}{2}}}{\hat{\sigma}^{z}_{j}(t)}{\sigma_{-\frac{L}{2}}\dots\sigma_{\frac{L}{2}}},$
(27)

where $\sigma_{j}=\pm 1$. We have therefore reduced the evaluation of a two-
time correlation function to the weighted sum of single-time functions using
product initial states, which can be realized in the experiment. In order to
evaluate Eq. (27) exactly, all $2^{L}$ product states in the z-basis would
have to be prepared, time evolved and measured to then obtain $C_{j}(t)$ from
the weighted sum. In practice, however, it suffices to sample a small number
$M\approx 10-100$ [29] of product states to evaluate Eq. (27) while replacing
$Z\rightarrow M$. In the main text, we considered the largest sector with
magnetization $\sim 0$. In Fig. 13 we compare the sampling of the largest
sector to the exact results of the thermal expectation value. The full trace
is evaluated using the typicality approach [43], which involves evolving Haar-
random states $\ket{\psi_{r}}$ in the entire Hilbert space and then averaging
over $5-10$ such states.

To further improve the convergence of the sampling of Eq. (27) we prepare the
conjugate configuration for each randomly drawn product state. In the
conjugate state, all spins are flipped except for the central one. This choice
of initial states reduces statistical noise by reproducing the exact property
$C_{j}(0)=\delta_{j,0}$ of the full trace for each pair of initial states.

## References

  * Gring  _et al._ [2012] M. Gring, M. Kuhnert, T. Langen, T. Kitagawa, B. Rauer, M. Schreitl, I. Mazets, D. A. Smith, E. Demler, and J. Schmiedmayer, Relaxation and prethermalization in an isolated quantum system, Science 337, 1318 (2012).
  * Bohnet  _et al._ [2016] J. G. Bohnet, B. C. Sawyer, J. W. Britton, M. L. Wall, A. M. Rey, M. Foss-Feig, and J. J. Bollinger, Quantum spin dynamics and entanglement generation with hundreds of trapped ions, Science 352, 1297 (2016).
  * Tang  _et al._ [2018] Y. Tang, W. Kao, K.-Y. Li, S. Seo, K. Mallayya, M. Rigol, S. Gopalakrishnan, and B. L. Lev, Thermalization near integrability in a dipolar quantum Newton’s cradle, Phys. Rev. X 8, 021030 (2018).
  * Bernien  _et al._ [2017] H. Bernien, S. Schwartz, A. Keesling, H. Levine, A. Omran, H. Pichler, S. Choi, A. S. Zibrov, M. Endres, M. Greiner, V. Vuletić, and M. D. Lukin, Probing many-body dynamics on a 51-atom quantum simulator, Nature 551, 579 (2017).
  * Guardado-Sanchez  _et al._ [2020] E. Guardado-Sanchez, A. Morningstar, B. M. Spar, P. T. Brown, D. A. Huse, and W. S. Bakr, Subdiffusion and heat transport in a tilted two-dimensional Fermi-Hubbard system, Phys. Rev. X 10, 011042 (2020).
  * [6] S. Scherg, T. Kohlert, P. Sala, F. Pollmann, B. H. M., I. Bloch, and M. Aidelsburger, Observing non-ergodicity due to kinetic constraints in tilted Fermi-Hubbard chains, arXiv:2010.12965 .
  * Schreiber  _et al._ [2015] M. Schreiber, S. S. Hodgman, P. Bordia, H. P. Lüschen, M. H. Fischer, R. Vosk, E. Altman, U. Schneider, and I. Bloch, Observation of many-body localization of interacting fermions in a quasirandom optical lattice, Science 349, 842 (2015).
  * Smith  _et al._ [2016] J. Smith, A. Lee, P. Richerme, B. Neyenhuis, P. W. Hess, P. Hauke, M. Heyl, D. A. Huse, and C. Monroe, Many-body localization in a quantum simulator with programmable random disorder, Nature Physics 12, 907 (2016).
  * Kaufman  _et al._ [2016] A. M. Kaufman, M. E. Tai, A. Lukin, M. Rispoli, R. Schittko, P. M. Preiss, and M. Greiner, Quantum thermalization through entanglement in an isolated many-body system, Science 353, 794 (2016).
  * Brydges  _et al._ [2019] T. Brydges, A. Elben, P. Jurcevic, B. Vermersch, C. Maier, B. P. Lanyon, P. Zoller, R. Blatt, and C. F. Roos, Probing Rényi entanglement entropy via randomized measurements, Science 364, 260 (2019).
  * Prüfer  _et al._ [2018] M. Prüfer, P. Kunkel, H. Strobel, S. Lannig, D. Linnemann, C.-M. Schmied, J. Berges, T. Gasenzer, and M. K. Oberthaler, Observation of universal dynamics in a spinor Bose gas far from equilibrium, Nature 563, 217 (2018).
  * Erne  _et al._ [2018] S. Erne, R. Bücker, T. Gasenzer, J. Berges, and J. Schmiedmayer, Universal dynamics in an isolated one-dimensional Bose gas far from equilibrium, Nature 563, 225 (2018).
  * Jurcevic  _et al._ [2017] P. Jurcevic, H. Shen, P. Hauke, C. Maier, T. Brydges, C. Hempel, B. P. Lanyon, M. Heyl, R. Blatt, and C. F. Roos, Direct observation of dynamical quantum phase transitions in an interacting many-body system, Phys. Rev. Lett. 119, 080501 (2017).
  * Zhang  _et al._ [2017] J. Zhang, G. Pagano, P. W. Hess, A. Kyprianidis, P. Becker, H. Kaplan, A. V. Gorshkov, Z.-X. Gong, and C. Monroe, Observation of a many-body dynamical phase transition with a 53-qubit quantum simulator, Nature 551, 601 (2017).
  * Lux  _et al._ [2014] J. Lux, J. Müller, A. Mitra, and A. Rosch, Hydrodynamic long-time tails after a quantum quench, Phys. Rev. A 89, 053608 (2014).
  * Bohrdt  _et al._ [2017] A. Bohrdt, C. B. Mendl, M. Endres, and M. Knap, Scrambling and thermalization in a diffusive quantum many-body system, New J. Phys. 19, 063001 (2017).
  * Cao  _et al._ [2011] C. Cao, E. Elliott, J. Joseph, H. Wu, J. Petricka, T. Schäfer, and J. E. Thomas, Universal quantum viscosity in a unitary Fermi gas, Science 331, 58 (2011).
  * Sommer  _et al._ [2011] A. Sommer, M. Ku, G. Roati, and M. W. Zwierlein, Universal spin transport in a strongly interacting Fermi gas, Nature 472, 201 (2011).
  * Schneider  _et al._ [2012] U. Schneider, L. Hackermüller, J. P. Ronzheimer, S. Will, S. Braun, T. Best, I. Bloch, E. Demler, S. Mandt, D. Rasch, and A. Rosch, Fermionic transport and out-of-equilibrium dynamics in a homogeneous Hubbard model with ultracold atoms, Nature Physics 8, 213 (2012).
  * Brown  _et al._ [2019] P. T. Brown, D. Mitra, E. Guardado-Sanchez, R. Nourafkan, A. Reymbaut, C.-D. Hébert, S. Bergeron, A.-M. S. Tremblay, J. Kokalj, D. A. Huse, P. Schauß, and W. S. Bakr, Bad metallic transport in a cold atom Fermi-Hubbard system, Science 363, 379 (2019).
  * Bandurin  _et al._ [2016] D. A. Bandurin, I. Torre, R. K. Kumar, M. Ben Shalom, A. Tomadin, A. Principi, G. H. Auton, E. Khestanova, K. S. Novoselov, I. V. Grigorieva, L. A. Ponomarenko, A. K. Geim, and M. Polini, Negative local resistance caused by viscous electron backflow in graphene, Science 351, 1055 (2016).
  * Crossno  _et al._ [2016] J. Crossno, J. K. Shi, K. Wang, X. Liu, A. Harzheim, A. Lucas, S. Sachdev, P. Kim, T. Taniguchi, K. Watanabe, T. A. Ohki, and K. C. Fong, Observation of the Dirac fluid and the breakdown of the Wiedemann-Franz law in graphene, Science 351, 1058 (2016).
  * Moll  _et al._ [2016] P. J. W. Moll, P. Kushwaha, N. Nandi, B. Schmidt, and A. P. Mackenzie, Evidence for hydrodynamic electron flow in PdCoO2, Science 351, 1061 (2016).
  * Zu  _et al._ [2021] C. Zu, F. Machado, B. Ye, S. Choi, B. Kobrin, T. Mittiga, S. Hsieh, P. Bhattacharyya, M. Markham, D. Twitchen, A. Jarmola, D. Budker, C. R. Laumann, J. E. Moore, and N. Y. Yao, Emergent hydrodynamics in a strongly interacting dipolar spin ensemble (2021), arXiv:2104.07678 .
  * Gromov  _et al._ [2020] A. Gromov, A. Lucas, and R. M. Nandkishore, Fracton hydrodynamics, Phys. Rev. Research 2, 033124 (2020).
  * Feldmeier  _et al._ [2020] J. Feldmeier, P. Sala, G. De Tomasi, F. Pollmann, and M. Knap, Anomalous diffusion in dipole- and higher-moment-conserving systems, Phys. Rev. Lett. 125, 245303 (2020).
  * Ljubotina  _et al._ [2017] M. Ljubotina, M. Žnidarič, and T. Prosen, Spin diffusion from an inhomogeneous quench in an integrable system, Nature Communications 8, 16117 (2017).
  * Bulchandani  _et al._ [2021] V. B. Bulchandani, S. Gopalakrishnan, and E. Ilievski, Superdiffusion in spin chains (2021), arXiv:2103.01976 .
  * Schuckert  _et al._ [2020] A. Schuckert, I. Lovas, and M. Knap, Nonlocal emergent hydrodynamics in a long-range quantum spin system, Phys. Rev. B 101, 020416 (2020).
  * Jurcevic  _et al._ [2014] P. Jurcevic, B. P. Lanyon, P. Hauke, C. Hempel, P. Zoller, R. Blatt, and C. F. Roos, Quasiparticle engineering and entanglement propagation in a quantum many-body system, Nature 511, 202 (2014).
  * Monroe  _et al._ [2021] C. Monroe, W. C. Campbell, L.-M. Duan, Z.-X. Gong, A. V. Gorshkov, P. W. Hess, R. Islam, K. Kim, N. M. Linke, G. Pagano, P. Richerme, C. Senko, and N. Y. Yao, Programmable quantum simulations of spin systems with trapped ions, Rev. Mod. Phys. 93, 025001 (2021).
  * [32] See Supplementary Material.
  * Zaburdaev  _et al._ [2015] V. Zaburdaev, S. Denisov, and J. Klafter, Lévy walks, Rev. Mod. Phys. 87, 483 (2015).
  * [34] S. Birnkammer, A. Bohrdt, F. Grusdt, and M. Knap, Characterizing topological excitations of a long-range Heisenberg model with trapped ions, arXiv:2012.09185 .
  * [35] D. Wei, A. Rubio-Abadal, B. Ye, F. Machado, J. Kemp, K. Srakaew, S. Hollerith, J. Rui, S. Gopalakrishnan, N. Y. Yao, I. Bloch, and J. Zeiher, Quantum gas microscopy of Kardar-Parisi-Zhang superdiffusion, to appear in the same arXiv posting .
  * Joshi  _et al._ [2020] M. K. Joshi, A. Fabre, C. Maier, T. Brydges, D. Kiesenhofer, H. Hainzer, R. Blatt, and C. F. Roos, Polarization-gradient cooling of 1d and 2d ion coulomb crystals, New Journal of Physics 22, 103013 (2020).
  * Torosov and Vitanov [2011] B. T. Torosov and N. V. Vitanov, Smooth composite pulses for high-fidelity quantum information processing, Phys. Rev. A 83, 053420 (2011).
  * Porras and Cirac [2004] D. Porras and J. I. Cirac, Effective quantum spin systems with trapped ions, Phys. Rev. Lett. 92, 207901 (2004).
  * Nevado and Porras [2016] P. Nevado and D. Porras, Hidden frustrated interactions and quantum annealing in trapped-ion spin-phonon chains, Phys. Rev. A 93, 013625 (2016).
  * Itano  _et al._ [1993] W. M. Itano, J. C. Bergquist, J. J. Bollinger, J. M. Gilligan, D. J. Heinzen, F. L. Moore, M. G. Raizen, and D. J. Wineland, Quantum projection noise: Population fluctuations in two-level systems, Phys. Rev. A 47, 3554 (1993).
  * Zarfaty  _et al._ [2019] L. Zarfaty, A. Peletskyi, E. Barkai, and S. Denisov, Infinite horizon billiards: Transport at the border between gauss and lévy universality classes, Phys. Rev. E 100, 042140 (2019).
  * Bachelard and Kastner [2013] R. Bachelard and M. Kastner, Universal threshold for the dynamical behavior of lattice systems with long-range interactions, Phys. Rev. Lett. 110, 170603 (2013).
  * Bartsch and Gemmer [2009] C. Bartsch and J. Gemmer, Dynamical typicality of quantum expectation values, Phys. Rev. Lett. 102, 110403 (2009).