Atomic three-body loss as a dynamical three-body interaction

Atomic three-body loss as a dynamical three-body interaction


We discuss how large three-body loss of atoms in an optical lattice can give rise to effective hard-core three-body interactions. For bosons, in addition to the usual atomic superfluid, a dimer superfluid can then be observed for attractive two-body interactions. The non-equilibrium dynamics of preparation and stability of these phases are studied in 1D by combining time-dependent Density Matrix Renormalisation Group techniques with a quantum trajectories method.

03.75.Lm, 42.50.-p

Cold atomic gases in optical lattices have proven a test bed for understanding novel quantum phases (1) and non-equilibrium many-body dynamics (2); (3). Recently, Syassen et al. (4) showed that a strong two-body loss process for molecules in an optical lattice (1) could produce an effective, elastic hard-core repulsion and thus a Tonks gas (4); (5). This is related to the quantum Zeno effect: a large loss dynamically suppresses processes creating two-body occupation on a particular site. Whilst elastic two-body interactions occur in many systems, regimes where elastic 3-body interactions dominate are rare in nature. Here we discuss how the ubiquitous, though normally undesirable three-body losses of atomic physics experiments can induce effective three-body interactions. These are associated with interesting quantum phases, including Pfaffian states (6), and could be used to stabilise three-component Fermi mixtures (7), assisting in the production of a colour superfluid state (8). We investigate Bosons in an optical lattice, where a three-body hard-core constraint stabilises the system with attractive two-body interactions, and a dimer superfluid phase emerges. We focus on the dynamics of this intrinsically time-dependent system, both testing the hard-core constraint for finite loss rates, and studying non-equilibrium properties including decay. In 1D, the exact evolution is computed by combining time-dependent density matrix renormalisation group methods (t-DMRG) (2); (3) with a quantum trajectories approach from quantum optics (9).

Figure 1: (a) Bosons in an optical lattice in the presence of three-body loss at a rate . (b) Example model parameters estimated for Cs at a magnetic field of 15 Gauss (where the scattering length is and the recombination length (10) where is the Bohr radius) as a function of lattice depth , showing (solid line), (dashed), and (dotted). Values of are obtained by integrating the measured three-body recombination rates in free space over a state with three particles in a single Wannier function. (c) The probability that at least one loss event has occurred at time , beginning with a single particle on each of 10 sites, and (solid line), 5 (dashed), and 10 (dotted), computed using t-DMRG methods (see text for details).

Three-body recombination (10) in an optical lattice corresponds to decay into the continuum of unbound states, and thus loss from the lattice. This can be described by a master equation in the Markov approximation (9), which for atoms in the lowest band of an optical lattice can be projected onto the corresponding basis of Wannier functions (1); (5), associated with bosonic annihilation operators on site . We can separate the master equation into terms which conserve particle number, corresponding to an effective Hamiltonian , and terms which remove three particles on a site:

where denotes the system density operator with atoms and . The dominant loss term is on-site three-body decay 1 and is the corresponding rate. The effective Hamiltonian is


with the nearest neighbor tunneling amplitude, the elastic two-body interaction, and the local potential. The Hamiltonian is valid in the limit where , with the band gap and the mean density. In an experiment these parameters, in particular , may be tuned whilst remains constant and large 2. In Fig. 1c we show example values of , , and using numbers for Caesium as a function of lattice depth.

If we begin in a pure state with particles, then loss processes lead to heating, in that they produce a mixed state of different particle numbers. Within a fixed particle number sector, the dynamics are described by . Three-body interactions emerge most clearly in the limit of rapid decay: . If we define the projector onto the subspace of states with at most two atoms per site and , then in second order perturbation theory we obtain the effective model


where , and denotes the set of nearest neighbours of site . The term describes the Hubbard dynamics, Eq. (2), supplemented by the hard-core constraint . Furthermore, the effective loss rates decrease as 3.

Thus, we see the clear emergence of a three-body hard-core constraint in the limit . We can study the physics of the projected model to obtain a qualitative understanding of the quantum phases associated with the projection. However, the residual loss processes make this system intrinsically time-dependent, and can give rise to heating. We therefore study the full non-equilibrium dynamics, by combining t-DMRG methods (2); (3) with an expression of the master equation as an average over quantum trajectories (9). Each stochastic trajectory begins from an initial pure state (sampled from the initial density matrix), and can be interpreted as describing a single experimental run, in which losses occurred at particular times and on sites . The evolution is described by the non-Hermitian , except for times , where losses (or quantum jumps) occur,


where the jump operator corresponds to three-body loss on site . In stochastic simulation of the master equation, the times are points where the norm of the state falls below a randomly chosen threshold. At these times, a random jump operator is selected according to the probabilities and applied. In this way we can both investigate individual trajectories and compute expectation values from the master equation. The latter is performed by stochastic average over both initial states and over jump events, which converges rapidly as the number of trajectories is increased.

The need to simulate many trajectories for convergence is offset by the efficiency of simulating states rather than density matrices. In contrast to directly applying t-DMRG to density matrices (11), we can simulate the master equation without squaring the local Hilbert space dimension, and by time-evolving states with fixed particle number, we also make use of existing optimisations for conserved quantities. Despite the application of local jump operators, we find the evolution quite efficient, especially for small numbers of jumps 4.

As an example of the suppression of loss, we consider preparing a homogenous initial state at unit filling in a deep optical lattice where . At time we suddenly ramp the lattice to a finite depth, and observe the probability that a single three-body loss event has occurred as a function of time. In Fig. 1b we plot this probability for different as a function of . We see a clear suppression of loss rates for large , and also a substantial decrease for larger , resulting from the decreased amplitude for doubly occupied sites.

In the limit of large , it is instructive to study the equilibrium phase diagram of the projected Hamiltonian . For , we observe the well-known Mott Insulator (MI) and atomic superfluid phases of the Bose-Hubbard model. However, the three-body hard-core condition will also stabilise the system for , where we find a dimer superfluid phase (see Fig. 2a). This is characterised by the vanishing of the order parameter signalling superfluidity of single atoms (ASF) (), while a dimer superfluidity (DSF) order parameter persists (). The superfluid regimes are connected via a quantum phase transition associated with the spontaneous breaking of a discrete symmetry, reminiscent of an Ising transition (12): the DSF order parameter transforms with the double phase compared to the ASF order parameter . Consequently, the symmetry exhibited by the DSF order parameter is broken when reaching the ASF phase.

We can obtain a qualitative mean-field picture using a homogeneous Gutzwiller ansatz wavefunction, given for the projected Hilbert space by , where . Normalisation implies , while the filling is . To examine the phases, we find the energy , , where and is the number of lattice sites. For any the energy is minimized for an integer multiple of . For (placing us in the ASF phase with ), the phase-locking expression contributes a “source” term linear in to the energy, and consequently the minimum of the energy cannot be located at . Thus, a finite atomic condensate always implies a dimer component . However, the reverse is not true and one may have pure dimer superfluidity without an atomic condensate.

At fixed , the energy is a function, e.g., of alone. For the second order transition found within our mean field theory, an instability for atomic superfluidity is indicated by its mass term crossing zero, . This leads to a critical interaction strength for the ASF-DSF transition, . Within the DSF phase the order parameter obeys independent of the interaction strength. For , we approach a MI state in a second order transition. At we find that the ASF-DSF transition takes place at the same coupling strength as the ASF-MI transition, but with the opposite sign. The complete mean field phase diagram in the plane of density and interaction strength is plotted in Fig. 2a.

From the last term in Eq. (3), we can estimate the initial loss rate from the ground state Gutzwiller wave function. We obtain the rate , which is zero in the MI, and for the DSF, . In the DSF phase, the critical temperature at low densities, and the energy density deposited by a single loss, . The number of independent loss events needed to melt the DSF is then proportional to , and the melting time strongly decreases for increasing density, proportional to .

These qualitative features are reproduced in 1D, as supported by numerical calculation of the ground state for . In Figs. 2b,c we show the characterisation of the crossover between the ASF and DSF regimes in 1D via the off-diagonal elements of the single particle density matrix, , and the dimer density matrix . In the MI regime, the off-diagonal elements of and decay exponentially. As we enter the superfluid regime, quasi-long range order is visible in the polynomial decay (linear on the logarithmic scale). As is made more negative, we see a return to exponential decay for the off-diagonal elements of , but the off-diagonal elements of still decay polynomially and, indeed, increase in magnitude. This characterises the DSF regime in 1D. Here, the transition to the DSF and MI regimes occurs at much smaller than in higher dimensions, but these two transitions again occur at similar for .

Figure 2: Equilibrium analysis of the projected Bose-Hubbard model . (a) Mean-field phase diagram as a function of and density, . (b,c) Magnitude of off-diagonal elements of (b) the single particle density matrix and (c) the dimer density matrix, , as a function of , for 10 (thin solid line), 5 (dotted), 0 (dot-dash), -5 (dashed), and -10 (thick solid line). These results are computed for 20 particles on 20 lattice sites in 1D, with box boundary conditions and , using imaginary time-evolution in t-DMRG, and plotted on a logarithmic scale.
Figure 3: Dynamics of adiabatic ramps into a dimer superfluid regime. (a) We begin with (i) a Mott-insulator state (ramping ), and (ii) a state with pre-prepared dimers in a superlattice (removing the superlattice). (b)-(c) The sum of kinetic () and interaction () energy and (inset) particle number as a function of time for two example trajectories, one with no loss events (dashed lines) and one with several loss events (solid lines). Here, (b) shows a ramp from to , with , with and ramp parameters, and (c) shows a ramp with a superlattice potential, , where , adjusted so that , with fixed . In each case, . For (b), we use 20 atoms on 20 lattice sites, for (c), 14 atoms on 23 lattice sites. (d) Plot showing the probability that no loss event has occurred after time for the ramps in (b) (dashed line) and (c) (solid).

A dimer superfluid phase could be prepared in an experiment beginning in states that have very small amplitude of three-body occupation due to two-body interactions. Hamiltonian parameters could then be ramped adiabatically across the phase boundary, making use of dynamical suppression of three-body occupation. We study two such scenarios as illustrated in Fig. 3a: (i) We begin from a Mott Insulator state, and then ramp from to in order to enter the dimer superfluid regime (which is intuitive, but associated with large probability of decay); or (ii) We apply a superlattice that raises the energy of certain sites and form a Mott Insulator with two particles per site in the lowest wells, then switch the interaction rapidly to on a timescale much faster than tunnelling between the lowest wells and ramp down the superlattice. In each case, we compute dynamics using t-DMRG methods with quantum trajectories, and consider an on-site three-body loss with . We ramp parameters on a timescale where without loss the ground state will be reached with minimal presence of excited states.

In Figs. 3b,c we show the time dependence of the sum of kinetic and interaction energies, and of the total particle number for example trajectories. For each ramp type, we choose a “lossless” trajectory, where the ground state is reached adiabatically, and a “lossy” trajectory, where three-body loss events lead to heating of the system (as holes are produced that correspond to excited states).

In Fig. 3d, we compare the probability for each type of ramp of producing the “lossless” trajectory. We see that for the ramp from the MI state, where we pass through a region of small , the probability of such a trajectory is essentially zero. For the superlattice ramp, on the other hand, it is much more likely that we will obtain the ground state from a randomly chosen trajectory. This is both because the superlattice allows us to use large , so that the two-body interaction reduces the probability of triply occupied sites, and because the superlattice facilitates the choice of a lower density in the system, which significantly reduces the effective loss rate.

In Fig. 4 we show the local density as a function of time and the final dimer density matrix for (a) the “lossless” and (b) the “lossy” trajectories of Fig. 3c. When a loss occurs it affects not just the density for the site on which it occurs, but also on neighbouring sites due to the knowledge that we obtain of the position of the remaining particles. We also see clearly the destruction of correlations in the region of the system where the loss occurs. Note, however, that off-diagonal correlations are still visible in parts of the system, and that a single loss event does not always destroy the properties of the final state entirely. As discussed above, the probability of an individual loss event increases with system size, but a single loss event will change the character of the final state less.

Figure 4: Comparison of (a) lossless and (b) lossy trajectories from Fig. 3c. We show the mean density as a function of position and time, and magnitude of the dimer correlation function () at the end of the ramp.

The three-body interactions discussed here could have applications to producing Pfaffian-like states, and also for fermions, especially stabilising three-component mixtures. The combination of quantum trajectories methods with t-DMRG has potential applications in the simulation of other classes of master equations. In addition, there are open questions regarding the nature of the ASF-DSF phase transition, which we will address in more detail within a quantum field theoretical treatment of the attractive Bose-Hubbard model with a three-body hard-core constraint (13).

We thank L.-M. Duan for discussions. Work in Innsbruck is supported by the Austrian Science Foundation (FWF) through SFB 15 and project I118_N16 (EuroQUAM_DQS), and by the DARPA OLE network. JMT is supported by a Pappalardo Fellowship at MIT.


  1. Off-site loss terms can arise from non-zero overlap of the corresponding Wannier functions, in analogy to off-site elastic interactions, which are discussed in L.-M. Duan, Europhys. Lett. 81, 20001 (2008). These are small for values used here.
  2. This is because non-universal behavior (where , with the scattering length) is observed away from a Feshbach resonance, where the rate of three-body recombination is large and weakly dependent on (10), while the two-body scattering length can be very small ().
  3. Note that off-site processes can place an upper bound on the useful value of , as increasing suppresses on-site loss, but increases the rate of loss from processes involving neighbouring lattice sites.
  4. We typically take the number of states retained in bipartite splittings (2), in the results presented here.


  1. I. Bloch, J. Dalibard, and W. Zwerger, Rev. Mod. Phys. 80, 885 (2008).
  2. G. Vidal, Phys. Rev. Lett. 93, 040502 (2004); F. Verstraete, V. Murg, and J. I. Cirac, Ann. Phys. 57, 143 (2008).
  3. A. J. Daley et al., J. Stat. Mech.: Theor. Exp. P04005 (2004); S.R. White and A.E. Feiguin, Phys. Rev. Lett. 93, 076401 (2004); F. Verstraete, J. J. Garcia-Ripoll, and J. I. Cirac, ibid. 93, 207204 (2004).
  4. N. Syassen et al., Science 320, 1329 (2008).
  5. J. J. Garcia-Ripoll et al., arXiv:0809.3679v1.
  6. B. Paredes, T. Keilmann, and J. I. Cirac, Phys. Rev. A 75, 053611 (2007).
  7. T. B. Ottenstein et al., arXiv:0806.0587v2.
  8. A. Rapp, G. Zarand, C. Honerkamp, and W. Hofstetter, Phys. Rev. Lett 98, 160405 (2007).
  9. H. J. Carmichael, An Open Systems Approach to Quantum Optics, (Springer, Berlin, 1993); C. W. Gardiner and P. Zoller, Quantum Noise, (Springer, Berlin, 2005).
  10. See, e.g., T. Kraemer et al., Nature 440, 315 (2006).
  11. M. Zwolak and G. Vidal, Phys. Rev. Lett. 93, 207205 (2004); F. Verstraete, J. J. Garcia-Ripoll, and J. I. Cirac, Phys. Rev. Lett. 93, 207204 (2004).
  12. M. W. J. Romans, R. A. Duine, Subir Sachdev, and H. T. C. Stoof, Phys. Rev. Lett. 93, 020405 (2004); L. Radzihovsky, J. I. Park, and P. B. Weichman, ibid. 92, 160402 (2004); L. Radzihovsky, P. B. Weichman, and J. I. Park, Ann. Phys. 323, 2376 (2008).
  13. S. Diehl et al., in preparation (2008).
This is a comment super asjknd jkasnjk adsnkj
The feedback cannot be empty
Comments 0
The feedback cannot be empty
Add comment

You’re adding your first comment!
How to quickly get a good reply:
  • Offer a constructive comment on the author work.
  • Add helpful links to code implementation or project page.