Free fermions at the edge of interacting systems
SSciPost Physics Submission
Free fermions at the edge of interacting systems
Jean-Marie St´ephan Univ Lyon, CNRS, Universit´e Claude Bernard Lyon 1, UMR5208, Institut CamilleJordan, F-69622 Villeurbanne, France* [email protected] 25, 2019
Abstract
We study the edge behavior of inhomogeneous one-dimensional quantum systems,such as Lieb-Liniger models in traps or spin chains in spatially varying fields. Forfree systems these fall into several universality classes, the most generic one beinggoverned by the Tracy-Widom distribution. We investigate in this paper the effectof interactions. Using semiclassical arguments, we show that since the densityvanishes to leading order, the strong interactions in the bulk are renormalized tozero at the edge, which simply explains the survival of Tracy-Widom scaling ingeneral. For integrable systems, it is possible to push this argument further, anddetermine exactly the remaining length scale which controls the variance of theedge distribution. This analytical prediction is checked numerically, with excellentagreement. We also study numerically the edge scaling at fronts generated byquantum quenches, which provide new universality classes awaiting theoreticalexplanation.
Contents a r X i v : . [ c ond - m a t . s t a t - m ec h ] A p r ciPost Physics SubmissionA Other universality classes 24 A.1 Tuning the dispersion relation 24A.2 Calogero-Sutherland models and β -matrix ensembles 24 References 25
The celebrated Tracy-Widom (T-W) distribution [1] was originally discovered while studyingthe largest eigenvalue of large random matrices. More precisely, it describes in this contextthe appropriately rescaled cumulative distribution function of the largest eigenvalue λ max N ofa random N by N gaussian hermitian matrix, in the limit N → ∞ : E ( s ) = lim N →∞ Proba (cid:32) λ max N − √ N − / N − / ≤ s (cid:33) . (1)The appearance of this distribution is not at all limited to random matrix theory. In fact,such a universal scaling occurs in edge problems as diverse as increasing subsequences of ran-dom permutations [2], growth models [3–5], dimer coverings on graphs [6], classical exclusionprocesses [7], or quantum quenches [8, 9], to name a few. In those problems, the T-W dis-tribution describes the edge properties of a macroscopic 2 d classical system at equilibrium,or the front of a 1 d system out of equilibrium. From a mathematical perspective T-W isbased on a determinantal point process (free fermions in physicist parlance), with a correla-tion kernel (propagator) known as the Airy kernel. While the diversity of problems where thisdistribution appears looks impressive, most of those are free fermions in disguise. A simplephysical picture was put forward in [3, 10] (see also [11] for an earlier related work). In such apicture the Airy kernel naturally emerges as a “filter” that projects onto the negative energyeigenstates of a free fermion model in a linear potential. Showing convergence to T-W inthose free problems then boils down to showing convergence of the correlation kernel to theAiry kernel, after appropriate edge rescaling.The aim of the present paper is to investigate several examples of physical 1 d interacting quantum mechanical models where the T-W distribution naturally appears in the groundstate. This will be done by combining heuristic semiclassical and thermodynamic BetheAnsatz arguments, supplemented by careful numerical checks. The main reason why this ispossible follows from a simple –but difficult to prove– renormalization argument: particles, sayin a trap, are typically diluted near the edge, so are less sensitive to the effects of interactionswhich might be otherwise very strong in the bulk. We will also investigate what happenswhen those interacting quantum systems are put out of equilibrium, which can lead to morecomplicated and much less understood universality classes.This long introduction is devoted to the free case, which helps put all the important con-cepts in place –once this is done treating interacting systems will prove no more complicated,since the edge will turn out to be free in the end. It is organized as follows. In section 1.1 weintroduce the free fermion model which has the Airy kernel as correlation kernel. We then2 ciPost Physics Submission present a derivation of the exact Fredholm determinant formula for the Tracy-Widom distri-bution (section 1.2), and briefly discuss various extensions. Finally, we explain on a simpleexample how T-W scaling occurs at the edge of a realistic fermion model (section 1.3). Themechanism for this is more important than the specific derivation, and follows from generalsemiclassical arguments. Let us stress that this introduction does not contain new results andfollows Ref. [10] to some extent; the only slight originality lies in the use of the language offield theory and Wick’s theorem. We consider the following second-quantized Hamiltonian on the real line H = (cid:90) R dx c † ( x ) (cid:18) − ∂ ∂x + x (cid:19) c ( x ) , (2)where the Dirac fields obey the anticommutation relations { c ( x ) , c † ( y ) } = δ ( x − y ), { c ( x ) , c ( y ) } =0 = { c † ( x ) , c † ( y ) } . This model is free, i.e. quadratic in the fermions operators, and can besolved exactly. Indeed, introducing the modes ψ † ( λ ) = (cid:90) R dx u ( λ, x ) c † ( x ) , ψ ( λ ) = (cid:16) ψ † ( λ ) (cid:17) † , (3)it is easy to show that [ H, ψ † ( λ )] = (cid:15) ( λ ) ψ † ( λ ), provided the single particle wave functions u ( λ, x ) satisfy the Schr¨odinger equation (cid:18) − ∂ ∂x + x (cid:19) u ( λ, x ) = (cid:15) ( λ ) u ( λ, x ) . (4)The solutions to this eigenvalue equation are well known to be Airy functions. Keeping onlythe eigenfunctions that decay to zero for | x | → ∞ : u ( λ, x ) = Ai( x + λ ) , Ai( x ) = (cid:90) R dq π e i ( qx + q / . (5)Those solutions are parametrized by a real number λ , the eigenenergies are given by (cid:15) ( λ ) = − λ . Hence the spectrum is continuous and unbounded. Due to the orthogonality relation (cid:82) R dx u ( λ, x ) u ( µ, x ) = δ ( λ − µ ), the modes satisfy the anticommutation relations { ψ ( λ ) , ψ † ( µ ) } = δ ( λ − µ ). The Hamiltonian becomes diagonal in terms of the new modes H = (cid:90) R dλ (cid:15) ( λ ) ψ † ( λ ) ψ ( λ ) , (cid:15) ( λ ) = − λ. (6)The ground state will play an important role in the following. It is a Dirac sea, obtained byfilling all the states with negative energies (corresponding to λ > (cid:104) ψ † ( λ ) ψ ( µ ) (cid:105) = δ ( λ − µ )Θ( µ ), where Θ is the Heavisidestep function. The propagator is given by G ( x, x (cid:48) ) = (cid:104) c † ( x ) c ( x (cid:48) ) (cid:105) (7)= (cid:90) ∞ dλ Ai( x + λ )Ai( x (cid:48) + λ ) . (8)3 ciPost Physics Submission This is known as the Airy kernel [1]. Of course, for free fermions problems the two point func-tion determines everything, more complicated observables reduce to determinants involvingthe propagator, by making use of Wick’s theorem [12]. The operator G Airy acting on functionsin L ( R ) as G Airy f ( x ) = (cid:82) R dyG ( x, y ) f ( y ) can be seen as a filter, that projects the function f ( x ) onto the subspace − d dx + x ≤ G Airy admits several generalizations, which we now briefly discuss. The first onecomes from introducing imaginary time operators c † ( x, τ ) = e τH c † ( x ) e − τH , with propagator G ( x, τ | x (cid:48) , τ (cid:48) ) = (cid:104) c † ( x, τ ) c ( x (cid:48) , τ (cid:48) ) (cid:105) (10)= (cid:90) ∞ dλ e − λ ( τ − τ (cid:48) ) Ai( x + λ )Ai( x (cid:48) + λ ) (11)for τ (cid:48) ≤ τ . This is known as the extended Airy kernel. The determinantal point process withcorrelation kernel G ( x, τ | x (cid:48) , τ (cid:48) ) is called the Airy process [3]. It is also possible to look atfinite temperature states, with averages taken as (cid:104) . (cid:105) β = Tr (cid:0) .e − βH (cid:1) , where the trace is takenover the underlying Fock space, and β is the inverse temperature. In that case the modeoccupation follows the Fermi-Dirac distribution, (cid:104) ψ † ( λ ) ψ ( µ ) (cid:105) β = δ ( λ − µ )1+ e − βµ , which leads to ageneralization (see e.g. [13–15]) that interpolates between the Airy kernel (zero temperature, β → ∞ ) and the Gumbel kernel (infinite temperature, β → Say we are interested in particle number fluctuations in an interval A = [ s, ∞ ) of R . Thenatural object to consider is the following generating functionΥ( α, s ) = (cid:68) e α (cid:82) ∞ s dx Q ( x ) (cid:69) , Q ( x ) = c † ( x ) c ( x ) , (12)which is known as full counting statistics [16] in condensed matter literature. A standardcomputation givesΥ( α, s ) = 1 + ∞ (cid:88) n =1 α n n ! (cid:90) ∞ s dx . . . (cid:90) ∞ s dx n (cid:104) Q ( x ) . . . Q ( x n ) (cid:105) (13)= 1 + ∞ (cid:88) n =1 ( e α − n n ! (cid:90) ∞ s dx . . . (cid:90) ∞ s dx n (cid:104) : Q ( x ) . . . Q ( x n ) : (cid:105) (14)= 1 + ∞ (cid:88) n =1 ( e α − n n ! (cid:90) ∞ s dx . . . (cid:90) ∞ s dx n det ≤ i,j ≤ n G ( x i , x j ) (15)= det s ( I + ( e α − G Airy ) . (16)Here : Q ( x ) . . . Q ( x n ): = c † ( x ) . . . c † ( x n ) c ( x ) . . . c ( x n ) denotes normal ordering, (14) followsfrom (13) by applying Wick’s theorem [12] and carefully rearranging the terms. In the last line,det s is the Fredholm determinant of an integral operator acting on L ([ s, ∞ )) with kernel ( e α − G ( x, y ). Eq. (15) can be taken as a definition of Eq. (16). The limit E ( s ) = lim α →−∞ Υ( α, s )4 ciPost Physics Submission is the probability that the interval A = [ s, ∞ ) contains no fermions. Obviously E ( −∞ ) = 0and E ( ∞ ) = 1. This emptiness formation probability (EFP) is given by the exact formula E ( s ) = det s ( I − G Airy ) . (17)The EFP is the cumulative distribution of what is known as the GUE Tracy-Widom distribu-tion , in particular it can be shown [1] to just equal (1). The corresponding probability densityfunction (pdf), p ( s ) = dE ( s ) ds , is illustrated in figure 1. It looks similar to a gaussian, but it is . . . . . − − − − s p ( s ) = dEds Figure 1: probability density function p ( s ) of the GUE Tracy-Widom distribution.slightly skewed; in fact, one can show that p ( s ) ∼ e − (4 / s / for s → ∞ and p ( s ) ∼ e − (1 / | s | for s → −∞ . It has mean (cid:104) s (cid:105) (cid:39) − . κ = (cid:104) ( s − (cid:104) s (cid:105) ) (cid:105) (cid:39) . κ / ( κ ) / (cid:39) . κ / ( κ ) (cid:39) . κ = (cid:104) ( s − (cid:104) s (cid:105) ) (cid:105) and κ = (cid:104) ( s − (cid:104) s (cid:105) ) (cid:105) − κ ) are thethird and fourth cumulants, respectively.It is not physically clear at this stage exactly of what this is a pdf in a realistic model. Toclarify this point, we discuss now a simple example where T-W emerges. The Hamiltonian (2) looks utterly unphysical at first sight: the potential is linear, and doesnot even confine particles to a given region of space. Another related complication lies in theDirac sea nature of the ground state, with infinite total particle number .The case of a harmonic potential is better behaved, and also of unquestionable experi-mental relevance, through its relation to the Tonks-Girardeau gas (see [17] for a review). Itturns out that the Airy Hamiltonian (2) describes the edge physics of the model in a harmonicpotential, through a mechanism that we discuss below. To be more concrete, we now considerthe Hamiltonian H = (cid:90) R dx c † ( x ) (cid:18) − ∂ ∂x + x − µ (cid:19) c ( x ) . (18)The parameter µ is a chemical potential, which allows to control the number of particles in theground state. This model can be solved in a similar way as the previous one, and the singleparticle wave functions may be expressed in terms of Hermite polynomials. The problem is, in One can show [11] that the particle number in the interval [ − a, ∞ ) diverges as π a / for a → ∞ . ciPost Physics Submission fact, formally identical to the well known quantum harmonic oscillator. Due to the confiningnature of the potential, the energy levels are now discrete. Using this approach, one can forexample show that the density of fermions in the ground state follows, when µ → ∞ , thecelebrated Wigner semi-circle law (cid:104) Q ( x ) (cid:105) = 1 π (cid:112) µ − x . (19) Bulk LDA—
It is enlightening to look at this problem using semiclassical analysis, some-times also known as local density approximation (LDA) in cold atom literature. The key as-sumption is separation of scales: we look at mesoscopic scales around some point x , namelywe look in an interval [ x − δx, x + δx ], where δx is much bigger that the mean distance be-tween particles, and much smaller than the system size (both to be determined at this stage).On such distances the system looks homogeneous, with a well defined effective chemical po-tential µ eff ( x ) = µ − x . The ground state propagator becomes the kernel of the projectiononto − d d ( δx ) ≤ µ − x , (20)which is easy to determine. Indeed, thinking in Fourier space, the above becomes k < k ,where k = (cid:112) µ − x , which defines a disk in phase space ( x, k ). Hence the desired projectoris given by (cid:104) c † ( x + δx ) c ( x + δy ) (cid:105) = (cid:90) k − k dk π e ik ( δx − δy ) (21)= sin k ( δx − δy ) π ( δx − δy ) , (22)consistent with the claimed density (19). The particle number is then determined self consis-tently as N = (cid:82) (cid:104) Q ( x ) (cid:105) = µ/
2, so the limit µ → ∞ , where LDA is expected to become exact,is the thermodynamic limit N → ∞ in the usual sense. The effective size of the system isthen L ∼ (2 N ) / , while the mean interparticle distance is of order a ∼ N − / . The result(22) is therefore valid in the limit N → ∞ , a (cid:28) δx, δy (cid:28) L , and k > Edge from LDA—
The behavior close to the edge is slightly more complicated, but canstill be obtained from semi-classical analysis (see Refs. [10,18] for discussions). To explore thisregime, we make the change of variable x = √ N + ˜ x , where the new variable is just assumedto be much smaller than system size for now. The propagator close to the edge becomes thekernel of the projection − d d ˜ x + (cid:16) √ N + ˜ x (cid:17) ≤ N. (23)Expanding the square, the term in ˜ x is subleading compared to 2 √ N ˜ x , so may be discarded.After a final change of variable ˜ x = (cid:96)u , (cid:96) = (8 N ) − / (24)the previous equation becomes − d du + u ≤ , (25)6 ciPost Physics Submission whose kernel is precisely the Airy kernel studied in section 1.1, see (9). Back in ˜ x coordinatesystem, this behavior occurs at scales of order (cid:96) = (8 N ) − / (cid:29) N − / , so does not contradictthe bulk LDA argument, even though we are in a different regime with lower density now .In this sense the limit is smooth, and LDA/semiclassics correctly predicts the edge behavioras a limit, since the result (25) can be proved by other means [1]. Semiclassically in phasespace, we go from a disk k + x ≤ µ for the bulk to a parabolic region q + u ≤ k, q are the momenta corresponding to x, u respectively. The scale (cid:96) — It is important to realize that the above derivation does not rely on the factthat the potential be harmonic. For any reasonably behaves potential, we expect the samescaling behavior, since any smooth potential can be linearized close to the edge, and this willprove the dominant contribution. For this reason, Airy scaling close to the edge is expectedto be quite generic. To the leading order, the only free parameter in such a mechanism isprecisely the scale (cid:96) we calculated above in a particular example.In terms of the bulk variables (recall a is interparticle distance, and L total system size)we have (cid:96)a ∝ (cid:0) La (cid:1) / , and for this reason the exponent 1 / Tracy-Widom—
T-W appears when looking at the distribution of the rightmost particle.It may be determined by looking at the emptiness formation probability, which is givenfor finite N by the Fredholm determinant E N ( x ) = det x ( I − K N ), where K N is the kernelassociated to the ground state propagator for the harmonic trap . As we have just established,this kernel scales to the Airy kernel in a suitable edge limit, which means the (suitably rescaled)distribution of the rightmost fermion, dE N /dx converges to the Tracy-Widom distribution.Physically, the scale (cid:96) also controls the standard deviation of the distribution of the rightmostparticle, which is given by (cid:112) (cid:104) x (cid:105) c = a(cid:96) + o ( (cid:96) ), where a (cid:39) . Relation to GUE—
The free fermion problem looked at in the previous subsection is infact formally identical to the random matrix problem where the T-W was originally discovered.Indeed, denote by | Ψ (cid:105) the N − particle ground state of (18). In first quantization language,the many-body wave function reads φ ( x , . . . , x N ) = (cid:104) c ( x ) . . . c ( x N ) | Ψ (cid:105) which is given by aSlater determinant. A direct calculation using properties of the Hermite polynomial and theVandermonde identity shows | φ ( x , . . . , x N ) | ∝ (cid:89) ≤ i The remainder of the paper is devotedto the effect of interactions. We study in section 2 interacting models in traps at equilibrium,which can be seen as generalizations of (18), and demonstrate that T-W scaling genericallysurvives at the edge (specific exceptions are discussed in appendix A). Section 3 tackles a morecomplicated quantum out of equilibrium problem, where the effects of interactions are subtle.In particular, we establish that the edge distribution has very long tail, in stark contrast withT-W. We conclude in section 4 and discuss some open problems. As already mentioned an obvious question, left unanswered in the introduction, lies in theeffect of interactions. We have discussed an explicit example, that has the Airy Hamiltonian(2) as effective edge Hamiltonian, but the free fermion structure was already built in, whichmeans the distribution of the last particle could always be expressed as a Fredholm deter-minant. Showing convergence to the T-W distribution, ignoring mathematical difficulties,amounts to showing convergence of the propagator to the Airy kernel.On the other hand, T-W is widely believed to be a universal distribution, and should alsoappear in problems where the free fermions structure is not already present in the microscopicmodel. For example, T-W scaling has been proved in the asymmetric exclusion processes(ASEP) for certain initial conditions, see e.g. [7]. The ASEP is related to the integrable XXZspin chain, but away from the free fermion point. Despite these notable exceptions, there arein general still very few rigorous or exact results in this direction.In the class of problems we look at, there is a simple argument explaining why T-W shouldappear at the edge of an interacting system (say) in a trap. Even if the underlying modelmay be extremely complicated, the edge is precisely the region where the density of particlesbecomes very low. In this region the quantum particles are diluted, and interactions withsufficiently fast decay, which might be very strong in the bulk, are expected to become weakerand weaker. Hence the particle become effectively free, and this makes generic T-W scalingbehavior quite plausible. This mechanism is not much different from the usual appearance ofa simpler effective field theory to describe the scaling limit of possibly extremely complicatedmicroscopic models.We discuss in this section an example where we are able to demonstrate this, and also,perhaps more importantly, are able to compute analytically the associated scale on which suchbehavior occurs. We do this using a combination of simple analytical arguments, backed byextensive numerical checks. Before doing that let us emphasize that the word generic in theprevious paragraph is important. In fact, two clear exceptions will be discussed in subsectionsA.1 and A.2. 8 ciPost Physics Submission The first example we look at is the Lieb-Liniger model in a harmonic trap, governed by thesecond-quantized Hamiltonian H = (cid:90) dx (cid:18) Ψ † ( x ) (cid:20) − (cid:126) m ∂ ∂x − µ + V ( x ) (cid:21) Ψ( x ) + g † ( x )Ψ ( x ) (cid:19) , (27)with g > x ) , Ψ † ( y )] = δ ( x − y ), [Ψ( x ) , Ψ( y )] = 0 = [Ψ † ( x ) , Ψ † ( y )]. This model is well-known to beintegrable in the absence of a trapping potential [24, 25]. The trap, however, typically breaksintegrability. In the following, we will consider the (integrability breaking) harmonic trapwhich is the most natural and experimentally relevant.Before proceeding any further, let us mention that this model is a natural generalizationof the Fermi gas looked at in section 1.3, in the following sense. In the limit of infinitelystrong repulsion, g → ∞ , the Tonks-Girardeau limit, the first quantized ground state bosonicwave function is given by [17] φ b ( x , . . . , x N ) = φ ( x , . . . , x N ) (cid:89) ≤ i Despite the fact that the system is not integrable, it is still possibleto rely on separation of scales. As before, we assume that the system is sufficiently uniformon mesoscopic scales, which means it looks, locally, identical to the ground state of the Lieb-Liniger model without an external potential. This observation allows us to use the groundstate thermodynamic properties of this Bethe-Ansatz integrable model. The thermodynamicBethe Ansatz (TBA) description of homogeneous ground state is well known, see e.g. [25],and has been used to predict density profiles [26] and more complicated correlation functions[27, 28] in the ground state.In the following we work in units where (cid:126) = m = 1. For a given chemical potential µ , theground state is parametrized by a set of rapidities, that satisfy Bethe equations [24]. In thethermodynamic limit the relevant quantity is the density of rapidities ρ ( k, µ ), which can beshown to satisfy the linear integral equation (LIE) ρ ( k, µ ) − π (cid:90) k F ( µ ) − k F ( µ ) K ( k, q ) ρ ( q, µ ) dq = 12 π , (29)with kernel K ( k, q ) = 2 g ( k − q ) + g . (30)Of great importance is also the energy of single particle excitations with quasimomentum k above the ground state. It can be shown to satisfy another LIE ε ( k, µ ) − π (cid:90) k F ( µ ) − k F ( µ ) K ( k, q ) ε ( q, µ ) dq = k − µ. (31) The first quantized form is H = − (cid:80) Nj =1 − (cid:126) m ∂∂x j − µ + V ( x j ) + g (cid:80) ≤ i From the previous argument, we have determined that the edge is simplylocated at x = ±√ µ , even though the full density profile can only be accessed in implicitform. The ground state is characterized in phase space by ε ( k, µ eff ( x )) ≤ , µ eff ( x ) = µ − x (33)where ε ( k, µ ) is given by (31). Now comes the following simple but crucial point: at theedge the contribution from the integral in (31) vanishes to the leading order, so introducing x = ˜ x + √ µ as before, we are left with the simple projection k + 2 √ µ ˜ x ≤ g/ √ µ , should notprevent the appearance of T-W scaling at the edge. We expect the subleading corrections dueto dressing to be no greater than those occuring already at the free fermion point. As before,the rightmost particle will be delocalized on scales of order (cid:96) = (4 µ ) − / . Since the densityprofile close to the edge follows from the same argument, the T-W scaling is tightly related,from a more pedestrian perspective, by the behavior of the density close to the edge, which isthe square-root scaling ∼ x →±√ µ π (cid:112) µ − x , which turns out to be independent on interactionshere. As a simple consequence, systems where the bulk density does not vanish as square-rootare unlikely to yield T-W scaling.It is also possible to interpret this result using field-theoretical language. An importantproperty of interacting inhomogeneous systems in the Luttinger liquid universality class is thatthe Luttinger parameter, which parametrizes the strength of interactions, depends on positionin the bulk [27, 28]. In such systems, the edge is precisely the place where it evaluates to one,the free fermion value (for inhomogeneous free fermions K = 1 throughout the system [29]).This argument should apply whenever the interaction between particles decays sufficientlyfast. To illustrate this last point we discuss in appendix A.2 an example with inverse squarelong-range interactions, for which the Luttinger parameter can take other values at the edge.10 ciPost Physics Submission We study here another similar but more general example, this time of discrete nature. TheHamiltonian we consider is that of the spin-1/2 XXZ chain on the infinite lattice H = (cid:88) x ∈ Z +1 / (cid:0) S x x S x x +1 + S y x S y x +1 + ∆( S z x S z x +1 − / − h ( x/R ) S z x (cid:1) , (35)where S αj = σ αj , and σ αj act as Pauli matrices on the j ’s copy of C and as identity on theothers (we take the Hilbert space ( C ) ⊗ L and implicitly assume L → ∞ ). Similar problemswith spatially varying magnetic fields have been considered in the literature [29–31]. Themagnetic field term depends on position, and plays a similar role as the trapping potentialbefore. Before investigating this, let us summarize known results in the case of a constantmagnetic field h . As is well known, the ground state has critical correlations for | h | < | h | > h ( u ), that also, for later convenience, satisfies h ( − u ) = h ( u ) and lim u →∞ h ( u ) = ∞ .The large parameter R in (35) defines an effective system size, set by the location where | h ( x/R ) | = 1 + ∆. Defining x e = Rh − (1 + ∆), inside the region [ − x e , x e ] the system isinhomogeneous with critical correlations, outside it is a fully polarized product state. Bulk and edge TBA— The TBA description of the ground state is also well known [25],and has a similar structure as the Lieb-Liniger one. It has also been checked numerically inRef. [31], on the example h ( u ) = u , that the LDA approach gives the correct density profiles.With this at hand it is straightforward to look at the edge behavior, the calculations areexactly the same as in the previous subsection. With x = x e + ˜ x , we find the edge behaviorin rapidity space k h (cid:48) (cid:16) x e R (cid:17) ˜ xR ≤ . (36)Assuming as before the emergence of Wick’s theorem at the edge means we get T-W scaling.Introducing the new scale (cid:96) ∆ = (cid:18) R h (cid:48) ( h − (1 + ∆)) (cid:19) / (37)and making the change of variables ˜ x = (cid:96) ∆ u , we recover the projector onto − d du + u ≤ 0. Thescale (cid:96) ∆ controls, as (cid:96) before, the standard deviation of the distribution of the last particle. Itis now of order R / , and depends explicitly on the interaction parameter ∆. This predictionis tested numerically in the next subsection. The analytical argument presented in the previous subsection is quite heuristic. Indeed, weassumed free fermion behavior at the edge, and determined the propagator (correlation kernel)by using a self-consistent TBA description. This makes a numerical confirmation necessary.Let us first note that numerical checks of Tracy-Widom scaling are notoriously difficult(see e.g. [32]). Since the associated scale is usually a power one third of the system size11 ciPost Physics Submission convergence is slow, even when reaching apparently very large system sizes. In classicalsetups Monte Carlo techniques are able to simulate large enough systems, however error barstend to blur the results, especially when trying to extrapolate the data. The situation inthe spin chain, we argue here, is slightly more favorable, which is one of the motivations forinvestigating T-W scaling in this quantum system. While the Hilbert space size naively growsexponentially fast, powerful variational techniques such as DMRG [33] are able to find theground state with very good accuracy for large enough R . Efficient DMRG libraries ableto implement continuous symmetries are now available in several programming languages(including Python [34] and C++ [35]), which simplifies our task considerably in the XXZ spinchain. The simulations shown below were performed using the C++ ITensor library [35].For the magnetic field we made the choice h ( u ) = u + au | u | , which satisfies the hypothesisexplained in the previous subsections. The term proportional to u | u | might seem artificial,however, its presence ensures that the length scale (37) associated to T-W depends on ∆ (forthe linear potential (cid:96) ∆ = ( R/ / unfortunately does not depend on ∆), and makes for astronger numerical test of our analytical argument. Ground state density profile— Let us first discuss the ground state magnetization profile (cid:104) S z x (cid:105) , which is shown in figure 2 for several values of ∆. The case ∆ < − > − 1. With our choice of magnetic field,an explicit computation solving a quadratic equation shows that the edge is located at x e = ± R (cid:112) a (1 + ∆) − a , (38)a prediction in very good agreement with numerics (note again the density profile for a = 0has already been checked in Ref. [31]). The whole profile is also invariant under reflectionsymmetry x → − x conjugated with up-down (particle-hole after a Jordan-Wigner transfor-mation) symmetry, due to the antisymmetry of the magnetic field we chose in (35). Suchprofiles are also related to equilibrium shapes of crystals in 2d, and have been investigatedmuch earlier in this context [36].We note in passing that another region develops in the middle of the chain for ∆ > > h , the homogeneous ground state is gapped with antiferromagnetic order.For | h | > ∆ − x e in (38) is not affected by this phenonenon, however, and this is what we focus on in thefollowing. Edge distribution— We now come to the actual check of our conjecture, which predictsT-W scaling with associated scale (cid:96) ∆ = (cid:32) R (cid:112) a (1 + ∆) (cid:33) / . (39)Accessing the edge distribution can be done in a straightforward way in DMRG. We studythe discrete analog of the “emptiness” formation probability, the probability that all spins at12 ciPost Physics Submission . . . . − − S z x + / x/R ∆ = − . 4∆ = 0 . 4∆ = 1∆ = 1 . edge x e Figure 2: Numerical density profile for R = 64, and a = 0 . 1, obtained using DMRG. Thedata is shown for four different values of ∆. In practice we use a system of total size L = 512,which is significantly larger than the effective size of the system 2 x e , outside of which the wavefunction is fully polarized. The central region with antiferromagnetic order is a specificity of∆ > 1, as mentionned in the text. In the following we are interested in the behavior at theedge x e , indicated with green arrows.position j ≥ x be up, E x = (cid:42) ∞ (cid:89) j = x (cid:18) S z j (cid:19)(cid:43) , (40)close the right edge , see figure 2. The discrete probability density function (dpdf) is thenreconstructed as p x = E x +1 − E x , and expected to converge to Tracy-Widom after properrescaling involving (cid:96) ∆ . This is shown in figure. 3 (left). As can be seen the agreementis excellent. Note however a slight shift along the horizontal axis. We interpret this as asubleading order one additive correction to the (cid:96) ∆ scaling, and checked that this is indeeda finite-size effect (not shown). The variance predicted by our analytical argument is alsoclearly confirmed by a finite-size scaling analysis, shown in figure. 3 (right). In this figure aswell as in later plots, the leading correction is expected to be of order R − / , and correspondsto the terms in ˜ x that were discarded around Eq. (23) or (36).To study more quantitatively the convergence to T-W, we also performed a finite-sizescaling analysis of the skewness and excess kurtosis, related to the third and fourth cumulants(for a gaussian all cumulants of order larger than two are zero). This is shown in figure 4,with very convincing agreement. Relative errors for the largest system sizes we could accessare typically 5% or less, depending on the value of ∆. After extrapolation this error falls The left edge analog would be the probability (cid:68)(cid:81) xj = −∞ (cid:16) − S z j (cid:17)(cid:69) that all spins are down, and leads tothe same results by symmetry. ciPost Physics Submission . . . . . − − − − ‘ ∆ p x ( x − x e ) /‘ ∆ TW ∆ = 0 . , R = 256∆ = 1 , R = 256∆ = 1 . , R = 256 0 . . . . . 95 0 0 . 02 0 . 04 0 . 06 0 . 08 0 . R − / ∆ = − . 4∆ = 0∆ = 0 . 4∆ = 1∆ = 1 . Figure 3: Left: rescaled dpdf for R = 256 and comparison with T-W (red line). Right:variance for a = 1 / 10 divided by the variance for a = 0, plotted as a function of R − / forseveral values of ∆. The analytical prediction, given by (1 + 4 a (1 + ∆)) − / is shown forcomparison in thick horizontal lines. The data is extrapolated by a straight line, shown indashed as a guide to the eye, with perfect agreement. The total chain length used for allsimulations is L = 8 R .well under a percent in all cases, which is remarkable given the numerical difficulties usuallyassociated to testing T-W. Of course, it is also possible to check higher order cumulants. . . . . . . . . . 18 0 0 . 01 0 . 02 0 . 03 0 . 04 0 . 05 0 . 06 0 . R − / ∆ = − . 4∆ = 0∆ = 0 . 4∆ = 1∆ = 1 . TW . . . . . . . . . 18 0 0 . 01 0 . 02 0 . 03 0 . 04 0 . 05 0 . 06 0 . R − / ∆ = − . 4∆ = 0∆ = 0 . 4∆ = 1∆ = 1 . TW Figure 4: Left: skewness divided by the T-W one ( (cid:39) . (cid:39) . R − / . As in the previousfigure, extrapolation to a straight line shows nearly perfect agreement.However, those probe finer and finer details of the distribution, which would not be visibleto the eye e.g. in figure 3. Since excess kurtosis shows larger errors than skewness, it isreasonable to expect finite-size effects to increase for higher order cumulants.Let us finally mention that it is possible to use the spin chain to simulate the Lieb-Linigermodel. This is done by considering the potential h ( u ) = u , and taking an appropriatelow density limit (see [38, 39]). In that case simulations are typically limited to less than ahundred particles, we also checked that for reasonable interactions strength the skewness iswithin ≤ 10% of T-W, with agreement improving for larger particle numbers.14 ciPost Physics Submission We have argued in previous sections that interactions renormalize to zero close to the edge.This implies the emergence of the fermionic Wick theorem, a key ingredient to get T-Wscaling. It is possible to check the fermionic Wick factorization property more explicitly,for example by looking at the entanglement entropy S ( x ) of an interval [ x, ∞ ) for x closeto the (right) edge. For generic interacting systems computing this exactly is extremelycomplicated, however, for free (Airy or not) fermions it may be simply determined from thepropagator [8, 40, 41], which leads us to conjecture that the formula S ( x e + (cid:96) ∆ s ) = − Tr s [ G Airy log G Airy + ( I − G Airy ) log( I − G Airy )] , (41)holds for any ∆ > − R → ∞ (the free case ∆ = 0 was previously established inRef. [42]). Here Tr s denotes trace on L [ s, ∞ ). Note once again that only the scale (cid:96) ∆ entersin the final result. Data for the rescaled entropy is shown in figure. 5. As can be seen theagreement is excellent and improves as R is increased. We observe slight deviations when s . . . . . . . . . − − − s TW R = 64 , ∆ = 0 . R = 256 , ∆ = 0 . Figure 5: Rescaled entanglement entropy close to the right edge, which is expected to convergeto (41) in the limit R → ∞ .becomes large negative. This is expected since the entropy still sees bulk effects for finite-size.We note that the bulk entanglement entropy is more complicated in inhomogeneous systems,with even the free case turning out to be nontrivial [43]. For interacting systems the fact thatthe Luttinger parameter depends on position makes a field-theoretic treatment more difficult(see [28] for a discussion for local operators). 15 ciPost Physics Submission We investigate in this section a different but related out-of-equilibrium setup, which showsinteresting edge behavior. We consider the infinite XXZ spin-1/2 model (35) in the absenceof a magnetic field. The system is initially prepared in the domain-wall state | Ψ (cid:105) = | . . . ↑↑↑↑↓↓↓↓ . . . (cid:105) , (42)and let evolve unitarily with the aforementioned Hamiltonian H (the wave function at time t isgiven by | Ψ( t ) (cid:105) = e − iHt | Ψ (cid:105) ). At long times, a non trivial magnetization profile develops, withdynamical edges that we wish to study. As we shall see, away from free fermions (∆ = 0) thiswill provide an example of a new universality class, beyond what is presently known. Beforeentering into specifics, let us remind once again that T-W is not the only known universalityclass even in equilibrium problems, even though it is probably the most frequent/natural. Toillustrate this, we discuss in appendix A two known exceptions to the scenario put forward inthe previous section. The example discussed here will be of different nature, however.The section is organized as follows. Several works have studied the spread of correlationsafter this quantum quench [44–55], we summarize the aspects that we need in section 3.1.Section 3.2 deals with previous results and claims for the edge behavior. We come in section3.3 to our new results. In particular, we numerically access the real-space distribution of therightmost up spin, the exact analog of what gave T-W in the previous section, or gives T-Wfor free fermions here. We show that this distribution is very delocalized, compared to otherclasses, and discuss in depth some of its properties. Finally, we summarize our findings insection 3.4. Despite the integrability of the XXZ chain and apparent simplicity of the quench protocol,exact computations of simple observables at finite time are extremely challenging, with onlythe return probability known in closed form [53]. A (generalized) hydrodynamic (GHD)description, able to tackle general such protocols, was put forward in Ref. [56, 57]. Thisapproach is expected to become exact for our quench in the limit x → ∞ , t → ∞ , x/t fixed,provided | ∆ | < 1, which we will assume from now on. It was used in Ref. [52] to compute thedensity profile analytically in that limit.For the convenience of the reader, numerical examples of such density profiles are shownin figure 6 for several values of ∆, and compared to the exact solution. The DMRG timeevolution is implemented using the method of [58] together with the higher order Trotterformulas of [59]. The GHD limit for this quench is quite peculiar, and the density profilein the bulk region turns out to be nowhere continuous as a function of ∆. This surprisingbehavior, reminiscent of Drude weight results [52, 60–66], which are believed to have also thisproperty, is ultimately related to the quantum group structure [67] underlying the XXZ chainat root of unity.We name the position x e where the GHD density profile vanishes the GHD edge . It isgiven by the simple formula x e /t = √ − ∆ [52, 53]. There can be subleading correctionsto this behavior. In fact, closer inspection of figure 6 (see in particular the inset) shows thatdensity decays slowly for x > x e before it hits another edge at x f = t [48]. For x > x f thedecay of the density appears to be super-exponential. Since the speed corresponding to x f canbe interpreted as the group velocity v f = 1 of a single magnon in a ferromagnetic background16 ciPost Physics Submission . . . . . . . . / h S z x ( t ) i + / x/t ∆ = 0∆ = 1 / 2∆ = √ / 2∆ = 0 . GHD free edge x f = tx f GHD edge x e = t √ − ∆ x e Figure 6: Numerical density profiles for x > 0, shown as a function of x/t at time t = 240. Thedata is compared to the GHD formula (cid:104) S x ( t ) (cid:105) = − q π arcsin (cid:16) sin πq sin γ xt (cid:17) for γ/π = p/q irrational, − x t sin γ otherwise. Inset: zoom on the region [ x e , x f ] = [ t √ − ∆ , t ] between the GHD edgeand free edges, where the GHD prediction vanishes but subleading corrections remain. For x > t , density vanishes super-exponentially fast for all values of ∆.and does not depend on interactions, we dub x f = t the free edge . It can also be seen as aLieb-Robinson-type bound in such a system. The fact that the GHD and free edge do notcoincide away from ∆ = 0 will play an important role in the following. T-W scaling for the edge front was established [8] at the free fermion point ∆ = 0 by anexact computation. However, such a scaling does not survive at the edge for ∆ (cid:54) = 0, aswas argued in Ref. [52], the simplest reason being the fact that the density profile is linearat the GHD edge, not square-root as in all the examples discussed in the present paper,see e.g. (19). Such a linear behavior was also observed numerically in more complicatedout-of-equilibrium setups [31]. An associated toy-model kernel [52], expected to qualitativelydescribe the GHD edge, was obtained from the exact computation of the density and currentprofiles. The calculation of those was formally identical to a different free fermion problemstudied in Ref. [48, 49], with time-dependent propagator C ( x, y | t ) = (cid:90) π − π dk π (cid:90) π − π dq π e it (cos k − cos q )+ ixk − iyq f ( k, q ) , (43)where f ( k, q ) = χ γ ( k, q )1 − e i ( k − q + i + ) + reg( k, q ) . (44)17 ciPost Physics Submission Here cos γ = ∆, and χ γ ( k, q ) = 1 if | k | , | q | ∈ { , γ } ∪ { π − γ, π } and zero otherwise. reg( k, q )denotes a function that is regular at k = q , but can have pointwise singularities, see [49] forexplicit expressions. The asymptotics are then studied using standard saddle point techniques,where the singular term dominates. The case γ = π/ x e = t is derived by cubic expansion around k, q = π/ γ (cid:54) = π/ k, q = γ yields C ( x, y | t ) = πγ √ t cos γ E ( x − t sin γ √ t cos γ , y − t sin γ √ t cos γ ), where E ( X, Y ) is the imaginary error kernel [52] E ( X, Y ) = (cid:90) ∞ dλE ( X + λ ) E ∗ ( Y + λ ) , E ( X ) = (cid:90) ∞ dQ π e iQX + iQ / . (45)The scaling behavior close to the edge is t / , instead of t / . In our language, this canbe naturally interpreted as the kernel of the projection − i ddX + X ≤ 0, consistent with thelinear behavior for the density profile and the edge free fermions assumption. This analytical . . . . . . . . . − − E ( X ) X t = 60 t = 120 t = 240 toy-model . . . . . . . . . − − E ( X ) X t = 60 t = 120 t = 240 toy-model Figure 7: Rescaled density profiles around the GHD edge x e = t √ − ∆ for increasing times t = 60 , , / / √ 2. In that case theagreement is only fair, the difference with the free fermion kernel does not seem to go awayin the limit t → ∞ .result in the toy-model is compared to numerical simulations in figure 7. As can be seen theagreement is decent for ∆ = 1 / 2, but gets worse for larger values of ∆, which means it isprobably not exact. The collapse as a function of √ t seems quite good, however, sufficient toconfidently exclude t / .What about the free edge, around which density is small but non-zero [48]? It was recentlystudied numerically in [55], where t / scaling close to x f was observed. The fact that a smallfraction of quasiparticles go faster than the TBA/GHD speed was interpreted as a consequenceof a slight order one excess in energy, due to the fact that (cid:104) Ψ | H | Ψ (cid:105) = − / 2, where GHDimplicitly assumes (cid:104) Ψ | H | Ψ (cid:105) = 0.We want to stress here that the observations made in Refs. [52, 55] are not incompatible,provided the results of [55] are interpreted carefully. First, the t / scaling is, in fact, alsopresent in the free fermion propagator (43). Indeed, the result (45) was obtained from (44) byneglecting the regular terms, which provide subleading contributions. However, in the region x e t < x < x f t this is not true anymore, since the indicator function χ γ in (44) vanishes inthe corresponding region of phase space. Close to x = x f = t the dominating saddle point is18 ciPost Physics Submission located at k, q = π/ 2, and yields a (subleading) product of two Airy functions, but not theAiry kernel. For x/t > 1, all correlations decay super-exponentially fast to zero.From the previous considerations, it is not clear how the distribution of the rightmostparticle would exactly look like, except for the fact that it should differ from T-W. This is thepurpose of the next subsection, where we study it numerically for the first time, and pointout an important analytical subtlety. As our previous analysis suggests, the t / contribution close to the second edge should onlyaccount for a small fraction of one real-space particle, since it is subleading compared to theAiry kernel (recall T-W accounts for exactly one particle). This means the distribution of thelast particle, the true analog of the Tracy-Widom distribution in our quench, should still bedominated by other effects, including diffusive effects in the neighborhood of the GHD edge x e . This can be checked by once again computing the EFP, and numerically reconstructingthe corresponding dpdf. The results are presented in figure 8 and show that the distributionis peaked around x e . The free edge x f is then simply the termination of the right tail of thedistribution. . . . . . − − − − . . x e x f Xt = 60 t = 120 t = 240 x/t t = 60 t = 120 t = 240 . . . . . − − − − . . x e x f Xt = 60 t = 120 t = 240 x/t t = 60 t = 120 t = 240 Figure 8: Rescaled distribution of the rightmost particle (rightmost up spin). As before theabscissa is X = x − t sin γ √ t cos γ , and data is shown for ∆ = 1 / / √ x e = t sin γ for the times we could access, even though its long right tail goes all theway to x free = t . The rescaled density profile for t = 240 is also shown in orange solid line forcomparison (it is appropriately normalized to allow for comparison with the distribution). Tohelp visualize the location of both edges (red x e and blue x f bullets), the same distribution isshown in the inset as a function of x/t .While the collapse as a function of √ t near the GHD edge seems fair, it is unlikely that thisfully describes the distribution of the last particle, due to the following analytical argument.Discarding the fact that the toy-model kernel (45) is unlikely to be exact for our quench, wefind that it behaves for large X, Y as E ( X, Y ) ∼ log X − log Y π ( X − Y ) , X (cid:54) = Y, (46) E ( X, X ) ∼ π X , (47)19 ciPost Physics Submission which is, importantly, not integrable for X → ∞ . This problem has to be cured by hand,introducing a hard cutoff at x = t to make the density profile consistent with Lieb-Robinsonbounds, but this would still mean that the figure above does not represent a true scalingfunction for the rescaled pdf. This suggests the possibility for logarithmic corrections infigure 8, which are hard to prove or disprove numerically.These corrections should affect transport properties also; for example the particle number N dil in the diluted region [ x e , ∞ ) was claimed to be of order one in Ref. [53], but if the truekernel decays as inverse distance as E does, then particle number should diverge logarithmi-cally with time. As shown in figure. 9, this looks plausible numerically, back in the interactingquench. . . . . . . . . . . t ∆ = 0∆ = 1 / 2∆ = 0 . 6∆ = 1 / √ 2∆ = 0 . Figure 9: Particle number N dil in the diluted region [ x e , ∞ ) on the right of the GHD edge, asa function of time (on a logarithmic scale). Data is shown for several values of ∆. The resultsare consistent with a slow logarithmic divergence of this particle number, except at the freefermion point, ∆ = 0, where it saturates very quickly.Pushing the numerics further than done here is unfortunately unlikely to pay huge divi-dends. Indeed, we observe that convergence is quite slow in general, worse than regular T-Wscaling encountered in this paper. In addition to the effects already mentionned, there areother competing terms, that are already present in the (probably simplified) free fermionmodel (43). In fact, we also checked that numerical convergence to the kernel (45) is alreadyvery slow in a discrete free fermion system modeling (45), even considering the very largetimes ( t > Let us summarize our main numerical observations for ∆ (cid:54) = 0. For most values of | ∆ | < 1, andall accessible times most of the probability distribution is concentrated near the GHD edge x e = t √ − ∆ . The distribution has an extremely long right tail, however, which extends allthe way to x f = t . In stark contrast, the free fermion T-W distribution is concentrated ona much smaller region of width t / near the free edge (which coincides with the GHD edge,since ∆ = 0 in that case).Motivated by the toy-model kernel of Eq. (45), which predicts a total particle number N dil ∝ (cid:82) x f x e + (cid:15) dxx − x e ∝ log( t/(cid:15) ) in the diluted region, we have observed numerically that particle20 ciPost Physics Submission number in the diluted region [ x e , ∞ ) grows with time, with a behavior consistent with alogarithmic divergence. This does suggest that the distributions shown in figure 8 might befar from converged, and might look different when the particle number becomes greater thanone (for | ∆ | ≤ . t = 10 , a time whichquickly increases as ∆ is decreased). We expect the distribution to shift to the right, possiblyeven move away from the GHD edge at extremely large times. Since (cid:82) x f x (cid:48) e dxx − x e does not divergeprovided √ − ∆ < x (cid:48) e /t < 1, we still expect the distribution to be at least supported on aninterval [ x (cid:48) e , x f ], with an extremely long right tail. In all cases the free edge will correspondto the termination of the tail, which means the distribution is very different from T-W.To further illustrate this last point, we have computed numerically the variance and skew-ness of the distribution, as a function of time. The results are shown in figure. 10. The . . . . . h x i / t t ∆ = 0∆ = 1 / 2∆ = 0 . 6∆ = 1 / √ 2∆ = 0 . . . . . . . . . s k e w n e ss t ∆ = 0∆ = 1 / 2∆ = 0 . 6∆ = 1 / √ 2∆ = 0 . T-W Figure 10: Left: variance divided by t , as a function of time t . As can be seen, the distributionis much more delocalized than expected from regular diffusion. The free fermion point is alsoshown for comparison, in that case we expect a decay as t − / . Right: skewness as a functionof time (on a logarithmic scale). It appears to grow slowly with time, possibly logarithmically.The skewness of T-W, which is much smaller, is also shown for comparison, and matches verywell the free fermion calculation.variance grows possibly as fast as t (or slightly slower), while skewness keeps on increasing:we find once again a behavior consistent with a logarithmic divergence, very different fromthe T-W finite value which is about 0 . local interactions.For similar reasons, it is not completely clear whether the final answer for correlations at theGHD will be free fermionic or not. Correlations near the free edge should be, however. In this paper, we have investigated a few simple inhomogeneous interacting quantum systemsin traps, and their edge properties. Our main result is extremely simple to formulate: at the21 ciPost Physics Submission edge the particle density goes to zero, so sufficiently local interactions are also renormalizedto zero. While this observation is well-known from standard TBA arguments, the fact that itholds at a subleading scale is perhaps underappreciated. This partly explains the universalityof such edge distributions, in particular T-W. In our case its appearance is ensured by thevalidity of the LDA (or semiclassical) hypothesis in the bulk, and then taking the edge limit,which is smooth. More importantly, the LDA/TBA approach also allowed us to computeexactly the length scale associated to T-W, essentially the only free parameter for such scaling.All those claims were carefully checked by large scale DMRG calculations in a spin chainmodel, that also admits Lieb-Liniger as a limit.It is of course difficult to prove our semiclassical treatment, since the system is non inte-grable, but already a proof for discrete inhomogeneous spin chains that map to free fermionswould be very interesting. Note also that the argument should carry over to inhomogeneousquantum systems whose homogeneous analogs are not integrable, but in that case we wouldnot be able to compute analytically the location of the edge and the scale associated to T-W,as we did in the present paper.There are several interesting directions for future investigations, let us mention someof those now. First, we only looked at ground states here, but it would be interesting toinvestigate finite-temperature effects, and see whether the Hamiltonian (2) still emerges atthe edge in the presence of interactions. Even though the edge effects are too small to beaccessible to current cold-atom experiments at small but finite temperature, such a resultwould nevertheless provide a clear experimental prediction.A better understanding of edge universality classes in out-of-equilibrium quantum prob-lems is obviously left as an important open problem. For the quench from a domain wall statethe edge distribution can in principle be computed exactly using the method put forwardin [53, 68–70], and applying it to the exact EFP in the six vertex model with domain wallboundary conditions, for which multiple integral exact formulas are available [71, 72]. Thismight provide a way to rigorously study those new edge universality classes for any valueof ∆, but technical difficulties, while we believe not insurmountable, remain formidable. Amore heuristic approach would be to better understand the corrections to GHD, which areless understood in the diluted regime we are interested in.The point ∆ = 1 is also of great interest, especially given the fact that the (sub-ballistic)transport properties of this point are still theoretically not so well understood for the purestates [50, 51, 53, 54, 73] encountered here. Studying the distribution of the rightmost particlein that case would be of great interest; since the signal still spreads ballistically, we expect aneven more spectacular long tail effect, possibly related to the difficulties in reliably extractinga (non-overestimated) transport exponent. The long tail effect should also be present for | ∆ | > 1, even though there is no transport in this quench.Let us finally emphasize that we have looked here at pure states, that have zero entropyin string-TBA language. Finite entropy states are more relevant when the system is preparedin a thermal density matrix. This means there is no direct connection between our edgebehaviors and the corrections to GHD studied in [74], or the transport studies [75, 76] atthe Heisenberg point. Finally, investigating edge distributions in chaotic systems would behighly desirable, also in relation to operator spreading. Looking at those problems with theperspective of the present paper should shed some light on these timely issues. Acknowledgments. I thank Filippo Colomo, Andrea De Luca, Jacopo De Nardis, J´erˆomeDubail, Alexandre Lazarescu, Gr´egoire Misguich, Herbert Spohn, Eric Vernier and Jacopo22 ciPost Physics Submission Viti for enlightening discussions. I am also grateful to J´erˆome Dubail and Herbert Spohnfor a careful reading of the manuscript. The DMRG calculations were performed using theITensor C++ library [35]. 23 ciPost Physics Submission A Other universality classes We have demonstrated in this paper how T-W naturally emerges at the edge of an inhomoge-neous interacting system. Our main motivation was to partially fill a gap in the literature, andfocus on interacting quantum system at equilibrium, which have not been much investigatedin this context. This does not mean that T-W scaling is systematic, as we briefly discuss here,however. In appendix A.1 we look at simple free fermions problems that do not exhibit T-Wbehavior, but are described by higher order free fermions kernels. Appendix A.2 deals withthe Calogero-Sutherland-Moser model, which belongs to the universality class of β -matrix en-sembles, which is not free fermions. An even more spectacular and less understood exceptionis discussed in section 3 in the main text. A.1 Tuning the dispersion relation Let us go back to the spin chain in a magnetic field studied in section 2.2. As is well known,the point ∆ = 0 can be mapped onto free fermions, upon performing a Jordan-Wigner trans-formation. In terms of lattice fermions { c i , c † j } = δ ij the Hamiltonian reads H = 12 (cid:88) x (cid:16) c † x +1 c x + c † x c x +1 − h ( x/R )(2 c † x c x − (cid:17) . (48)The homogeneous case (constant h ) can be solved by going to Fourier space. The dispersionrelation reads in that case ε ( k ) = cos k − h . For a varying magnetic field, LDA tells us theground state propagator is the kernel of the projection cos k − h ( x/R ) < 0. Near the edge x e = ± Rh − (1), the cosine may be expanded to second order at k = 0 , π , and we recoverT-W scaling.It is of course possible to consider different dispersion relations, which correspond to addingnext nearest neighbors hoppings. For example the choice ε ( k ) = cos k − cos 2 k is quarticaround k = 0, (cid:15) ( k ) = 3 / − k / O ( k ). This means the corresponding edge behavior willbe governed by the kernel of the projection18 d dx + xR ≤ , (49)which implies R / scaling at the edge, instead of R / . The distribution of the rightmostparticle will then be given by a different distribution, built with a kernel constructed fromfunctions A ( u ) = (cid:82) R dq π e iqu + iq / , instead of Airy functions. This kernel has been studied ina slightly different free fermions context in [77].Several other examples have been found in statistical mechanical literature, in particularin relation to limit shapes. Those include the Pearcey kernel [78] for quartic singularities(Airy is cubic), or the tacnode kernel [79] (which includes, roughly speaking, quadratic bandtouching). We refer to [80] for a review of these free fermionic universality classes. A.2 Calogero-Sutherland models and β -matrix ensembles Another exception to our previous discussion is provided by the Calogero-Sutherland-Mosermodel [81] in a harmonic trap, with first quantized Hamiltonian H = N (cid:88) j =1 (cid:32) − ∂ ∂x j + x j (cid:33) + (cid:88) i (cid:54) = j β ( β/ − x i − x j ) . (50)24 ciPost Physics Submission This is a long range interacting system for β (cid:54) = 2, to which our previous renormalizationargument does not apply. Contrary to models with short-range interactions such as Lieb-Liniger, the diluted particles close to the edge might still interact strongly with their bulkcounterparts, so we do not necessarily expect free fermions factorization.It can be shown analytically that this is precisely what happens. For inverse squareinteractions the wave function can be obtained exactly, and its modulus square given by | ψ cs ( x , . . . , x N ) | ∝ (cid:89) ≤ i Level-spacing distributions and the airy kernel , Communi-cations in Mathematical Physics (1), 151 (1994), doi:10.1007/BF02100489.[2] J. Baik, P. Deift and K. Johansson, On the distribution of the length of the longestincreasing subsequence of random permutations , J. Amer. Math. Soc. (4), 1119 (1999),doi:10.1090/S0894-0347-99-00307-0.[3] M. Pr¨ahofer and H. Spohn, Scale invariance of the png droplet and the airy process ,Journal of Statistical Physics (5), 1071 (2002), doi:10.1023/A:1019791415147.[4] M. Pr¨ahofer and H. Spohn, Universal distributions for growth processes in dimensions and random matrices , Phys. Rev. Lett. , 4882 (2000),doi:10.1103/PhysRevLett.84.4882.[5] K. Johansson, Shape fluctuations and random matrices , Communications in Mathemat-ical Physics (2), 437 (2000), doi:10.1007/s002200050027.[6] K. Johansson, The arctic circle boundary and the airy process , Ann. Probab. (1), 1(2005), doi:10.1214/009117904000000937.[7] C. A. Tracy and H. Widom, Asymptotics in asep with step initial condition , Communi-cations in Mathematical Physics (1), 129 (2009), doi:10.1007/s00220-009-0761-0.25 ciPost Physics Submission [8] V. Eisler and Z. R´acz, Full counting statistics in a propagating quantumfront and random matrix spectra , Phys. Rev. Lett. , 060602 (2013),doi:10.1103/PhysRevLett.110.060602.[9] N. Allegra, J. Dubail, J.-M. St´ephan and J. Viti, Inhomogeneous field theory inside thearctic circle , Journal of Statistical Mechanics: Theory and Experiment (5), 053108(2016).[10] H. Spohn, Exact solutions for kpz-type growth processes, random matrices, and equilib-rium shapes of crystals , Physica A: Statistical Mechanics and its Applications (1), 71(2006), doi:10.1016/j.physa.2006.04.006, Fundamental Problems in Statistical Physics.[11] V. L. Pokrovsky and A. L. Talapov, Ground state, spectrum, and phase dia-gram of two-dimensional incommensurate crystals , Phys. Rev. Lett. , 65 (1979),doi:10.1103/PhysRevLett.42.65.[12] G. C. Wick, The evaluation of the collision matrix , Phys. Rev. , 268 (1950),doi:10.1103/PhysRev.80.268.[13] K. Johansson, From gumbel to tracy-widom , Probability Theory and Related Fields (1), 75 (2007), doi:10.1007/s00440-006-0012-7.[14] D. S. Dean, P. Le Doussal, S. N. Majumdar and G. Schehr, Noninteracting fermionsat finite temperature in a d -dimensional trap: Universal correlations , Phys. Rev. A ,063622 (2016), doi:10.1103/PhysRevA.94.063622.[15] D. Betea and J. Bouttier, The periodic Schur process and free fermions at finite temper-ature , arXiv e-prints 1807.09022 (2018), https://arxiv.org/abs/1807.09022 .[16] L. S. Levitov and G. B. Lesovik, Charge distribution in quantum shot noise , SovietJournal of Experimental and Theoretical Physics Letters , 230 (1993).[17] M. A. Cazalilla, R. Citro, T. Giamarchi, E. Orignac and M. Rigol, One dimensionalbosons: From condensed matter systems to ultracold gases , Rev. Mod. Phys. , 1405(2011), doi:10.1103/RevModPhys.83.1405.[18] T. Tao, Topics in Random Matrix Theory , American Mathematical Society, ISBN9780821874301, doi:10.1365/s13291-013-0059-1 (2012).[19] P. J. Forrester, Log-gases and random matrices. , London Mathematical Society Mono-graphs. Princeton University Press, ISBN 978-0-691-12829-0 (2010).[20] V. Eisler, Universality in the full counting statistics of trapped fermions , Phys. Rev.Lett. , 080402 (2013), doi:10.1103/PhysRevLett.111.080402.[21] R. Marino, S. N. Majumdar, G. Schehr and P. Vivo, Phase transitions and edge scalingof number variance in gaussian random matrices , Phys. Rev. Lett. , 254101 (2014),doi:10.1103/PhysRevLett.112.254101.[22] P. Calabrese, P. Le Doussal and S. N. Majumdar, Random matrices and en-tanglement entropy of trapped fermi gases , Phys. Rev. A , 012303 (2015),doi:10.1103/PhysRevA.91.012303. 26 ciPost Physics Submission [23] D. S. Dean, P. Le Doussal, S. N. Majumdar and G. Schehr, Noninteracting fermions ina trap and random matrix theory , arXiv e-prints 1810.12583 (2018), https://arxiv.org/abs/1810.12583 .[24] M. Gaudin, The Bethe Wavefunction , Cambridge University Press,doi:10.1017/CBO9781107053885 (2014).[25] V. E. Korepin, N. M. Bogoliubov and A. G. Izergin, Quantum Inverse Scattering Methodand Correlation Functions , Cambridge Monographs on Mathematical Physics. Cam-bridge University Press, doi:10.1017/CBO9780511628832 (1993).[26] V. Dunjko, V. Lorent and M. Olshanii, Bosons in cigar-shaped traps: Thomas-fermiregime, tonks-girardeau regime, and in between , Phys. Rev. Lett. , 5413 (2001),doi:10.1103/PhysRevLett.86.5413.[27] D. M. Gangardt and G. V. Shlyapnikov, Local correlations in a strongly interactingone-dimensional bose gas , New Journal of Physics (1), 79 (2003), doi:10.1088/1367-2630/5/1/379.[28] Y. Brun and J. Dubail, The inhomogeneous gaussian free field, with applicationto ground state correlations of trapped 1d bose gases , SciPost Phys. , 37 (2018),doi:10.21468/SciPostPhys.4.6.037.[29] J. Dubail, J.-M. St´ephan and P. Calabrese, Emergence of curved light-conesin a class of inhomogeneous luttinger liquids , SciPost Phys. , 019 (2017),doi:10.21468/SciPostPhys.3.3.019.[30] A. G. Volosniev, D. Petrosyan, M. Valiente, D. V. Fedorov, A. S. Jensen and N. T.Zinner, Engineering the dynamics of effective spin-chain models for strongly interactingatomic gases , Phys. Rev. A , 023620 (2015), doi:10.1103/PhysRevA.91.023620.[31] V. Eisler and D. Bauernfeind, Front dynamics and entanglement in the xxz chain with agradient , Phys. Rev. B , 174301 (2017), doi:10.1103/PhysRevB.96.174301.[32] C. B. Mendl and H. Spohn, Searching for the tracy-widom distribution in nonequilibriumprocesses , Phys. Rev. E , 060101 (2016), doi:10.1103/PhysRevE.93.060101.[33] S. R. White, Density matrix formulation for quantum renormalization groups , Phys.Rev. Lett. , 2863 (1992), doi:10.1103/PhysRevLett.69.2863.[34] J. Hauschild and F. Pollmann, Efficient numerical simulations with tensor net-works: Tensor network python (tenpy) , SciPost Phys. Lect. Notes p. 5 (2018),doi:10.21468/SciPostPhysLectNotes.5.[35] Intelligent tensor library (itensor) , http://itensor.org/ .[36] I. M. Nolden and H. van Beijeren, Equilibrium shape of bcc crystals: Thermal evolutionof the facets , Phys. Rev. B , 17224 (1994), doi:10.1103/PhysRevB.49.17224.[37] S. Chhita and K. Johansson, Domino statistics of the two-periodic aztec diamond , Ad-vances in Mathematics , 37 (2016), doi:10.1016/j.aim.2016.02.025.27 ciPost Physics Submission [38] B. Golzer and A. Holz, The nonlinear schrodinger model as a special continuum limitof the anisotropic heisenberg model , Journal of Physics A: Mathematical and General (11), 3327 (1987), doi:10.1088/0305-4470/20/11/035.[39] B. Schmidt and M. Fleischhauer, Exact numerical simulations of a one-dimensionaltrapped bose gas , Phys. Rev. A , 021601 (2007), doi:10.1103/PhysRevA.75.021601.[40] I. Peschel, Calculation of reduced density matrices from correlation functions , Jour-nal of Physics A: Mathematical and General (14), L205 (2003), doi:10.1088/0305-4470/36/14/101.[41] M.-C. Chung and I. Peschel, Density-matrix spectra of solvable fermionic systems , Phys.Rev. B , 064412 (2001), doi:10.1103/PhysRevB.64.064412.[42] V. Eisler and I. Peschel, Surface and bulk entanglement in free-fermion chains ,Journal of Statistical Mechanics: Theory and Experiment (4), P04005 (2014),doi:10.1088/1742-5468/2014/04/p04005.[43] J. Dubail, J.-M. St´ephan, J. Viti and P. Calabrese, Conformal field theory for inhomo-geneous one-dimensional quantum systems: the example of non-interacting fermi gases ,SciPost Phys. , 002 (2017), doi:10.21468/SciPostPhys.2.1.002.[44] T. Antal, Z. R´acz, A. R´akos and G. Sch¨utz, Transport in the xx chain at zero tem-perature: Emergence of flat magnetization profiles , Phys. Rev. E (5), 4912 (1999),doi:10.1103/PhysRevE.59.4912.[45] D. Gobert, C. Kollath, U. Schollw¨ock and G. Sch¨utz, Real-time dynamics in spin- chains with adaptive time-dependent density matrix renormalization group , Phys. Rev.E , 036102 (2005), doi:10.1103/PhysRevE.71.036102.[46] T. Antal, P. L. Krapivsky and A. R´akos, Logarithmic current fluctuationsin nonequilibrium quantum spin chains , Phys. Rev. E , 061115 (2008),doi:10.1103/PhysRevE.78.061115.[47] J. Mossel and J.-S. Caux, Relaxation dynamics in the gapped xxz spin-1/2 chain , NewJournal of Physics (5), 055028 (2010), doi:10.1088/1367-2630/12/5/055028.[48] T. Sabetta and G. Misguich, Nonequilibrium steady states in the quantum xxz spin chain ,Phys. Rev. B , 245114 (2013), doi:10.1103/PhysRevB.88.245114.[49] J. Viti, J.-M. St´ephan, J. Dubail and M. Haque, Inhomogeneous quenches in a freefermionic chain: Exact results , EPL (Europhysics Letters) (4), 40011 (2016),doi:10.1209/0295-5075/115/40011.[50] M. Ljubotina, M. Znidaric and T. Prosen, Spin diffusion from an inhomoge-neous quench in an integrable system , Nature Communications , 16117 (2017),doi:10.1038/ncomms16117.[51] M. Ljubotina, M. ˇZnidariˇc and T. Prosen, A class of states supporting diffusive spindynamics in the isotropic heisenberg model , Journal of Physics A: Mathematical andTheoretical (47), 475002 (2017), doi:10.1088/1751-8121/aa8bdc.28 ciPost Physics Submission [52] M. Collura, A. De Luca and J. Viti, Analytic solution of the domain-wall nonequilibriumstationary state , Phys. Rev. B , 081111 (2018), doi:10.1103/PhysRevB.97.081111.[53] J.-M. St´ephan, Return probability after a quench from a domain wall initial state in thespin-1/2 xxz chain , Journal of Statistical Mechanics: Theory and Experiment (10),103108 (2017), doi:10.1088/1742-5468/aa8c19.[54] G. Misguich, K. Mallick and P. L. Krapivsky, Dynamics of the spin- heisen-berg chain initialized in a domain-wall state , Phys. Rev. B , 195151 (2017),doi:10.1103/PhysRevB.96.195151.[55] V. B. Bulchandani and C. Karrasch, Subdiffusive front scaling in interacting integrablemodels , arXiv e-prints 1810.08227 (2018), https://arxiv.org/abs/1810.08227 .[56] O. A. Castro-Alvaredo, B. Doyon and T. Yoshimura, Emergent hydrodynamics inintegrable quantum systems out of equilibrium , Phys. Rev. X , 041065 (2016),doi:10.1103/PhysRevX.6.041065.[57] B. Bertini, M. Collura, J. De Nardis and M. Fagotti, Transport in out-of-equilibrium xxz chains: Exact profiles of charges and currents , Phys. Rev. Lett. , 207201 (2016),doi:10.1103/PhysRevLett.117.207201.[58] M. P. Zaletel, R. S. K. Mong, C. Karrasch, J. E. Moore and F. Pollmann, Time-evolvinga matrix product state with long-ranged interactions , Phys. Rev. B , 165112 (2015),doi:10.1103/PhysRevB.91.165112.[59] K. Bidzhiev and G. Misguich, Out-of-equilibrium dynamics in a quantum impurity model:Numerics for particle transport and entanglement entropy , Phys. Rev. B , 195117(2017), doi:10.1103/PhysRevB.96.195117.[60] X. Zotos, Finite temperature drude weight of the one-dimensional spin-1/2 heisenbergmodel , Phys. Rev. Lett. , 1764 (1999), doi:10.1103/PhysRevLett.82.1764.[61] A. Kl¨umper and K. Sakai, The thermal conductivity of the spin-1/2 xxz chain at arbi-trary temperature , Journal of Physics A: Mathematical and General (9), 2173 (2002),doi:10.1088/0305-4470/35/9/307.[62] T. Prosen, Open xxz spin chain: Nonequilibrium steady state and a strict bound on ballis-tic transport , Phys. Rev. Lett. , 217206 (2011), doi:10.1103/PhysRevLett.106.217206.[63] T. Prosen and E. Ilievski, Families of quasilocal conservation laws and quantum spintransport , Phys. Rev. Lett. , 057203 (2013), doi:10.1103/PhysRevLett.111.057203.[64] E. Ilievski and J. De Nardis, Microscopic origin of ideal conductivity in integrable quan-tum models , Phys. Rev. Lett. , 020602 (2017), doi:10.1103/PhysRevLett.119.020602.[65] V. B. Bulchandani, R. Vasseur, C. Karrasch and J. E. Moore, Bethe-boltzmann hy-drodynamics and spin transport in the xxz chain , Phys. Rev. B , 045407 (2018),doi:10.1103/PhysRevB.97.045407.[66] A. Urichuk, Y. Oez, A. Kl¨umper and J. Sirker, The spin Drude weight of the XXZ chainand generalized hydrodynamics , arXiv e-prints 1808.09033 (2018), https://arxiv.org/abs/1808.09033 . 29 ciPost Physics Submission [67] C. G´omez, M. Ruiz-Altaba and G. Sierra, Quantum Groups in Two-DimensionalPhysics , Cambridge Monographs on Mathematical Physics. Cambridge University Press,doi:10.1017/CBO9780511628825 (1996).[68] B. Pozsgay, The dynamical free energy and the loschmidt echo for a class of quantumquenches in the heisenberg spin chain , Journal of Statistical Mechanics: Theory andExperiment (10), P10028 (2013), doi:10.1088/1742-5468/2013/10/P10028.[69] L. Piroli, B. Pozsgay and E. Vernier, From the quantum transfer matrix to the quenchaction: the loschmidt echo in xxz heisenberg spin chains , Journal of Statistical Mechanics:Theory and Experiment (2), 023106 (2017), doi:10.1088/1742-5468/aa5d1e.[70] L. Piroli, B. Pozsgay and E. Vernier, Non-analytic behavior of the loschmidtecho in xxz spin chains: Exact results , Nuclear Physics B , 454 (2018),doi:10.1016/j.nuclphysb.2018.06.015.[71] F. Colomo and A. Pronko, Emptiness formation probability in the domain-wall six-vertexmodel , Nuclear Physics B (3), 340 (2008), doi:10.1016/j.nuclphysb.2007.12.016.[72] F. Colomo and A. G. Pronko, The arctic curve of the domain-wall six-vertex model ,Journal of Statistical Physics (4), 662 (2010), doi:10.1007/s10955-009-9902-2.[73] O. Gamayun, Y. Miao and E. Ilievski, Domain-wall dynamics in the landau-lifshitzmagnet and the classical-quantum correspondence for spin transport , Phys. Rev. B ,140301 (2019), doi:10.1103/PhysRevB.99.140301.[74] J. De Nardis, D. Bernard and B. Doyon, Hydrodynamic diffusion in integrable systems ,Phys. Rev. Lett. , 160603 (2018), doi:10.1103/PhysRevLett.121.160603.[75] E. Ilievski, J. De Nardis, M. Medenjak and T. Prosen, Superdiffusion inone-dimensional quantum lattice models , Phys. Rev. Lett. , 230602 (2018),doi:10.1103/PhysRevLett.121.230602.[76] S. Gopalakrishnan and R. Vasseur, Kinetic theory of spin diffusion and superdiffusionin XXZ spin chains , arXiv e-prints 1812.02701 (2018), https://arxiv.org/abs/1812.02701 .[77] P. Le Doussal, S. N. Majumdar and G. Schehr, Multicritical edge statistics for themomenta of fermions in nonharmonic traps , Phys. Rev. Lett. , 030603 (2018),doi:10.1103/PhysRevLett.121.030603.[78] C. A. Tracy and H. Widom, The pearcey process , Communications in MathematicalPhysics (2), 381 (2006), doi:10.1007/s00220-005-1506-3.[79] M. Adler, P. L. Ferrari and P. van Moerbeke, Nonintersecting random walks in theneighborhood of a symmetric tacnode , Ann. Probab. (4), 2599 (2013), doi:10.1214/11-AOP726.[80] K. Johansson, Edge fluctuations of limit shapes , arXiv e-prints 1704.06035 (2017), https://arxiv.org/abs/1704.06035 .[81] B. Sutherland, Beautiful Models , World Scientific, doi:10.1142/5552 (2004).30 ciPost Physics Submission [82] J. A. Ramirez, B. Rider and B. Vir´ag, Beta ensembles, stochastic airy spectrum, and adiffusion , J. Amer. Math. Soc. , 919 (2011), doi:10.1090/S0894-0347-2011-00703-0.[83] G. Borot, B. Eynard, S. N. Majumdar and C. Nadal, Large deviations of the maximaleigenvalue of random matrices , Journal of Statistical Mechanics: Theory and Experiment (11), P11024 (2011), doi:10.1088/1742-5468/2011/11/P11024.[84] L. Dumaz and B. Vir´ag, The right tail exponent of the tracy–widom β distribution , Ann.Inst. H. Poincar´e Probab. Statist.49