Persistent oscillations after quantum quenches in $d$ dimensions

We obtain analytical results for the time evolution of local observables in systems undergoing quantum quenches in $d$ spatial dimensions. For homogeneous systems we show that oscillations undamped in time occur when the state produced by the quench includes single-quasiparticle modes and the observable couples to those modes. In particular, a quench of the transverse field within the ferromagnetic phase of the Ising model produces undamped oscillations of the order parameter when $d>1$. For the more general case in which the quench is performed only in a subregion of the whole $d$-dimensional space occupied by the system, the time evolution occurs inside a light cone spreading away from the boundary of the quenched region as time increases. The additional condition for undamped oscillations is that the volume of the quenched region is extensive in all dimensions.


Introduction
A main question for the nonequilibrium dynamics of extended quantum systems is whether time evolution eventually leads to some form of relaxation [1] or can produce a different behavior. The difficulty of the question calls for analytical study in the basic nonequilibrium setting, which corresponds to quantum quenches 1 . Here an isolated extended system is in the ground state |0 of its Hamiltonian H 0 until the time t = 0, when the instantaneous change of an interaction parameter leads to the new Hamiltonian H that rules the unitary evolution for t > 0. The procedure dynamically generates the post-quench state and leaves no ambiguities in the initialization of the nonequilibrium evolution.
The way to analytically study the quench dynamics on general grounds was found in [4], where the case of translation invariant one-dimensional systems was considered. The theory is perturbative in the quench size λ, and crucially incorporates the fact that normally the quasiparticle modes interact (even at equilibrium, namely at λ = 0). In particular, it was found that, under conditions determined by the theory and including quasiparticle interaction, local observables (e.g. the order parameter) exhibit oscillations undamped in time. It was argued in [5] that these oscillations appearing at leading order in λ are not washed out by higher orders and stay undamped also when λ is not small. Undamped oscillations were then numerically observed over hundreds of periods within the time scale accessible to the simulation of a large quench of the longitudinal magnetic field in the Ising chain [6], a basic candidate discussed in [4,7].
In this paper we consider the case of d spatial dimensions. We show that for a translation invariant system with post-quench Hamiltonian H = H 0 + λ dx Ψ(x), the large time limit of the one-point function of an operator Φ takes the form where Φ eq λ is the equilibrium expectation value in the theory with the Hamiltonian H, M > 0 the quasiparticle mass, F O 1 the matrix element of O between |0 and the singlequasiparticle state, and α ≥ d/2 .
Throughout the paper we adopt natural units in which the maximal velocity of the quasiparticles is v max = 1. In case of systems possessing several quasiparticle species the term in the square bracket is summed over the species. (1) shows, in particular, the presence of undamped oscillations under the same conditions determined in [4] for the case d = 1 (for which α is generically 3/2), namely when the excitations produced by the quench include a single-quasiparticle mode 2 (F Ψ 1 = 0) and the observable couples to this mode (F Φ 1 = 0). In turn, (1) is a particular case (D = R d ) of the more general situation that we are going to consider, namely that of a quench performed only in a subregion D of the full space R d occupied by the system. For such an inhomogeneous quench, corresponding to the post-quench Hamiltonian we show that the large time limit (1) generalizes to where α is the same entering (1), while rules the large time behavior t −β D of the integral term. If L d , with L → ∞, is the volume of R d , there are undamped oscillations (β D = 0) when the volume of D, vol(D), is of order L d , a condition weaker than D = R d . On the other hand, β D = d/2 when vol(D) is finite. Intermediate values of β D occur depending on the number of dimensions in which vol(D) is extensive. We also show that the effect of the quench is appreciable only inside a light cone that at time t contains D and the external region within distance t from the boundary ∂D. These analytical results -whose degree of generality is uncommon in the framework of nonequilibrium quantum dynamics -lead to the following physical deductions. The remainder of (1) (terms of order λ 2 and higher) has one of the following behaviors for t → ∞: indefinite growth in modulus, approach to a constant, or undamped oscillations. Hence, if we exclude that Φ(x, t) indefinitely grows in modulus -something that is not expected on physical grounds -(1) leads to undamped oscillations also when λ is not small (provided that F Ψ 1 F Φ 1 = 0). The same applies to (4) when β D = 0. To understand the role of D in this respect, we observe that the oscillations propagate inside the light cone, and are sustained by the energy produced in the quench and carried by the quasiparticles. In order the oscillations to stay undamped, a nonzero energy density inside the light cone is needed at late times (besides F Ψ 1 F Φ 1 = 0). Since the energy produced by the quench and conserved by the time evolution is proportional to vol(D), and since the volume enclosed by the light cone becomes L d as t → ∞, a nonzero energy density at late times requires that vol(D) is of order L d , and this is the case corresponding to β D = 0.
The paper is organized as follows. The main theoretical results are derived in the next section, while in section 3 we illustrate the features of the time evolution also with a number of examples of quenched regions D. Section 4 illustrates the role played by internal symmetries through the basic example of the quantum Ising model, and the last section contains some concluding remarks.

Post-quench state and one-point functions
We consider a d-dimensional system occupying the whole space R d . Before the quench the system is translation invariant and in the ground state |0 of the Hamiltonian H 0 . We perform the theoretical analysis exploiting the complete basis of asymptotic quasiparticle states |p 1 , . . . , p n of the pre-quench theory, with p i denoting the d-dimensional momenta of the quasiparticles. The asymptotic states are eigenstates of H 0 with eigenvalues equal to the sum of the quasiparticle energies E p i = M 2 + p 2 i . M > 0 is the quasiparticle mass and measures the distance from a quantum critical point. In order to simplify the notation we refer to the case of a single quasiparticle species; generalizations are straightforward and will be discussed when relevant.
The quench at t = 0 is performed changing the Hamiltonian to (3), and for D = R d breaks translation invariance. Since the quench excites quasiparticle modes, the prequench state |0 evolves into the state |ψ 0 = S λ |0 , where (T denotes chronological ordering) is the operator whose matrix elements n|S λ |0 give the probability amplitude that the quench induces the transition from |0 to |n . Here we adopt the compact notation |n = |p 1 , . . . , p n . To first order in the quench parameter λ we have where we defined introduced the notation with P the momentum operator and O a generic local operator. As usual, an infinitesimal imaginary part is given to the energy to make the time integral in (6) convergent 3 . The result (7) shows that the quench produces excitation modes with any number of quasiparticles and all possible momenta. Only when the quasiparticles do not interact, so that H 0 and H are quadratic in the quasiparticle modes and F Ψ n ∝ δ n,2 , a post-quench state with quasiparticles organized in pairs is obtained.
The one-point function of a local observable Φ is given by the expectation value Φ(x, t) on the post-quench state (7). In the formalism of asymptotic states the spacetime dependence is carried by the operator and is extracted exploiting (12). The variation δ Φ(x, t) of the one-point function of a hermitian observable with respect to the prequench value is given, at first order in λ, by where we took into account that normalizing by ψ 0 |ψ 0 = 1 + O(λ 2 ) is immaterial at first order, and added the term to ensure continuity at t = 0, namely the condition δ Φ(x, 0) = 0, which has no reason to be automatically satisfied.

Light cone from the quenched domain and behavior at large times
We can use (9) to rewrite (13) as For t large the rapid oscillation of the exponential suppresses the integrals over momenta unless the phase is stationary, namely unless where we introduced the quasiparticle velocities Since |v i | < 1, the stationarity condition (15) is satisfied when This means that, for any point y ∈ D, the effect of the quench is appreciable only within a distance t from y, namely the maximal distance that the quasiparticles excited by the quench at the point y could reach at time t. Hence, the time evolution takes place inside a light cone 4 containing D and the external region within distance t from ∂D (figure 1). For D finite, x fixed and t large enough, the stationarity condition v i = (y−x)/t implies that the time dependence in (13) receives a significant contribution only when all momenta p i are small, and then can be evaluated with g D → vol(D) and E p i → M + p 2 i /2M . The matrix elements (10) do not diverge at small momenta (see e.g. [9]), and in such a limit [F Ψ n ] * F Φ n will behave as momentum to a power 2α n ≥ 0. In particular α 1 = 0, since F O 1 is a real constant for the scalar hermitian operators that we consider. With this information, it is easy to rescale the momenta and see that the n-quasiparticle contribution in (13) behaves at large times as t −(nd/2+αn) .
More generally, suppose that D goes from −∞ to +∞ in k of the d spatial dimensions (say x 1 , x 2 , . . . , x k ). Now g D ∝ δ(P 1 ) · · · δ(P k ), and this gives an extra contribution −k to the counting of the powers of momentum, so that the large time behavior is modified to t −[(nd−k)/2+αn] . As long as F Ψ 1 F Φ 1 = 0 the leading contribution comes from the singlequasiparticle mode n = 1 and goes as t −(d−k)/2 . For an homogeneous quench (D = R d , i.e. k = d) we have undamped oscillations and the result in the square bracket of (1), with the lower bound (2) coming from n = 2. A generic k gives the time dependence in the square bracket of (4), behaving as t −β D , with β D = (d − k)/2 that satisfies (5).
On the other hand, suppose that D differs from R d for the subtraction of a finite region, namely D ∪D = R d withD finite. In this case g D = g R d − gD, and the dependence at large times corresponds to the difference of the previous cases k = d and k = 0. Hence, we have β D = 0 and undamped oscillations at sufficiently large times produced by the integral in (4). The conclusion that in general there are undamped oscillations when is nonzero corresponds to the physical picture anticipated in the introduction: the energy produced by the quench is proportional to vol(D) and spreads in time inside the light cone, so that the energy density will be asymptotically proportional to ρ D ; a nonzero ρ D is able to keep the oscillations undamped. The amplitude of the undamped oscillations goes to zero if ρ D goes to zero. The condition ρ D > 0 amounts to vol(D) extensive in all dimensions.
We finally explain the origin of the time-independent term in (1) and (4). For this purpose observe that in the equilibrium theory with Hamiltonian (3) the first order contribution in λ to the one-point function Φ(x) eq λ is 5 where we used (12), expanded over asymptotic states, and finally compared with (14).
Hence, recalling what we just concluded about the time-dependent part, we have when β D = 0. When β D = 0, namely when vol(D) is extensive in all dimensions, the r.h.s. of (20) is the value around which the undamped oscillations take place.
3 Following the time evolution

General features
We rewrite (4) in the form where f (x, t) = M 2 (2π) d dp 1 5 The subscript c indicates the connected part of the two-point function.
is a dimensionless function. For F Ψ 1 F Φ 1 = 0, f (x, t) determines the large time behavior of the one-point function for small quenches, until a time scale t λ that goes to infinity as λ is reduced 6 . On the other hand, the analysis of section 2.2 does not allow to neglect the terms with n > 1 in (13) when t is not large. However, it is known that in d = 1 [7,5] the contribution of the n-quasiparticle state is normally rapidly suppressed as n increases, so that f (x, t) provides a good approximation also for short times. If this is true also for d > 1 we should expect, recalling (19) and (14), that the first order contribution to the equilibrium expectation value is well approximated as We will soon explicitly illustrate that (23) indeed yields the result expected for δ Φ(x) eq λ , namely a function that is essentially constant for x ∈ D and zero otherwise. With this anticipation, we see that (21) (without the term O(t −α−β D )) can be used not only for large times but, with good approximation, also for short times, meaning that the function (22) yields a global view of the time evolution for small quenches. We now consider this function for a number of quenching domains D.

Rotationally invariant quenched domains
If D is the d-dimensional sphere of radius b centered in the origin, (9) and (22) yield and respectively, with J α (z) the Bessel function and r = |x|. The function (25) is plotted at different times in d = 2 (D = disk) in figure 2. The result for t = 0 confirms what we anticipated about (23) and the ability of f (x, t) to give an accurate view of the time evolution also for short times. As t increases the figure clearly shows the spreading of the light cone located at distance t from the boundary of the disk. The boundary modes propagate also inside the disk, but for t < b they leave unaffected the central region with r < b − t. Here the function essentially behaves as for a homogeneous quench, namely is spatially constant with undamped oscillations in time. The presence of the light cone is also illustrated in figure 3. Figure 4 shows the damping of the oscillations at large times, with the suppression t −d/2 expected for the case D finite. A straightforward generalization is that of a quench with D the d-dimensional spherical shell b 1 < r < b 2 , which yields The time evolution is illustrated in figures 5 and 6 for b 2 finite, and in figure 7 for b 2 infinite.  Figure 8 shows the corresponding function f (x, t) at different times in d = 2 for a square domain. Again, the image at t = 0, with the function essentially constant inside the square and vanishing outside, illustrates (through (23)) that the time evolution is described with good approximation also at short times. For t > 0, the general result that the light cone is located at distance t from the boundary of D leads to a rounding in correspondence of the corners of the square. This example gives an idea of the patterns that can be expected when ∂D increasingly deviates from a smooth surface.

Role of internal symmetries: Ising model
The model-dependent information for the use of (1) and (4) is whether the matrix elements F Ψ 1 and F Φ 1 vanish or not. They are generically nonzero, unless a symmetry forces them to vanish. In this section we illustrate the role of symmetries through the basic example   We will consider quenches in which the system, which for t < 0 is in the ground state |0 of the Hamiltonian H 0 given by (28) with h x = 0, evolves for t > 0 with the Hamiltonian (3) with 7 λ Ψ(x) equal to δh z σ z (x) (δh z ≪ h z ) or to h x σ x (x). Depending on the fact that we start from the paramagnetic or ferromagnetic phase, we have the four quenches depicted in figure 9.
The specialization of (1) and (4) to the different quenches proceeds through symmetry considerations in the pre-quench theory. In the paramagnetic phase, the fundamental quasiparticle excitation is created by the order parameter operator σ x (x) and then is Z 2odd; it follows that F σ z 1 = 0 and F σ x 1 = 0. In the ferromagnetic phase the symmetry is spontaneously broken and both F σ z 1 and F σ x 1 are nonzero in d > 1. The case d = 1 is special because the excitations in the ferromagnetic phase have a topological nature (they are kinks, see [10] for a review) and can couple to σ x and σ z only in topologically neutral pairs, with the consequence that F σ z 1 = F σ x 1 = 0. This information about F σ z 1 and F σ x 1 determines the time evolution of σ x and σ z through (1) when the quenches 8 I-IV are performed in the whole space R d , and through (4) when they are performed only in a subregion D. In particular we saw that, if vol(D) is extensive in all dimensions, undamped oscillations at large time occur when F Ψ 1 F Φ 1 = 0, and we explicitly provide this information in table 1 for the case Φ = σ x . The peculiarity of the cases III and IV in d = 1 follows from the previous observation about kinks.
The topological nature of the excitations of the ferromagnetic phase in d = 1 gives rise to an additional caveat about quench III in the Ising chain. While (1) and (4) hold 7 Our theory is formulated in the continuum, which in the present case is accessed working not too far from the critical point. 8 For quench I the symmetry implies σ x = 0 at all orders in λ.  figure 9. The information about F Ψ 1 F σ x 1 determines the time evolution of the order parameter σ x for homogeneous quenches through (1), and for quenches in the subregion D through (4). The peculiarities of the case d = 1 are discussed in the text. also in this case and do not yield undamped oscillations at first order in λ = h x , the longitudinal field makes kinks unstable and confines them into topologically neutral pairs [11,12,10], thus generating nonperturbatively the single-quasiparticle modes able to produce undamped oscillations of the order parameter on a time scale that becomes accessible for h x not too small 9 . These oscillations have indeed been observed numerically in [16]. Hence, we see that for quench III kink confinement in d = 1 produces for h x large enough the undamped oscillations that in d > 1 are already present at first order in h x . For quench IV, instead, the stability of the kinks precludes 10 in d = 1 the undamped oscillations that arise in d > 1. It is worth recalling that the presence of kinks in the spontaneously broken phase of systems with discrete symmetry, as well as their confinement under explicit symmetry breaking, are generic in d = 1 (see [19,20]), so that considerations analogous to those we made for Ising apply more generally.

Quench
In d = 1 the theory (28) possesses a single species of quasiparticles at h x = 0, so that (1) and (4) hold as they are. In d > 1, on the other hand, more species may be present, in which case the term where a labels the different species; a similar generalization occurs in (4). For the ferromagnetic phase in d = 2 there is numerical consensus [21,22,23,24,25,26] about the existence, 9 As recalled above, the results at first order in λ quantitatively hold until a time scale that goes to infinity as λ goes to zero. Ref. [7] contains a more detailed discussion of Ising quenches in d = 1, including the agreement with other analytical results [13,14] available for quenches I and IV (noninteracting fermions) when D = R. More generally, for d = 1 one can exploit the exact knowledge of all matrix elements (10) when the pre-quench theory is integrable [4]. Following the suggestion of an anonymous referee, we also mention that it has been argued in [15] that the dynamics following a specific initial condition with exponential structure -not originated by a quench -can lead to oscillation damping. 10 Another manifestation of the different nature of the excitations in the spontaneously broken phase of the Ising model in d = 1 and d > 1 arises in the physics of interfaces [17,18].
besides the lightest quasiparticle with mass M 1 , of a second stable quasiparticle with mass M 2 ≈ 1.8M 1 . A spectral analysis of the undamped oscillations expected for quench IV should provide an alternative, possibly more accurate, way of determining this mass ratio. It was also observed in [5] that in presence of several quasiparticle species oscillation frequencies M a − M b can arise at order λ 2 . In addition, the comparison performed in [7] with the result of [14] for noninteracting fermions indicates that corrections of order λ 2 and higher lead to the replacement of the pre-quench masses with the post-quench ones, and numerical evidence in this sense for interacting quasiparticles was given in [6,27]. By post-quench masses we mean the masses of the equilibrium theory with the post-quench values of the couplings; their difference from the pre-quench masses is of order λ. The mass ratio M 2 /M 1 to be observed in quench IV in d = 2 does not depend on λ in the scaling region.
The time evolution of the order parameter σ x (in our notations) following a homogeneous quench in the Ising ferromagnetic phase has been numerically investigated in d = 2 in [28,29], for a pre-quench value h i z = 0 that maximizes the distance from the continuum limit of our analytical study. Still, confirming a robustness of our results, the plots show undamped oscillations, at least for post-quench values h f z not too close to h c z . When h f z approaches h c z the correlation length becomes large, the long cylinder with a six-site circumference used in the simulation no longer approximates the plane, and a crossover to one-dimensional behavior (with damping of the oscillations) can be expected 11 . In perspective, it would be very interesting to have numerical results for small quenches to perform the first quantitative comparison with analytical results in d = 2.

Conclusion
We studied quantum quenches of systems in d spatial dimensions that are initially in the ground state of a spatially homogeneous Hamiltonian H 0 . The quench is performed instantaneously changing an interaction parameter inside a spatial region D. The analytical results that we derived provide unique benchmarking for numerical and experimental methods, as well as a general picture of the time evolution of local observables. In the first place the evolution takes place inside a light cone that originates from the boundary of D at the moment of the quench and spreads outwards as time increases. Inside the light cone the observable undergoes oscillations with frequency equal to the quasiparticle mass, which persist undamped at late times under two types of conditions. The first type, of dynamical nature, requires that the state produced by the quench includes single-quasiparticle modes and that the observable couples to them. These requirements are generically fulfilled when the quasiparticles interact, unless internal symmetries of the system cause the vanishing of some matrix elements, a mechanism that we illustrated through the paradigmatic example of the quantum Ising model. The second condition, of geometrical type, is that the energy density does not tend to zero at large times, and is fulfilled when the volume of D is extensive in all dimensions. The wavefronts spreading from the boundary of D are increasingly structured as the boundary deviates from a smooth surface. Our formulae apply to any choice of D, and we provided explicit illustrations of these features through some examples.