# Quantum phase-space analysis of population equilibration in multi-well ultracold atomic systems

###### Abstract

We examine the medium time quantum dynamics and population equilibration of two, three and four-well Bose-Hubbard models using stochastic integration in the truncated Wigner phase-space representation. We find that all three systems will enter at least a temporary state of equilibrium, with the details depending on both the classical initial conditions and the initial quantum statistics. We find that classical integrability is not necessarily a good guide as to whether equilibration will occur. We construct an effective single-particle reduced density matrix for each of the systems, using the expectation values of operator moments, and use this to calculate an effective entropy. Knowing the expected maximum values of this entropy for each system, we are able to quantify the different approaches to equilibrium.

###### pacs:

03.75.Lm, 03.75.Kk, 67.25.du, 03.75.Gg## I Introduction

The relaxation to equilibrium of closed quantum systems and their temporal evolution are important areas of study, as seen in, for example thermal1 (); thermal2 (); thermal3 (), with a beautiful experiment by Kinoshita et al. Kinoshita () having shown that relaxation to equilibrium does not happen in a trapped one dimensional Bose gas with contact interactions. This was not unexpected for a one dimensional untrapped Bose gas with point interactions, which is known to be an integrable system, but it had been thought that practical features such as the harmonic trap and imperfectly point-like interactions would compromise the integrability and the system would relax. On the other hand, there are closed quantum systems which are known to relax to a generalised equilibrium, without any interactions with a thermal cloud or other reservoir thermal2 (); zhang (); lfs (). To the best of our knowledge, there is as yet no established consensus on the mechanism by which this happens.

For the relaxation of closed quantum systems to a generalised equilibrium, one proposal is the eigenstate thermalisation hypothesis (ETH), in which every eigenstate of the Hamiltonian implicitly contains a thermal state srednicki (); thermal3 (), with the initial coherence between the eigenstates being lost by dephasing during the dynamics. Srednicki, when introducing this hypothesis, claimed that a necessary condition was the validity of Berry’s conjecture that the energy eigenfunctions behave like Gaussian random variables Berry (), which is expected to hold for systems which exhibit classical chaos in at least a large majority of the classical phase space.

In this work we study the relaxation or lack thereof in bosonic ultra-cold atomic systems held in two Joel (), three Nemoto (), and four-well Anglin (); nosso () tunnel-coupled potentials. For these sytems, the only constants of the motion are the total number of atoms and the total energy, so that in the classical sense we would expect only the twin well model to be integrable. Following the general wisdom, we would therefore not expect the two well system to relax to equilibrium, whereas we have previously shown that the four-well system, at least in one particular configuration, will demonstrate this feature over medium times nosso ().

We use stochastic integration in the truncated Wigner representation Wminus1 (); Steel () to perform phase-space analyses of the quantum dynamics and calculate the expectation values of the well populations and other operator products which allow us to construct effective reduced single-particle density matrices. From these matrices we are able to calculate a pseudoentropy which we previously found useful for the four well system nosso (). For the two-site model we are also able to calculate exact quantum results using a matrix equation for the number state coefficients, finding excellent agreement with the stochastic results over most of the time domain. We also calculate a type of Lyupanov exponent classically, using coupled Gross-Pitaevskii equations, and investigate the extent to which this allows us to predict which systems will relax. We find that all the systems we investigate will reach a generalised, but not necessarily static, equilibrium, but not for arbitrary initial quantum states. We find that the differences can be predicted qualitatively via the off-diagonal elements of the reduced density matrices at the beginning of the time evolution.

## Ii Physical model, Hamiltonian and equations of motion

We consider three different physical arrangements, firstly with two wells of equal depth, secondly with three equal wells at the vertices of an equilateral triangle, and thirdly with four equal wells in a square arrangement. We will follow a procedure which is effectively equivalent to changing the Hamiltonians at , in that all our initial conditions will have at least one well unoccupied. We outline here the the standard procedure for developing an effective Hamiltonian for two wells Joel (), where we consider separate wells with an independent condensate in each of the at the beginning of our investigations. The procedures for three or four wells are a simple extension of this nosso (). The Hamiltonian for a condensate in an external trapping potential, , may be written as

(1) |

where is the field operator for the condensate, and the non-linear interaction parameter is , where is the s-wave scattering length describing two-body collisions within the condensate, and is the atomic mass. In the case where the external potential provides a two well confinement for the condensate, we may simplify the above Hamiltonian by making use of the two-mode approximation. At zero temperature all atoms in the system are condensed and if the ground state energies of the condensate in either well are sufficiently separated from the energies of the condensate in all other excited single particle states, transitions to or from the modes of interest and these higher lying states can be neglected. We may then expand the field operator as

(2) |

where the are bosonic annihilation operators in each of the wells, and the are the ground state spatial wave functions of the condensate in each of the wells.

Using this in Eq. (1), we find an effective Hamiltonian

(3) | |||||

where we have neglected the spatial overlap of the different well densities. The single well bound state energies, , are

(4) |

, the tunnel coupling, is

(5) |

and effective non-linear interaction term is

(6) |

We set the single well bound state energies equal because we will consider only symmetric potentials where we can set .

This process simplifies the analyses by approximating the atoms in each separate well as being in a single mode, meaning that our equations are much easier to solve than would otherwise be the case. Generalising to the three and four site models, we will use as bosonic annihilation operators in each of the two or three wells, and and , with , in each side of the four-well system. We parametrise time by setting , so that dimensionless time as displayed in the results is labelled as . A schematic of the three well configuration is shown in Fig. 1, with the four well system being treated as shown in Fig. 2.

We may now also write the effective Hamiltonians for the three and four site systems, generalising from that for twin wells. For the three-well system we find the effective Hamiltonian,

(7) | |||||

Finally, for the the four-well system we have,

(8) | |||||

In order to calculate the quantum dynamics of these three systems, we will make use of the truncated Wigner representation Wminus1 (); Steel (), which gives results indistinguishable from those of the full quantum matrix equations for the number state coefficients DanChris () up until at least the relaxation time in the case of two and three wells, and is much easier to use for four-wells, where the full matrix equations become very difficult to use. Following the standard methods QNoise (), we find generalised Fokker-Planck equations for the Wigner pseudoprobability functions. As these contain third order derivatives, they cannot be mapped onto stochastic differential equations. Although methods have been developed to allow a mapping onto stochastic difference equations nossoEPL (), the numerical integration of these equations is even more unstable than that of the positive-P representation equations Steel (); P+ (), so we will not follow this route here. Instead, we truncate the third order terms in the generalised Fokker-Planck equations and make a mapping to coupled equations for the Wigner variables corresponding to the operators found in the Hamiltonians. We note here that classical averages of the Wigner variables correspond to symmetrically ordered operator expectation values, so that the necessary reordering must be undertaken before we arrive at solutions for physical quantities, for which normal ordering is more appropriate.

Having done this, we find the sets of coupled equations given immediately below. For the two-well system, the truncated Wigner equations of motion are,

(9) |

for the three-well system they are

(10) |

and for the four-wells they are

(11) |

Although the three sets of equations above might look classical, the Wigner variables themselves are drawn from appropriate distributions for the desired initial quantum states, with the stochasticity coming from the initial conditions for the chosen quantum states. The truncated Wigner equations are solved numerically by taking averages over a large number of stochastic trajectories, with initial conditions sampled probabilistically nosso (); AxMuzza ().

## Iii Classical integrability and stability

Before we investigate the quantum dynamics further we will examine the classical systems, which are described by coupled Gross-Pitaevskii equations which have the same form as above, but have deterministic initial conditions. Berry’s conjecture Berry () states that a quantum system can thermalise if the related classical system is unstable or chaotic in at least a large majority of the classical phase space. To examine the applicability of this conjecture, we numerically calculate approximate Lyupanov exponents for the three and four mode systems Eckmann (), which we define for the dynamics of one of the coupled wells as

(12) |

where

(13) |

where is an initial condition slightly perturbed from . Note that we have used numbers here rather than complex amplitudes, as it is the numbers which are directly observable. In practice, we obviously cannot integrate the equations for infinite time, so we integrate the coupled GPE type equations over a reasonably long time and look at the development of and hence . We also make the caveat here that because the system is both bounded and periodic, these effective Lyupanov exponents are not as definitive as in a non-periodic system.

In these systems as we define them, the only constants of the motion are the total energy and the total number of atoms. The obvious classical degrees of freedom are the number of atoms in each well, which, given the constraint on total number, would seem to be one less than the number of wells. This approach suggests that the three-well system has the same number of classical degrees of freedom as it does constants of the motion, and therefore should be classically integrable. The four-well system, as we have previously shown, is not stable, at least when two different tunneling rates are present nosso ().

In Fig. 3\subreffig:Ly3 we show the exponents for the three-well system, beginning in a situation far from the expected equilibrium, which would have one third of the atoms in each well. The classical solutions are again regular and periodic, with full oscillations between and atoms in well three, with the other two oscillating between and . The perturbed classical solutions, while still regular and periodic, never see more than or less than atom in wells one and two, while well three oscillates between and . Fig. 3\subreffig:Ly4 shows that the four-well system also behaves in a similar manner.

## Iv Relaxation

In a state of relaxed equilibrium, we expect that the entropy of the system will be maximised. In a quantum system, the entropy is normally defined using the density matrix for the system. While this works very well for small systems, such as those that are of interest in discrete variable quantum information applications, it becomes difficult to calculate the full density matrix for a many-body interacting quantum system, which will have a much larger Hilbert space. As in our previous work nosso (), we can define something which behaves as a reduced single-particle density matrix, which then allows us to calculate an effective entropy. Importantly for possible experimental investigations, all the quantities needed are in principle measurable using techniques developed by Ferris et al. Andyhomo (). For the four-well system we define,

(14) |

with the matrices for two and three wells being obvious truncations of the above. It is then an easy matter to calculate a single-particle pseudoentropy from this matrix,

(15) |

which will have a maximum value of when the atoms are equally distributed throughout the wells, which is statistically the most probable situation.

In Fig. 4 we show the results for the two-well system, beginning with all atoms initially in one of the wells. From Fig. 4\subreffig:number2 we see that the atoms equalise between the two wells, with the initial quantum states having almost no discernible effect on the process. This is also clear when we look at Fig. 4\subreffig:entropy2, where the pseudoentropy rises to a value very close to , which is the maximum single-particle value expected of this system in thermal equilibrium. We therefore see that, once quantum effects are included, this system appears to relax to something close to its equilibrium state at zero temperature without any contact with a reservoir, despite the fact that it is classically integrable. We have also integrated the exact quantum equations for the number state coefficients using a matrix method for initial Fock states DanChris (); Tania () and see a partial revival of the oscillations centred around , as shown in Fig. 5\subreffig:revive2. We used this method out to , finding one more very partial revival at around . In Fig. 5\subreffig:prob2, we show the at for well one of this system. The true equilibrium ground state for repulsive interatomic interactions would be a binomial distribution between the two wells, with the distribution in each well closely approximating a Gaussian centred on . We see that this is not the case here, with the highest probably actually being for all the atoms in well one. This shows that, despite the pseudoentropy being close to its maximum and the mean populations being equal at this time, the system is not in a true equilibrium state. We will return to this issue below.

When we look at the triple well system, we find differences which depend entirely on the quantum statistics, in a similar manner to previous investigations of quantum superchemistry RCsuper (); Wigsuper (). As shown in Fig. 6, when we start with a total of atoms, half in one well and half in another, with the third well unoccupied, the subsequent behaviours for initial coherent and Fock states are qualitatively different. The initial Fock states equilibriate to an equal number in each well and the pseudoentropy climbs to very close to , but initial coherent states evolve to a different final distribution. We find approximately atoms in each the wells which were originally occupied, while the initially unoccupied well holds a final number of . The final value of the pseudoentropy is well below the possible maximum, as the state of the three wells at the end of our integration time is by no means the statistically most probable distribution. A clue to the explanation of this behaviour can be found in the reduced single-particle density matrices at the final time of . For the initial Fock states we find,

(16) |

while for initial coherent states the final result is

(17) |

from which we see that the off-diagonal elements are much larger for the initial coherent states. As these elements represent coherences between the modes, they contain information and therefore, the larger they are, the further the system is from the maximum entropy equilibrium state.

For the four-well system, as shown in Fig. 7, we see that the maximum pseudoentropy () is closely approached for both initial Fock and coherent states, but also that the population dynamics are quite different. The small oscillations about the equilibrium value in well die down much slower for initial coherent states than for initial Fock states. The off-diagonal elements of the matrix at the end of the integration time were also larger for initial coherent states, although the number distribution is obviously tending towards a quarter of the total in each well.

We note here that for initial Fock states in all these systems, the off-diagonal elements of the reduced density matrix are all zero at , while this is not the case for initial coherent states. This can be easy seen by considering some of the typical initial matrix elements for each system. For the two well system with an initial coherent state in one side and nothing in the other, we find

(18) |

while for a Fock state in one side and nothing in the other, we have

(19) |

so that all the off-diagonal elements are originally zero. For the three-well system with our coherent state initial conditions, we find two of the off-diagonal elements are non-zero,

(20) |

with expectation values of any operator products containing or being zero. For the initial Fock states in wells one and two, we find that all the off-diagonal elements are again initially zero. With the four-well system and the two Fock states, all off-diagonal elements are initially zero, whereas for the two coherent state occupations, we find the non-zero elements,

(21) |

We therefore see that, where the initial reduced density matrices are the same for each choice of quantum states, as in the two-well example, we see almost indistinguishable evolution of both the populations and the pseudoentropy. Where they exhibit the most difference, as in the three-well example, the subsequent dynamics are very different. The four-well case, where the initial matrices are less different, lies between these two extremes.

## V Conclusions and Discussion

Our stochastic phase-space investigations of the quantum dynamics of these three different Bose-Hubbard models shows that all three may equilibrate for given initial conditions. We have also shown that an effective reduced single-particle density matrix may be calculated using phase-space methods. While this matrix does not have all the information contained in the full quantum density matrices, it does allow for the calculation of a pseudoentropy and an investigation of the equilibriation process. Perhaps more importantly, all the measurements necessary to construct this density matrix are possible in principle. We have shown that classical integrability is not necessarily a good guide as to whether a given system will relax or not, as seen by the two-well system, which does relax to a close to equilibrium state, with only minor revivals over the time we investigated. We have also shown that the presence of initial coherences between the different modes can affect the equilibriation process, as shown by the presence of the off-diagonal elements of the density matrix. We have also demonstrated that the stochastic phase-space methods developed for quantum optics can be useful for the theoretical study of statistical mechanical processes in interacting atomic systems where the Hilbert space becomes so large as to make density matrix methods extremely difficult to use.

What we cannot show, due to the limitations of our numerical methods, is whether the oscillations in these systems ever revive fully in finite time, after their initial collapse. We suffer from two restrictions here. The first is that our exact method, although applicable to the two and three well models, can still not be used for arbitrary times. The second is that the truncated Wigner method gives us accurate information on the collapse of the oscillations, but not on the revivals. However, taken together, these two methods allow us to investigate the medium time dynamics and also allow us to determine the effects of different initial quantum states on the mean-field dynamics. In none of the systems we have examined here do the solutions of the mean-field classical equations even approximate the true dynamics of the mean-fields. In that sense these are truly quantum systems. On a final note, it is interesting to consider the systems we have examined as having non-Markovian reservoirs connected to them at , when the unoccupied wells become accessible to the atoms. Whether the relaxation we see can be described accurately as a smaller subsytem coming into some type of equilibrium with this bath will be investigated in later work.

## Acknowledgments

This research was supported by the Australian Research Council under the Centres of Excellence and Future Fellowships Programs. The authors thank Matthew Davis for stimulating discussions.

## References

- (1) J.M. Deutsch, Phys. Rev. A43, 2046 (1991).
- (2) M. Rigol, V. Dunjko, V. Yurovsky, and M. Olshanii, Phys. Rev. Lett. 98, 050405 (2007).
- (3) M. Rigol, V. Dunjko and M. Olshanii, Nature 452, 854 (2008).
- (4) T. Kinoshita, T. Wenger, and D.S. Weiss, Nature 440, 900 (2006).
- (5) J.M. Zhang, C. Shen, and W.M. Liu, arXiv:1102.2469v1 (2011).
- (6) L.F. Santos, A. Polkovnikov, and M. Rigol, Phys. Rev. Lett. 107, 040601 (2011).
- (7) M. Srednicki, Phys. Rev. E50, 888 (1994).
- (8) M.V. Berry, J. Phys. A 10 (1977).
- (9) G.J. Milburn, J.F. Corney, E.M. Wright, and D.F. Walls, Phys. Rev. A55, 4318 (1997).
- (10) K. Nemoto, C.A. Holmes, G.J. Milburn, and W.J. Munro, Phys. Rev. A63, 013604 (2000).
- (11) M.P. Strzys and J.R. Anglin, Phys. Rev. A81, 043616 (2010).
- (12) C.V. Chianca and M.K. Olsen, Phys. Rev. A83, 043607 (2011).
- (13) R. Graham, Springer Tracts in Modern Physics, 66, 1 (1973).
- (14) M.J. Steel, M.K. Olsen, L.I. Plimak, P.D. Drummond, S.M. Tan, M.J. Collett, D.F. Walls, and R. Graham, Phys. Rev. A58, 4824 (1998).
- (15) D.F. Walls and C.T. Tindle, J. Phys. A 4, 534 (1972).
- (16) C.W. Gardiner, Quantum Noise, (Springer-Verlag, Berlin, 1991).
- (17) L.I. Plimak, M.K. Olsen, M. Fleischhauer, and M.J. Collett, Europhys. Lett. 56, 372 (2001).
- (18) P.D. Drummond and C.W. Gardiner, J. Phys. A 13, 2353 (1980).
- (19) M.K. Olsen and A.S. Bradley, Opt. Commun. 282, 3924 (2009).
- (20) J.-P. Eckmann and D. Ruelle, Rev. Mod. Phys. 57, 617 (1985).
- (21) A.J. Ferris, M.K. Olsen, E.G. Cavalcanti, and M.J. Davis, Phys. Rev. A78, 060104(R) (2008).
- (22) T.J. Haigh, A.J. Ferris, and M.K. Olsen, Opt. Commun. 283, 3540 (2010).
- (23) M.K. Olsen and L.I. Plimak, Phys. Rev. A68, 031603 (2003).
- (24) M.K. Olsen, Phys. Rev. A69, 013601 (2004).