New Journal of Physics

New Journal of Physics

Quick search Find article
Quick search
Find article
New J. Phys. 9 (2007) 221
doi:10.1088/1367-2630/9/7/221
PII: S1367-2630(07)44932-1

Fast initialization of a high-fidelity quantum register using optical superlattices

B Vaucher1, S R Clark, U Dorner and D Jaksch

Clarendon Laboratory, University of Oxford, Parks Road, OX1 3PU, United Kingdom

1 Author to whom any correspondence should be addressed.

Email: benoit.vaucher@merton.ox.ac.uk

Received 27 February 2007
Published 10 July 2007

Abstract. We propose a method for the fast generation of a quantum register of addressable qubits consisting of ultracold atoms stored in an optical lattice. Starting with a half filled lattice we remove every second lattice barrier by adiabatically switching on a superlattice potential which leads to a long wavelength lattice in the Mott insulator state with unit filling. The larger periodicity of the resulting lattice could make individual addressing of the atoms via an external laser feasible. We develop a Bose-Hubbard-like model for describing the dynamics of cold atoms in a lattice when doubling the lattice periodicity via the addition of a superlattice potential. The dynamics of the transition from a half filled to a commensurately filled lattice is analysed numerically with the help of the time evolving block decimation algorithm and analytically using the Kibble–Zurek theory. We show that the timescale for the whole process, i.e. creating the half filled lattice and subsequent doubling of the lattice periodicity, is significantly faster than adiabatic direct quantum-freezing of a superfluid into a Mott insulator for large lattice periods. Our method therefore provides a high-fidelity quantum register of addressable qubits on a fast timescale.

Contents

1. Introduction

Systems of cold atoms trapped in optical lattices provide the unique opportunity to coherently manipulate a large number of atoms [1]–[3]. The remarkable degree of experimental control offered by these systems, as well as the possibility to use the internal hyperfine states of the atoms to encode qubits, make them particularly suited for quantum information processing (QIP). In this context, optical lattices in a Mott-insulating (MI) state with unit filling can be viewed as the realization of a quantum register, and it is possible to collectively manipulate the qubits stored in such a register experimentally [4, 5]. However, in many quantum computing schemes based on neutral atoms stored in optical lattices the application of single qubit gates [3] or single qubit measurements [6] requires the ability to address single atoms with a focused laser beam. This remains experimentally challenging since these operations have to be performed without perturbing the state of other atoms in their vicinity.

A number of strategies have been proposed to circumvent this problem by using global operations [7]–[9], e.g. via `marker atoms' which are moved to a particular lattice site and interact with the corresponding register qubit such that an external laser affects only that qubit [10]. Another way is simply to use a quantum register in which the atoms are distant enough such that they can be addressed individually by a laser. This method requires an optical lattice with filling factor n  =  1 in the MI state with a sufficiently large distance between the atomsNote2 . The initialization time of a MI state is proportional to the tunnelling time of the atoms between neighbouring sites. Therefore, by using the conventional method of quantum-freezing a superfluid (SF) state [4], [12]–[14], it scales exponentially with the lattice spacing [15]–[17].

In this paper, we propose an alternative method to generate a long wavelength lattice with one atom per site in the MI state. Starting with a one-dimensional (1D) lattice with a short period—and hence a short initialization time—and filling factor n  =  1/2 we remove every second potential barrier by adiabatically turning on a superlattice. This superlattice potential has already been experimentally realized [18, 19]. This procedure eventually leads to a long wavelength lattice in which the periodicity has been doubled and where the atoms are in a MI state with n  =  1 (see figure 1). This scheme does not require changing the angles of the intersecting laser beams. Furthermore, we show that our method allows for the initialization of a MI state with unit filling factor on a timescale which, although scaling exponentially with the final lattice spacing, is approximately one order of magnitude smaller than the direct quantum-freezing method. Although we only consider the case of a homogeneous lattice, the results presented in this paper extend to the case of weak harmonic confinements, quartic [20] and box traps [21].

Figure 1

Figure 1. (a) The initial profile of the lattice is associated with the value of s  =  1. The periodicity of the lattice is then progressively doubled (b) until it reaches its final profile (c) associated with the parameter s  =  0. The number of sites M corresponds to the number of unit cells in the large lattice limit. By starting with a filling factor of n  =  1/2, this procedure leads to a lattice with filling factor of n  =  1.

This paper is structured as follows. In section 2, we introduce the model used to describe the system dynamics. In section 3.1, we discuss ground state properties, particularly two-site correlation functions and quasi-momentum distributions. Also, we present and discuss numerical results for the probability of staying in the ground state during the transition depending on the speed of the ramping. In section 3.2, we apply the analytical Kibble–Zurek (KZ) theory and compare it with our numerical results. Finally, we summarize and conclude in section 4.

2. Model

We consider a gas of interacting ultracold bosonic atoms loaded into a 3D optical lattice. The lattice is formed by pairwise orthogonal standing wave laser fields and its optical potential is given by [22]

Equation (1)

Here VT is the depth of the potential in the y- and z-directions created by pairs of lasers with wavenumber k  =  2π/λ, wave length λ and period aT  =  λ/2. The extension to the case of different optical potentials in the y- and z-directions is straightforward. In the x-direction two pairs of laser beams with a long wavelength λL and short wavelength λS  =  λL/2 are applied. The potential in the x-direction is thus given by

Equation (2)

with V0 the depth of the lattice, kL  =  2π/λL and kS  =  2π/λS. The depths of the potentials will be expressed in units of the recoil energy ER  =  kL2/2m with m the mass of the atoms (taking hslash  =  1 throughout). The parameter s in [0,1] is determined by the relative intensities of the two pairs of lasers. By changing s from 1 to 0 the lattice profile is continuously transformed from a sinusoidal potential with a small period aS  =  λS/2 to one with a long period a  =  λL/2, thus halving the number of lattice sites per unit length along the x-direction (see figure 1). The lattice constant a corresponds to the size of a unit cell for s  <  1.Note3  We refer to the lattice profile with parameters s  =  1 and s  =  0 as to the small lattice limit and the large lattice limit, respectively.

The Hamiltonian of the system in second quantization reads

Equation (3)

where \skew2\hat{h}_0(\bar{\bf r})=-(1/{2 m }) \nabla^2 + V_{\rm OL}(\bar{\bf r}) is the one-particle Hamiltonian. The symbol \bar{\bf r}=({\bf r},s) represents the position variable {\bf r} and lattice parameter s. The interaction between the atoms is modelled by s-wave scattering with g  =  4πas/m where as is the s-wave scattering length. The bosonic field operators obey the usual commutation relations [\psi(\bar{\bf r}),\psi^{\dagger}({\bar{\bf r}^{\prime}})]=\delta({\bf r}-{\bf r^{\prime}}) with δ denoting the Dirac delta function.

We restrict our considerations to the case where VT is sufficiently large so that motion along the y- and z-directions is frozen. The dynamics of the system is then effectively 1D along the x-direction. As shown in figure 2(a) the two lowest bands of the Hamiltonian \hat{{h}}_{0,x}=-(1/2 m)\left({{\rm d}}/{{\rm d}x}\right)^2 + V_{\rm S}(x,s) are separated in energy by much less than the typical motional excitation energy E_{\rm ex} = \sqrt{4 V_0 E_{\rm R}}  [22] for values s ≈ 1.Note4  Therefore, despite assuming that atoms loaded into the lattice are ultracold, we have to consider the two lowest Bloch bands in x-direction to obtain an accurate description of the atomic dynamics. However, excitations to higher bands in the y- and z-directions are neglected in our investigations since we assume that the temperature of the atomic cloud is kBT  ll  Eex. We then expand the bosonic field operator as

Equation (4)

where M is the number of lattice sites and \hat{{\alpha}}^\dagger_i (α  =  a, b) creates a particle in the mode associated with the localized function \phi_{\alpha,i}(\bar{\bf r}) centred at site i. The mode functions \phi_{\alpha,i}(\bar{\bf r})=w_{\alpha,i}(x,s) W_{i,0}(y)W_{i,0}(z) are factorized into a product of well-localized Wannier functions (WFs) Wi,0 of the lowest Bloch band in the y- and z-directions [23] and mode functions wα,i(x,s) in x-direction.

Figure 2

Figure 2. Band structure along the x-direction in (a) the small and (b) large lattice limit for V0  =  10ER. The points represents the values of q for q  =  0, π/a, π/2a. The value of epsilonα,q correspond to the energy of single-particles with momentum q in the α-th Bloch band. In the small lattice limit, the two first Bloch bands are connected.

The aim of the next section is to describe the single particle dynamics in the tight binding (TB) approximation. If we were to use WFs for wα,i(x,s) this approximation would restrict our model to sinusoidal Bloch bands [24]. Because of the deviation of the lowest two bands from a sinusoidal dispersion relation (see figure 2(a)) when s ≈ 1 we instead use generalized WF (GWFs) for wα,i(x,s) (see appendix A for a detailed definition and a description of their properties) [25]. By exploiting the optimization procedure described in appendix A, we calculate GWFs wα,i(x,s) which are well localized at lattice sites i. Typical shapes of GWFs and the effects of optimizing their localization are shown in figure 3. We note that these GWFs are in general composed of superpositions of Bloch orbitals of both bands and are not related to WFs by a local transformation. Only when s ≈ 0 the wα,i(x,s) are equivalent to the WFs of the first and the second Bloch band, respectively. Finally, they are always (anti-)symmetric with respect to the centre of the lattice site i for α  =  a (α  =  b).

Figure 3

Figure 3. (a) Optimized (thick line) and non-optimized (dashed line) GWFs associated with the first mode for V0  =  30ER and the lattice profiles s  =  1. The area between the optimized and non-optimized mode functions is shaded in order to illustrate how the localization procedure reduces the spread of the optimized function. (b) The square of the mode functions for V0  =  30ER and s  =  1. By combining the mode functions shown in (b), we can construct two new mode functions corresponding to particles localized in either the left (c) or the right (d) well of a given site.

2.1. Single-particle Hamiltonian

Inserting the approximate field operator equation (4) into the first term of equation (3) yields the single-particle part of the Hamiltonian in terms of \hat{{a}}^\dagger_i and \hat{{b}}^\dagger_i. Applying the TB approximation, which amounts to keeping only nearest-neighbour hopping terms, the single-particle Hamiltonian can be approximated by

Equation (5)

where

Equation (6)

is the hopping matrix element between neighbouring sites along the x-axis and

Equation (7)

is the local on-site energy of a particle in mode α. Note that hopping between modes a and b within one site is not allowed by the symmetry properties of the GWFs. However, the inclusion of nonzero hopping matrix elements Jab and Jba is essential to accurately reproduce the single particle behaviour of the full Hamiltonian (3). The symmetry properties of WFs would not allow the inclusion of these terms [24].

For periodic boundary conditions, the parameters Vα and Jα,β can be found from the band structure without explicit calculation of the mode functions (see appendix B). The numerical values of Vα and Jα,β for V0  =  30ER are shown in figure 4(a). Using these parameters, we find that \hat{{H}}_0(s) very accurately reproduces the band structure of the exact Hamiltonian for all values of s, thus justifying the utilization of the TB approximation and the corresponding GWFs.

Figure 4

Figure 4. Parameters of the effective Hamiltonian \hat{{H}}_{\rm eff} as functions of s for a lattice depth of V0  =  30ER and VT  =  60ER with \tilde{a}=a\, a_{\rm T}^2. The hopping matrix elements (a) have been calculated using the method described in appendix B, while the on-site interaction energies (b) are calculated using optimized GWFs.

2.2. Interaction Hamiltonian

To calculate the interaction matrix elements of \hat{{H}} the explicit form of the localized GWFs is needed. We find that for V0  >  10ER, off-site interaction terms are at least two orders of magnitude smaller than on-site interactions. We therefore only keep the dominant on-site terms and find the interaction Hamiltonian \hat{{H}}_{\rm I}

Equation (8)

where \hat{{n}}^a_i=\hat{{a}}^\dagger_i \hat{{a}}_{i} and \hat{{n}}^b_i=\hat{{b}}^\dagger_i \hat{{b}}_{i} are the site-occupation number operators. The on-site interaction matrix elements are given by

Equation (9)

The numerical values of Uα,β as a function of the lattice profile s are shown in figure 4(b). Figure 4(b) shows that their values become equal as s → 1 for sufficiently large V0.

Combining the single- and two-particle contributions the effective Hamiltonian describing the system dynamics is given by

Equation (10)

By using \hat{{H}}_{\rm eff}(s) for s varying in time, we implicitly assume that the system adiabatically follows changes in the mode functions wα,i. For all dynamical calculations carried out in this work we have carefully chosen the time dependence of s so that such non-adiabatic contributions can safely be neglected.

2.3. Limiting cases

In the small lattice limit (s  =  1), the superpositions \varphi_{{\rm L},i}=(w_{a,i}-w_{b,i})/\sqrt{2} and \varphi_{{\rm R},i}=(w_{a,i}+w_{b,i})/\sqrt{2} correspond to mode functions localized in the left and in the right well of site i respectively (see figures 3(c) and (d)). The associated bosonic operators are defined by

Equation (11)

Given that the new mode functions are sufficiently localized within each well, the parameters of \hat{{H}}_{\rm eff} for s  =  1 can be written as

Equation (12)

where J=\int{\rm d}x \,\varphi_{{\rm R},i}^* \hat{{h}}_{0,x} \varphi_{{\rm L},i+1}, E=\int{\rm d}x \,\varphi_{{\rm L},i}^* \hat{{h}}_{0,x} \varphi_{{\rm L},i} and U=g \int{\rm d}x \,|\varphi_{{\rm R},i}|^4. Notice that the parameters shown in figure 4 are consistent with these equations. Expressing the Hamiltonian (10) in terms of the operators \hat{{c}}^\dagger_{\alpha,i} (α  =  L, R) and using the parameters (12) we find that

Equation (13)

where

Equation (14)

Therefore, as expected, \hat{{H}}_{\rm eff}(s=1) is equivalent to the one-band Bose–Hubbard model (BHM) [22] with 2M sites.

The mode functions keep their symmetry with the two peaks moving towards the centre of the cell when s is decreased. When s ≈ 0 is reached wa,i(x, s) and wb,i(x, s) are equivalent to the WFs of the first and the second Bloch band, respectively. In this limit we obtain a standard two band BHM for M sites.

3. Timescale for the preparation of the quantum register

In this section, we present numerical as well as analytical results characterizing the ground state properties of the system and the timescale necessary to initialize the quantum register.

3.1. Numerical results

The numerical calculations have been carried out using both the exact matrix representation of the Hamiltonian \hat{{H}}_{\rm eff} and the time-evolving Block decimation (TEBD) algorithm (see appendix C for details).

We have evaluated the ground state and dynamical properties of our system for two different values of g corresponding to different interaction regimes. These values have been chosen such that in the large lattice limit (with filling factor n  =  1) the ground state is a Mott-insulator state for g1 and g2 with Jaa/Uaa  <  0.3 [26]. In the small lattice limit, g1 and g2 produce ground states corresponding to a strongly interacting Tonks–Girardeau (TG) gas (U/J  =  215) and a SF (U/J  =  7), respectively.

3.1.1. Ground states properties In the large lattice limit, the ground state | {\psi_{0,{\rm L}}}\rangle of the system is populated exclusively by particles in the lowest Bloch band, i.e. a-mode particles. Therefore, in this limit, we only consider the one-particle density matrix given by

Equation (15)

where |ψrangle is the state of the system. The quasi-momentum distribution of particles in the a-mode is given by [27]

Equation (16)

As expected, in this limit and for commensurate filling n  =  1 the ground state is a Mott-insulator for both values of g (see figures 5(a)–(d)). The quasi-momentum distributions show that in this limit particles are uniformly distributed in the first band (see figures 5(a) and (c)).

Figure 5

Figure 5. Ground state one-particle density matrix and quasi-momentum distribution for M  =  12 lattice sites and the parameters of panels parameters. The panels (a)–(d) correspond to the large lattice limit and the panels (e)–(h) to the small lattice limit. The quasi-momentum distribution (a), (e) and one-particle density matrix (b), (f) are for g  =  g1. The quasi-momentum distribution (c), (g) and one-particle density matrix (d), (h) are for g  =  g2.

In the small lattice limit, the one-particle density matrix is given by

Equation (17)

where the operators \hat{{c}}^\dagger_{i^{\prime}} are constructed from the \hat{{a}}^\dagger_i and \hat{{b}}^\dagger_i operators using the transformations (11). The quasi-momentum distribution is given by

Equation (18)

In this limit (with filling factor n  =  1/2), the characteristics of the ground state | {\psi_{0,S}}\rangle depend on the value of g. For g  =  g1, the ratio U/J  =  215 and the ground state's correlations as well as the quasi-momentum distributions (see figures 5(e) and (f)) are characteristic of a TG gas (see e.g. [28] and references therein). The value g  =  g2 yields the ratio U/J  =  7 and the system exhibits the behaviour of a SF (see figures 5(g) and (h)), with particles occupying mainly the q  =  0 momentum state. Notice that 1D systems described by the BHM with a fixed filling factor n  =  1 cross the MI–SF phase boundary at the critical point (J/U)c ≈ 0.3 [26]. Thus, in our case the SF behaviour of the system in the small lattice limit is due to the fractional filling of the lattice.

3.1.2. Simulation of the dynamics Starting from a system with half-filling and a lattice profile s  =  1, we investigate the timescale required to obtain a nearly perfect MI state—or quantum register—with filling factor n  =  1 by ramping the lattice profile down to s  =  0. The quality of the register is determined by the fidelity

Equation (19)

defined as the overlap between the state of the system | {\psi}\rangle at the end of the ramping process and the ground state in the large lattice limit | {\psi_{0,{\rm L}}}\rangle. Furthermore, we calculate the fluctuations of the number of particles in the a-mode at site i

Equation (20)

where \left\langle \circ \right\rangle =\langle {\psi}|\circ | {\psi}\rangle. Since particle-number fluctuations are suppressed in a MI state, nonzero fluctuations indicate the presence of excitations, such as double occupancies or particles in the second band, in the final state.

We test different ramps by simulating the system dynamics between ti  =  0 to t_{\rm f}=\tau_{\rm Q} where τQ is the ramping time (the time required to complete the ramping process from s  =  1 to s  =  0). Here, each ramp σ corresponds to a function s  =  sσ(t) varying from sσ(0)  =  1 to sσ(τQ)  =  0.

For linear ramping we use s_{\rm lin}(t)=(\tau_{\rm Q}-t/E_{\rm R})/\tau_{\rm Q} . The fidelity and particle-number fluctuations obtained using this strategy for different quench times and g  =  g1 and g  =  g2 are shown in figures 6(a)–(d). The linear ramp is shown in figure 7(a).

Figure 6

Figure 6. Dynamical simulation of \hat{{H}}_{\rm eff} using the parameters shown in figure 4 for two different values of g. The ramping strategies sgap and sman are shown in figure 7(a). Simulation results for the fidelity and particle-number fluctuations using: (a) and (b) a linear ramp and g  =  g1; (c) and (d) a linear ramp and g  =  g2; (e) and (f) the ramp sgap and g  =  g1; (g) and (h) the ramp sman and g  =  g2. The vertical lines indicate the value of the particle tunneling time ( t_{\rm tun}\sim 1/J_{aa} ) in the large lattice limit.

Figure 7

Figure 7. (a) Different ramps used in our numerical calculations. (b) Transition probability between the ground and the first (labelled by e) and second (labelled by 2e) excited states as a function of time for the ramps slin and sman. The transition probabilities have been calculated via the exact diagonalization of \hat{{H}}_{\rm eff} for M  =  4 and a quench time of τQ  =  20/ER.

Another ramping strategy we use consists of adapting the velocity of the ramp proportionally to the energy gap between the ground and the first excited state. This ramp is denoted by sgap(t). We evaluate the ramp function sgap(t) numerically (see appendix D) for a system with M  =  4 sites and V0  =  30ER. The ramp sgap(t) for g  =  g1 is shown in figure 7(a). We expect this ramping strategy to be more efficient than the linear one, since accelerating the ramp when the gap is large while slowing it down when the gap is small should suppress transitions of particles to excited levels. Our numerical calculations have shown that when g  =  g1, the utilization of this ramp does indeed reduce the quench time needed to obtain a nearly perfect fidelity to a half of the tunnelling time in the large lattice limit (see figures 6(e) and (f)). For g  =  g2, we find that compared to slin(t), this strategy only marginally improves the fidelity.

In order to reduce the time required to obtain a given fidelity for g  =  g2, a more sophisticated ramping strategy is needed. In the following, we provide a simple method to estimate the efficiency of different ramps without running a complete dynamical simulation of the system. The transition probability between the ground and some excited state | {k}\rangle at time t for a ramp σ is approximately given by [29]

Equation (21)

where | {k}\rangle =| {k(\sigma,t)}\rangle is the kth instantaneous eigenstateNote5  of \hat{{H}}_{\rm eff}(s_\sigma(t)) and \omega_{0k} is the transition frequency between the ground state and |krangle. Therefore, the assessment of the efficiency of a ramp can be done by evaluating the functional

Equation (22)

where the index k runs over all the values associated with levels connected to the ground state. The functional A(\sigma,\tau_{\rm Q}) corresponds to the average transition probability per unit time for a given strategy σ and quench time τQ. We calculate the value of the functional A(σ,τQ) numerically via exact diagonalization of \hat{{H}}_{\rm eff} for a small system. This method allows to optimize ramps by minimizing the value of A(σ, τQ). The optimized ramp for a small system is then used in the simulation of larger systems. For instance, the strategy sman(t) shown in figure 7(a) was designed and optimized manually using this method. For g  =  g2, we find that A({\rm lin},\tau_{\rm Q})/A({\rm man},\tau_{\rm Q})\approx 2.3 for \tau_{\rm Q}=20 /E_{\rm R} (see figure 7(b)), thus showing the better efficiency of the strategy sman(t) compared to slin(t). As shown in figure 6, dynamical simulations of the system with g  =  g2 confirm that this strategy reduces the time required to obtain a given fidelity.

For systems initially in the SF regime (g  =  g2), the fidelity curves exhibit small oscillations (see figure 6(e)). These can be understood from time-dependent perturbation theory as oscillations occurring when some of the frequencies involved in the Fourier decomposition of the perturbation Hamiltonian enter into resonance with system frequencies. This is expected since SF have a dense spectrum at low energies and are therefore likely to enter into resonance with one of the frequencies of the perturbation Hamiltonian [30]. Hence, the amplitude of the oscillations in the fidelity curve associated with a ramping strategy s1(t) should be smaller than those associated with a ramping strategy s2(t) if A(1,\tau_{\rm Q}) <  A(2,\tau_{\rm Q}) for all τQ. This is what is observed from our numerical simulations (see figures 6(f)–(h)).

The quasi-momentum distribution of the particles in the a-mode at the end of the different ramping processes for a system with g  =  g1 are shown in figure 8. For the linear ramp, the quasi-momentum distribution of particles shown in figure 8(a) does not correspond to that of a MI state for the quench times considered. In contrast, for the ramp sgap the quasi-momentum distribution becomes approximately flat for quench times of τQ  >  200/ER. Even for the fastest ramps, we find that the occupation of the b-mode is less than 2% of the total number of particles. Thus, the experimental measure of the register fidelity can be made by comparing the quasi-momentum distribution of particles in the final state with, e.g. that shown in figure 5(c).

Figure 8

Figure 8. The quasi-momentum distributions at the end of a ramp for a system of M  =  12 sites using the parameters shown in figure 4 for g  =  g1. (a) slin; (b) sgap.

3.1.3. Discussion of the numerical results In the BHM with U/J  ll  1, the tunnelling time 1/J determines the adiabatic timescale of the system. However, as soon as many-body interactions are sufficiently large, this timescale often becomes very non-adiabatic [31]. The main observation that can be drawn from the numerical calculations presented in the last section is that by preparing the system in a TG state (g  =  g1), and using an efficient ramping strategy, it is possible to initialize a very deep MI state on a timescale equivalent to half the tunnelling time in the large lattice limit (see figure 6(b)). The time required to initialize a MI state with unit filling as well as a TG state with half filling is approximately ten times the tunnelling time of the final system [15, 16, 31, 32]. Since the tunnelling time in the large lattice limit is two orders of magnitude larger than in the small lattice limit, the total time required to initialize a MI state using our procedure is an order of magnitude faster than the direct quantum-freezing method. In this estimation we assume that the initial BEC has zero temperature, i.e. we do not take the effect of defects present in the initial state into account.

The experimental realization of initial states with g  =  g1 and g  =  g2 can be achieved using Feshbach resonances. For a magnetic Feshbach resonance fluctuations in the magnetic field result in fluctuations of the size of the gap between the ground and the first excited state which will affect the performance of our scheme. For instance, in the case g  =  g1 magnetic field fluctuations of 10 mG will change the gap by approximately 0.5% for 85Rb or 133Cs atoms [33, 34]. We assume adiabatic evolution of the system and thus these fluctuations will have negligible repercussions on the fidelity of the final state. To realize the SF regime with g  =  g2 very stable magnetic fields are required. Magnetic field fluctuations of e.g. 1 mG will lead to gap fluctuations of approximately 1% in 23Na and 85Rb. We finally remark that our scheme could also be used without employing a Feshbach resonance. This case corresponds to an intermediate value of g between g1 and g2. While a detailed analysis of the intermediate regime is beyond the scope of the present work, we do not expect qualitative differences compared to the interaction strengths considered here.

3.2. Analytical results

In this section we derive an approximate expression for the quench time required to obtain a given fidelity in the case of a linear ramp.

The energy spectrum of systems in the TG and SF regime is gaplessNote6  (see e.g. [35]), while in the MI regime, the gap between the ground and first excited state is proportional to the on-site interaction energy [30]. Since the relaxation time τ(t) of the system—the time required by the system to adjust to a change of parameters at time t—is inversely proportional to the gap between the ground and the first excited state, the relaxation time in the small lattice limit is large, while it is small and finite in the large lattice limit. This observation suggests that the adiabatic-impulse (AI) assumption from KZ theory can be used to evaluate the adiabaticity of a ramp with respect to the quench time [36]–[38].

The AI approximation is based on the following considerations: (i) when the gap between the ground and the first excited state is large, the relaxation time of the system is short and thus a system starting its evolution in the ground state remains in the ground state, i.e. its evolution is adiabatic. (ii) when the gap between the ground and the first excited state is small, the system's relaxation time is large and the system no longer adapts to changes of the Hamiltonian's parameters and its state becomes effectively frozen. The system is then in the impulse regime. The instant \hat{t} at which the system passes from the impulse to the adiabatic regime, and inversely, is defined by the equation [36, 37, 39]

Equation (23)

where \alpha=\mathcal{O}(1) is a constant. Note that \hat{t} is a time and not an operator. In the AI approximation, the time-evolution of the system is either adiabatic or impulse. Thus, the density of defects {\mathcal D}, which corresponds to the density of excitations caused by a change of parameters which drive the system from the impulse to the adiabatic regime can be approximated by [37]

Equation (24)

where | {\Psi_g(0)}\rangle and | {\Psi_e(\hat t)}\rangle are the ground and first excited states at the initial time t  =  0 and at time t=\hat t, respectively. Hence, without solving the time-dependent Schrödinger equation it is possible to make predictions for the density of defects (24) resulting from a given dynamical process.

In order to apply the KZ theory to our problem, we develop an effective model describing the system dynamics. A similar model was recently examined by Cucchietti et al  [38]. We find from numerical calculations (see figure 7(b)) that most of the excitations created in the system are caused by transitions from the ground to the first excited state. Furthermore, we examined the form of the eigenvectors of the Hamiltonian (10) in both limits for different system sizes via exact diagonalization. This revealed that both the ground state and the first accessible excited state can be approximated by an expansion in only two basis states. The elements of this reduced basis set are given by

Equation (25)

Equation (26)

where, e.g. |{2;0;1;\cdots;1}\rangle =| {2}\rangle _1\,\otimes\,|{0}\rangle_2\,\otimes\,| {1}\rangle _3\,\otimes\,\cdots\,\otimes\,| {1}\rangle _M with | {n}\rangle _i=(1/\sqrt{n!})(\hat{{a}}^\dagger_i)^n | {{\rm vac}}\rangle. The basis state | {2}\rangle corresponds to a superposition of all possible states of a system of M particles in the a-modes with (M–2) singly occupied sites, one doubly occupied site, and one empty site. In the limit M → ∞, the matrix representation \hat{{H}}_{\rm R} of Hamiltonian (10) in the reduced basis \left\{ | {1}\rangle ,| {2}\rangle \right\} reads, up to a constant energy Va

Equation (27)

The instantaneous eigenstates of equation (27) associated with the energies of the ground and first excited levels are given by | {g(t)}\rangle =-\sin (\theta (t)/2) | {1}\rangle + \cos (\theta (t)/2)| {2}\rangle and | {e(t)}\rangle =\cos (\theta (t)/2) | {1}\rangle + \sin(\theta (t)/2)| {2}\rangle, respectively, with \cos (\theta (t))=-U_{aa}(t)/[U_{aa}(t)^2+8J_{aa}(t)^2]^{{1}/{2}}, θ in [0,π]. Furthermore, we approximate the parameters Jaa and Uaa as linear functions of time

Equation (28)

where \Delta J_{aa}=|J_{aa}(s=1)-J_{aa}(s=0)|, \Delta U_{aa}=|U_{\rm init}-U_{aa}(s=0)| and U_{{\rm init}}={\rm min}[U_{aa}].

We further simplify the Hamiltonian \skew3\hat{H}_{\rm R} by replacing the hopping term Jaa(t) by its time average. Setting J_{aa}(t)=\skew6\bar{J}_{aa} with \skew6\bar{J}_{aa}=(1/\tau_Q)\int_0^{\tau_Q} {\rm d}t\,J_{aa}(t) and rescaling the time as t \rightarrow t^{\prime} -(U_{{\rm init}} \tau_Q/\Delta U_{aa}) yields the transformation

Equation (29)

which turns \skew3\hat{H}_{\rm R} into the Landau–Zener form \skew3\hat{H}_{\rm R}=\skew3\hat{A}t^{\prime} +\hat B, where \skew3\hat{A} and \skew3\hat{B} are Hermitian matrices and \skew3\hat{A} is diagonal [40]–[42]. The energy spectrum of the Hamiltonian (29) reproduces approximately the features of the spectrum of \skew3\hat{H}_{\rm eff} except that the gap is overestimated in the small lattice limit.

Starting at t  =  0 in the small lattice limit where the system is impulse, we use the AI approximation to derive the density of defects at the end of a linear ramp which drives the system to the large lattice limit, where the system is adiabatic. Defining the relaxation time as the inverse of the energy gap 1/ΔE between the levels of \skew3\hat{H}_{\rm R}, with \Delta E= [(\Delta U_{aa}t/\tau_{\rm Q})^2+ (2\sqrt{2} \skew6\bar{J}_{aa})^2]^{{1}/{2}}, equation (23) can be solved analytically and the instant \hat{t} at which our system exits the impulse regime is given by

Equation (30)

where \eta = 4\skew6\bar{J}_{aa}^2/ \Delta U_{aa}. We evaluate the density of defects {\mathcal D} using equation (24) with | {\Psi_g(t_i)}\rangle =|{g(0)}\rangle and | {\Psi_e(\hat{t})}\rangle =|{e(\hat{t})}\rangle. In order to simplify the expression of {\mathcal D}, we have set θ(0)  =  π/2, which is a greater value than the one we would obtain using the real parameters of the system. This has no significant physical consequences in our case as it is actually equivalent to considering a smaller energy gap in the small lattice limit. The density of defects is then given by

Equation (31)

where \hat{\theta} =\theta(\hat{t}). Inverting equation (31), we obtain

Equation (32)

which gives an approximate analytical expression of the quench time τQ required to obtain a density of defects {\mathcal D}.

In our system, defects correspond mainly to doubly occupied sites. Thus, the number of defects is approximately measured by the operator

Equation (33)

Hence, the density of defects is given by {{\mathcal D}}=\langle \hat{{K}} \rangle /M=\langle(\hat{{n}}_i^a)^2 \rangle- \langle \hat{{n}}_i^a\rangle. Numerical calculations show that the density of defects {{\mathcal D}} has the same scaling behaviour as Δna,i2. A comparison between equation (31) and the numerical values of the density of defects in the large lattice limit for different values of the ramping time τQ is shown in figure 9. For M  =  12 particles and g  =  g2, we find that the analytical formula for {{\mathcal D}} fits the numerical data well for α  =  1.41. For g  =  g1 the fit is less accurate. This is expected since for this value of g, the system has a less distinct separation between the adiabatic and impulse regime than for g  =  g2. For the number of particles we have been able to simulate, the fit improves as we increase the number of particles for both values of g. In addition to this, we see from equation (32) that the time required to initialize a register with a given fidelity—and thus the adiabatic time for small {\mathcal D}—scales with the ratio \Delta U_{aa}/\skew6\bar{J}_{aa}^2.

Figure 9

Figure 9. Comparison between the analytical formula for the density of defects (solid line) and the results of the TEBD calculations for the density of kinks (dots) {\mathcal D}=\langle \hat{{K}} \rangle /M, for a system of M  =  12 sites with (a) g  =  g2 and \alpha=\sqrt{2} and (b) g  =  g1 and α  =  154. Using the parameters shown in figure 4, we have Uinit/g2  =  6.8, ΔUaa/g2  =  16.4, \skew6\bar{J}_{aa}/E_{\rm R}=0.027.

4. Conclusion

We have shown that the dynamics of an optical lattice whose periodicity is doubled via superlattice potentials is very well described by a two-mode Hubbard-like Hamiltonian. The parameters of this Hamiltonian have been evaluated in the TB approximation using optimally localized GWFs. The doubling of the period removes half of the lattice sites and doubles the filling factor. We have shown that this doubling can be used for the fast initialization of a quantum register. By starting from a half filled lattice in the small lattice limit filled by either a TG (g  =  g1) gas or a SF (g  =  g2), a commensurate MI state corresponding to an atomic quantum register can be obtained on timescales shorter than those achieved by direct quantum freezing of a SF with same lattice spacing. Furthermore, we derived an analytical expression for the density of defects as a function of the quench time for linear ramping of the superlattice. We found that the time required to achieve a given density of defects is proportional to the ratio \Delta U_{aa}/\skew6\bar{J}_{aa}^2.

Our numerical calculations of ground state properties suggest that doubling the lattice period drives the system through a quantum phase transition for large lattices M → ∞. The eventual abrupt change in the ground state properties might be observable by time-of-flight measurements. An investigation of whether such a quantum phase transition indeed exists is beyond the scope of the current work but will be investigated in future work.

In this work we concentrated on the transition from filling factor n  =  1/2 to n  =  1. We finally note that the idea developed in this paper may be extended by considering lattices with an initial filling factor of n  =  1/2ℓ (where ℓ is an integer). Subsequently removing every second barrier will create a lattice with period 2a and filling factor 1/2ℓ–1. This procedure could be repeated ℓ times providing a lattice with filling factor n  =  1 and large lattice spacing 2ℓa.

Acknowledgments

This work was supported by the EU through the STREP project OLAQUI and a Marie Curie Intra-European Fellowship within the 6th European Community Framework Programme. The research was also supported by the EPSRC (UK) through the QIP IRC (GR/S82716/01) and project EP/C51933/01. DJ thanks the Beijing International Center for Mathematical Research at Peking University for hospitality while carrying out parts of this work.

Appendix A. Definition and localization properties of GWFs

In the TB limit, the effective single particle Hamiltonian of our system in momentum space (in the basis | {\alpha}\rangle _q=[(1/\sqrt{M})\sum_{i=1}^M\, {\rm e}^{{\rm i}{q a i}} \hat{{\alpha}}^\dagger_i]|{0}\rangle with α  =  a, b) can be written as

Equation (A.1)

The elements of the matrices \hat{{\varepsilon}}(\mu) are given by [43]

Equation (A.2)

For μ  =  0 and 1, the elements of the matrices \hat{{\varepsilon}}(\mu) correspond to the local site energy and hopping matrix elements between neighbouring sites, respectively. For the TB approximation to be accurate, we need the eigenvalues Eα,q of \skew3\hat{H}_{0,q} to reproduce very closely the band structure of the exact single particle Hamiltonian \skew2\hat{h}_{0,x} for all values of the lattice profile s. If we were to use WFs as mode functions wα,i, the matrices \hat{{\varepsilon}} would be diagonal [24] and the dispersion relations of the two Bloch bands sinusoidal. In order to obtain a more accurate description, we use GWFs as mode functions. The definition of GWFs is given by

Equation (A.3)

where \tilde{\Psi}_{\alpha,q}=\sum_{\beta=a,b}U^{(q)*}_{\beta,\alpha}\,\Psi_{\beta,q} is called the generalized-Bloch orbital with Ψα,q the Bloch function associated with the band α and Ri is the centre of site i [25, 44]. The rows of the 2 × 2 matrix U(q) contain the (real) normalized eigenvectors of \skew3\hat{H}_{0,q} associated with the eigenvalues Eα,q, that is \sum_{\mu=a,b}(\skew3\hat{H}_{0,q})_{n,\mu} U^{(q)}_{\alpha,\mu}=E_{\alpha,q}U^{(k)}_{\alpha,n}  [43]. Inserting GWFs in equation (A.2), we recover the elements epsilonαβ(μ) and, hence, GWFs correspond to the mode functions associated with the effective Hamiltonian \skew3\hat{H}_{0,q} [43]. Notice that the definition of GWFs reduces to that of WFs for U^{(q)}={\leavevmode\hbox{\small1\kern-3.5pt\normalsize1}}.

A.1. Localization properties of GWFs

Given a valid set of GWFs, another equally valid set of GWFs can be obtained by applying the following transformation on the U(q) matrices

Equation (A.4)

where phiα(q) are (real) functions of q which can be chosen freely as long as they do not introduce discontinuities in the generalized Bloch function [45]. The gauge transformation (A.4) is equivalent to re-phasing each Bloch function as \Psi_{\alpha,q}\,\rightarrow\,{\rm e}^{{\rm i}{\phi_\alpha(q)}}\Psi_{\alpha,q}. Notice that gauge transformations do not affect the value of the parameters Vα and Jα,β calculated using the relations (7) and (6), respectively, but they alter the localization properties—the spread—of the GWFs. Following the convention suggested by Blount [45], we set the phase functions phiα(q) in a manner that leads to maximally localized GWFs. That is, we choose the phase functions such that the resulting GWFs minimize the spread functional

Equation (A.5)

where in our case α  =  a,b while \left < x^2 \right > _\alpha=\langle {w_{\alpha,i}}|x^2 |{w_{\alpha,i}}\rangle and \left < x \right > _\alpha=\langle {w_{\alpha,i}}|x | {w_{\alpha,i}}\rangle correspond to the centre of a GWF and its second moment, respectively.

Expanding \tilde{u}_{\alpha,q}(x,s)={\rm e}^{-{\rm i}{q x}} \tilde{\Psi}_{\alpha,q}(x,s) into plane waves yields

Equation (A.6)

where Kj  =  2πj/L with L  =  M. Invoking the translational symmetry of the lattice and the convolution theorem [45], the value of the functional Ω is minimized when the expansion coefficients Gα,j(q, s) are chosen real, which is always possible when the lattice possesses mirror symmetry [43]. This is equivalent to the choice of purely real GWFs for even generalized-Bloch functions and purely imaginary ones for odd generalized-Bloch functions. Notice that this conclusion is in agreement with a conjecture of Marzari et al  [25] on the real nature (up to a global phase) of maximally localized WFs.

We have numerically evaluated the phases phiα(q) using the algorithm described in [25] for the special case of 1D WFs. This procedure minimized the functional Ω in the limit of very fine sampling of the q-space. The effect of this localization procedure is illustrated in figure 3(a).

Appendix B. Parameters of the single-particle Hamiltonian

We diagonalize \skew3\hat{H}_0 in momentum space for periodic boundary conditions. Using the Fourier transformations \skew3\hat{a}_i=(1/\sqrt{M})\sum_{q} \, {\rm e}^{{\rm i}{q a i}} \skew3\hat{a}_q and \skew3\hat{b}_i=(1/\sqrt{M}) \sum_{q} \, {\rm e}^{{\rm i}{q a i}} \skew3\hat{b}_q, \skew3\hat{H}_0 becomes block diagonal

Equation (B.1)

where q  =  2πν/Ma, ν=1  ... M. In the basis \{| {a}\rangle_q,| {b}\rangle _q\}, the operator \skew3\hat{H}_{0,q} reads

Equation (B.2)

with | {a}\rangle_q=\skew3\hat{a}^\dagger_q| {0}\rangle and | {b}\rangle_q=\skew3\hat{b}^\dagger_q| {0}\rangle. For simplicity, the explicit dependence of the parameters on the lattice profile has been dropped.

Due to the periodicity of the lattice, the eigenvalues of \skew3\hat{H}_{0,q} exhibit a band structure. By choosing the points q=0,\,\pi/a,\,\pi/2a in the Brillouin zone, we derive and solve a set of equations for the parameters Vα and Jα,β as functions of the eigenvalues Eα,q of \skew3\hat{H}_{0,q}

Equation (B.3)

For the eigenvalues of \skew3\hat{H}_{0,q} to reproduce the band structure of the exact single-particle Hamiltonian along the x-direction \skew2\hat{h}_{0,x}, we evaluate the parameters Vα and Jα,β using the eigenvalues epsilonα, q of \skew2\hat{h}_{0,x} obtained via exact numerical diagonalization for the same points in the Brillouin zone (see figure 2). The numerical values of the parameters obtained via this procedure are shown in figure 4(a).

The accuracy of the TB approximation can be tested by evaluating the standard deviation between the exact and approximated band structure. That is, taking Nq different points qi on each band, we define the standard deviation between the exact and approximate band structure for a given lattice profile by \sigma^2_s=(1/2N_q)\,\sum_{i=1}^{N_q} ( \Delta\varepsilon_{a,q_i}^2+\Delta\varepsilon_{b,q_i}^2), with \Delta\varepsilon_{\alpha,q_i}=\varepsilon_{\alpha,q_i}-E_{\alpha,q_i}. Averaging over Ns different lattice profiles si, we obtain \left < \sigma\right > _s=[(1/N_s)\sum_{i=1}^{N_s}\sigma^2_{s_i}]^{{1}/{2}}=3.4\times10^{-2} E_{\rm R} for a lattice depth of V0  =  10 ER. This excellent agreement improves further as we increase the value of V0, and hence fully justifies the TB approximation.

Appendix C. Dynamical and ground state calculations using the TEBD algorithm

The TEBD algorithm is based on directly manipulating a matrix product representation of the many-body wavefunction. Here, we shall briefly describe the key aspects of this algorithm and refer the reader to some of the recent literature [46]–[49] for more detail.

An arbitrary state of a 1D quantum lattice system composed of M sites can be written as

Equation (C.1)

where c_{{j_1 \cdots j_M}} is a set of {\ddo}^M complex amplitudes and |{j_m}\rangle is a basis spanning the local {\ddo}-dimensional Hilbert space of site m. Within time-dependent density matrix renormalization group (DMRG) the amplitudes c_{j_1 \cdots j_M} are constructed from a product of tensors

Equation (C.2)

where \{\alpha\}=\{\alpha_1,\ldots,\alpha_{M-1}\}, \{\chi\} =\{\chi_1,\ldots,\chi_{M-1}\} and with Γ and λ tensors chosen to be constructed from the set of M-1 Schmidt decompositions for contiguous partitions of the system. Specifically, the elements of \lambda^{\lbrack m\rbrack}_{\alpha} are taken to be Schmidt coefficients of the bipartite splitting after site m in | {\psi}\rangle=\sum_{\alpha=1}^{\chi_m}\lambda^{\lbrack m\rbrack}_{\alpha}|{L^{m}_{\alpha}}\rangle |{R^{m}_{\alpha}}\rangle with Schmidt rank χm. The Schmidt states |{L^{m}_\alpha}\rangle and |{R^{m}_\alpha}\rangle spanning the left \{1,\ldots,m\} and right \{m+1,\ldots,M\} subsystems of sites respectively are then specified by the corresponding sums remaining in equation (C.2).

The usefulness of this representation is based on the observation that for 1D systems with a Hamiltonian composed of nearest neighbour terms the ground state and low-lying excited states have Schmidt coefficients \lambda^{\lbrack m\rbrack}_{\alpha} which rapidly decay with α when arranged in descending order. Consequently, rather than allowing the Schmidt ranks χm to grow to their maximum permissible value a much smaller fixed upper-limit χ can be imposed truncating the representation while still providing a near unit overlap with the exact state | {\psi}\rangle. Fixing the Schmidt ranks results in the number of parameters scaling as O({\ddo}\chi^2M) and so curtails the possible exponential growth with M seen for general coefficients c_{j_1 \cdots j_M}.

The matrix product representation also permits the efficient update of the state after the action of a unitary operator on any two neighbouring lattice sites. This proceeds by modifying the Γ tensors associated to the sites and the λ tensor linking them and requiring a number of operations which scales as O({\ddo}^4\chi^3). The resulting tensors are then systematically truncated back to a maximum rank of χ.

Dynamical simulations can be performed by decomposing the time evolution operator \exp(-{\rm i}H\delta t), for small time step \delta t, into a sequence of pairwise unitaries via a Suzuki-Trotter expansion. Given the properties outlined such a calculation is likely to be accurate for a practical value of χ if both the initial state and the states generated by the dynamics remain in the low-energy manifold of the system. To determine the appropriate χ calculations are repeated with increasing values of χ until the final result converges and are unaffected by further increases. For practical purposes the convergence is usually quantified by the robustness of the expectation values calculated. The accuracy of a calculation is also gauged by the sum of the discarded Schmidt coefficients at each time step—a quantity which should necessarily be small—and the deviation of normalization of the final state from unity which indicates the accumulated effect of truncation.

Finally, initial states are typically taken to be the ground state of the system which are found either by applying the DMRG procedure or, as in this work, by simulating imaginary time evolution through the repeated application of exp(–Hδt) and subsequent renormalization of the state. In our simulations, we have used χ  =  40.

Appendix D. The ramp sgap

The ramp sgap can be evaluated as follows. The energy gap between the ground and the first excited state associated with a given lattice profile s is given by {\rm gap}(s)=\omega_{0e} where \omega_{0k}=(\varepsilon_{{\rm eff},k}-\varepsilon_{{\rm eff},0}) with \varepsilon_{{\rm eff},\alpha} the αth instantaneous eigenvalue of \skew3\hat{H}_{\rm eff}(s). The gap is evaluated numerically via the direct diagonalization of \skew3\hat{H}_{\rm eff} for a small number of sites. The function sgap(t) is defined by the equation {\rm d}s_{\rm gap}(t)/{\rm d}t={\rm gap}(s_{\rm gap}(t))/K where K^{-1}=\int_{0}^{\tau_Q}{\rm d}t\, {\rm gap}(s_{\rm gap}(t)) is a normalization constant.

References

[1]
Jaksch D and Zoller P 2005 The cold atom Hubbard toolbox Ann. Phys. 315 52 
CrossRef
[2]
Lewenstein M, Sanpera A, Ahufinger V, Damski B, Sen De A and Sen U 2006 Ultracold atomic gases in optical lattices: mimicking condensed matter physics and beyondPreprint cond-mat/0606771
Preprint
[3]
Cirac J I and Zoller P 2004 New frontiers in quantum information with atoms and ions Phys. Today 57 38 
CrossRef
[4]
Jaksch D, Briegel H-J, Cirac J I, Gardiner C W and Zoller P 1999 Entanglement of atoms via cold controlled collisions Phys. Rev. Lett. 82 1975 
CrossRef
[5]
Mandel O, Greiner M, Widera A, Rom T, Hänsch T W and Bloch I 2003 Coherent transport of neutral atoms in spin-dependent optical lattice potentials Phys. Rev. Lett. 91 010407 
CrossRefPubMed
[6]
Raussendorf R and Briegel H-J 2001 A one-way quantum computer Phys. Rev. Lett. 86 5188 
CrossRefPubMed
[7]
Vollbrecht K G H, Solano E and Cirac J I 2004 Ensemble quantum computation with atoms in periodic potentials Phys. Rev. Lett. 93 220502–4 
CrossRefPubMed
[8]
Zhou Z-W, Han Y-J and Guo G-C 2006 Local operations in qubit arrays via global periodic manipulation Phys. Rev. A 74 052334–5 
CrossRef
[9]
Kay A and Pachos J K 2004 Quantum computation in optical lattices via global laser addressing New J. Phys. 6 126 
IOPscience
[10]
Calarco T, Dorner U, Julienne P S, Williams C J and Zoller P 2004 Quantum computations with atoms in optical lattices: marker qubits and molecular interactions Phys. Rev. A 70 012306 
CrossRef
[11]
Häffner H et al  2005 Nature 438 643 
CrossRefPubMed
[12]
Greiner M, Mandel O, Esslinger T, Hänsch T W and Bloch T W 2002 Quantum phase transition from a superfluid to a Mott insulator in a gas of ultracold atoms Nature 415 39 
CrossRefPubMed
[13]
Julienne P S, Williams C J, Band Y B and Trippenbach M 2005 Loading Bose–Einstein-condensed atoms into the ground state of an optical lattice Phys. Rev. A 72 053615–7 
CrossRef
[14]
Sklarz S E, Friedler I, Tannor D J, Band Y B and Williams C J 2002 Flat-phase loading of a Bose–Einstein condensate into an optical lattice Phys. Rev. A 66 053620 
CrossRef
[15]
Clark S R and Jaksch D 2004 Dynamics of the superfluid to Mott-insulator transition in one dimension Phys. Rev. A 70 43612 
CrossRef
[16]
Blakie P B and Porto J V 2004 Adiabatic loading into optical lattices Phys. Rev. A 69 013603 
CrossRef
[17]
Rey M 2004 A Ultra-cold bosonic atoms in optical lattices PhD Thesis  (University of Maryland at College Park) 
[18]
Peil S, Porto J V, Laburthe T B, Obrecht J M, King B E, Subbotin M, Rolston S L and Phillips W D 2003 Patterned loading of a Bose–Einstein condensate into an optical lattice Phys. Rev. A 67 051603 
CrossRef
[19]
Sebby-Strabley J, Anderlini M, Jessen P S and Porto J V 2006 Lattice of double wells for manipulating pairs of cold atoms Phys. Rev. A 73 033605 
CrossRef
[20]
Gygi O, Katzgraber H G, Troyer M, Wessel S and George B G 2006 Simulations of ultracold bosonic atoms in optical lattices with anharmonic traps Phys. Rev. A 73 063606 
CrossRef
[21]
Meyrath T P, Schreck F, Hanssen J L, Chuu C-S and Raizen M G 2005 Bose–einstein condensate in a box Phys. Rev. A 71 041604 
CrossRef
[22]
Jaksch D, Bruder C, Cirac J I, Gardiner C W and Zoller P 1998 Cold bosonic atoms in optical lattices Phys. Rev. Lett. 81 3108 
CrossRef
[23]
Ashcroft N W and Mermin N D 1976 Solid State Physics  (Philadelphia: Brooks/Cole Thomson Learning) 
[24]
Schnell I, Czycholl G and Albers R C 2002 Hubbard-u calculations for cu from first-principle Wannier functions Phys. Rev. B 65 075103 
CrossRef
[25]
Marzari N and Vanderbilt D 1997 Maximally localised generalised Wannier functions for composite energy bands Phys. Rev. B 56 12847 
CrossRef
[26]
Mühner T D, White S R and Monien H 2000 One-dimensional Bose–Hubbard model with nearest-neighbor interaction Phys. Rev. B 61 12474 
CrossRef
[27]
Roth R and Burnett K 2003 Phase diagram of bosonic atoms in two-color superlattices Phys. Rev. A 68 023604 
CrossRef
[28]
Paredes B, Widera A, Murg V, Mandel O, Folling S, Cirac I, Shlyapnikov G V, Hansch T W and Bloch I 2004 Tonks–Girardeau gas of ultracold atoms in an optical lattice Nature 429 277–81 
CrossRefPubMed
[29]
Schiff L I 1968 Quantum Mechanics 3rd edn  (New York: McGraw-Hill) 
[30]
Clark S R and Jaksch D 2006 Signature of the superfluid to mi-insulator transition in the excitation spectrum of cold atoms New J. Phys. 8 160 
IOPscience
[31]
Hecker D J, Simsarian J E, Haeffner H, McKenzie C, Browaeys A, Cho D, Helmerson K, Rolston S L and Phillips W D 2002 A Bose–Einstein condensate in an optical lattice J. Phys. B: At. Mol. Opt. Phys. 35 3095–110 
IOPscience
[32]
Griessner A, Daley A J, Clark S R, Jaksch D and Zoller P 2007 Dissipative dynamics of atomic hubbard models coupled to a phonon bath: dark state cooling of atoms within a bloch band of optical lattices New J. Phys. 9 44 
IOPscience
[33]
Kohler T, Goral K and Julienne P S 2006 Production of cold molecules via magnetically tunable feshbach resonances Rev. Mod. Phys. 78 1311 
CrossRef
[34]
Weber T, Herbig J, Mark M, Nagerl H-C and Grimm R 2003 Bose–Einstein condensation of cesium Science 299 232–5 
CrossRefPubMed
[35]
Pethick C J and Smith H 2002 Bose–Einstein Condensation in Dilute Gases  (Cambridge: Cambridge University Press) 
[36]
Damski B 2005 The simplest quantum model supporting the Kibble–Zurek mechanism of topological defect production: Landau–Zener transitions from a new perspective Phys. Rev. Lett. 95 035701 
CrossRefPubMed
[37]
Damski B and Zurek W H 2006 Adiabatic-impulse approximation for avoided level crossings: from phase transition dynamics to Landau–Zener evolutions and back again Phys. Rev. A 73 063405 
CrossRef
[38]
Cucchietti F M, Damski B, Dziarmaga J and Zurek W H 2007 Dynamics of the Bose–Hubbard model: transition from Mott insulator to superfluid Phys. Rev. A 75 023603 
CrossRef
[39]
Zurek W H 1985 Cosmological experiments in superfluid helium Nature 317 505 
CrossRef
[40]
Zener C 1932 Non-adiabatic crossing of energy levels Proc. R. Soc. Lond. Ser. A 137 696 
CrossRef
[41]
Sinitsyn N 2004 Generalisation of the Landau–Zener theory in the physics of nanoscale systems PhD Thesis  (Texas AM University) 
[42]
Rubbmark J R, Kash M M, Littman M G and Kleppner D 1981 Dynamical effects at avoided level crossings: a study of the Landau–Zener effect using Rydberg atoms Phys. Rev. A 23 3107 
CrossRef
[43]
Sporkmann B and Bross H 1994 Calculation of Wannier functions for fcc transition metals by Fourier transformation of Bloch functions Phys. Rev. B 49 10869 
CrossRef
[44]
Bross H 1971 A general interpolation scheme for the one-electron energies of a solid Z. Phys. 243 311–25 
CrossRef
[45]
Blount E I 1962 Formalisms of band theory Solid State Physics: Advances in Research and Applications  (New York: Academic) p 305 
[46]
Vidal G 2003 Efficient classical simulation of slightly entangled quantum computations Phys. Rev. Lett. 91 147902 
CrossRefPubMed
[47]
Vidal G 2004 Efficient simulation of one-dimensional quantum many-body systems Phys. Rev. Lett. 93 40502 
CrossRef
[48]
Vidal G 2007 Classical simulation of infinite-size quantum lattice systems in one spatial dimension Phys. Rev. Lett. 98 070201 
CrossRefPubMed
[49]
Daley A J, Kollath C, Schollwöck U and Vidal G 2004 Time-dependent density-matrix renormalization-group using adaptive effective hilbert spaces J Stat. Mech.: Theory Exp. 04 P04005 
IOPscience

Notes

Note2  For example, in ion trap quantum computing, distances in the m regime are sufficient to address individual atoms see [11].

Note3  To avoid any discontinuity, we work with a lattice periodicity of a even in the case s  =  1.

Note4  The two lowest bands form two segments of the lowest Bloch band if a cell size of a/2 is used in the case s  =  1.

Note5  Note that neglecting non-adiabatic changes of the GWFs does not correspond to P0kσ(t)  =  0 at all times.

Note6  In finite size systems, the spectrum is not gapless, only very dense.



Please login to access our web services, or create an account if you don't yet have one.

You must have cookies enabled in your web browser to be able to login.

Username
Password

Forgotten your password? Get a new one here.