## Abstract

We present a semiclassical approach to many-body quantum propagation in terms of coherent sums over quantum amplitudes associated with the solutions of corresponding classical nonlinear wave equations. This approach adequately describes interference effects in the many-body space of interacting bosonic systems. The main quantity of interest, the transition amplitude between Fock states when the dynamics is driven by both single-particle contributions and many-body interactions of similar magnitude, is non-perturbatively constructed in the spirit of Gutzwiller’s derivation of the van Vleck propagator from the path integral representation of the time evolution operator, but lifted to the space of symmetrized many-body states. Effects beyond mean-field, here representing the classical limit of the theory, are semiclassically described by means of interfering amplitudes where the action and stability of the classical solutions enter. In this way, a genuinely many-body echo phenomenon, coherent backscattering in Fock space, is presented arising due to coherent quantum interference between classical solutions related by time reversal.

## 1. Introduction

Interference effects, as a hallmark of quantum behaviour, are so essential for the description of the fundamental physical processes around us that the founding fathers of quantum mechanics included interference already at the kinematic level of modern quantum theory in the form of a central postulate: the state-space structure of quantum mechanics is that of a vector space, and, therefore, the sum of two possible states represents also a valid state. This kinematic concept extends its implications into the domain of unitary dynamics in isolated quantum systems, which is postulated to be represented by linear operators acting on state vectors and, therefore, allowing for the coherent superposition of different histories when describing a single physical process. Given the fragile nature of quantum coherent effects, much has been learned from the formal correspondence between classical field theory, in particular from Maxwell’s electrodynamics, and the solutions of both the stationary and time-dependent single-particle Schrödinger equation describing (predominantly) the orbital degrees of freedom of a particle in real space. Many of the benchmark coherent effects, now routinely studied in the framework of quantum transport, were, for instance, first observed within a purely classical wave context. In this respect, echo phenomena were very similar: there is a perfect correspondence between classical wave echoes and quantum wave mechanical signatures of coherent superpositions that reveal the information about the initial state.

In exact analogy with classical wave systems where ray dynamics serves as a backbone supporting wave propagation in the limit of short wavelengths, wave mechanics of single-particle quantum systems can be described in the limit of large classical actions (compared with ) by an asymptotic version of quantum mechanics where only classical information is required (apart from ). This semiclassical programme, initiated by Bohr and cast in its elaborate, modern form by Gutzwiller [1], accomplishes this idea to its ultimate limits: not only the average incoherent aspects of a quantum process are given by their classical expectations (in accordance with the correspondence principle) but also the quantum fluctuations around them can be unambiguously calculated using only classical information. The seminal contributions of van Vleck, Feynman and Gutzwiller consist in the explicit construction of the machinery transforming classical information into quantum mechanical predictions. The ultimate success of the theory, based in particular on the underlying propagator, is reflected in the physical understanding and quantitative analysis of quantum interference phenomena in terms of classical solutions, their actions and their stabilities.

The fruitful analogy between quantum amplitudes and classical fields cannot be simply lifted to the many-body level where an effective single-particle description ceases to be valid. For interacting many-body systems the quantum state is by definition a high-dimensional object: many-body interference takes place in the space of quantum states without a real-space analogue. Therefore, the obvious question arises whether it is possible to export the highly successful semiclassical techniques into the realm of interacting quantum systems. More precise questions concern the very existence of a semiclassical limit and the technical issue of the implementation of a Gutzwiller approach in many-body space.

These sorts of questions might have remained academic if the experimental studies of interference effects, where semiclassics has proved its value, had not advanced so rapidly. The successful realization of genuine many-body coherent effects such as Rabi oscillations and Anderson localization in cold atom quantum gases and the first steps towards Fock space tomography of cold atoms in optical lattices [2,3] triggered theoretical methods able to capture many-body interference phenomena in their characteristic non-perturbative form.

As the classical limit of systems of interacting bosonic cold atoms, representing discrete bosonic quantum fields, has been known for some time in the form of a nonlinear classical field equation, one might think that a semiclassical approach could have been easily implemented. However, in the context of interacting bosons this Gross–Pitaevskii equation is derived from a variational approach describing the ground state of the system [4]; although within this regime it has enjoyed remarkable success, it has been implicitly assumed that it cannot by construction deal with full-quantum effects and that it is restricted to ground-state properties. Here, we show that the same mean-field equation emerges as the rigorous limit of a quantum theory of interacting bosons when the latter is analysed semiclassically within a field-theoretical formalism. In this approach, one finds naturally that the classical limit corresponds to , where *N* is the total number of particles, thus providing the connection to the standard derivations of the mean-field equations. Raising the Gross–Pitaevskii equation from a mean-field picture to the notion of a classical limit makes it amenable to semiclassical approaches and thereby allows classical information encoded in the mean-field equations to be used in an alternative way to describe both many-body interference and quantum effects of excited states beyond the mean-field description. In particular, the whole quantum chaos machinery that enables us to study the quantum manifestations of classical Hamiltonian structures such as tori, stability islands and bifurcations, and the emergence of universal features in the quantum description due to instability of the classical motion [1], can now, in principle, be applied to the many-body domain.

Here we provide a detailed account of such a semiclassical theory of interacting bosonic fields, where quantum indistinguishability and inter-particle interactions are accounted for. This approach, which was achieved only recently [5], follows and extends Gutzwiller’s work [1] based on the asymptotic analysis of the exact path integral representation of the quantum propagator within a suitable basis. Equipped with a van Vleck–Gutzwiller-type propagator in Fock space, we will apply this approach to the calculation of an echo effect due to the coherent interference between quantum amplitudes associated with different classical mean-field solutions related by time-reversal symmetry. We dubbed this novel many-body echo effect the *coherent backscattering (CBS) in Fock space*.

Echoes serve as prominent probes for unitary quantum evolution and in particular deviations from it, with spin [6] and polarization [7] echoes in nuclear magnetic resonance as ground-breaking examples. The generalization of these concepts, originally devoted to spin-like properties, to the quantum time evolution and in particular the study of its sensitivity with respect to perturbations has led to the closely related notions of fidelity [8] and the Loschmidt echo [9] and a whole corresponding research field (for reviews see [10–12]). This broad interest in echo phenomena was partly triggered by the pioneering work of Jalabert & Pastawski [13], who considered the Loschmidt echo to be related to (the overlap of) the spatial degrees of freedom of wave functions in wave-chaotic single-particle systems employing a semiclassical approximation of the propagator. Here we raise such a semiclassical method to the level of propagation in many-body Fock space; however, we do not consider an echo upon time inversion as in the case above, but instead a many-body generalization of dynamical echoes [14,15] such as CBS.

To this end, we provide a comprehensive and detailed derivation (only briefly outlined in [5]) of the semiclassical propagator in Fock space that describes, for example, bosonic systems on a lattice where the number of single-particle states is finite. The paper is structured as follows: in §2, we introduce so-called quadrature eigenstates, which are used to first derive a path integral and then a semiclassical approximation for the propagator. After that, we demonstrate how to perform a basis change from quadrature to Fock states in §3. Finally, the resulting semiclassical propagator is used to derive CBS in Fock space in §4, before we conclude in §5.

## 2. The propagator in quadrature representation

### (a) Quadrature states

In this section, we will introduce the quadrature eigenstates and their most important properties without any proof. For more details we refer the reader to, for example, [16].

Quadrature eigenstates are states in Fock space of a many-body system that is built from a given set of single-particle states labelled by (a finite number of) integers *l*=1,2,3,…,*L*, which can be, for example, eigenstates of the corresponding single-particle problem or individual sites (wells) of an optical lattice. The *l*th creation operator is then defined to create a particle in state , while the *l*th annihilation operator annihilates it again, such that a Fock state , with *n*_{l}=0,1,2,… being the occupation number of state , satisfies
2.1*a*and
2.1*b*

Equation (2.1) implies the commutation relations
2.2*a*and
2.2*b*as well as
2.3and
2.4where is the vacuum, i.e. the state of the empty system.

The quadrature operators and are finally deduced from the creation and annihilation operators according to
2.5where *b* is for now an arbitrary, but real, parameter. In the following their eigenstates
2.6will be referred to as quadrature (eigen)states. These quadrature eigenstates can be expanded in Fock states by using their overlaps
2.7*a*and
2.7*b*where *H*_{n} denotes the *n*th Hermite polynomial. For large occupation numbers, *n*_{l}≫1, these overlaps can in the oscillatory region, , be approximated by the Wentzel-Kramers-Brillouin (WKB) approximation [17],
2.8where *F*(*q*,*n*) is the generating function for the canonical transformation from (*q*,*p*) to (*n*,*θ*) with and and is given by
2.9It is worth noting that the generating function can also be written as an integral,
2.10

where . In fact, this is the same generating function as the one for the canonical transformation from position and momentum to action angle variables, showing that the quadrature states are formally equivalent to the position and to the momentum of a particle in a harmonic potential.

As a generating function of a canonical transformation, the derivatives of *F*(*q*,*n*) satisfy
2.11*a*and
2.11*b*Note that, by this definition of *θ*, the phase is only defined between 0 and *π*, and therefore |*p*|=*p*.

The overlaps of different quadrature eigenstates are given by
2.12*a*
2.12*b*
2.12*c*Finally, the resolution of unity in terms of quadrature eigenstates is given by
2.13Thus, the quadrature states are a proper choice for the construction of a path integral, because the resolution of unity is given by an *integral* rather than a sum. Moreover, for the semiclassical approximation following this construction, they seem to be a better choice than coherent states, because their eigenvalues are real and therefore the complexification necessary for the coherent state propagator can be avoided.

### (b) The path integral

In quadrature representation, the exact quantum mechanical propagator is given by the matrix element
2.14with
2.15being the quantum mechanical evolution operator for a many-body Hamiltonian , which may be, for example, a typical Bose–Hubbard Hamiltonian of the form
2.16Here, the first term can be considered as the single-particle part of the Hamiltonian while the second term is the most general form of possible two-body interactions. The matrix *h* occurring in the single-particle term, in general, has the energies of the ground-lying single-particle states *ϵ*_{l}=*h*_{ll}, which later on will also be called on-site energies, on its diagonal and describes transitions between single-particle states (later the term *hopping* will be used) via its off-diagonal entries *h*_{l,l′≠l}.

The path integral for equation (2.15) is constructed by first using Trotter’s formula [18],
2.17where *τ*=(*t*_{f}−*t*_{i})/*M*, and inserting one unit operator in terms of quadratures *q* and *p*, equation (2.13), for each time step *τ*. Following this procedure, however, as already noted in the first quantized case for charged particles in a magnetic field by Feynman in his original derivation of the path integral [19] and reviewed in detail in [20], one has to be careful with the terms of the Hamiltonian containing products of and operators when letting them act on a or state, respectively (see also appendix A). After performing all these steps, as is done in appendix A, in the continuous time limit the path integral for the propagator in quadrature representation is given by [21,22]
2.18with ** ψ**(

*t*)=(

**q**(

*t*)+

*i*

**p**(

*t*))/(2

*b*) and the boundary conditions

**q**(

*t*

_{i})=

**q**

^{(i)}as well as

**q**(

*t*

_{f})=

**q**

^{(f)}. Here, the classical Hamiltonian is defined by 2.19which for a Bose–Hubbard Hamiltonian of the form (2.16) can be evaluated as 2.20with 2.21Thus, the classical Hamiltonian equation (2.20) is determined by the quantum to classical replacement rule 2.22Note that the classical Hamiltonian also contains constant terms, which give rise to a constant phase in the propagator. Such phases have been found for the semiclassical coherent state propagator [23–30], and have been denoted as the Solari–Kochetov extra-phase. It should also be noted that the classical Hamiltonian depends on the chosen ordering of the propagator.

### (c) Semiclassical approximation

The semiclassical approximation to the quantum propagator in quadrature representation is obtained by applying the stationary phase approximation to equation (2.18). This is to expand the exponent around the stationary points up to second order in **q**^{(m)} and **p**^{(m)}, and then perform the resulting Gaussian integral.

The stationarity conditions in **q**^{(m)} and **p**^{(m)} are equivalent to those in *ψ*^{(m)} and *ψ*^{(m)}*. In the time-continuous limit , these are given by the equations of motion
2.23These equations of motion are the time-dependent mean-field equations corresponding to the full-quantum Hamiltonian, which for a Bose–Hubbard Hamiltonian, equation (2.16), are determined by the mean-field Hamiltonian equation (2.19). The apparent complex character of these equations must be understood as a convenient, compact way of representing the set of real Hamilton equations
2.24*a*and
2.24*b*From the point of view of dynamical systems, equations (2.23) and (2.24) admit the whole zoo of effects characteristic of Hamiltonian systems, from integrability to chaos passing through mixed-phase space [31–33]. It is easy to check that, if the quantum Hamiltonian commutes with the number operator,
2.25Equation (2.23) preserves the (classical) number of particles
2.26and therefore the system is not only Hamiltonian but also automatically bounded. The symmetry of the classical Hamiltonian corresponding to this conserved quantity is the *U*(1) gauge symmetry, i.e. the invariance of the Hamiltonian when the global phase of the mean-field wave function ** ψ**(

*t*) is changed, 2.27Moreover, these equations of motion, equation (2.23), are generally the ingredients for the so-called truncated Wigner method [34–41].

Finally, in the case of constant on-site interactions, *U*_{ll′l′′l′′′}=*Uδ*_{ll′}*δ*_{ll′′}*δ*_{ll′′′}, the equations of motion are (up to renormalization of the single-particle energies due to ordering) given by the discrete Gross–Pitaevskii equations [4].

Thus, the semiclassical propagator in quadrature representation will be given by a sum over mean-field wave functions, whose dynamics is described by the time-dependent mean-field equation (2.23), and is subject to the following boundary conditions on the real parts of ** ψ**(

*t*) in time: 2.28The thus obtained shooting problem has, in general,

*several*solutions, which will be labelled by an index

*γ*and will—in accordance with the notation for the van Vleck propagator in configuration space—be denoted as

*trajectories*, although they are not trajectories in the conventional sense. Still thinking in terms of ‘usual trajectories’ may help when evaluating expressions containing the semiclassical propagator.

We note that, contrary to the semiclassical propagators derived for coherent states [23], where the boundary conditions for the classical trajectories fix both the real and imaginary parts of the mean-field wave function at initial and final time, and thus the shooting problem is overdetermined, no complexification of ** ψ** is necessary in the approach presented here, i.e.

*** is still the complex conjugate of**

*ψ***.**

*ψ*Evaluating the exponent along the stationary points will yield the phase that the trajectory contributes to the propagator and is given by the classical action of *γ*,
2.29

where ** ψ**(

*t*)=(

**q**(

*t*)+

*i*

**p**(

*t*))/(2

*b*) is determined by the solution of the equations of motion, (2.23), under the boundary conditions (2.28), such that the semiclassical propagator is given by 2.30The semiclassical amplitude finally is determined by the integration over quantum fluctuations in the path integral around the stationary (classical) trajectories. The evaluation of these integrals, which is performed in appendix B, is very similar to those presented for the semiclassical coherent state propagator in [23,24] and yields (see appendix B) 2.31and thus the semiclassical propagator in quadrature representation is given by 2.32

It is important to note that here **q**(*t*) and **p**(*t*) for the classical trajectories—in contrast with the coherent state propagators derived in [23–26,42]—are still real variables, because only the initial and final *real parts* are fixed by the boundary conditions, i.e. in the approach presented here, no complexification is necessary.

Finally, we state the derivatives of the action with respect to the boundary conditions, which can be calculated by using the equations of motion and are given by
2.33Furthermore for a time-independent Hamiltonian, the time derivative of the action is given by
2.34where *E*_{γ}=*H*^{(cl)}(** ψ***(

*t*

_{i}),

**(**

*ψ**t*

_{i})) is the conserved energy of the trajectory.

## 3. Basis change to Fock states

Starting from the semiclassical propagator in quadrature representation, we can now also derive a semiclassical expression for the propagator in Fock representation according to
3.1by using the asymptotic form equation (2.8) for the overlap of a quadrature with a Fock state, neglecting the contributions to the integral from ,^{1} and evaluating the integrals in stationary phase approximation.

Using equations (2.11a) and (2.33), it is straightforward to find that the stationarity conditions are given by
3.2*a*and
3.2*b*Note that the two possible signs of **p**(*t*_{i}) and **p**(*t*_{f}), which are not encoded in the generating function itself, are obtained from the two possible signs of the generating function when expressing the cosine in equation (2.8) by exponentials.

However, it turns out that the boundary conditions (3.2) give rise to continuous families of trajectories originating from the *U*(1)-gauge invariance, i.e. from the fact that if ** ψ**(

*t*) solves the equations of motion under the boundary conditions (3.2) then

**(**

*ψ**t*)e

^{iθ}also does. Therefore, one of the 2

*L*integrations, e.g. the one over , cannot be performed in stationary phase approximation.

In fact, it turns out to be beneficial to transform the integrations over into integrations over *θ*^{(i/f)}_{l} defined by
3.3*a*
3.3*b*
3.3*c*Then one can perform all the integrations except the one over *θ*^{(i)}_{1} in stationary phase approximation, yielding again the stationarity conditions (3.2) with
3.4being automatically fulfilled due to conservation of particles. Thus, these integrations select trajectories *γ* satisfying these boundary conditions and additionally
3.5

The evaluation of the exponent at the stationary points yields the classical action of the resulting classical trajectories
3.6where *θ*_{l}(*t*) and *n*_{l}(*t*) are the real variables defined by
3.7It is easy to see that—again due to conservation of particles—the classical action (3.6) is independent of the initial global phase .

Furthermore, the semiclassical prefactor, which is derived in appendix C, turns out to be independent of , too, and can be written as 3.8

where the prime at the determinant indicates that the determinant is taken on the (*L*−1)×(*L*−1) sub-matrix obtained by removing the first line and first row of the full matrix.

The integration over *θ*^{(i)}_{1} is therefore trivial and the final semiclassical expression for the propagator in the Fock basis is given by
3.9where the sum runs now over the *continuous families* of trajectories joining the initial and final occupations.

## 4. Coherent backscattering

A natural application of the semiclassical propagator is the transition probability from an initial Fock state **n**^{(i)} to a final one **n**^{(f)},
4.1i.e. the probability of measuring the Fock state **n**^{(f)} at time *t*_{f}, if the system was prepared in the Fock state **n**^{(i)} at the initial time *t*_{i}.

Inserting the semiclassical propagator (3.9) yields for the transition probability a double sum over the trajectories
4.2Under the disorder average, i.e. an average over the on-site energies *ϵ*_{l}=*h*_{ll}, due to the scaling of the action with the total number of particles, the action difference in the exponential gives rise to strong oscillations. Most contributions to the double sum cancel out on the disorder average, unless the two actions for trajectories *γ* and *γ*′ are correlated.

### (a) Diagonal approximation

The most obvious pairs of trajectories, which are correlated, are the identical ones, i.e. *γ*=*γ*′ (figure 1). For those pairs, the action difference is obviously zero and the transition probability in diagonal approximation reduces to
4.3where ** θ**(

*t*

_{i}) is the vector containing the initial phases for the trajectory

*γ*.

By inserting an integration over the initial phases together with a Dirac delta, ensuring that only those phases contribute which match the ones given by the trajectories and using the properties of the delta function, one can prove the sum rule
4.4with ** ψ**(

**n**

^{(i)};

*θ*^{(i)};

*t*

_{f}) being the classical time evolution of the initial state

*ψ*^{(i)}with components 4.5

Therefore, the diagonal approximation yields the classical transition probability 4.6which in chaotic systems can be expected to be for large enough propagation times roughly constant as a function of the final occupations, as long as their interaction energies, 4.7are roughly the same. That this is in fact the case has also been observed numerically in [5].

### (b) Coherent backscattering contribution

A second class of correlated pairs of trajectories are those shown in figure 2, namely , where denotes time reversal. For Bose–Hubbard systems, time reversal is achieved by complex conjugation, such that the time reversal of a solution ** ϕ**(

*t*+

*t*

_{i}) of the equations of motion is

**(**

*ϕ**t*

_{f}−

*t*)*. Of course, these pairs exist only if the Hamiltonian itself is time-reversal symmetric, which for a Bose–Hubbard system corresponds to a real Hamiltonian. It turns out that a Bose–Hubbard chain with nearest neighbour hopping only is always time-reversal symmetric, because it can be transformed into a real one by multiplying each component of the wave function by a certain time-independent phase. In the same way, for certain combinations of the phases of the hopping parameter in a Bose–Hubbard ring, the Hamiltonian can be transformed into a real one.

Having identified the time reverse of a trajectory *γ*, it is straightforward to check that the actions of both trajectories are the same,
4.8However, in (4.2) pairing *γ* with its time reverse is possible only if **n**^{(f)}=**n**^{(i)}. Therefore, one first has to replace *γ* and *γ*′ by nearby trajectories and , which both start and end at **n**^{(0)}=(**n**^{(i)}+**n**^{(f)})/2. In the course of this replacement, one also has to expand the actions in **n**^{(f)} and **n**^{(i)} up to linear order around **n**^{(0)},
4.9In the prefactor one can simply replace the second derivative of the action of the original trajectory by the corresponding one of the new trajectory.

If is the time reverse of , the initial and final phases of are related to those of by
4.10and
4.11which, after applying the sum rule (4.4), yields for the CBS contribution to the transition probability
4.12where *δ*_{TRI} is one if the system is time-reversal invariant and zero otherwise.

For a classically chaotic system, within a disorder average one can replace the final phases ** θ**(

**n**

^{(i)},

*θ*^{(i)};

*t*

_{f}) by identically independent in [0,2

*π*] uniformly distributed random variables, such that for

*P*

^{(cbs)}one can perform an average over the final phases yielding a factor . Thus, the averaged CBS contribution to the transition probability is given by 4.13

As long as time-reversal symmetry is the only discrete symmetry of the system, pairing a trajectory with itself or its time reverse are the only two generic possibilities to get a contribution with vanishing action difference. All further possible pairs of trajectories give contributions which are suppressed by a factor of 1/*N* compared with the latter two. In the semiclassical limit, the averaged transition probability is then given by
4.14i.e. the quantum transition probability is enhanced for **n**^{(f)}=**n**^{(i)} compared with the classical one by a factor of two due to the constructive interference between time-reversed paths.

However, only paths which are not self-retracing may contribute to the CBS peak, because the time reverse of a self-retracing path is again the path itself and therefore this pairing would already be included in the diagonal approximation. For short times, one can expect that all trajectories coming back to the initial occupations are all self-retracing and therefore a certain minimal time is required to observe CBS (see supplementary material of [5]). This behaviour as well as the existence of CBS itself has been shown numerically in [5].

We remark that, by its very definition, CBS in Fock space requires a dephasing mechanism (in our case, the combined effect of dynamical instability and average over disorder realizations), and it is therefore a robust effect against external perturbations as long as they do not produce decoherence. The existence of coherent but robust interference effects in many-body systems is all the more remarkable considering the extreme delicacy of the associated quantum coherences. Also, the prediction and numerical confirmation of this dynamical echo effect has deep implications for the very active field of thermalization of closed, interacting quantum systems (see [43] for a recent review). There, the main question is whether local observables of many-body systems reach an equilibrium value in agreement with the expectations from quantum statistical physics. Our findings show that, even under strong ensemble average (which is expected to speed up thermalization) if the system preserves time-reversal symmetry, the memory of the many-body initial state remains for arbitrarily long times. Although this observation is not in contrast with thermalization as understood in [43], as the effect of CBS on the expectation values of *local* observables is a 1/*N* effect, it indicates that time reversal may play an unexpected role in the thermal equilibration of many-body systems.

## 5. Conclusion

In this paper, we have rigorously constructed an approximation for the quantum mechanical transition amplitude connecting two many-body states of a lattice system of interacting bosons occupying a finite set of single-particle orbitals. We applied this method to predict an echo-type backscattering effect due to quantum interference in Fock space. Our approximation is semiclassical in spirit, in the sense that it relies on a stationary phase analysis of the exact path integral representation of the quantum mechanical propagator, but corresponds to an asymptotic, classical limit that is different from the usual limit . Rather, we show that, when written in the language of fields, an interacting bosonic system has as its classical limit a nonlinear classical field theory, namely the mean-field Gross–Pitaevskii equation. The semiclassical method then unambiguously prescribes how to use the solutions of this equation to describe quantum processes in terms of the actions and stabilities of the classical solutions, and in this way we find that the mean-field limit, commonly used only in the context of ground-state properties due to its usual derivation from a variational approach, actually contains all the information one needs to understand quantum mechanical effects for excited states, including quantum fluctuations and interference. This conceptual step is one of the main results of this paper.

Once the semiclassical propagator is derived and its domain of validity is clarified, we turned to its application to give a semiclassical picture of many-body interference in Fock space. Coherent sums over semiclassical amplitudes then lead to constructive interference of the transition probability that is not captured in any mean-field level of approximation. Following the steps learned from the study of quantum interference in mesoscopic systems, we then calculated the Fock-space version of the mesoscopic dynamical echo, and found a universal enhancement of the return probability due to the interplay between interference and interactions. This prediction, successfully confirmed in numerical calculations presented in [5], shows the power of the semiclassical thinking in many-body systems, a road that has been shown to be useful for other set-ups like the 1-site Bose–Hubbard system [44], the WKB approach for the 2-site case of [45–49], bosonic transport in optical lattices [50], connecting soliton-like solutions of discrete nonlinear equations with properties of the quantum spectra [51], investigating spectral statistics of fully chaotic systems [52–55] and going beyond the classical truncated Wigner method to describe dynamical processes in many-body systems [56]. While the extension of our methods to non-relativistic Fermionic fields on the lattice (partially initiated in [57]) is the subject of work in progress, the ultimate goal of addressing fully fledged quantum fields in the continuum, or equivalently an infinite number of possible single-particle orbitals, will face fundamental aspects of such systems in the presence of interactions^{2} like renormalization issues, well beyond the present tools of semiclassical analysis. However, non-relativistic bounded one-dimensional systems such as interacting bosons with delta interactions in a ring (the Lieb–Liniger model) are of high interest because of their experimental realization in cold atom systems; exploring the formal construction of the semiclassical approximation there is a realistic goal.

## Authors' contributions

T.E. derived the path integral as well as the semiclassical approximation for the propagator and drafted the manuscript. J.D.U. computed the transition probability using the semiclassical propagator in Fock representation. J.D.U. and K.R. together conceived the study, supervised and guided the work. All authors worked on the manuscript.

## Competing interests

The authors declare that they have no competing interests.

## Funding

This work was financially supported by the Deutsche Forschungsgemeinschaft within the research unit FOR 760.

## Acknowledgements

The authors particularly thank Peter Schlagheck for attendance, cooperation and numerous useful conversations and discussions during the whole course of this project.

## Appendix A. The path integral in quadrature representation

In this appendix, the derivation of the path integral, equation (2.18), is presented starting from the Trotter formula, equation (2.17). The first step is to introduce resolutions of unities in terms of quadrature states, equation (2.13), after each factor in the product of equation (2.17), A 1

In order to evaluate the matrix elements of the right-hand side of equation (.1), a Hubbard–Stratonovich transformation has to be performed,
A 2where the integration runs over Hermitian matrices *σ* and (*U*^{−1})_{klmn} is chosen such that . This transformation can be proved by shifting and using the relations and *U*_{klmn}=*U*_{lkmn}=*U*_{klnm} induced by the Hermiticity of the Hamiltonian (2.16).

Applying equation (.2) to the matrix elements in equation (.1) for a Hamiltonian of the form equation (2.16) yields
A 3with *h*_{σ} being the effective single-particle Hamiltonian matrix with entries given by
A 4Next, after replacing and by and , one has to replace the exponential containing and operators by a product of exponentials, which contain only or only operators. The product has to be arranged such that the operators can be applied directly to the quadrature states at the right or at the left. To this end recognize that
A 5Note here that is antisymmetric and thus its diagonal is zero.

The last term in equation (.5) still contains terms . Thus before splitting up the exponentials, we need to get rid of these terms. This is possible by means of a method known as ‘uncompleting the square’ or ‘Gaussian trick’,
A 6When splitting up the exponential, one should recognize that not only the last term in equation (.6) is proportional to , but due to the Gaussian weight of the integration . Thus by applying
A 7terms proportional to *τ* are correctly taken into account and yields for the matrix elements in equation (.3) after also introducing a unit operator in terms of
A 8

Now, the operators can simply be replaced by and the operators in the most left exponential by , while those in the most right exponential can be replaced by . After that, first the integration over **u** and finally the one over *σ* can be carried out, thus undoing the previously performed transformations, yielding
A 9

where *ψ*^{(m)}=(**q**^{(m)}+*i***p**^{(m)})/(2*b*), and
A 10

Note that in the continuous limit, *τ*→0, and therefore
A 11where *H*^{(cl)} is given by equation (2.20).

## Appendix B. The semiclassical prefactor in quadrature representation

In order to derive the semiclassical prefactor, equation (2.31), according to the stationary phase approximation, the second variation of the exponent around the classical trajectory is kept for the integration yielding Gaussian integrals,
B 1which after substituting **x**^{(m)}=2*b***x**^{(m)}, **y**^{(m)}=2*b***y**^{(m)} and introducing the short-hand notation reads
B 2where is the *L*×*L* unit matrix and the result has to be evaluated at **x**^{(M)}=0.

For the evaluation of these integrals we define the matrices *V* ^{(m)},*W*^{(m)} and *Y* ^{(m)}, such that, after integrating **x**^{(1)},**y**^{(1)},…**x**^{(m′−1)},**y**^{(m′−1)}, the semiclassical amplitude is given by
B 3Obviously, *V* ^{(m)}, *W*^{(m)} and *Y* ^{(m)} fulfil the initial conditions
B 4Performing the integrals over **x**^{(m′)} and **y**^{(m′)} in equation (B 3) and using the inverse of a matrix in block form,
B 5then shows that also the conditions
B 6hold with .

Defining and keeping only terms up to linear order in *τ*, these conditions can be written as
B 7with the initial condition
B 8In the continuous limit, *τ*→0, the recursion relation for *X* gets transformed into a first-order differential equation,
B 9with initial condition *X*(*t*_{i})=0 and is solved by
B 10Therefore, the integral over the fluctuations yields
B 11

Using the differential equation for ∂**q**/∂**p**(*t*_{i}), the exponent can be written as
B 12where in the second steps the equations of motion for ∂**p**(*t*)/∂**p**(*t*_{i}) have been used. Thus, the semiclassical prefactor of the bosonic propagator in quadrature representation is given by the van Vleck–Gutzwiller determinant
B 13

## Appendix C. The semiclassical prefactor in Fock representation

Starting from equation (3.1) and using equations (2.8), (2.11a), (2.11b), (2.32), the semiclassical propagator in Fock state representation can be written as
C 1As the asymptotic formula for the overlap is valid only for the oscillating region , one can substitute and for *l*=2,…,*L*, as well as for *l*=1,…,*L*. In the following, for reasons that will become obvious after the stationary phase conditions, we will refer to the *θ*^{(i/f)}_{l}’s as ‘phases’. Note that these definitions correspond to those in equation (2.11b). Together with the sums over **s**^{(i)} and **s**^{(f)}, the integrals over *θ*^{(i/f)}_{l} then run from −*π* to *π* and the phase *θ*^{(i)}_{1} does not enter in the Hamiltonian, which is part of the action.

The second derivative of the exponent in (.1) with respect to the final phases at the stationary point then yields

Therefore, together with the prefactor given by the second derivative of the generating function, the semiclassical prefactor after performing the integration over the final phases is given by C 2where we defined evaluated at the stationary point.

Noting that
where **1** is the vector, for which every entry is equal to 1, one can write (.2) also as
C 3Substituting and in the kinetic part of the action, as required by the stationary phase condition equation (3.2), in consideration of the preservation of the total number of particles by the classical equations of motion yields for the classical action in Fock state representation
C 4Moreover, we have
C 5and therefore
C 6It is important to note that the two matrices in the prefactor have different dimensions. The first one is an (*L*−1)×(*L*−1)-matrix, while the second one is *L*×*L* dimensional. However, as *θ*^{(f)} does not depend on *θ*^{(i)}_{1}, the second one has the form
C 7where is the vector resulting from *θ*^{(i)} by skipping its first entry. Moreover, we can match the dimensions of the two matrices by using
C 8Then,
C 9The determinant of this matrix however can be further simplified by using Sylvester’s determinant theorem, which states that for an *n*×*m* matrix *A* and an *m*×*n* matrix *B*
C 10where is the *n*×*n* unit matrix and applying it to the right-hand side of
C 11With this, one can show that
C 12Finally, we see that the integrand in (.6) is independent of , such that the final result reads
C 13

## Footnotes

One contribution of 12 to a theme issue ‘Loschmidt echo and time reversal in complex systems’.

↵1 The overlap decays exponentially outside this region.

↵2 In the (integrable) non-interacting case, an exact analytical formula for the semiclassical propagator can be found, such that eventually the limit can be taken.

- Accepted February 5, 2016.

- © 2016 The Author(s)

Published by the Royal Society. All rights reserved.