Non-equilibrium dynamics of isolated quantum systems

The non-equilibrium dynamics of isolated quantum systems represent a theoretical and experimental challenge raising many fundamental questions with applications to different fields of modern physics. In these proceedings, we briefly review some of the recent findings on the subject, with particular emphasis to the existence of stationary expectation values of local observables and to their statistical mechanics description. It turns out that the appropriate statistical ensemble describing these asymptotic values depends on whether the Hamiltonian governing the time evolution is integrable or not.


Introduction
Recent years witnessed an increasing theoretical and experimental effort in condensed matter and statistical physics in order to describe and understand the fate of isolated many-body quantum systems which are prepared in a non-equilibrium initial state (see e.g.[1,2] for reviews).When one of these many-body quantum systems is prepared in a non-equilibrium state |Ψ 0 and it is let evolve with a Hamiltonian H, the time evolved state is (setting = 1) A number of theoretical and experimental investigations have unambiguously shown that for large times and in the thermodynamic limit, the expectation values of the some observables relax to stationary values, although the dynamics governing the evolution is unitary and the initial state is pure.Under some specific conditions, these stationary values are the same as those obtained in a thermal ensemble in which the (effective) temperature is fixed by the expectation value of H in the initial state.This is in fact the problem of thermalisation of an isolated quantum system, that was initiated in 1929 by Von Neumann [3], but that for decades was considered only an academic question, because for all realistic condensed matter systems the coupling to the environment is unavoidable on the time scale of observation of the system.The situation completely changed after some recent and pioneering coldatomic experiments [4][5][6][7][8][9], as we are going to report soon.Furthermore the general findings of these studies have the potential to shed new light on long standing problems in apparently disconnected fields such as cosmology, nuclear and particle physics.For example, in high energy collisions of particles, as well as nuclei, there is a long standing "paradox" (see e.g.[10][11][12] and references therein), supported by a huge experimental evidence, that few-body a e-mail: calabrese@sissa.itobservables, in the energy range of the soft scale of QCD (O(GeV)), are well described by thermal mixed distribution, although the initial state is pure and the dynamics is unitary.It is reasonable to believe that the explanation of this apparent paradox is exactly the same for statistical systems which we are going to review in the following.
In the spirit of statistical physics, one of the main question emerging from these studies is to understand the most appropriate statistical ensemble describing the stationary behaviour without having to solve, on a case by case basis, the complicated non-equilibrium dynamics.
These proceeding are organised as follows.In Sec. 2 we briefly review two cold-atomic experiments in order to show the change of paradigm that they brought into statistical physics and condensed matter.In Sec. 3 we introduce the concept of reduced density matrix and we explain in which sense an isolated many-body quantum system can attain a steady state.In Sec. 4 some simple results for free field theories are explicitly derived, while in Sec. 5 more complicated results for interacting integrable models are only reviewed without any derivation.Finally, in Sec.6, we briefly mention some interesting topics of the nonequilibrium dynamics of quantum systems which have not been covered here for lack of space.

A brief report on cold atoms experiments
Since the pioneering work by Greiner et al. [4] showing the possibility of engineering and measuring interacting cold atomic gases under unitary non-equilibrium dynamics, many experiments (see e.g.[5][6][7][8][9]) have been performed, providing a large number of new results and constantly driving the theoretical research.Among all these experiments, we decided to limit our discussion to two particular ones which highlight the main features that we will describe theoretically in the following.The first one is the famous quantum Newton cradle [5].In this experiment, a cloud of a few hundreds dilute 87 Rb atoms is initially prepared in a harmonic trapping potential at a very low temperature (which can be considered practically zero for what follows).The density of the gas is such that the interaction can be taken to be point-like (i.e. a δ-function two body interaction), but it is strong enough and non negligible.The cloud is split by a laser beam in two counter-propagating clouds with opposite momentum.The two clouds then climb the harmonic potential up to the maximum value allowed by energy conservation, subsequently move back toward the centre of the trap where they interact; after this interaction the clouds climb again the potential and the process is repeated many times until the system becomes stationary (when it does).The system is recorded for many of these oscillations (obviously, since the measures are destructive, in the experiment this procedure is repeated many times, see the original reference for all details and for a suggestive graphical representation of the experiment).It has also been argued that the time evolution is essentially unitary during the whole probed time window.The results of this experiment are considered milestones.Indeed, it has been shown that the momentum distribution function attains for large time a stationary distribution for arbitrary space dimensionality.The details of the stationary values do depend on spatial dimensionality.In one dimension, the system relaxes slowly in time to a non-thermal distribution, while in two and three dimensions, systems relax very quickly and thermalise.It has been suggested that the one-dimensional case is special because the system is almost integrable, as we will discuss in the following.
A second cold-atom experiment [6] explores the existence of a stationary state in a one-dimensional bosonic lattice system starting from a non-equilibrium state.The Hamiltonian modelling the evolution -the Bose-Hubbard Hamiltonian-is known to be non integrable.The experiment shows that some observable approach stationary values that are compatible with the thermal ones, suggesting that integrability more than the space dimensionality is the crucial ingredient to establish asymptotic thermal behaviour.Indeed nowadays the common belief is that generic systems thermalise [13] while integrable ones attain stationary values described by a generalised Gibbs ensemble [14] in which the nontrivial integrals of motion have to be taken into account.

Stationary state and reduced density matrix
We still need to clarify in which sense some observables in the these non-equilibrium protocols can be described by a mixed state such as the thermal one.The crucial concept to solve this apparent paradox is the reduced density matrix.Let us imagine to take an arbitrary extended quantum system (in arbitrary dimension) and ideally divide into two complementary spatial parts (i.e.take a bipartition) denoted as A and A respectively.Under the unitary time evolution, the entire system A ∪ A will always be in the pure state |Ψ(t) given by Eq. ( 1), but this is not the case if we limit to consider only the subsystem A. Indeed starting from the density matrix of the entire system we can define the reduced density matrix of the subsystem A by tracing over the degrees of freedom in A, i.e. considering the so called reduced density matrix From the knowledge of ρ A (t) all the correlation functions local within A can be obtained.Thus as long as we are interested in local observable we do not have to retain information about the whole system, but we can limit to consider ρ A .The importance of ρ A stems from the fact that it is the quantity which generically displays a stationary behaviour described by some statistical ensemble, while the full density matrix ρ(t) always corresponds to a pure state with zero entropy.Indeed the von Neumann entropy of ρ A (t) has a non-vanishing extensive part (which is often called entanglement entropy, see e.g.[15]).Following Refs.[16][17][18][19], it is usually said that a system reaches a stationary state if a long time limit of the reduced density matrix exists, i.e. if the limit exists for any finite subsystem A. This limit can exist only in the case when we have already taken the thermodynamic limit for the entire system, otherwise phenomena like quantum revivals and recurrence would prevent strict relaxation (but it is still possible that time average quantities attain thermal values).We will say that the stationary state ρ A (∞) is described by a given statistical (mixed state) ensemble with full density matrix ρ E , if its reduced density matrix restricted to A equals ρ A (∞), i.e. if defined ρ A,E = Tr A (ρ E ), for any finite subsystem A, it holds In particular, this implies that arbitrary local multi-point correlation functions within subsystem A can be evaluated as averages with the density matrix ρ E .By no means this implies that ρ E equals the full density matrix of the system which is clearly impossible being the former a mixed state and the latter a pure one.When a system thermalises, ρ E is the standard Gibbs distribution ρ E ∝ e −βH and this is expected to be the case when the model is non-integrable.Here the inverse temperature β is not a free parameter but it is fixed by the requirement that the energy, which is a constant of motion, must have the same value in the initial state and in the asymptotic one, i.e.Tr[Hρ E ] = Ψ 0 |H|Ψ 0 .It is usually suggestively stated that the infinite part A of the system "acts as an heat bath for A".This expectation (i.e. the thermalisation of a non-integrable system) is supported by a few theoretical arguments such as the eigenstate thermalisation hypothesis [13,20] and it is compatible with a large number of simulations obtained with the most advanced numerical techniques, see e.g.[21][22][23][24] and references therein.
Conversely, for an integrable model, the proper statistical ensemble describing the system for long time is a generalised Gibbs ensemble (GGE) rather than a thermal one.The density matrix of the GGE is defined as [14] where I n is set of commuting integrals of motion, i.e.
[I n , I m ] = 0, and Z is a normalisation constant Z = Tr e − n λ n I n .Even in this case, the Lagrange multipliers λ n are not free parameters, but are fixed by the requirement that the integrals of motion assume the same value in the stationary state and in the initial one, i.e.
Tr[I n ρ GGE ] = Ψ 0 |I n |Ψ 0 for all I n .In this case, the infinite part A acts as a very peculiar bath on A since infinite information about the initial state is retained.
It is important for a proper definition of the GGE to specify which conserved charges enter in the GGE density matrix above.Indeed any quantum model has too many integrals of motion, regardless of its integrability.For example, all the projectors on the eigenstates O n = |E n E n |, are conserved for all Hamiltonians, but these cannot contribute to the stationary state in the most general case because otherwise (among the many strange consequences) no system will ever thermalise.Furthermore, taking all these projector in to account would be equivalent to timeaverage the complicated many-body dynamics, but this will result in a trivial prediction that has no connection with the economy of an ensemble description of statistical physics.Indeed such an ensemble would retain all of the information about the initial state, rather than information about only a minimal set of integrals of motion.To clarify this point, it has been understood recently [16,19,25] that only local integrals of motion should be used in Eq. ( 6) as long as we are interested in the expectation values of local observables such as the reduced density matrix (an integral of motion is said to be local if it can be written as an integral-sum in the case of a lattice model-of a given local current, in the spirit of Neither theorem).Thus, from now on, we will always refer to the GGE in Eq. ( 6) in which only local integrals of motion are included.
Finally, it is worth spending few words on what happens for finite systems.While the existence of the stationary state strictly relies on having taken the thermodynamic limit (as it is the case for everything in statistical physics), it is definitely important to understand how these concepts can be applied to cold-atomic experiments (or also nuclei) with only a few hundreds constituents.Denoting by v the typical velocity for the spreading of correlations through the system, it is evident that for times such that vt L (where L is the linear size of the finite system) the boundary conditions cannot affect any expectation value.Thus, for all practical purposes, the system can be considered as infinite.Consequently, as long as L is large enough to guarantee the existence of a time window such that the relaxation time of a given observable τ O is τ O L/v, then this observable will attain an almost stationary value also in a finite system.

Some simple examples: Quantum quenches in free field theories
Let us start by considering the simplest possible instance of a non-equilibrium quantum system which is an isolated harmonic oscillator, with Hamiltonian (we set m = 1) We consider the case in which the system is prepared for t < 0 in the ground state of a harmonic oscillators with a given frequency ω 0 and at time t = 0 the frequency is suddenly quenched to another value ω (this is commonly referred to as a quantum quench [26], because a Hamiltonian parameter has been suddenly changed).An elementary quantum mechanics calculation (either solving the Heisenberg equations of motion or playing with the ladder operators) leads to the time dependent mean square displacement Not surprisingly the harmonic oscillator oscillates forever fully deserving its name.However, from this example and without making any further calculation, we can construct a non-trivial case in which the local observables of a many-body system approach stationary values [26,27].Indeed it is enough to couple many harmonic oscillators with a nearest neighbour coupling to form a harmonic chain with Hamiltonian and with lattice spacing a.The variables π i and φ i satisfy standard commutation relations [φ i , φ j ] = [π i , π j ] = 0 and [φ i , π j ] = iδ i j .This is nothing but a lattice discretisation of a scalar free field theory with mass m described by Klein-Gordon equation.The Hamiltonian can be straightforwardly diagonalised in momentum space where we dropped all uninteresting additive constants.In this quadratic form, each mode is a harmonic oscillator with momentum dependent frequency and annihilation operator Let us now consider the case in which the harmonic chain is prepared in the ground state of the Hamiltonian (9) for a given values of the mass m 0 which at time t = 0 is quenched to another value m.Each momentum mode evolves according to Eq. ( 8) with frequencies given by Ω k and Ω 0k (i.e.Eq. ( 10) defined with m and m 0 respectively) and so it will never attain a stationary value.Conversely, the correlation function in real space (that for fixed distance r is a local observable) can be calculated by Fourier transform, obtaining φ r (t)φ 0 (t) − φ r (0)φ 0 (0) where we also subtracted the initial value only to have a more compact form.In the one-dimensional chain the subscript BZ just refer to the trivial Brillouin zone [−π/a, π, a], but the above formula is valid for a lattice in arbitrary dimension when the appropriate Brillouin zone is considered and most of the following considerations indeed remain valid in arbitrary dimension (see [27] for details).Eq. ( 12) is extremely instructive.By a simple stationary phase argument it is straightforward to see that the correlation function has a large time limit given by showing indeed that a stationary value is attained.But not only, it is also very simple to show that this stationary value is indeed the same as the one predicted by a generalised Gibbs ensemble ρ GGE ∝ e − k λ k n k constructed by using as integrals of motion the mode occupation number n k = a † k a k .For the elementary derivation of the GGE, we refer the interested reader to Ref. [27].Here we limit to stress that there is no contradiction between this result and what stated in the previous section due to the fact that the mode occupation operators are non-local.Indeed, one can prove [16,19] that despite being non local, the mode occupations are linear combinations of some local charges I m , thus it is always possible to rewrite k λ k n k = m γ m I m and, consequently, the GGEs built with n k and with I m are equivalent.
At this point it is natural to wonder how robust these prediction are under a change the initial condition, which in general we would like to be an arbitrary eigenstate of some Hamiltonian, even corresponding to a non-integrable model.Indeed, it is not difficult to think about very specific initial states for which a stationary state cannot exist.For example a linear superposition of a finite number of eigenstates of the post-quench Hamiltonian will provide persistent oscillation in time.Can we find some conditions for the initial state/Hamiltonian guaranteeing the existence of steady state (for local operators) which is described by a GGE?A simple general requirement has been recently shown for free theories [47]: a steady state exists and it is described by a GGE if the initial state satisfy the cluster decomposition property lim (14) Although this condition has not been generalised to interacting models (both integrable and not), it is reasonable to expect that cluster decomposition property must play an important role in the general non-equilibrium evolution of an isolated quantum system.

Some results about interacting integrable systems
The exact treatment of the non-equilibrium dynamics of interacting integrable models is immensely more complicated than the free theories, despite of the exact solvability of the models (i.e. an exact knowledge of the spectrum and of the eigenstates, but not often in a very practical form) and although a few methods have been developed specifically to tackle these out of equilibrium problems [48][49][50][51][52].
Due to these difficulties, many studies attempted the explicit exact construction of the GGE [53][54][55] starting from specific initial conditions.These GGE predictions could be (and have been in some cases) successively tested in numerical simulations or even in actual experiments.Before reviewing some recent results, let us briefly recall what are the needed building blocks to study the quench dynamics in a generic situation in such a way to grasp the technical difficulties one can encounter.The first problem one faces is to write the initial state |Ψ 0 in terms of the eigenstates of the Hamiltonian H governing the time evolution.Let us generically denote the normalised eigenstates as |n , then that the initial state can be written as where a n are the overlaps a n ≡ n|Ψ 0 between the initial state and the eigenstates.Consequently, the time evolved state is where E n is the energy of the state |n .This provides the time dependent expectation value of an arbitrary observable O, in terms of the form factors n|O|m , as Summing up, in order to characterise the non-equilibrium dynamics, the needed ingredients are (i) a complete characterisation of all eigenstates |n of a Hamiltonian and their energies; (ii) the norms of the eigenstates and the form factors of relevant operators in this basis; (iii) the overlaps between the initial states and the eigenstates.

EPJ Web of Conferences
For integrable models, the Bethe Ansatz is a very efficient tool to obtain most of these ingredients since it provides a full set of eigenstates (called Bethe states) with their energies, the norms and the form factors of the most relevant local operators [56].What is not (yet) known in general is how to obtain the overlaps between Bethe states and generic initial states.Up to now, only very few exact results exist for these overlaps with very specific initial states [57][58][59][60][61][62][63].Clearly, finding compact and tractable expressions for the overlaps between Bethe states and more generic initial states would allow exact calculations for a variety of potentially interesting situations and experiments.
However, the ingredients listed above are only the starting point for the description of the quench dynamics, because the sum ( 17) is still to be performed and this is a very difficult step.Fortunately, a recently-proposed method (termed either representative state approach or quench action formalism) gives an exact analytical description of the post-quench steady state in the thermodynamic limit [52].The essential building blocks of this method are once again the overlaps between the initial states and Bethe states.
There are two paradigmatic integrable models that have been (mostly) considered in the literature.The first one is the so called Lieb-Liniger model, a onedimensional Bose gas with pairwise delta interaction on a ring of circumference L with periodic boundary conditions (PBC), i.e. with second quantised Hamiltonian [64] where φ(x) is a canonical bosonic field, c the coupling constant and we set = 2m = 1.The second model is the anisotropic Heisenberg spin chain (also called XXZ spinchain) with Hamiltonian where σ α are the Pauli matrices at the site and periodic boundary conditions σ α L+1 = σ α 1 are always assumed.The two models present pros and cons.Indeed, the Lieb-Liniger model has the advantage to display a simpler Bethe ansatz solution, but the XXZ spin-chain can be straightforwardly simulated on a computer because it is defined on the lattice.
The non-equilibrium dynamics of the Lieb-Liniger model starting from the non interacting ground state (i.e.c = 0) has been firstly studied in Ref. [53] where the GGE was approximately constructed with an ingenious integrable lattice regularisation of the model.Subsequently the overlaps have been analytically derived [63], allowing (thanks to the quench action method of Ref. [52]) the explicit evaluation of the long time limit of many observables [63] whose values agree with the previously constructed GGE.It is worth mentioning that the quench from c = 0 to c = ∞ has been treated with more elementary techniques [46], whose results are reproduced by the proper limit of the general construction [63].A typical example of the comparison between GGE and numerical simulations.In the figure the data for the correlation σ x j σ x j+k at distance k = 1, 2, 3 are reported as function of time for a quench starting in the Neel state and evolving with the XXZ Hamiltonian (19) with Δ = 4.The numerical data (full lines) clearly attain stationary values which are compatible with GGE (dashed lines), but not with the thermal ones (dotted lines).The figure is extracted from Ref. [55] where many more initial states (dimer, tilted ferromagnets, tilted Neel and a few more) and final Hamiltonians for different values of Δ > 1 have been studied, finding similar results.
The anisotropic Heisenberg spin chain has a similar story, but with a different end.Indeed explicit GGE predictions (based on all the known local charges) have been derived in Refs.[54,55] for quenches starting from many initial product states such as Neel states, dimer states, tilt ferromagnetic state and a few more.For some technical reason, it has been possible to obtain GGE predictions only in the gapped antiferromagnetic phase with Δ > 1.In Ref. [55], these analytic predictions have been compared with numerical simulations (based on density matrix renormalisation group and time evolving block decimation algorithms) finding a very good agreement in all cases.In Fig. 1 we report a typical example of comparison between the GGE predictions and the time dependent numerical simulations (obtained by means tensor network algorithms, see [55] for details).
However, this is not the end of the story, because the overlaps of the Neel and dimer states with the general Bethe states were finally derived in Ref. [59] (in a form useful for the computation because they were implicitly known already from Ref. [58]).From these overlaps, the infinite time limit of the evolution has been explicitly calculated for the Neel state in Ref. [65,66] and for the dimer state in [67], again exploiting the representative state approach of Ref. [52].It has been very surprising to discover that these exact results differ from the previously known GGE of Ref. [55] (reported in Fig. 1).However, in most XLIV International Symposium on Multiparticle Dynamics (ISMD 2014) 08001-p.5 cases, the difference between exact and GGE results is extremely small and it is impossible to observe it in any numerical simulation, as for example in the case of a quench from the Neel state.Just to quote an example the nearestneighbour correlation function (i.e. the one reported in Fig. 1) in the GGE and for infinite time can be expanded close to the Δ = ∞ point obtaining [66] 21) with the first three orders being identical.The difference is more pronounced for the quench starting from the dimer state, although it remains relatively small.However, in this case it is possible to observe the difference between the two predictions [67], and the numerical simulation confirm indeed that the GGE built is Ref. [55] yields wrong prediction, while the representative state approach correctly reproduce the numerical results.
At this point we should critically ask what is going wrong in the GGE construction (several checks have been performed excluding categorically computational mistakes).This issue has been considered already in a few manuscripts appeared in the last months [69][70][71][72][73], but a conclusive answer has not yet been provided.The most probable scenario concerns the presence of additional charges which are not yet explicitly known and that do not appear in the standard GGE construction [55].For example, it is well known that in the gapless regime (|Δ| < 1) and in the presence of a current, there are additional quasi-local (i.e. with a kernel decaying with the distance and not with support on a finite interval as the standard ones) which must be added to the GGE [71].These charges however do not explain the GGE failure in the gapped regime.It has been only recently pointed out that in the continuum quantum field theory there are additional charges which are not ultra-local [73], but that anyhow are local enough to be included in the construction of the GGE.The lattice counterpart of these integrals of motion are still unknown.

Conclusions
We have briefly reviewed here some of the recent results about the non-equilibrium dynamics of isolated systems showing that they currently represent a theoretical and experimental challenge raising many fundamental questions in many-body quantum mechanics.We have discussed mainly the issue of existence and description of a stationary state starting from a given initial state.However, this is only a small part of a very broad research on non equilibrium dynamics of isolated quantum systems.
A very important issue which has not been covered here is how these stationary values are approached.In Refs.[26,27] a light cone spreading has been proposed according to which correlations are generated by quasiparticles, entangled over regions of the order of the correlation length in the initial state, which after the quench propagate classically through the system.This idea has been confirmed by a recent experimental investigation [7].Another timely issue concern what happen in a system when integrability is slightly broken (as in the case of the famous quantum Newton cradle, cf.Sec. 2).It has been shown that for very long time scales (depending on the integrability breaking parameter), the system is trapped in a metastable state which can be effectively described by a deformed GGE constructed with the perturbed conserved charges of the integrable model [74].Only for longer time the system is expected to display the asymptotic thermal behaviour, but how this crossover takes place is not yet clear.Oppositely, when integrability is strongly broken the approach to the thermal values appear to be very fast.

Figure 1 .
Figure1.A typical example of the comparison between GGE and numerical simulations.In the figure the data for the correlation σ x j σ x j+k at distance k = 1, 2, 3 are reported as function of time for a quench starting in the Neel state and evolving with the XXZ Hamiltonian(19) with Δ = 4.The numerical data (full lines) clearly attain stationary values which are compatible with GGE (dashed lines), but not with the thermal ones (dotted lines).The figure is extracted from Ref.[55] where many more initial states (dimer, tilted ferromagnets, tilted Neel and a few more) and final Hamiltonians for different values of Δ > 1 have been studied, finding similar results.