# Quasi-adiabatic dynamics of ultracold bosonic atoms in a one-dimensional optical superlattice

## Abstract

We study the quasi-adiabatic dynamics for a one-dimensional system of ultracold bosonic atoms loaded in an optical superlattice. Focusing on a slow linear variation in time of the superlattice potential, the system is driven from a conventional Mott insulator phase, to a superlattice-induced Mott insulator, crossing in-between a gapless critical superfluid region. Due to the presence of a gapless region, a number of defects depending on the velocity of the quench appear. Our findings suggest a power-law dependence similar to the Kibble-Zurek mechanism for intermediate values of the quench rate. For the temporal ranges of the quench dynamics that we considered, the scaling of defects depends non trivially on the width of the superfluid region.

###### pacs:

03.75.Nt, 05.10.Cc, 05.30.Jp, 73.43Nq## I Introduction

Ultracold atomic gases in optical lattices provide a unique platform to probe a wide range of quantum phenomena with a high degree of controllability. The success of the Bose-Hubbard (BH) Hamiltonian (1); (2) in elucidating the behavior of ultracold bosons in optical lattices has stimulated a great deal of interest, from theoretical and experimental points of view. In particular, the seminal work on the Mott insulator (MI) to superfluid (SF) quantum phase transition (QPT) (3) paved the way for a number of studies that led to the observation of various exotic quantum phases (4); (5).

Non-equilibrium quantum effects in such systems can be investigated by varying in time parameters such as the optical lattice depth or the magnetic field close to a Feshbach resonance. Such experimental possibilities have spurred a renewed interest in the study of both sudden and quasi-adiabatic quenches (6); (7). The latter would provide important insights into non-equilibrium quantum phase transitions. In the presence of a ground-state energy gap always being finite during a very slow time evolution, the adiabatic theorem guarantees that the system will remain in the instantaneous ground state of the corresponding time-dependent Hamiltonian. However if a gapless region is crossed, the system will be unable to stay in its equilibrium ground state, regardless of how slowly it is quenched. The non-adiabatic evolution inevitably excites the system and a number of defects will appear in the evolved state. The mechanism of such defects formation had been first addressed by Kibble and Zurek (KZ) in the context of classical phase transitions in the early universe (8); (9), and more recently extended to the quantum regime for the case of adiabatic quenches across a single quantum critical point (10); (11).

The possibility to apply this kind of quenches has led to a number of theoretical studies addressing different types of many-body systems, including spin chains and ultracold quantum gases (see, e.g., Refs. (12); (13); (14); (15); (16); (17); (18); (19); (20); (21); (22); (23); (24); (25); (26); (27); (28); (29); (30); (31); (32); (33)). Despite the large body of work in this field, several aspects involving the response of such systems to slow quenches have not been completely understood and deserve further investigation. It is believed that, in presence of non-isolated critical points or of extended critical regions, the validity of the KZ mechanism is a priori not obvious, even if in some cases it is still possible to predict the defect density by identifying a dominant critical point, or by using scaling arguments (20); (22); (23).

The dynamics of ultracold bosons in an optical lattice subjected to a quasi-adiabatic quench has been theoretically analyzed for the MI-SF as well as the reverse transition (12); (14). The SF-MI transition for a slow quench has been also studied taking the effects of the parabolic trapping potential into account (30). These results highlighted the emergence of a scaling behavior for the characteristic length scale as a function of the quench rate, which is well approximated by a power law. However it has been later shown that, for a phase transition of the Kosterlitz-Thouless type (as is the case for the MI-SF transition in the one-dimensional BH model), the exponents depend on that rate, and are generally different from the KZ prediction, based on the critical exponents that are relevant for asymptotically long quench times (33). Experimental evidence in support of the growth of the condensate excitations with a power-law dependence on the quench rate has been observed for the MI-SF transition for ultracold bosons in an optical lattice (31). Furthermore, for the same transition in a similar system, the observation of the emergence of coherence and a power-law dependence of the correlation length on the quench rate for intermediate quenches has recently been reported (34).

The feasibility of superposing different optical lattices with distinct frequencies (35) also made the study of local relaxation dynamics possible in such superlattice setups (36); (37). An interesting property of these composite structures when they are loaded with bosonic atoms is that they can facilitate the generation of multiple lattice-modulated MI phases, which can isolate SF regions in the parameter space of the system (38); (39). To the best of our knowledge, slow quenches for QPTs in optical superlattices starting and ending in insulating phases, and crossing a superfluid region in between, have not been addressed so far, and this is the focus of our present work.

Here we consider a one-dimensional (1D) system of ultracold bosonic atoms loaded in an optical superlattice, formed by two superimposed optical standing waves with different frequencies. At zero temperature, this system exhibits different quantum phases: MI, SF and a superlattice-induced MI (SLMI) with periodically modulated onsite occupation (38); (39). The superlattice potential is constructed as to vary linearly in time, and is chosen in such a way that it crosses a gapless region in between two insulating regions. We consider the formation of defects in the final state after the quench, and demonstrate a non-trivial scaling of the excess energy as a function of the quench rate. We tackle this problem by means of the time-dependent density matrix renormalization group (DMRG) method, in the formalism of matrix-product-state ansatz (40).

The paper is organized as follows. We start introducing our model and discussing the static properties of its ground-state phase diagram that are relevant to the ongoing discussion (Sec. II). In Sec. III we define our dynamical protocol and discuss the formation of defects and the behavior of two-point correlation functions at the end of the protocol. Finally, in Sec. IV we draw our conclusions.

## Ii Model and phase diagram

The model is described by the following Hamiltonian:

(1) |

where denote the creation and annihilation operators on site satisfying the usual bosonic commutation rules, with being the corresponding number operator. The parameter denotes the hopping amplitude, is the on-site repulsive interaction strength and quantifies the superlattice potential depth. For the period-two optical superlattice that we have considered, has a finite value of for odd sites, and it is zero for even sites. Hereafter we work in units of and set as the energy scale. The critical value for the MI-SF transition is located at for (41), as shown in Fig. 1.

The ground-state phase diagram of model (1) has been studied by means of mean field theory (38), quantum Monte Carlo techniques (42) and the DMRG method (39). In 1D and at integer filling , this is given in the - plane as in Fig. 1. Here we identify the various phases by analyzing the behavior of the ground-state energy gap as a function of the system size in the following way. First we observe that, for a MI, this is finite and coincides with the charge gap , that is the difference between the energy cost to add () and to remove () a particle from the system. On a chain of finite length , the numerical evaluation of the Mott gap has been thus obtained by performing three DMRG iterations, with projections on different number sectors . The corresponding ground states respectively give the desired energies , . In the SF region, this gap vanishes as the inverse of the system size . The critical points have been extracted as those in which the product for the smallest () and the largest () considered lengths differed more than four percent (see the insets of Fig. 1).

The nature of the insulating phase (MI or SLMI) depends on the relative strength of and . A SF region is present in between these two insulators, and its width decreases with increasing . The intervening SF phase arises because of the competition between the superlattice potential and the on-site two-body interaction . For , there is only a transition from a gapless SF to a gapped SLMI at a critical value of . But for , there are three possible scenarios. If is much smaller than , the system remains in the MI phase. When it becomes comparable to , the system makes a transition to the SF phase, and for large values of it enters the SLMI phase.

## Iii Quasi-adiabatic quench dynamics

In order to probe the slow quench dynamics of ultracold atoms in an optical superlattice, Eq. (1), it is necessary to analyze the excitations that are generated when the gapless SF phase is crossed. In view of the specific features of the phase diagram, it is convenient to fix a value of and increase the parameter , so as to drive the system across the MI-to-SF and then the SF-to-SLMI phase transitions. As can be seen in Fig. 1, the width of the critical region changes with . This reflects into a non trivial dependence of the rate of defects generation with , as discussed later.

We adopt a linear variation in time of the superlattice potential , which is given by

(2) |

Here denotes the time of the quench, while and are respectively the initial and the final values of the superlattice potential. After fixing the value of , we choose and such that the system starts from a MI and ends in a SLMI phase (except for the cases with , where there is no MI phase). In between the initial and the final insulating phases, there is a SF region, whose width depends on . Due to the presence of a gapless region at the thermodynamic limit , a certain number of defects in the final state after the evolution will appear, no matter how slowly the quenching is performed (43). Below we shed light on these defects.

The system wavefunction evolves according to the time-dependent Schrödinger equation. We computed at the final time after the quench (2), using a time-evolving block-decimation strategy (40); (44). We simulated systems up to sites with no more than bosons per site, and used open boundary conditions. The time interval has been discretized into many slices of time-length , where is assumed to be constant. The corresponding time evolution operator has been expanded by means of a sixth-order Suzuki-Trotter decomposition (45). We have been able to consider a time-step , and reach a threshold for the discarded states , by using a bond-link dimension for all our simulations (46).

### iii.1 Excess energy

To quantify the defects generation due to the non-adiabatic crossing of the SF region during the time evolution, we focus on the residual energy , defined as the excess energy above the ground state:

(3) |

where denotes the energy of the system at time , while is the instantaneous ground-state energy for Hamiltonian at time . This quantity serves as analog of the defect density originally considered by Kibble and Zurek (see, e.g., Refs. (18); (19); (20)). Let us now discuss its behavior after a time modulation of the superlattice depth , from to , as dictated by Eq. (2). In particular we focus on the scaling of with for different values of interaction .

The typical scenario is depicted in Fig. 2a, where we are able to distinguish three distinct behaviors as a function of the quench time. For , the dynamics is ruled by the adiabatic regime typical of slow quenches: the time-evolved wavefunction remains very close to the instantaneous ground state of . The residual energy follows a power-law behavior

(4) |

with . This exponent can be obtained within the Landau-Zener formalism (47), where the quantum evolution is studied by means of an effective two-level approximation with an avoided level crossing. We point out that the adiabatic regime can occur only in the presence of an instantaneous ground-state energy gap which remains finite along the sweeping (2) (i.e., for quenches much slower than the inverse square of the minimum crossed gap, ). Therefore this is a behavior related to a finite-size effect, which disappears at the thermodynamic limit where the gap in the superfluid region is rigorously zero. In the opposite regime of fast quenches (), the dynamics is strongly non adiabatic and the initial state is essentially frozen during the evolution. The excess energy thus saturates with . The intermediate regime in between is the most interesting one, since it is crucially affected by the critical properties of the region crossed by the system.

In the intermediate regime, our data display a power-law scaling of the type in Eq. (4). This fairly agrees with the general behavior predicted by KZ and verified in many cases, when the system is adiabatically driven across isolated quantum critical points (6). The KZ mechanism roughly identifies two types of evolution, either adiabatic or impulsive, according to the distance from the critical point. The time (distance from the critical point) at which the system switches to the impulsive regime depends on the quench velocity. This simple argument indeed predicts a power-law scaling form for as a function of , with a rate expressed in terms of the critical exponents dictating the phase transition. However, for the crossing of continuous phase transitions with extended critical regions, as it is in our case, the KZ scaling may still give insightful information, but cannot be regarded as ultimately predictive. In the specific case of a Kosterlitz-Thouless transition analogous to that occurring in the BH model at , it has been shown that the exponentially slow gap closure induces a power-law scaling which generally depends on the quench rate. The exponent differs from that obtained with the usual KZ mechanism using the critical exponents of the transition (33). In our specific situation, the system is quenched from a gapped to a gapless phase, and then to another gapped phase, thus crossing two QPTs and an extended critical (SF) region (Fig. 1, for fixed ). This is an even more complex scenario, where it is impossible to grasp the quantitative power-law scaling predicted by KZ, and hence we expect the emergence of a more complex and inhomogeneous behavior, in terms of the size of the crossed critical region.

In the case of we clearly identify an intermediate scaling region with , as extracted from fits of the numerical data (see the solid line in Fig. 2a, denoting the best fit of the data series at the size ). Our results do not display a significant dependence of on . Notice also that the excess energy per unit length is universal in the scaling region, whose width increases with , due to the gap closure in the SF phase (inset of Fig. 2a). A qualitatively similar behavior is observed for a quench from SF to SLMI phase, at fixed —Fig. 2b, blue circles, shows an example with . The value in the intermediate region obtained for that case is considerably larger than that for .

These observations support the evidence that any appropriate scaling analysis should depend non trivially on , while a simple KZ argument cannot predict this complex behavior. We point out that it is also not guaranteed that does not change with the quench rate, as theoretically predicted for the MI-SF transition at (33). In the range of we were able to address, we did not observe such dependence. However, as reported in the inset of Fig. 2b, for quenches across the phases MI SF SLMI with , we found a rather complex dependence of on , and hence on the width of the intermediate SF region. In particular the defects production rate decreases monotonically as a function of the time during which the system is crossing the gapless region.

We also checked the dependence of on the starting and ending points inside the insulator (we varied and , for fixed ). Results in Fig. 3 indicate a tendency toward a slight decrease of , if the gapped region crossed by the quench increases.

### iii.2 Correlation functions

Finally we examined the behavior of the two-point correlation function at the end of the quasi-adiabatic dynamics, and observed an exponential decay with the distance , as shown in Fig. 4. The two points and have been chosen in a symmetric way with respect to the center of the chain in order to minimize boundary effects, such that , for odd , and , for even (for instance, for sites, corresponds to , ; corresponds to , ; corresponds to , , and so on).

In the inset we plot the correlation function

(5) |

This clearly exhibits a non-monotonic behavior as a function of the quench rate (48). In particular we notice that is increasing initially in the intermediate scaling region. This can be attributed to the persistence of quasi-long range order which the system gained while quenching through the SF phase. But after some critical value of , we observe that starts to decrease. Such a behavior is ascribed to the onset of adiabatic regime. Since the system ends up in an insulating phase (SLMI), where the correlation function decays exponentially, we expect in the adiabatic regime the correlation length to be small. In particular, the value of at which the transition from the intermediate regime to the adiabatic regime takes place, obtained from the analysis of the residual energy (Fig. 2) coincides with that seen in the correlation length (Fig. 4).

## Iv Summary

We have theoretically analyzed the slow quench dynamics of ultracold bosonic atoms in a one-dimensional optical superlattice. By considering a linear time dependence of the superlattice potential, we showed that, when crossing a gapless superfluid region, the system has the tendency to generate defects. This fact is due to the adiabaticity loss during the time evolution, which occurs even when the system is quenched very slowly. Our results show a complex dependence of the rate of defect generation on the quench velocity, which cannot be understood in terms of the Kibble-Zurek physics underlying the crossing of a single critical point.

From an experimental point of view, the behavior of the excess energy could be verified by means of time-of-flight measurements of the correlation length. This, in turn, may reveal itself as a simple indicator of the presence or absence of the power-law scaling regime for the defects production as a function of the quench velocity. Trapping ultracold bosonic atoms in optical standing waves is probably the most favourable setup to probe this kind of physics. Recent experiments have already verified the power-law behavior (31); (34) in the framework of the Bose-Hubbard model. Moreover, the capability that has been demonstrated by a variety of out-of-equilibrium experiments with great accuracy and for large coherence times, ranging from sudden quenches to adiabatic variation of the system’s parameters, puts our results arising from the quenching of the superlattice potential in a favourable position for verification in the laboratory.

###### Acknowledgements.

We thank R. Fazio for useful discussions. AD and BPD thank D. Sen and A. Dutta for stimulating discussions. AD acknowledges the Academy of Finland through its Center of Excellence Programme (2012-2017) and under Project No. 13251748. DR acknowledges support from the Italian MIUR through FIRB Project RBFR12NLNA. The computations were performed on the Intel Cluster at the Indian Institute of Astrophysics.### References

- M. P. A. Fisher, P. B. Weichman, G. Grinstein, and D. S. Fisher, Phys. Rev. B 40, 546 (1989).
- D. Jaksch, C. Bruder, J. I. Cirac, C. W. Gardiner, and P. Zoller, Phys. Rev. Lett. 81, 3108 (1998).
- M. Greiner, O. Mandel, T. Esslinger, T. W. Hänsch, and I. Bloch, Nature 415, 39 (2002); ibid. 419, 51 (2002).
- M. Lewenstein, A. Sanpera, V. Ahufinger, B. Damski, A. Sen(De), and U. Sen, Adv. Phys. 56, 243 (2007).
- I. Bloch, J. Dalibard, and W. Zwerger, Rev. Mod. Phys. 80, 885 (2008).
- J. Dziarmaga, Adv. Phys. 59, 1063 (2010).
- A. Polkovnikov, K. Sengupta, A. Silva, and M. Vengalattore, Rev. Mod. Phys. 83, 863 (2011).
- T. W. B. Kibble, J. Phys. A: Math. Gen. 9, 1387 (1976).
- W. H. Zurek, Nature 317, 505 (1985).
- W. H. Zurek, U. Dorner, and P. Zoller, Phys. Rev. Lett. 95, 105701 (2005).
- A. Polkovnikov, Phys. Rev. B 72, 161201(R) (2005).
- S. R. Clark and D. Jaksch, Phys. Rev. A 70, 043612 (2004).
- J. Dziarmaga, Phys. Rev. Lett. 95, 245701 (2005).
- R. Schützhold, M. Uhlmann, Y. Xu, and U. R. Fischer, Phys. Rev. Lett. 97, 200601 (2006).
- R. W. Cherng and L. S. Levitov, Phys. Rev. A 73, 043614 (2006).
- F. M. Cucchietti, B. Damski, J. Dziarmaga and W. H. Zurek, Phys. Rev. A 75, 023603 (2007).
- L. Cincio, J. Dziarmaga, M. M. Rams, and W. H. Zurek, Phys. Rev. A 75, 052321 (2007).
- T. Caneva, R. Fazio and G. E. Santoro, Phys. Rev. B 76, 144427 (2007).
- T. Caneva, R. Fazio and G. E. Santoro, Phys. Rev. B 78, 104426 (2008).
- F. Pellegrini, S. Montangero, G. E. Santoro, and R. Fazio, Phys. Rev. B 77, 140404(R) (2008).
- K. Sengupta, D. Sen, and S. Mondal, Phys. Rev. Lett. 100, 077204 (2008).
- S. Deng, G. Ortiz, and L. Viola, Europhys. Lett. 84, 67008 (2008).
- U. Divakaran, A. Dutta, and D. Sen, Phys. Rev. B 78, 144301 (2008).
- A. Polkovnikov and V. Gritsev, Nat. Phys. 4, 477 (2008).
- C. De Grandi, R. A. Barankov and A. Polkovnikov, Phys. Rev. Lett. 101, 230402 (2008).
- E. Canovi, D. Rossini, R. Fazio, and G. E. Santoro, J. Stat. Mech. P03038 (2009).
- V. Mukherjee and A. Dutta, Europhys. Lett. 92, 37004 (2010).
- M. Haque and F. E. Zimmer, Phys. Rev. A 87, 033613 (2013).
- C. Trefzger and K. Sengupta, Phys. Rev. Lett. 106, 095702 (2011).
- J.-S. Bernier, G. Roux, and C. Kollath, Phys. Rev. Lett. 106, 200601 (2011); J.-S. Bernier, D. Poletti, P. Barmettler, G. Roux, and C. Kollath, Phys. Rev. A 85, 033641 (2012).
- D. Chen, M. White, C. Borries, and B. DeMarco, Phys. Rev. Lett. 106, 235304 (2011).
- D. Pertot, A. Sheikhan, E. Cocchi, L. A. Miller, J. E. Bohn, M. Koschorreck, M. Köhl, and C. Kollath, Phys. Rev. Lett. 113, 170403 (2014).
- J. Dziarmaga and W. H. Zurek, Sci. Rep. 4, 5950 (2014).
- S. Braun, M. Friesdorf, S. S. Hodgman, M. Schreiber, J. P. Ronzheimer, A. Riera, M. del Rey, I. Bloch, J. Eisert, and U. Schneider, Proc. Natl. Acad. Sci. 112, 3641 (2015).
- S. Fölling, S. Trotzky, P. Cheinet, M. Feld, R. Saers, A. Widera, T. Müller, and I. Bloch, Nature 448, 1029 (2007).
- M. Cramer, A. Flesch, I. P. McCulloch, U. Schollwöck, and J. Eisert, Phys. Rev. Lett. 101, 063001 (2008); A. Flesch, M. Cramer, I. P. McCulloch, U. Schollwöck, and J. Eisert, Phys. Rev. A 78, 033608 (2008).
- S. Trotzky, Y.-A. Chen, A. Flesch, I. P. McCulloch, U. Schollwöck, J. Eisert, and I. Bloch, Nat. Phys. 8, 325 (2012).
- A. Dhar, M. Singh, R. V. Pai, and B. P. Das, Phys. Rev. A 84, 033631 (2011).
- A. Dhar, T. Mishra, R. V. Pai, and B. P. Das, Phys. Rev. A 83, 053621 (2011).
- U. Schollwöck, Ann. Phys. 326, 96 (2011).
- T. D. Kühner, S. R. White, and H. Monien, Phys. Rev. B 61, 12474 (2000).
- V. G. Rousseau, D. P. Arovas, M. Rigol, F. Hebert, G. G. Batrouni, and R. T. Scalettar, Phys. Rev. B 73, 174516 (2006).
- Note that, by first taking the thermodynamic limit , we are guaranteed that the system in the SF phase will be gapless, therefore the evolution will never be adiabatic and defects will necessarily be formed at any finite value of . On the contrary, by first taking the asymptotically slow quench limit , the evolution will stay adiabatic at any size , and defects will not form.
- G. Vidal, Phys. Rev. Lett. 93, 040502 (2004).
- M. Suzuki, Commun. Math. Phys. 51, 183 (1976); H. F. Trotter, Proc. Am. Math 10, 545 (1959).
- We have checked that, on the scale of all the figures presented here, our results do not present significant deviations, if we take a time step corresponding to half of the value used throughout the simulations shown in this paper. This produces a discrepancy in the decay rates which turns out to be less than .
- C. Zener, Proc. R. Soc. A 137, 696 (1932); L. D. Landau and E. M. Lifshitz, Quantum Mechanics, Vol. 3 (Pergamon Press, 1958).
- The two summations in Eq. (5) have been performed neglecting ten lattice sites per each border, so to reduce their influence (i.e., we considered ).