# Dynamical phase diagram of Gaussian wave packets in optical lattices

###### Abstract

We study the dynamics of self-trapping in Bose-Einstein condensates (BECs) loaded in deep optical lattices with Gaussian initial conditions, when the dynamics is well described by the Discrete Nonlinear Schrödinger Equation (DNLS). In the literature an approximate dynamical phase diagram based on a variational approach was introduced to distinguish different dynamical regimes: diffusion, self-trapping and moving breathers. However, we find that the actual DNLS dynamics shows a completely different diagram than the variational prediction. We numerically calculate a detailed dynamical phase diagram accurately describing the different dynamical regimes. It exhibits a complex structure which can readily be tested in current experiments in BECs in optical lattices and in optical waveguide arrays. Moreover, we derive an explicit theoretical estimate for the transition to self-trapping in excellent agreement with our numerical findings, which may be a valuable guide as well for future studies on a quantum dynamical phase diagram based on the Bose-Hubbard Hamiltonian.

###### pacs:

03.75.Gg, 03.75.Lm, 67.85.HjBose-Einstein condensates (BECs) trapped in periodic optical potentials have proven to be an invaluable tool to study fundamental and applied aspects of quantum optics, quantum computing and solid state physics Bloch (2005); Bakr et al. (2009, 2010); Simon et al. (2011); Morsch and Oberthaler (2006). In the limit of large atom numbers per well, the dynamics can be well described by a mean-field approximation which leads to a lattice version of the Gross-Pitaevskii equation, the discrete nonlinear Schrödinger equation (DNLS) Pethick and Smith (2008); Morsch and Oberthaler (2006). One of the most intriguing features of the dynamics of nonlinear lattices is that excitations can spontaneously stably localize even for repulsive nonlinearities. This phenomenon of discrete self-trapping, also referred to as the formation of discrete breathers (DB), is a milestone discovery in nonlinear science that has sparked many studies (for reviews see Flach and Gorbach (2008); Campbell et al. (2004); Hennig and Fleischmann (2013)). DBs have been observed experimentally in various physical systems such as arrays of nonlinear waveguides Eisenberg et al. (1998); Morandotti et al. (1999) and Josephson junctions Trías et al. (2000); Ustinov (2003), spins in antiferromagnetic solids Schwarz et al. (1999); Sato and Sievers (2004), and BECs in optical lattices Eiermann et al. (2004). Self-trapping has also been shown to exist in the dynamics of the Bose-Hubbard Hamiltonian Hennig et al. (2012) and in calculations beyond the Bose-Hubbard model, which include higher-lying states in the individual wells Sakmann et al. (2009); Trujillo-Martinez et al. (2009).

Under which conditions will a Gaussian distributed initial condition in an optical lattice become diffusive, or self-trapped or a moving breather after sufficient propagation time? This question was addressed in a seminal work Trombettoni and Smerzi (2001a) that has become a standard reference in both experimental and theoretical studies involving self-trapping in optical lattices including BECs and optical waveguide arrays in the past decade Morsch and Oberthaler (2006); Christodoulides et al. (2003); Fleischer et al. (2003a); Albiez et al. (2005); Cataliotti et al. (2001); Fleischer et al. (2003b); Eiermann et al. (2004); Morsch et al. (2001); Burger et al. (1999); Smerzi et al. (2002); Anker et al. (2005); Iwanow et al. (2004); Eiermann et al. (2003); Reinhard et al. (2013); Bersch et al. (2012); Dienst et al. (2013); Dong et al. (2013). By means of a variational approach that approximates the DNLS dynamics it was shown that the dynamical phase diagram is divided in in different regimes (diffusion, self-trapping and moving breathers) Trombettoni and Smerzi (2001a). The ‘dynamical phase diagram’ distinguishes between qualitatively different steady state solutions. However, a recent study Franzosi et al. (2011) indicated that numerical simulations based on the DNLS show strong deviations from the variational dynamics. Moreover, as we show in Fig. 2a, the actual parameter regions in which the different dynamical phases can be observed are completely different from those predicted in ref. Trombettoni and Smerzi (2001a). Therefore, a new theory of the self-trapping of Gaussian wave packets is needed that can predict the dynamical regimes.

Here, we study Gaussian wave packets of BECs in deep optical lattices in the mean-field limit (described by the DNLS). We both analytically and numerically present a detailed and accurate dynamical phase diagram that separates the different dynamical regimes (diffusion, self- trapping, moving breathers and multi-breathers).

Although our focus is on BECs in optical lattices, where our dynamical phase diagram can be readily probed with single-site addressability Gericke et al. (2008); Albiez et al. (2005); Morsch and Oberthaler (2006), our results apply as well other systems described by the DNLS including optical waveguide arrays Christodoulides et al. (2003); Christodoulides and Eugenieva (2001). Our results may as well be a valuable guide for studies that aim to understand the deviations of the correlated self-trapping dynamics as described by the Bose-Hubbard Hamiltonian from the mean field dynamics. Examples of such deviations have been observed both experimentally Reinhard et al. (2013) and theoretically Hennig et al. (2012).

Model.– In the limit of large atom numbers per well, the dynamics of dilute Bose-Einstein condensates trapped in deep optical potentials are well described by the mean-field Bose-Hubbard Hamiltonian Trombettoni and Smerzi (2001b); Pethick and Smith (2008); Buonsante and Penna (2008)

(1) |

where is the lattice size, is the norm (number of atoms at site ), denotes the onsite interaction (between two atoms at a single lattice site), is the onsite chemical potential and is the strength of the tunnel coupling between adjacent sites. The corresponding dynamical equation is the DNLS which read in its dimensionless form

(2) |

for , and . In the numerics we present we use periodic boundary conditions (), however, we checked that our results equally hold for closed boundary conditions. The DNLS describes a high dimensional chaotic dynamical system. Its dynamics, however, in general is far from being ergodic and e.g. shows localization in the form of stationary DB (localized excitations pinned to the lattice) and so called moving breathers traversing the lattice Flach and Gorbach (2008) (see Fig. 1).

The different characteristic types of dynamics are exemplified in Fig. 1: An initial condition that rapidly disperses, i.e. shows diffusive behavior, can be seen in Fig. 1a. Moving breathers are strictly speaking not supported by the DNLS Oxtoby and Barashenkov (2007), they radiate norm and will eventually be pinned to the lattice or diffusively spread. However, as in Fig. 1b, these losses of norm can become infinitesimally small so that these solutions remain traveling for extremely long times and can therefore be regarded as actual moving breather solutions for practical use. Finally an example of a stationary DB emerging from a Gaussian initial condition is shown in Fig. 1c.

Exact DB solutions were calculated analytically in Hennig et al. (2010); Flach and Gorbach (2008) for a system of three sites (trimer). For larger lattices the DB solutions can be numerically calculated using the anti-continuous method Flach and Gorbach (2008); Aubry (1997); Marin and Aubry (1996) or a Newton method Proville and Aubry (1999). To prepare an experimental system in exact breather states is usually not feasible in practice. In contrast, Gaussian initial conditions – for which we derive a dynamical phase diagram – can be well controlled in experiments on BECs with single-site addressability Gericke et al. (2008); Albiez et al. (2005); Morsch and Oberthaler (2006).

Variational approach.– In Trombettoni and Smerzi (2001a, b) the dynamics of Gaussian wave packets defined by

(3) |

was studied using a variational collective coordinate approach (a technique very successfully applied in a variety of fields Campbell et al. (1983); Campbell (1980)). Here, and are the center and the width of the Gaussian distribution and and their associated momenta. The variational approach leads to approximate equations of motion for the conjugate variables and Trombettoni and Smerzi (2001a, b). It is important to note that the variational approach assumes that the excitation is well approximated by a Gaussian at all times. The effective (approximate) Hamiltonian is Trombettoni and Smerzi (2001b)

(4) |

with depending on the initial values of the center and width parameter as well as and their conjugate momenta.

In ref. Trombettoni and Smerzi (2001a) a dynamical phase diagram was derived based on the variational approach Eq. 4. These theoretical predictions however disagree with our numerical simulations of the actual DNLS dynamics, because the final dynamical state will in most cases be highly non-Gaussian. Hence, the variational approach breaks down, e.g., Fig. 1a shows diffusion, while the variational approach predicts a self-trapped state. In Fig. 1b we find a moving breather while the variational approach yields diffusive behavior. Moreover, the entire phase diagram Fig. 2b(top) for the DNLS is completely different from the prediction of the variational approach (shown as dashed and dotted lines in Fig. 2a).

Numerical construction of the phase diagram.– In order to construct a dynamical phase diagram based on the DNLS criteria to distinguish the different dynamical regimes are needed. This separation is non-trivial, as part of the atom cloud can, e.g., remain trapped in one region while another part diffuses in the remainder of the lattice. We can separate the regimes by defining a local norm and local average current

(5) |

and two order parameters by

(6) |

where is the central site at which assumes its maximum. Both quantities were evaluated at time and averaged over the last 10% of the time.

In the diffusive regime both order parameters are small. Moving breathers are characterized by large and large , while in the self-trapping regime one finds large but nearly vanishing currents due to the stationarity of the self-trapping solutions. Fig. 2 shows and as a function of the initial phase difference and the nonlinearity for (see Supplemental Material for ). The examples of Fig. 1 are marked in Fig. 2 by different symbols to demonstrate how the above criteria can be used to distinguish the dynamical regimes (a cut along can be found in the Supplemental Material).

We construct the dynamical phase diagram which assigns every initial condition to a specific dynamical regime by defining suitable thresholds for the order parameters and . Figure 2a(top) shows remarkably sharp transitions in the maximum local probability current. We use the thresholds to separate moving breathers from self-trapping and diffusion and to delimit diffusion from self-trapping. By setting a threshold for the probability density function of to be observed in a diffusive state (see Supplemental Material for details), we find for and for .

With these thresholds we obtain the dynamical phase diagram shown in Fig. 2b. Our results are not sensitive to the exact value of the parameter . For comparison we show a phase diagram obtained with in the Supplemental Material which differs only minimally from Fig. 2b. To identify higher order self-trapped excitations of more complex structure (see Fig. 1d) in our dynamical phase diagram, we use a third order parameter, the local energy in the vicinity of measured relative to the fixed point energy (which will be introduced below), i.e. . For higher order self-trapping we found in contrast to for discrete breathers centered around a single site. We chose and marked this regime in black in Fig. 2b.

Theory.– In the following, we derive an analytical estimate for the self-trapping transition separating the different phases in the dynamical phase diagram. The exact DB solution forming a single peak centered at a single lattice site (denoted cp for central peak) is a local maximum in the energy and can be numerically constructed using the methods described in Refs. Proville and Aubry (1999); Flach and Gorbach (2008); Aubry (1997); Marin and Aubry (1996). The cp-solutions correspond to elliptic fixed points in phase space Rumpf (2004) which are separated from moving and chaotic solutions by saddle points in the energy landscape. These saddle points correspond to another kind of stationary solutions centered between two lattice sites called central bond (cb) states. They are unstable Rumpf (2004) and have lower energy than cp-solutions. Due to continuity the cb-solutions are considered intermediate states for an excitation to hop between lattice sites.

A necessary condition for self-trapping is , where is the energy of the initial Gaussian and is the energy of the cb-solution. To derive an expression for , note that the energy of an arbitrary state of norm compared to the state of the same shape but of unit norm scales as . Thus reads

(7) |

with . We approximate the energy of the initial Gaussian by with (see Eq. 4), which we found agrees very well with the energy of the DNLS Hamiltonian (Eq. (1)), although the dynamics strongly differs.

We make the following Ansatz to estimate the critical nonlinearity at the transition to self trapping

(8) |

where and . The critical nonlinearity at the transition to self-trapping reads

(9) | |||||

Our analytical estimate of the self-trapping transition Eq. (9) is shown in Fig. 2a (thick black line) and Fig. 2b (thick red line) and is in excellent agreement with the simulations of the DNLS dynamics. By fixing we can determine the constant numerically. We find and . We chose for the estimate since the phase difference between neighboring sites for an exact solution of a discrete breather is . Below, we will present an independent way to obtain by calculating directly.

How can we interpret the rescaled nonlinearity in Eq. 8? The ability of localized excitations to move inside a lattice can be analyzed using the concept of the Peierls-Nabarro (PN) energy barrier which we shortly review: A crucial condition for a localized excitation to move across the lattice is that the initial energy is lower than the energy of a cb-solution such that this intermediate state cannot act as a barrier. Self-trapping on the other hand requires the initial energy to be between the maximum energy state cp and cb which can thus act as a barrier in phase space and inhibits migration. This energy gap between cp and cb is referred to as the PN barrier Kivshar and Campbell (1993); Hennig et al. (2010); Rumpf (2004).

We report the PN barrier in Fig. 3. Self-trapping requires a non-zero PN barrier, which we find for nonlinearities larger than in Fig. 3. Since the transition to self-trapping is reflected by the onset of a non-zero PN barrier, we identify . With we find in excellent agreement with the numerical value . Note that only a fraction of the norm of the initial Gaussian will actually be trapped when a breather is created (and this norm is smaller the more the phase differences of the Gaussian deviates from the phase difference of the breather fixed point, i.e. ). Our theoretical analysis focuses on the self-trapped fraction of norm and thus the existence of a non-zero PN barrier corresponding to the respective fraction of norm . This of course leaves the background unattended yet proves to be a valid approximation.

In conclusion, we calculated a dynamical phase diagram for the different dynamical regimes (diffusion, moving breathers and self-trapping) of Gaussian initial conditions in periodic optical potentials. By defining two order parameters – the maximum local atomic density and the maximum local current – we numerically construct the dynamical phase diagram. We derived an explicit expression (Eq. 9) for the nonlinear interaction that separates the dynamical regimes, in very good agreement with the numerical results. Gaussian initial conditions can be well controlled in experiments on BECs with single-site adressability Gericke et al. (2008); Albiez et al. (2005); Morsch and Oberthaler (2006) where our predictions can be readily tested experimentally. We hope our study further stimulates experiments and theoretical work on the transition between self-trapping, diffusion and solitons in the quantum regime, e.g., based on the Bose-Hubbard Hamiltonian.

###### Acknowledgements.

HH acknowledges financial support by the German Research Foundation (DFG), grant no. HE 6312/1-1. We thank Jérôme Dorignac and David K. Campbell for useful discussions.## References

- Bloch (2005) I. Bloch, Nat. Phys. 1, 23 (2005).
- Bakr et al. (2009) W. S. Bakr, J. I. Gillen, A. Peng, S. Fölling, and M. Greiner, Nature 461, 74 (2009).
- Bakr et al. (2010) W. S. Bakr, A. Peng, M. E. Tai, R. Ma, J. Simon, J. I. Gillen, S. Fölling, L. Pollet, and M. Greiner, Science 329, 547 (2010).
- Simon et al. (2011) J. Simon, W. S. Bakr, R. Ma, M. E. Tai, P. M. Preiss, and M. Greiner, Nature 472, 307 (2011).
- Morsch and Oberthaler (2006) O. Morsch and M. Oberthaler, Rev. Mod. Phys. 78, 179 (2006).
- Pethick and Smith (2008) C. J. Pethick and H. Smith, Bose-Einstein Condensation in Dilute Gases (Cambridge University Press, Cambridge, UK, 2008).
- Flach and Gorbach (2008) S. Flach and A. V. Gorbach, Phys. Rep. 467, 1 (2008).
- Campbell et al. (2004) D. K. Campbell, S. Flach, and Y. S. Kivshar, Phys. Today 57, 43 (2004).
- Hennig and Fleischmann (2013) H. Hennig and R. Fleischmann, Phys. Rev. A 87, 033605 (2013).
- Eisenberg et al. (1998) H. S. Eisenberg, Y. Silberberg, R. Morandotti, A. R. Boyd, and J. S. Aitchison, Phys. Rev. Lett. 81, 3383 (1998).
- Morandotti et al. (1999) R. Morandotti, U. Peschel, J. S. Aitchison, H. S. Eisenberg, and Y. Silberberg, Phys. Rev. Lett. 83, 2726 (1999).
- Trías et al. (2000) E. Trías, J. J. Mazo, and T. P. Orlando, Phys. Rev. Lett. 84, 741 (2000).
- Ustinov (2003) A. V. Ustinov, Chaos 13, 716 (2003).
- Schwarz et al. (1999) U. T. Schwarz, L. Q. English, and A. J. Sievers, Phys. Rev. Lett. 83, 223 (1999).
- Sato and Sievers (2004) M. Sato and A. J. Sievers, Nature 432, 486 (2004).
- Eiermann et al. (2004) B. Eiermann, T. Anker, M. Albiez, M. Taglieber, P. Treutlein, K.-P. Marzlin, and M. K. Oberthaler, Phys. Rev. Lett. 92, 230401 (2004).
- Hennig et al. (2012) H. Hennig, D. Witthaut, and D. Campbell, Phys. Rev. A 86, 051604 (2012).
- Sakmann et al. (2009) K. Sakmann, A. I. Streltsov, O. Alon, and L. S. Cederbaum, Phys. Rev. Lett. 103, 220601 (2009).
- Trujillo-Martinez et al. (2009) M. Trujillo-Martinez, A. Posazhennikova, and J. Kroha, Phys. Rev. Lett. 103, 105302 (2009).
- Trombettoni and Smerzi (2001a) A. Trombettoni and A. Smerzi, Phys. Rev. Lett. 86, 2353 (2001a).
- Christodoulides et al. (2003) D. N. Christodoulides, F. Lederer, and Y. Silberberg, Nature 424, 817 (2003).
- Fleischer et al. (2003a) J. W. Fleischer, M. Segev, N. K. Efremidis, and D. N. Christodoulides, Nature 422, 147 (2003a).
- Albiez et al. (2005) M. Albiez, R. Gati, J. Fölling, S. Hunsmann, M. Cristiani, and M. Oberthaler, Phys. Rev. Lett. 95, 010402 (2005).
- Cataliotti et al. (2001) F. S. Cataliotti, S. Burger, C. Fort, P. Maddaloni, F. Minardi, A. Trombettoni, A. Smerzi, and M. Inguscio, Science 293, 843 (2001).
- Fleischer et al. (2003b) J. Fleischer, J. Fleischer, T. Carmon, T. Carmon, M. Segev, M. Segev, N. Efremidis, N. Efremidis, D. Christodoulides, and D. Christodoulides, Phys. Rev. Lett. 90, 023902 (2003b).
- Morsch et al. (2001) O. Morsch, J. Müller, M. Cristiani, D. Ciampini, and E. Arimondo, Phys. Rev. Lett. 87, 140402 (2001).
- Burger et al. (1999) S. Burger, K. Bongs, S. Dettmer, W. Ertmer, K. Sengstock, A. Sanpera, G. V. Shlyapnikov, and M. Lewenstein, Phys. Rev. Lett. 83, 5198 (1999).
- Smerzi et al. (2002) A. Smerzi, A. Trombettoni, P. Kevrekidis, and A. Bishop, Phys. Rev. Lett. 89, 170402 (2002).
- Anker et al. (2005) T. Anker, M. Albiez, R. Gati, S. Hunsmann, B. Eiermann, A. Trombettoni, and M. Oberthaler, Phys. Rev. Lett. 94, 020403 (2005).
- Iwanow et al. (2004) R. Iwanow, R. Schiek, G. Stegeman, T. Pertsch, F. Lederer, Y. Min, and W. Sohler, Phys. Rev. Lett. 93, 113902 (2004).
- Eiermann et al. (2003) B. Eiermann, P. Treutlein, T. Anker, M. Albiez, M. Taglieber, K.-P. Marzlin, and M. Oberthaler, Phys. Rev. Lett. 91, 060402 (2003).
- Reinhard et al. (2013) A. Reinhard, J.-F. Riou, L. Zundel, D. Weiss, S. Li, A. Rey, and R. Hipolito, Phys. Rev. Lett. 110, 033001 (2013).
- Bersch et al. (2012) C. Bersch, G. Onishchukov, and U. Peschel, Phys. Rev. Lett. 109, 093903 (2012).
- Dienst et al. (2013) A. Dienst, E. Casandruc, D. Fausti, L. Zhang, M. Eckstein, M. Hoffmann, V. Khanna, N. Dean, M. Gensch, and S. Winnerl, nature materials 12, 535 (2013).
- Dong et al. (2013) G. Dong, J. Zhu, W. Zhang, and B. A. Malomed, Phys. Rev. Lett. 110, 250401 (2013).
- Franzosi et al. (2011) R. Franzosi, R. Livi, G. Oppo, and A. Politi, Nonlinearity 24, R89 (2011).
- Gericke et al. (2008) T. Gericke, P. Würtz, D. Reitz, T. Langen, and H. Ott, Nat. Phys. 4, 949 (2008).
- Christodoulides and Eugenieva (2001) D. N. Christodoulides and E. D. Eugenieva, Phys. Rev. Lett. 87, 233901 (2001).
- Trombettoni and Smerzi (2001b) A. Trombettoni and A. Smerzi, J. Phys. B 34, 4711 (2001b).
- Buonsante and Penna (2008) P. Buonsante and V. Penna, J. Phys. A 41, 175301 (2008).
- Oxtoby and Barashenkov (2007) O. Oxtoby and I. Barashenkov, Phys. Rev. E 76, 036603 (2007).
- Hennig et al. (2010) H. Hennig, J. Dorignac, and D. Campbell, Phys. Rev. A 82, 053604 (2010).
- Aubry (1997) S. Aubry, Physica D 103, 201 (1997).
- Marin and Aubry (1996) J. Marin and S. Aubry, Nonlinearity 9, 1501 (1996).
- Proville and Aubry (1999) L. Proville and S. Aubry, Eur. Phys. J. B 11, 41 (1999).
- Campbell et al. (1983) D. Campbell, J. Schonfeld, and C. Wingate, Physica D 9, 1 (1983).
- Campbell (1980) D. K. Campbell, Annals of Physics 129, 249 (1980).
- Rumpf (2004) B. Rumpf, Phys. Rev. E 70, 016609 (2004).
- Kivshar and Campbell (1993) Y. S. Kivshar and D. K. Campbell, Phys. Rev. E 48, 3077 (1993).
- Ng et al. (2009) G. S. Ng, H. Hennig, R. Fleischmann, T. Kottos, and T. Geisel, New J. Phys. 11, 073045 (2009).

## Appendix A Supplemental Material

### a.1 A cut through the dynamical phase diagram

The dynamical phase diagram of Fig. 2 shows a remarkable transition from diffusive to moving breather behavior back to diffusive motion and than to self trapping for . In Fig. 4 we shows curves of the order parameters as a function of the nonlinearity for , illustrating the sharp features in the order parameters at the transitions between the different regimes.

### a.2 Definition of the threshold

In order to find a suitable upper threshold of the maximal local norm in the diffusive regime we calculate the cumulative probability density function (pdf) of in the diffusive case where the single site norms are known to be exponentially distributed Ng et al. (2009). Assuming statistical independents in this regime gives the pdf of the maximum single site norm as . A DB however is localized over several sites around the maximum. Therefore we need to calculate the PDF of the local norm within a range of sites on either side of the maximum. The PDF of the sum of two single-site norms is given by the convolution of the two PDFs. For sites with exponentially distributed norms we can thus calculate the PDF of norms iteratively by , where denotes convolution. The PDF of the local norm in the range of sites around the maximum can consequently be expressed as . Let us examine an initial condition which leads to vanishing current . We will consider it to be self-trapped if its evolution leads to a maximum local norm that is sufficiently unlikely to be found in the diffusive state. We define as the threshold the norm for which the PDF falls below . We find for and for . However, our results are not sensitive to the exact value of this parameter . Assuming yields for and for which result in the dynamical phase diagrams shown in Fig. 5, which differ only minimally from those of Fig. 2b.