On unitary time evolution out of equilibrium

We consider $d$-dimensional quantum systems which for positive times evolve with a time-independent Hamiltonian in a nonequilibrium state that we keep generic in order to account for arbitrary evolution at negative times. We show how the one-point functions of local operators depend on the coefficients of the expansion of the nonequilibrium state on the basis of energy eigenstates. We express in this way the asymptotic offset and show under which conditions oscillations around this value stay undamped at large times. We also show how, in the case of small quenches, the structure of the general results simplifies and reproduces that known perturbatively.


Introduction
The problem of the fate at large times of an extended and isolated quantum system out of equilibrium was addressed analytically in [1] for the case of noninteracting fermions in one dimension.Understanding the case of interacting quasiparticles turned out to be a difficult problem that could be faced through the perturbative approach introduced in [2] and further developed in [3,4,5,6].The perturbative analysis applies to the basic way of dynamically generating nonequilibrium evolution, which has been called "quantum quench" [7,8] in analogy with thermal quenches of classical statistical systems: the system is in a stationary state until the sudden change of an interaction parameter leads to the new Hamiltonian that rules the unitary time evolution thereafter.The theory is perturbative in the quench size (see Eq. ( 35)), while the interaction among the quasiparticles can be arbitrarily strong.It allows to actually determine the nonequilibrium state generated by the quench and to follow analytically the time evolution.One result is that, in any spatial dimension, interaction leads to persistent oscillations of one-point functions (e.g. the order parameter) when the state produced by a homogeneous quench includes a one-quasiparticle mode and no internal symmetry prevents the observable to couple to that mode [2].It was observed in [4] that, if undamped oscillations are present at first order in the quench size, they will remain as a feature of the nonperturbative result, a prediction that found a remarkable confirmation in [9], where no decay of the oscillations was observed in a simulation of the Ising chain reaching times several orders of magnitude larger than the perturbative timescale.The undamped oscillations of [2] have also been observed in simulations performed over shorter time intervals, see [10,11,12,13,14] for a nonexhaustive list and [4,5] for a discussion of the different instances.
It is the purpose of this paper to further investigate the properties of unitary nonequilibrium evolution at large times.To this end, we consider a quantum system in d spatial dimensions with Hamiltonian      where H does not depend on time and is also translation invariant in space.We are interested in the time evolution of the system for t > 0, and expand the nonequilibrium state on the basis of the quasiparticle states |p 1 , . . ., p n of the Hamiltonian H, with coefficient functions f n (p 1 , . . ., p n ) which know about the evolution of the system since t = −∞.The main question we want to answer is how the one-point function of a local operator Φ depends on the coefficients f n .Hence, in order to study the features of the generic case, these coefficients are not specified and we perform a nonperturbative analysis relying on the structural properties of unitary time evolution of quasiparticle modes.While we analyze the problem in a more general way, we anticipate here the result for the fully translation invariant case (no spatial inhomogeneities are inherited from t < 0) and large times; we obtain where F Φ 1 is the one-quasiparticle form factor of the operator Φ and M is the quasiparticle mass.We show that all the f n 's enter the determination of h 1 and of the asymptotic offset A Φ , in a way that we find out.On the other hand, for dynamically generated nonequilibrium evolution, the condition h 1 F Φ 1 = 0 for the presence of the undamped oscillations in (2) amounts to f 0 f 1 F Φ 1 = 0 on symmetry grounds.For small quenches, which are a particular case of (1), we show how the structure of the present general results simplifies and reduces to that known from perturbation theory.While the result (2) is written for the case of a single quasiparticle species, we will also discuss the generalization to several species.
The paper is organized as follows.In the next section we consider the problem of one-point functions in the general case of the Hamiltonian (1), while in section 3 we analyze the fully translation invariant case.The way the general results reduce to those known perturbatively for small quenches is shown in section 4 before providing some final remarks in the last section.An appendix reviews some historical developments.

General setting
We consider a d-dimensional quantum system with the Hamiltonian (1), where H does not depend on time and is also translation invariant in space.We are interested in the time evolution for positive times, which we investigate considering the expectation values Φ(x, t) of local, scalar, Hermitian operators Φ(x, t) on the quantum state of the system.We consider that at t = −∞ the system is in the ground state |Ω of its Hamiltonian H 0 (x, −∞), with the normalization Ω|Ω = 1.Then, if I denotes the identity operator, the condition holds at t = −∞ and is henceforth preserved by the time evolution.We have where P denotes the momentum operator.The state |ψ of the system can be generally expanded on the basis of asymptotic quasiparticle 1 states |p 1 , . . ., p n of the theory with Hamiltonian H, which are eigenstates of energy and momentum with eigenvalues respectively.Energy and momentum of the quasiparticles are related as where M > 0 is the quasiparticle mass and measures the distance from a quantum critical point 2 ; we also adopt the state normalization The expansion of the state |ψ on the basis of asymptotic quasiparticle states takes the form 1 Quasiparticles are collective excitation modes commonly exhibited by statistical systems.For instance, a very general case in which they arise is the proximity of critical points (e.g. in spin models or quantum gases). 2 In general, statistical systems allow for critical points in their parameter space, and relativistic invariance emerges in their proximity.At the same time the formalism includes the nonrelativistic case as a low energy limit.For the sake of notational simplicity we refer to the case of a single quasiparticle species; generalizations will be discussed when relevant.

Translation invariant case
So far we allowed for t > 0 the presence of spatial inhomogeneities inherited from the time evolution at negative times.We now consider the case in which such inhomogeneities are absent and the state of the system is translation invariant.This is a particular case in which the coefficient functions in (8) take the form5 so that we have and the expression (11) for the expectation value becomes Use of ( 15) now leads to an expansion in terms of connected matrix elements which for the subtracted expectation value (20) reads Notice that the terms coming from the connected part of the original matrix elements (13) contain delta functions that do not mix the momenta {p i } and {q j }, while the terms coming from the disconnected parts involve delta functions of differences of these momenta.Hence, the expansion in terms of the connected matrix elements can be written in the form with coefficient functions h m,n (q 1 , ..., q m , p 1 , ..., p n ) and hm,n (q 1 , ..., q m , p 1 , ..., p n ) related to the coefficients fn of (24) as and 6hm,n (q 1 , ..., q m , p 1 , ..., where For t → ∞ the integrand of ( 26) rapidly oscillates because of the exponential factor e i( Ẽ−E)t and suppresses the integrals unless the phase is stationary, namely unless the momenta are small.The coefficient functions and the matrix elements in (26) generically go to constants in this limit 7 .The behavior for t large enough of the contribution to (26) with (m, n) quasiparticles is then obtained using the nonrelativistic expression of the energies and rescaling the momentum components by √ t; this gives where B Φ m,n and BΦ m,n are constants and we took into account that the term (m, n) = (0, 0) is not present in (26).We see that the leading time dependence comes from (m, n) equal (1, 0) and (0, 1) and corresponds to undamped oscillations.Notice that in absence of the delta function in (23) (i.e. in the generic inhomogeneous case of previous section) the oscillations coming from the (1, 0) and (0, 1) contributions would be damped as t −d/2 .
Since the only relativistic invariant that can be formed from the energy and momentum of a single particle is a constant, F Φ m,n with m + n = 1 is a constant.Besides (17), the matrix elements F Φ m,n with m and n interchanged are related by crossing symmetry [15], which amounts to analytic continuation in the momenta; this leads to the real constant Putting all together, the large time limit of the one-point function ( 21) is given by (2) with and In the physical dynamical problems we consider, the fn 's in (23) are nonzero unless an internal symmetry 8 forces some of them to vanish.In the current case of a single particle species, a Z 2 symmetry may lead to the vanishing of the fn 's with n even or of those with n odd.Since h 1 is a sum of terms containing f * n fn+1 , we have h 1 = 0 unless f0 f1 = 0. Hence, the condition for the presence of undamped oscillations in (2) is f0 f1 F Φ 1 = 0. Notice that the asymptotic offset (33) differs from the constant Φ(0) − G s Φ (0) entering (21) for the subtraction of the term (m, n) = (1, 1) in the sum; the reason is that in (21) this term is canceled by the (1, 1) contribution to (26), which has Ẽ = E and is time-independent.Equation (30) shows that the first subleading contributions to (2) come from (m, n) equal (0, 2), (2, 0), (1, 2) and (2, 1) and correspond to damped oscillations.
The analysis we performed above is easily generalized along the same lines to the case of several quasiparticle species a = 1, 2, ..., k with masses M a .In particular, the oscillations that remain undamped at large times take the form where the first sum generalizes the term present in (2), with F Φ 1a and h 1a corresponding to (31) and (32) with the specification of the species of the particle.The second sum, on the other hand, is a contribution arising from the fact that the term (m, n) = (1, 1) is no longer time-independent when the two particles have different masses; h 1b,1a corresponds to h 1,1 of ( 27) with the specification of the quasiparticle species 9 .We also see that the condition for the presence of this second type of undamped oscillations is f1b f1a F Φ 1b,1a = 0. Once again this condition involves one-quasiparticle states and is satisfied unless an internal symmetry causes the vanishing of one of the three factors.This clarifies the role of symmetries for undamped oscillations, a role that had been debated in the literature (see [17] and references therein).In the perturbative theory of quantum quenches, the undamped oscillations with frequencies M b − M a arise at second order in the quench size [4], as will also be seen in the next section.

Comparison with perturbative results
The Hamiltonian (1) includes as a particular case that in which the negative and positive time Hamiltonians differ for the change of an interaction parameter, namely the homogeneous quench 9 For a = b, F Φ 1b,1a = F Φ,c 1b,1a follows from the fact that particles of different species cannot annihilate each other.The contribution multiplying h1b,1a is damped at large times.
We then refer to Ψ(x) as the quench operator and to λ as the quench size.A general perturbative analysis can be performed in the quench size [2], in any dimension d [5] and for arbitrarily strong interaction among the quasiparticles.When the system is in the ground state of H 0 for negative times, the post-quench state is given by ( 23) with10 [2, 5] It then follows from (32), ( 27), ( 28) and ( 31) that and from (2) that undamped oscillations show up already at leading order in the quench size, as originally shown in [2].Notice also that ( 27) leads to so that in perturbation theory the undamped oscillations with frequences M a − M b in (34) arise at second order in the quench size, as observed in [4].While the expressions (39) and (40) coincide with those of the perturbative calculations of [2,4], a subtlety has to be pointed out: those perturbative calculations were performed in the basis of the quasiparticle states of the pre-quench theory (i.e. the unperturbed theory λ = 0), while the basis we use in this paper is that of the t > 0 theory.The point, however, is that the difference between the two bases can be ignored when working at leading order in perturbation theory 11 .
We see from (36), (37), ( 27) and ( 28) that h m,n are of order λ 2 unless m or n vanish, and that hm,n are in any case of order λ 2 .It follows that the asymptotic offset (33) takes the form where Φ(0) is now the expectation value on the ground state of the pre-quench theory, and It was shown in [3] that (41), (42) lead to where Φ λ is the expectation value on the ground state of the post-quench theory.The nonperturbative expression (33) suggests that in general there is no simple way of expressing the offset.
In our analysis of the large time behavior of one-point functions in the previous section we used the fact that the matrix elements ( 14) generically go to some finite constant when the momenta tend to zero.In d = 1, however, the quasiparticles often possess fermionic statistics 12 , and for the matrix elements (14) this leads to a zero when q i = q j or p i = p j , and to a pole when q i = p j .The poles are accompanied by an iǫ prescription [19] which anyway displaces them from the integration path over momenta 13 .The effect of the zeros in the matrix elements (14) and in the coefficient functions which multiply them in the expressions such as (26) is to enhance the time decay in (30), without affecting the undamped oscillations in (2) which are generally derived in any dimension.The perturbative results recalled above are of course consistent with this fact, and indicate that fermionic statitistics in d = 1 enhances the decay of the remainder in (2) (t −3/2 instead of t −1/2 ) [2,5].If the quasiparticles have fermionic statistics the suitable sign factors have to be introduced in (15) and will affect the combinatorial prefactor in (28).

Conclusion
In this paper we studied the nonequilibrium dynamics of quantum statistical systems in d spatial dimensions which for positive times evolve with a Hamiltonian H which is timeindependent and translation invariant in space.The nonequilibrium state was expanded on the basis of energy eigenstates (asymptotic quasiparticle states) of H, with coefficient functions f n which were left generic in order to account for arbitrary evolution for t < 0 under some Hamiltonian H 0 (x, t).We then showed how the evolution for positive times of the one-point functions of local operators (e.g. the order parameter) depends on the f n 's.
While the theory shows that the large time dynamics is determined by low-energy modes, our framework ensures that the results hold also in the vicinity of quantum critical points.It also allows to appreciate the role played by the connectedness structure of matrix elements, a circumstance noted since [20], where this structure was shown to account for the light cone spreading of correlations in two-point functions.
In the generic case (1), in which translation invariance is absent, the theory leads to oscillations of the one-point functions that normally decay as t → ∞.This is the case, in particular, when the system is confined in a finite region of space before this spatial constraint is removed for t > 0 (release from a trap 14 ).In such a situation, the energy density carried by the quasiparticle excitations goes to zero as t → ∞ (local dissipation) and is insufficient to sustain the oscillations at large enough times.A first illustration of this phenomenon was given perturbatively in [4,5] in the framework of inhomogeneous quantum quenches.
On the other hand, when the analysis is specialized to the case in which no spatial inhomogeneity is inherited from negative times, the theory shows that one-point functions exhibit undamped oscillations when no internal symmetry prevents a one-quasiparticle contribution to the nonequilibrium state or the coupling of the operator to this contribution.This result confirms the one obtained perturbatively since [2,4] for the case of the instantaneous change of an interaction parameter.
We also obtained the expression (33) of the asymptotic offset of one-point functions in terms of the matrix elements of the operator and of the coefficients specifying the nonequilibrium state.We showed how the structure of this result simplifies in the particular case of a small quench from the ground state and allows, up to higher order corrections in the quench size, the resummation originally shown in [3].While no similar resummations seem likely for the full result (33), it is known that expansions over quasiparticles modes often converge rapidly providing very good approximations from the first few terms 15 .It will be interesting to investigate to which extent this happens in the present case.

A Appendix
Since energy is conserved, nothing can force collective oscillation modes of an isolated homogeneous statistical system to always decay.In this respect, the undamped oscillating term (39) obtained in [2] for one-point functions after a quench of size λ provided the first analytical result.The fact that F Ψ 1 vanishes for noninteracting quasiparticles naturally accounted for the circumstance that the undamped oscillations were not found in [1] nor in the many studies devoted in more recent years to the transverse field Ising chain (free fermions).The result (39) also explained that interaction alone is not sufficient to produce such oscillations; it is also necessary that no internal symmetry causes the vanishing of the one-quasiparticle matrix elements F Ψ 1 and F Φ 1 .It was then pointed out in [2] that the simplest system where to look for those oscillations is the Ising chain with longitudinal field.Indeed, the longitudinal field introduces interaction among the quasiparticles and leaves no internal symmetry.The author of [2] did not know that "strong numerical evidence" of undamped oscillations had already been observed precisely in this model in [10].
The analysis of [2] was perturbative in the quench size λ, and it was pointed out in the same paper that the results of finite order perturbation theory can be quantitatively accurate up to the timescale t λ ∼ 1/λ 1/(d+1−X Ψ ) , where X Ψ is the scaling dimension of the quench operator.While t λ can be made arbitrarily large reducing the size of the quench, it is finite for fixed λ and the question of the fate of the oscillations beyond this timescale has to be posed.It was then observed in [4] that, for F Ψ 1 F Φ 1 = 0, the first order result (39) actually implies undamped oscillations in the full result (first order plus higher orders) for the one-point functions of usual physical interest.To see this, let us focus on the remainder of (39), namely the resummation of all terms beyond first order.Mathematically, this will be a function that for t → ∞ can either diverge, or approach a constant value, or itself exhibit undamped oscillations.On the other hand, since the contribution of order λ is bounded, the first possibility cannot occur for ordinary physical observables such as a local magnetization: at a given point of space this can grow in time at most to the limit of maximal ordering, but cannot diverge.Having discarded the possibility of divergence, the other two possibilities lead to undamped oscillations for the complete result (first order plus remainder) as long as these are present at first order16 [4].
What we just recalled clearly indicated that it should be possible to derive the presence of undamped oscillations of one-point functions in isolated homogeneous systems with interacting quasiparticles without reference to perturbation theory, and this is what we have done in this paper.The result (2) is mathematically clear: for t → ∞ the remainder vanishes at least as t −d/2 and there will be undamped oscillations as long as h 1 F Φ 1 = 0.While in this respect this is the end of the story, it may be worth spending some words on a point that might be confusing.Before [2] there was no way to analytically determine the state produced by the quench (35) in the generic case of interacting quasiparticles and some studies based on trial states were proposed.In particular, with the idea that this could be relevant for models with interacting quasiparticles and integrable equilibrium dynamics in d = 1, it was proposed (see e.g.[22]) to consider states of the kind arising in a different physical problem, namely that of perfectly reflecting boundary conditions in integrable models [23].These states are made of pairs of quasiparticles with opposite momenta and this "pair structure" exponentiates.Hence, we refer to them as exponential states; sometimes they are also called "integrable states" for the analogy we recalled.Their relevance for the quench problem (35) with interacting quasiparticles could not be shown and, in the words of [24], this issue was "simply set aside".On the other hand, the results of [2] recalled in section 4 (remember footnote 10) showed that such a peculiar feature as the pair structure -not to mention its exponentiation -does not arise in the quench (35) if the quasiparticles interact, and that in such a case the full complexity of the state (23) needs to be faced 17 .
Hence, we see that in the case of interacting quasiparticles these exponential states can only be interpreted as initial conditions in which at t = 0 the state is fine tuned by hand on the exponentiated pair structure; they do not arise in our problem (1) in which nonequilibrium is generated dynamically by the time evolution since t = −∞ and no fine tuning is possible.This essential physical difference completely changes also the technical nature of the problem with respect to our case.Indeed, since an exponential state does not allow the notions of negative time and of pre-quench state18 , one-point functions are normalized dividing by the norm of the exponential state itself.Divergences specifically due to the fine tuning on the pair structure appear both in the numerator and denominator, and a "linked cluster expansion" (in the amplitude of the quasiparticle pairs assumed small) as well as specific regularization steps are devised in order to produce a finite result (see [24,25,26]).This recipe yields some positive powers of time retaining memory of the assumed exponential structure of the state, and these powers are conjectured to resum into an exponential.In this way an exponential damping of one-quasiparticle oscillation modes was proposed in [25].However, we see that this exponential damping of oscillations is a conjecture relying on a chain of assumptions, and that already the first of these assumptions (the pair structure) is ruled out in the problem of dynamically generated nonequilibrium with interacting quasiparticles considered in the present paper.
It may also be worth pointing out that this example of exponential states illustrates something more general.The derivation of the present paper shows that the possibility of undamped collective oscillation modes follows from general properties of unitary time evolution in isolated homogeneous quantum systems in any dimension, as required by the fact that nothing can force oscillation modes to always decay in presence of energy conservation.
No fine tuning of the coefficient functions f n (p 1 , . . ., p n ) in the state (8) to some special form can change the conclusion valid for the generic case.In addition, no such fine tuning can be physically justified for the class (1) -or its subclass (35) -of dynamical problems we considered.
The irrelevance of fine tuning arguments for the problem at hand is by now well known and they are not mentioned in [27], a recent paper considering persistent oscillations.These authors, clearly unaware that the observation of [4] recalled earlier in this appendix rules out this possibility, asked whether oscillations undamped at first order in the quench size can altogether disappear in the full result.It is then not surprising that their treatment by a truncated BBGKY hierarchy of a dimerized XXZ chain in a staggered magnetic field was inconclusive: in their words, "it is not impossible that higher-order corrections would cause the oscillations to remain at late times" 19 .Of course, from the general point of view, the possibility of undamped oscillations is not questioned in [27] and realizations are listed 20 .In the present paper we showed nonperturbatively how and under which conditions they emerge.
Since undamped oscillation modes are clearly possible in energy-conserving homogeneous systems and are theoretically understood, it is not surprising that they are observed in numerical simulations.As a matter of fact, "strong numerical evidence that such states do not relax even at very long times" was obtained already several years ago [10] for the Ising chain with longitudinal field 21 .For its simplicity, this model minimizes the numerical efforts 19 It is worth stressing that the perturbative parameter involved in the considerations of [27] is not the quench size of [2] but the interaction strength.This is problematic, since the oscillations in their model are attributed to a bound state.The latter is nonperturbative in the interaction strength and a consistent perturbative study of the oscillations starting from the noninteracting theory is awkward in principle before than in practice.This problem is absent in the perturbation theory in the quench size, where the quasiparticle spectrum is already that of the interacting theory and can include bound states, topological excitations, etc. [2,3,5]. 20See the published version. 21From the point of view of the historical development of numerical studies, it is interesting to recall that the perturbative results of [2] (amplitude and frequency of the oscillations) were first numerically confirmed in [28] for a quench of the longitudinal field in the Ising chain.A decrease of the amplitude in time was observed in the same paper for a double quench (i.e. both in the longitudinal and transverse field), for which no perturbative formula is available.Since only the first three oscillations were observed, this initial decrease and leaves no margins of interpretation about the timescales.Outstanding evidence on the absence of damping was then given in [9], where the time evolution was followed over hundreds of oscillations and times exceeding by three orders of magnitude the two timescales (t λ and the inverse mass gap) allowed by the analytical formulation of the quench 22 .Lack of relaxation of one-point functions has by now been observed in several other numerical works in d = 1 (see e.g.[11,12,13,14]), and its search has started in d = 2 (see [31] and references therein) despite the size limitations that still affect simulations in this case.