Kardar-Parisi-Zhang physics in the quantum Heisenberg magnet

TL;DR

This paper demonstrates that the KPZ universality class describes magnetization dynamics in the quantum Heisenberg spin chain, revealing a new quantum analog of classical fluctuation phenomena.

Contribution

It provides the first evidence of KPZ universality in a quantum many-body system, using novel theoretical and numerical methods.

Findings

KPZ scaling function describes magnetization dynamics in the quantum Heisenberg chain

Energy conservation does not influence KPZ physics in this quantum system

Introduces new tools for analyzing quantum nonequilibrium correlations

Abstract

Equilibrium spatio-temporal correlation functions are central to understanding weak nonequilibrium physics. In certain local one-dimensional classical systems with three conservation laws they show universal features. Namely, fluctuations around ballistically propagating sound modes can be described by the celebrated Kardar-Parisi-Zhang (KPZ) universality class. Can such universality class be found also in quantum systems? By unambiguously demonstrating that the KPZ scaling function describes magnetization dynamics in the SU(2) symmetric Heisenberg spin chain we show, for the first time, that this is so. We achieve that by introducing new theoretical and numerical tools, and make a puzzling observation that the conservation of energy does not seem to matter for the KPZ physics.

Equations

$$h_{r,r+1}=J\left(s_r^{\mathrm{x}}{s}_{r+1}^{\mathrm{x}}+s_r^{\mathrm{y}}{s}_{r+1}^{\mathrm{y}}+s_r^{\mathrm{z}}{s}_{r+1}^{\mathrm{z}}\right),$$

$$U^{†}{M}_{2r}U-M_{2r}=j_{2r-1}^{\mathrm{o}}-j_{2r+1}^{\mathrm{o}},$$

$$U^{†}{M}_{2r-1}U-M_{2r-1}=j_{2r-2}^{\mathrm{e}}-j_{2r}^{\mathrm{e}},$$

$$j_{2r-1}^{\mathrm{o}}=2\sin (J)j_{2r-1}-\frac{1}{2}{\sin }^2(J/2)(s_{2r}^{\mathrm{z}}-s_{2r-1}^{\mathrm{z}}),$$

$$\rho (t=0)$$

$$⟨s_r^{\mathrm{z}}(t){⟩}_{\mu }=-\mu \sum _{r^{\prime }}{\theta }_{r^{\prime }}⟨s_r^{\mathrm{z}}(t)s_{r^{\prime }}^{\mathrm{z}}⟩+\mathcal{O}({\mu }^2),$$

$$⟨s_{r-1}^{\mathrm{z}}(t){⟩}_{\mu }-⟨s_r^{\mathrm{z}}(t){⟩}_{\mu }\approx \mu ⟨s_r^{\mathrm{z}}(t)\sum _{r^{\prime }}{\theta }_{r^{\prime }}(s_{r^{\prime }}^{\mathrm{z}}-s_{r^{\prime }+1}^{\mathrm{z}})⟩$$

$$=2\mu ⟨s_r^{\mathrm{z}}(t)s_0^{\mathrm{z}}⟩-2\mu ⟨s_r^{\mathrm{z}}(t)s_{-L/2}^{\mathrm{z}}⟩.$$

$$⟨s_0^{\mathrm{z}}(0)s_r^{\mathrm{z}}(t)⟩=\underset{\mu \to 0}{\lim }\frac{⟨s_{r-1}^{\mathrm{z}}(t){⟩}_{\mu }-⟨s_r^{\mathrm{z}}(t){⟩}_{\mu }}{2\mu }.$$

$${\partial }_{t}h=\frac{1}{2}\lambda {\left({\partial }_{r}h\right)}^2+\nu {\partial }_r^{2}h+\sqrt{\Gamma }\zeta ,$$

$$g(\phi )$$

$$f(\phi )$$

$${\partial }_{r}z(r,t)=\frac{a\mu }{t^{2/3}}f\left(\frac{br}{t^{2/3}}\right),$$

$$j(r,t)$$

$$h(\phi )$$

Full text

Introduction.– Universality – where different systems can be described by the same underlying mathematical structure – is at the core of our understanding of nature. For instance, the properties of any thermalizing system can be described by the same equilibrium ensembles of statistical physics. Out of equilibrium less is known in general, in a way justifiably so, simply because the world of nonequilibrium is much richer. One of the more famous universality classes that can (among other) describe various nonequilibrium phenomena livibook is that of the Kardar-Parisi-Zhang (KPZ) equation. The KPZ equation was originally introduced to describe stochastic growth of surfaces kardar86 , and is a diffusion equation with the simplest possible nonlinearity (relevant at large scales) and an additional white noise term (equivalently, the surface’s slope is described by the stochastic Burgers equation). Besides describing surface dynamics it can be found in various contexts, ranging from exclusion processes to random matrix theory, for review see kpzrev1 . The KPZ equation itself harbors rich mathematical problems kpzrev2 .

Nonequilibrium physics is one of the more propulsive areas of today’s theoretical physics. Close to equilibrium one can use Green-Kubo formulas and express nonequilibrium properties in terms of equilibrium correlation functions kubobook . A downside to such an approach is that the calculation of spatio-temporal correlation functions is often very complicated. Any possible universality in their long-time behavior would therefore be highly appreciated. For classical fluids in one dimension such a picture has in fact been put forward beijeren12 ; spohn14 in a form of nonlinear fluctuating hydrodynamics spohnrev that describes (anomalous) fluctuations around sound peaks due to nonlinearity in one-dimensional systems that have 3 conservation laws (momentum, energy and mass), and are in general nonintegrable. That fluctuations are indeed described by the KPZ scaling function spohn04 has been verified in a number of classical systems kulkarni13 ; mendl14 ; das14 ; das19 ; unpub . So-far there has been no observation of the KPZ universality class scaling function in quantum systems.

In this Letter we observe the KPZ scaling functions in an integrable quantum model that does not have any ballistic component. Namely, we show with an unprecedented accuracy (an order of magnitude larger than in simulations of classical systems) that an infinite temperature spin-spin correlation function in a paradigmatic SU(2) symmetric quantum Heisenberg chain has a KPZ form. Such accuracy is a result of two novelties: (i) using linear response formulation we show that one can calculate the equilibrium correlation function as an expectation value in a particular nonequilibrium state whose time evolution is easier to calculate, (ii) we directly treat an ensemble evolution, avoiding statistical averaging (as done in classical simulations), and which is, even more importantly, structurally stable. In addition, to discern the role played by conserved quantities, we show that in an integrable trotterized Floquet generalization vanicat18 of the model, that does not conserve the energy, the same KPZ scaling is observed. We note that the KPZ scaling exponents have been observed in various stochastic quantum settings, like random quantum circuits nahum17 ; nahum18 or noisy evolution lamacraft18 .

The model.– In classical systems the KPZ scaling function describes fluctuations around a sound mode, whose width scales as $\sim t^{1/z}$ with a dynamical exponent $z=\frac{3}{2}$. Therefore, to observe it one has to move to a ballistically moving reference frame, which, if the velocity is not known analytically, can introduce numerical inaccuracies. We are therefore going to look for KPZ physics at infinite temperature in the one-dimensional Heisenberg spin-$\frac{1}{2}$ chain at zero magnetization (half-filling) where the ballistic contribution is zero due to the spin-flip (particle-hole) symmetry and the spin transport shows a KPZ dynamical scaling exponent $z=\frac{3}{2}$. Namely, such superdiffusive magnetization transport has been observed in a nonequilibrium steady state where the current scales as $j\sim 1/L^{z-1}$ prl11 as well as in the spreading of an inhomogeneous initial state where the width scales as $\sigma\sim t^{1/z}$ ljubotina17 . The local Hamiltonian density is

[TABLE]

where $s^{\alpha}_{r}=\frac{1}{2}\sigma^{\alpha}_{r}$, $\alpha\in\{{\rm x},{\rm y},{\rm z}\}$, are spin operators (Pauli matrices) at site $r\in\{-\frac{L}{2},\ldots,\frac{L}{2}-1\}$. Theoretical explanation of the scaling exponent $z=\frac{3}{2}$ is still lacking, but consistent derivations within assumptions of generalized hydrodynamics were recently given vasseur18 . In particular, it is possible to estimate the diffusion constant markomprl ; denardis18 ; vasseur18 and prove its divergence, i.e. $z<2$ enejprl .

Here, in order to observe precise spatio-temporal profiles of spin and current densities, we will consider two dynamical setings: continuous time evolution $U^{t}=e^{-iHt}$ generated by $H=\sum_{r=-L/2}^{L/2-2}h_{r,r+1}$ (where we set $J=1$) or discrete time evolution with one step propagator ${U}={U}_{\rm e}{U}_{\rm o}$, with ${U}_{\rm o}=e^{-i\sum_{r}{h}_{2r-1,2r}}$ and ${U}_{\rm e}=e^{-i\sum_{r}{h}_{2r,2r+1}}$ (where we use $J=\frac{\pi}{2}$, and where one also observes the superdiffusive scaling $z=\frac{3}{2}$ ljubotina19 ). Both settings are characterized by both a global $SU(2)$ symmetry and integrability.

In order to study transport we must derive the expressions for the local spin current density operators for both the continuous-time and discrete-time models. The former is the standard spin current in the Heisenberg model ${j}_{r}={s}_{r}^{\rm x}{s}_{r+1}^{\rm y}-{s}_{r}^{\rm y}{s}_{r+1}^{\rm x}$ which fulfills the continuity equation $\frac{{\rm d}{s}^{\rm z}_{r}}{{\rm d}t}={j}_{r-1}-{j}_{r}$. The current in the discrete-time model turns out to be slightly more complicated, with the operator being different on odd and even sites due to the staggered nature of the propagator $U$. The two currents densities satisfy a pair of continuity equations

[TABLE]

where ${M}_{r}={s}_{r}^{\rm z}+{s}_{r+1}^{\rm z}$. The simpler odd current can then be seen to take the form

[TABLE]

whereas the even current is simply the odd current propagated by half a time step ${j}_{2r}^{\rm e}={U}^{\dagger}_{\rm e}{j}_{2r}^{\rm o}{U}_{\rm e}$ and acts on 4 adjacent sites.

We begin by preparing our system in a weakly polarized domain-wall mixed initial state ljubotina17

[TABLE]

An example of time evolution for both models is shown in Fig. 1, using the scaling variable $\xi=\frac{r}{t^{1/z}},z=\frac{3}{2}$. While this choice of the initial state provides a numerically stable and efficient way to study spin transport ljubotina17 , we emphasize that for our purposes it provides us with an efficient way to study the infinite-temperature spin-spin correlation function $\langle{s}^{\rm z}_{0}{s}^{\rm z}_{r}(t)\rangle$, where $A(t)\equiv U^{-t}AU^{t}$ and $\langle\,\boldsymbol{\cdot}\,\rangle\equiv 2^{-L}{\rm tr}(\cdot)$ denotes the infinite-temperature expectation value. We explain that in the following section.

Linear response.– We start by expanding the initial state (4) to linear order in $\mu$, evolving it in time, and writing down the expectation value for a single spin,

[TABLE]

where we introduced $\langle\,\boldsymbol{\cdot}\,\rangle_{\mu}={\rm tr}[(\cdot)\rho_{\mu}]/{\rm tr}\rho_{\mu}$ as the expectation value in the weak domain-wall initial state (4) and $\theta_{r}\equiv 1(-1)$ for $r\geq 0(<0)$. Accounting for the translational invariance of the infinite-temperature expectation value we obtain

[TABLE]

In the thermodynamic limit $L\to\infty$ the second term vanishes as there are no correlations across infinite distances, and using the cyclic property of the trace we get

[TABLE]

This is our first main result.

It shows that a weak domain wall initial state can be seen as a trick that allows us to calculate the infinite-temperature spin-spin correlation. We next recall spohnrev why the LHS of Eq.(7) is in certain classical systems described by the KPZ scaling function.

Kardar-Parisi-Zhang equation.– The KPZ stochastic partial differential equation was initially suggested to model the growth of surface $h(r,t)$ through random deposition kardar86

[TABLE]

where $\zeta(r,t)$ is a space-time uncorrelated noise.

Of particular interest to us will be the correlation function $C(r,t)=\langle\left[h(r,t)-h(0,0)-t\langle\partial_{t}h\rangle\right]^{2}\rangle$ – representing the fluctuations of the height around the expected value – and its second derivative $\frac{1}{2}\partial_{r}^{2}C(r,t)=\langle\partial_{r}h(0,0)\partial_{r}h(r,t)\rangle$ – describing the slope correlations (here brackets denote noise averaging). In terms of scaling functions $g(\varphi)$ and $f(\varphi)$ one has

[TABLE]

These can be obtained from the exact solution of the polynuclear growth model spohn04 (a model in the KPZ universality class), and have been tabulated with high precision in Ref. prahoferwww . Nonlinear fluctuating hydrodynamics predicts that the correlation function of a conserved quantity, in our case $\langle{s}_{0}^{\rm z}(0){s}_{r}^{\rm z}(t)\rangle$, should be given by the so-called KPZ scaling function $f(\varphi)$.

Using Eq.(7) this correlation function is equal to the magnetization difference on consecutive sites in the state $\rho(t)\propto U^{t}\rho_{\mu}U^{-t}$ (see Fig. 1). However, taking the discrete derivative increases numerical errors, so alternatively, one can also look at the scaling form of the current $j(r,t)=\langle{j}_{r}(t)\rangle_{\mu}$. In a diffusive process, the scaling forms of both the current as well as of the magnetization difference are Gaussian. Relation between the two in a general non-diffusive situation can be derived from the continuity equation.

Defining a shorthand notation $z(r,t)=\langle{s}_{r}^{\rm z}(t)\rangle_{\mu}$, and $\varphi=b\xi$, we write an ansatz

[TABLE]

where we introduced two system-dependent parameters $a$ and $b$, and use continuum notation for the magnetization difference. Taking into account the continuity equation $\partial_{t}z=-\partial_{r}j$, one may obtain the shape of the spin current profile. Expressing everything in terms of $g(\varphi)$ (using per-partes integration and Eq. (9)) we get

[TABLE]

The form of $j(r,t)$, i.e. the function $h(\varphi)$, is therefore uniquely determined by the form of $\partial_{r}z(r,t)$, i.e., the KPZ function $f(\varphi)$.

We employ extensive numerical simulations foot3 using the time-evolving block decimation algorithm vidal03 ; vidal04 ; schollwock11 for matrix-product density operator in order to study the time evolution of a domain-wall like initial state in both the continuous and discrete time Heisenberg models. This allows us to compute the infinite-temperature spin-spin correlations (7) in a numerically stable way with manageable bond dimensions $\chi$. Fig. 2 shows the results and the best-fitting KPZ profile for both the spin and spin current. Due to higher numerical accuracy we only fit the data for the current, obtaining $a$ and $b$ (11), which automatically fixes the spin difference profiles (10). In order to avoid even-odd staggering in the discrete-time model we take the difference of two consecutive pairs of spins, rather than a difference of two spins, and appropriately scale the continuity equation. For comparison we also show best-fitting Gaussians. Because the KPZ scaling functions $f(\varphi)$ and $h(\varphi)$ are rather close to Gaussians for not too large arguments, one in fact needs at least two decades of accuracy to be able to distinguish the two. With our numerics we have accuracy over about three decades in the continuous model and about four in the discrete one. We can clearly confirm that the KPZ scaling functions emerge at sufficiently long times.

Free parameters $a$ and $b$ are found to be $a=b\approx 0.67$, conjectured to be $\frac{2}{3J^{2/3}}$, for the continuous-time model. Similarly, for the discrete-time model we find $a=b\approx 0.43$, data for other values of $J$ are well described by the formula $a=b\approx\frac{2^{1/3}}{3|\tan(J/2)^{2/3}|}$.

Because the KPZ $f(\varphi)$ and $h(\varphi)$ are not Gaussian, their ratio $h/f\equiv w(br/t^{2/3})$ which appears in a relation $j(r,t)=[2t^{1/3}/(3b^{2})]w(br/t^{2/3})\,\partial_{r}z(r,t)$ is not a constant. Therefore, Fick’s law, even with a time-dependent diffusion constant foot2 , is violated (Fig. 3).

We also show the dependence of current profiles on the bond dimension $\chi$ used in simulations, Fig. 4. In the discrete-time case we use slightly smaller $\chi$, however the acquired times are larger (Fig. 2), as well as the sizes ($L=7200$ vs. $L=400$). As a net result the wall-times of discrete model simulations are about half as long as for a continuous one despite about a decade better accuracy (Fig. 4). We stress that in the best classical simulations (hard-point gas mendl14 ) slightly less than two decades of agreement with KPZ are achieved. What distinguishes our quantum simulations is that we directly work with an ensemble, encoded in the many-body density matrix $\rho(t)$, so no averaging is needed. It is an interesting open problem how to do such efficient ensemble simulations for classical many-body models, in particular since for continuous variables the local function spaces are infinitely dimensional.

Lastly, we note that taking a slightly larger domain-wall step $\mu=0.02$ we are even able to observe (data not shown) second order $\mu^{2}$ corrections to the dynamics in the form of a small ballistically spreading front, traveling away from the site of the quench.

Discussion.– We have shown that the infinite temperature spin-spin correlation function in the isotropic Heisenberg spin-$\frac{1}{2}$ model obeys the Kardar-Parisi-Zhang scaling. This is the first such observation in a deterministic quantum model. We stress that in order to reliably show the KPZ physics one has to look at the full distribution function of fluctuations and not e.g. just the dynamical scaling exponent being $z=\frac{3}{2}$. For instance, a related spreading exponent $\frac{1}{3}$ generically appears in free or dilute models, see e.g. racz04 ; vir18 .

High accuracy of over four decades was achieved by using a trick where we simulate the melting of a slightly polarized domain wall by directly evolving the density operator, which is, through linear response, equivalent to studying the equilibrium spin-spin correlation function.

Besides providing a method to efficiently probe spatio-temporal correlation functions in quantum models, several new directions are opened. The most important is the question of universality. Namely, in nonlinear fluctuating hydrodynamics, so-far verified only in classical models, the KPZ universality is associated to the existence of 3 conservation laws. It is not clear which 3 conserved quantities (if at all) are responsible for the observed behavior. By studying a kicked Floquet generalization of the isotropic Heisenberg model which does not conserve the energy, but nevertheless shows the KPZ physics, we show that the energy is not one of them. It remains to be seen if the observed behavior is in any way related to integrability and the SU(2) symmetry of the model.

Acknowledgements.– We acknowledge useful related discussions with J. De Nardis, E. Ilievski, M. Medenjak, and H. Spohn. The authors acknowledge support by the European Research Council (ERC) through the advanced grant 694544 – OMNES and the grants P1-0402 and J1-7279 of the Slovenian Research Agency (ARRS).