text:
# Effective mass analysis of Bose-Einstein condensates in optical lattices — stabilization and levitation

## I Introduction

The dynamics of classical and quantum systems in periodic potentials is a central paradigm of physics, finding applications from condensed matter physics to optics, and from nonlinear dynamics to atomic and to plasma physics. In recent years, the experimental and theoretical study of quantum-degenerate atomic systems in periodic potentials has opened up new avenues of investigation. In particular, Bose-Einstein condensates are a macroscopic quantum system that is amenable to exquisite experimental control. As a result, many phenomena studied in solid state systems can be re-examined in a more direct and dramatic fashion. Even more importantly perhaps, it is now possible to realize experimentally model systems that had previously been the object of considerable theoretical studies, but were all but impossible to test experimentally. One particularly beautiful example is the experimental realization of the Hubbard model leading to the demonstration of the superfluid to Mott insulator transition in a Bose condensate of 87Rb atoms mott .

The first experiments involving the dynamics of Bose-Einstein condensates (BEC) in periodic potentials were carried out by Anderson and Kasevich, who used this approach to demonstrate a mode-locked atom laser modelock  and observe atomic Josephson oscillations modelock ; josephson . The list of phenomena that were subsequently experimentally demonstrated and/or theoretically investigated includes the generation of atomic number squeezing squeez , the observation of the superfluid-Mott insulator phase transition, the generation of discrete dis  or gap solitons gap , the prediction of modulational instabilities stability  and superfluid flow flow , the observation of Bloch oscillations blochosc , the analysis and observation of coherent acceleration acceleration , studies of magnetism mag , etc.

It is well known that a particle confined to an infinite periodic potential and acted upon by an external force behaves as if possessing an effective mass that can be substantially different from its true mass, and may even take negative values. In particular, it is this property that is at the core of proposals to generate bright matter-wave solitons in BECs with repulsive interactions gap . The purpose of the present paper is to extend these studies by analyzing the temporal evolution of a condensate in the periodic potential provided by an optical lattice. The initial state of the condensate is chosen to be a (approximate) Bloch state modulated by a slow-varying Gaussian envelope, and we compare two situations where the Bloch state is either associated with a positive or a negative effective mass. We show that in situations where the dynamics of the condensate is well approximated by a negative effective mass, the periodic potential can lead to the stabilization of an otherwise unstable condensate. We also demonstrate theoretically the levitation of condensates of negative effective masses. These studies further allow us to determine the impact of self-interactions and of finite condensate widths on the general usefulness of the effective mass concept.

The paper is organized as follows. Section II gives a brief review of the linear problem of a particle inside a periodic sinusoidal potential and introduces the concept of effective mass. Section III applies and extends these ideas to the dynamics of a BEC in an optical lattice. We discuss the equation of motion for the slowly varying condensate envelope function, from which we can gain useful physical insights into the dynamics of the system. A full numerical solution is presented in Sec. IV, which demonstrates in particular the stabilization and the levitation of a condensate of negative effective mass, in agreement with the prediction of Sec. III. Finally, Section V presents concluding remarks on the usefulness of the effective mass concept and an outlook.

## II Linear properties of infinite periodic potential

In this section, we briefly review important aspects of the linear problem of a particle of mass $`m`$ inside a one-dimensional infinite periodic potential of the form

$$V(x)=V_0\mathrm{cos}^2(k_0x),$$
(1)
with corresponding time-independent Schrödinger equation

$$\left[\frac{\mathrm{}^2}{2m}\frac{d^2}{dx^2}+V_0\mathrm{cos}^2(k_0x)\right]\varphi (x)=E\varphi (x).$$
(2)
We proceed by introducing the dimensionless quantities

$$\theta =k_0x;b=\frac{2mE}{\mathrm{}^2k_0^2};h^2=\frac{2mV_0}{\mathrm{}^2k_0^2},$$
(3)
in terms of which Eq. (2) becomes

$$\frac{d^2\varphi }{d\theta ^2}+\left(bh^2\mathrm{cos}^2\theta \right)\varphi =0.$$
(4)
Eq. (4) has the form of Mathieu’s equation, whose solutions are well known mathieu . As a warm-up we sketch out the main features of its solutions which, according to Bloch’s theorem, must be of the form

$$\varphi _{n,s}(\theta )=e^{is\theta }F_n(\theta ),$$

where $`s`$ is real and arbitrary (in standard textbook language, $`s=k/k_0`$ where $`\mathrm{}k`$ is the quasi-momentum). In our scaled units, the periodic functions $`F_n(\theta )`$ have a period of $`\pi `$. The energy spectrum associated with the periodic potential $`V(\theta )`$ exhibits a band structure familiar from solid state physics. Each value of $`s`$ gives a discrete spectrum whose structure is periodic with respect to $`s`$. This property allows one to restrict the discussion to the first Brillouin zone, $`1<s1`$. As $`s`$ is increased within that zone, the energy levels trace out curves that are restricted to a small energy band. We restrict our discussion to the first two bands, as depicted in Fig. 1(a). Note that although the probability distribution $`|\varphi |^2`$ is a periodic function of $`\theta `$, $`\varphi `$ itself is not periodic unless $`s`$ takes integer values — e.g., at the center or the edges of the Brillouin zone. For a detailed discussion of the solutions to Eq. (4), see e.g. Ref. mathieu .

For $`s=0`$ or 1, the eigenstates of the Mathieu’s equation are periodic and take the form of standing waves. In particular, at the edge $`s=1`$ of the first Brillouin zone, the wave functions $`\varphi _{n,s}`$ for the first two energy bands, $`n=1,2`$, can be expressed as the Fourier series

$$\varphi _{1,s=1}(\theta )=\underset{n=0}{\overset{\mathrm{}}{}}A_{2n+1}\mathrm{sin}(2n+1)\theta ,$$
(5)
and

$$\varphi _{2,s=1}(\theta )=\underset{n=0}{\overset{\mathrm{}}{}}B_{2n+1}\mathrm{cos}(2n+1)\theta ,$$
(6)
where

$$G_{2n+1}=k_{2n1}\frac{1}{G_{2n1}},n>1,$$

$$k_n=(4b2h^24n^2)/h^2,$$

$$G_{2n+1}=A_{2n+1}/A_{2n1},\mathrm{or}B_{2n+1}/B_{2n1},$$

and $`A_3=k_1+1`$, $`B_3=k_11`$.

The effective mass, $`m^{}`$, which characterizes the response of the particle to external perturbations, is defined as

$$m^{}=\frac{\mathrm{}^2}{^2E/k^2}=\frac{2m}{^2b/s^2}.$$
(7)
Fig. 1(b) shows the variation of the effective mass of the first two bands with the quasi-momentum. One observes that it is near infinite (i.e., $`1/m^{}=0`$) for certain values of the quasi-momentum $`s`$, and even becomes negative in certain regions of the zone. In more than one dimension, $`m^{}`$ acquires a tensorial character, with elements given by

$$m_{\mu \nu }^{}=\frac{\mathrm{}^2}{^2E/(k_\mu k_\nu )}.$$

In the next section, we apply the effective mass formalism to the case of a self-interacting Schrödinger field described by the Gross-Pitaevskii equation. Specifically, we show that a BEC prepared in a Bloch state with negative effective mass behaves as if the signs of its self-interaction and of the external confining potential had been reversed.

## III Effective mass equation

We consider a BEC in a sinusoidal optical lattice potential. Near zero temperature, the system is described to an excellent degree of approximation by a time-dependent Gross-Pitaevskii equation for the normalized condensate wave function $`\psi (𝐫,t)`$, generalized to include the effects of three-body collisions,

$$i\mathrm{}\frac{\psi (𝐫,t)}{t}=\left[\frac{\mathrm{}^2}{2m}^2+V(𝐫)+U(𝐫)+Ng|\psi |^2iN^2K_3|\psi |^4\right]\psi (𝐫,t).$$
(8)
Here, $`V(𝐫)`$ is the periodic lattice potential, $`U(𝐫)`$ is an additional external potential that is taken to be slowly varying on the scale of the lattice period (e.g., a confining potential), $`N`$ is the total number of atoms, $`g`$ is the strength of the nonlinear inter-atomic interaction, and $`K_3`$ gives the rate of three-body recombination loss loss .

We proceed by expanding the condensate wave function $`\psi (𝐫,t)`$ on the complete set of Bloch functions $`\varphi _{n𝐤}(𝐫)`$ that are the stationary solutions of the associated linear Schrödinger equation

$$\left[\frac{\mathrm{}^2}{2m}^2+V(𝐫)\right]\varphi _{n𝐤}(𝐫)=E_{n𝐤}\varphi _{n𝐤}(𝐫).$$
(9)
The subscripts $`n`$ and $`𝐤`$ represent the band index and quasi-momentum, respectively. The Bloch functions satisfy the orthonormality condition:

$$_{\mathrm{cell}}𝑑𝐫\varphi _{m𝐤}^{}(𝐫)\varphi _{n𝐤^{}}(𝐫)=\frac{\mathrm{\Omega }}{(2\pi )^3}\delta _{mn}\delta _{\mathrm{𝐤𝐤}^{}},$$

where the integral is taken over a single period of the lattice of volume $`\mathrm{\Omega }`$. Our use of the resulting expansion,

$$\psi (𝐫,t)=\underset{n,𝐤}{}A_{n𝐤}(𝐫,t)\varphi _{n𝐤}(𝐫)e^{iE_{n𝐤}t/\mathrm{}},$$
(10)
is motivated by the the fact that the Bloch functions can capture the rapid oscillations of the condensate wave function, while the slowly varying envelope functions $`A_{n𝐤}`$ will describe the slow center-of-mass motion of the condensate.

The equations of motion for the amplitudes $`A_{n𝐤}`$ are found by inserting the expansion (10) into the generalized Gross-Pitaevskii equation (8). The resulting equations are however of little use if a very large number of $`A_{n𝐤}`$ is required to accurately describe the dynamics of the condensate. Here, we consider a more specialized situation where the matter-wave field is characterized by a central wave vector $`𝐤_0`$ corresponding to the mean velocity of the condensate (this would correspond to the carrier wave in conventional optics). Hence, we expand $`\psi (𝐫,t)`$ in a way reminiscent of the slowly varying envelope approximation of quantum optics as

$$\psi (𝐫,t)=\underset{n}{}f_n(𝐫,t)\varphi _{n𝐤_0}(𝐫)e^{iE_{n𝐤_0}t/\mathrm{}}.$$
(11)
Inserting Eq. (11) into Eq. (8) and applying the effective mass book  or multiple scales sipe  methods yields a set of equations governing the time evolution of the condensate envelope functions $`f_n(𝐫,t)`$ steel ,

$$i\mathrm{}\left(\frac{f_n}{t}+𝐯_gf_n\right)=\left[\frac{\mathrm{}^2}{2m_{\mu \nu }^{}}\frac{^2}{x_\mu x_\nu }+U(𝐫)+Ng^{}|f_n|^2iN^2K_3^{}|f_n|^4\right]f_n.$$
(12)
The derivation of Eq. (12), the so-called effective mass equation, is given in Appendix A. Here $`x_\mu `$ and $`x_\nu `$ are cartesian components of $`𝐫`$ and we use an implicit summation over repeated indices. The velocity

$$𝐯_g=\frac{1}{m}\varphi _{n𝐤_0}|\widehat{𝐩}|\varphi _{n𝐤_0},$$
(13)
is the drift velocity of the $`n`$-th band contribution to the condensate envelope. As shown in Appendix A, it is equal to the gradient of the energy $`E_{n𝐤}`$ with respect to the quasi-momentum $`𝐤`$, evaluated at $`𝐤_0`$. We note that $`𝐯_g`$ vanishes at the extreme points of the band, in particular at the band edges. The coefficients $`g^{}`$ and $`K^{}`$ give the renormalized nonlinear interaction strength and three-body loss rate, respectively. They are given by

$$g^{}=\frac{(2\pi )^3}{\mathrm{\Omega }}g_{\mathrm{cell}}𝑑𝐫|\varphi _{n𝐤_0}(𝐫)|^4,$$
(14)
and

$$K_3^{}=\frac{(2\pi )^3}{\mathrm{\Omega }}K_3_{\mathrm{cell}}𝑑𝐫|\varphi _{n𝐤_0}(𝐫)|^6$$
(15)
The lattice potential $`V(𝐫)`$ does not appear in these expressions, its effects being already incorporated into the effective mass $`m^{}`$ which appears in the first term on the right-hand side of Eq. (12) as well as in the drift term on the left-hand side of that equation.

The effects of a negative mass, $`m^{}<0`$, can be immediately inferred from the complex conjugate of Eq. (12),

$$i\mathrm{}\left(\frac{f_n^{}}{t}+𝐯_gf_n^{}\right)=\left[\frac{\mathrm{}^2}{2|m_{\mu \nu }^{}|}\frac{^2}{x_\mu x_\nu }U(𝐫)Ng^{}|f_n^{}|^2iN^2K_3^{}|f_n^{}|^4\right]f_n^{}.$$
(16)
This equation shows that as compared to the case of a positive effective mass, the condensate density envelope profile $`|f_n(t)|^2`$ evolves now under the influence of an inverted external potential $`U(𝐫)`$ and an inverted nonlinearity $`g^{}`$. As expected, though, the three-body recombination rate does not change sign and still represents a loss term. The change in sign of $`g^{}`$ and $`U`$ leads to a number of consequences. In particular, past work has shown how this property can be exploited to launch bright matter-wave solitons in condensates with positive scattering length gap  and induce modulational instability stability . The following section explores further consequences of this property.

## IV Effects of Negative effective mass

### IV.1 Inverted nonlinearity

We first investigate the effect of an inverted nonlinearity. As is well known, this will change the two-body interaction from being repulsive to attractive, and vice versa. Since condensates with attractive interactions are normally unstable against collapse, the sign change of $`m^{}`$ is clearly expected to have a significant effect.

To reduce the scale of the computation, we restrict ourselves to a two-dimensional system consisting of a pancake-shaped condensate whose axial motion is frozen to the trap ground state. (We note that condensates with reduced dimensions have already been realized in the laboratory reduce .) Furthermore, to focus on the effects of the nonlinearity, we assume that except for the lattice potential

$$V(x,y)=V_0[\mathrm{cos}^2(k_0x)+\mathrm{cos}^2(k_0y)],$$
(17)
there is no additional external potential in the transverse dimensions.

The evolution of the transverse component $`\psi _{}`$ of the condensate wave function is then given by the generalized Gross-Pitaevskii equation

$`i\mathrm{}{\displaystyle \frac{\psi _{}}{t}}`$
$`=`$ $`\left[{\displaystyle \frac{\mathrm{}^2}{2m}}_{}^2+V(x,y)+Ng|\psi _{}|^2iN^2K_3|\psi _{}|^4\right]\psi _{},`$
where

$$_{}^2=^2/x^2+^2/y^2.$$

For the initial state of the condensate, we assume a gaussian transverse wave function of the form

$$\psi _{}(x,y,t=0)=𝒩f(x,y)\varphi _{n𝐤}(x,y),$$
(18)
where $`f(x,y)=\mathrm{exp}[(x^2+y^2)/w^2]`$ is a slowly varying Gaussian envelope function of width $`w`$, $`\varphi _{n𝐤}`$ is a Bloch wave function of the linear Hamiltonian, and the constant $`𝒩`$ normalizes $`\psi _{}`$ to unity. Figs. 2 and 3 illustrate the dynamics of the condensate for two cases, one corresponding to a negative effective mass and the other to a positive one. We choose the Bloch wave function $`\varphi _{n𝐤}`$ associated with a negative effective mass $`m^{}`$ as

$$\varphi _{n𝐤}\mathrm{sin}(k_0x)\mathrm{sin}(k_0y).$$
(19)
This corresponds approximately to the Bloch state at the top of the first energy band, whose full expression is given by $`\varphi _{1,s=1}`$ of Eq. (5). The approximate form (19) retains only the dominant term in the Fourier series (5). Its overlap with $`\varphi _{1,s=1}`$ is larger than 99% for the parameters of our simulations. Similarly, for the state of positive effective mass $`m^{}`$ we use

$$\varphi _{n𝐤}\mathrm{cos}(k_0x)\mathrm{cos}(k_0y),$$
(20)
which corresponds roughly to the Bloch state $`\varphi _{2,s=1}`$ at the bottom of the second band, see Eq. (6). Both states (19) and (20) are located at the edge of the first Brillouin zone, with corresponding drift velocities $`𝐯_g=0`$. Hence the center of the condensate stays at $`x=y=0`$ during the course of its time evolution.

The two columns in Fig. 2 show snapshots of the density distribution of the condensate for positive (left column) and negative (right column) effective masses $`m^{}`$. In both cases the nonlinear interaction strength $`g`$ is negative, hence in the absence of the lattice the condensate would be unstable and subject to collapse. This however can be changed if the effective mass of the BEC becomes negative, as demonstrated in the right column. In that case, the condensate behaves as if experiencing a repulsive two-body interaction, its width increasing in time. This is in sharp contrast to the situation for a positive effective mass, in which case the condensate width quickly implodes note . This implosion is characterized by a complex dynamical behavior, such as the formation of atomic “bursts,” observed at later times. (We remark that while an attractive condensate can expand if it possesses a sufficiently large kinetic energy, this is not the case of Fig. 2, where the two initial states have nearly identical kinetic energies.) The collapse dynamics depicted here is rather similar to that of an attractive condensate in a harmonic potential, a system which has gained much interest recently as a result of the capability of tuning the scattering length from positive to negative via Feshbach resonance jila ; saito ; santos . We note however that the use of the generalized Gross-Pitaevskii equation becomes questionable in the regime of extreme collapse. This point though it is important, is not central to our results.

We conclude this section by observing that the number of condensate atoms decreases in time, due to the presence of three-body loss. This is particularly evident for the state with positive effective mass, as shown in Fig. 3. This is as expected, since the atomic density is higher in that case. Finally, we note that as a result of the tensorial character of the effective mass, it can readily be positive along some dimensions and negative along others. For example, in two dimensions a condensate wave function of the form

$$\varphi _{n𝐤}\mathrm{sin}(k_0x)\mathrm{cos}(k_0y),$$

has an effective mass that is negative along $`x`$ and positive along $`y`$-direction. This leads to a situation where the condensate will expand along $`x`$ and contract along $`y`$, as has been numerically confirmed.

### IV.2 Inverted external potential

We mentioned that if a condensate of negative effective mass is subject to an external potential $`U(𝐫)`$ that varies slowly over the lattice period, then this potential will appear inverted to the condensate. For example, under the effect of gravity, the center of mass of the condensate will climb up, rather than fall down. To demonstrate this effect, we solve the one-dimensional Gross-Pitaevskii equation

$$i\mathrm{}\frac{\psi }{t}=\left[\frac{\mathrm{}^2}{2m}\frac{^2}{z^2}+V(z)+U(z)+Ng|\psi |^2\right]\psi ,$$
(21)
where $`V(z)=V_0\mathrm{cos}^2(k_0z)`$ is the lattice potential and

$$U(z)=Gz,$$

represents the gravitational potential, $`G`$ being a constant. We neglect the three-body loss term, which is not essential in the present discussion.

Again, we assume that the initial condensate wave function is a broad Gaussian envelope function centered at $`z=0`$, modulated by a sinusoidal function. Taking that modulation of the form $`\mathrm{sin}(k_0z)`$ yields a state of negative effective mass $`m^{}`$, while $`\mathrm{cos}(k_0z)`$ gives a state of positive effective mass. Figure 4(a) illustrates the evolution of the center of mass of the condensate for these two cases. As expected from the previous discussion, we see that it initially climbs up for negative $`m^{}`$ and falls down for positive $`m^{}`$.

As the center of mass of the condensate of negative effective mass climbs up against the gravitational field, its gravitational potential energy increases. This increase must be offset by a decrease in kinetic energy (assuming for the sake of argument that the nonlinear interaction energy is negligible). Since the initial condensate possesses a finite kinetic energy, roughly given by $`\mathrm{}^2k_0^2/(2m)`$, conservation of energy therefore sets an upper limit on how high its center of mass can climb. This is illustrated in Fig. 4(a). A good estimate for the maximum height can be obtained from

$$\frac{\mathrm{}^2k_0^2}{4m}mGz_{\mathrm{max}},$$

where the left hand side is about half of the initial kinetic energy and the right hand side is the maximum gain in gravitational energy. For lithium atoms in an optical lattice formed by laser light of wavelength of 0.5 $`\mu `$m, one finds that $`z_{\mathrm{max}}`$ is about 0.4 mm, a rather significant distance in the context of integrated atom optics iao . Figure 4(b) compares the condensate density distribution at the maximum height to its initial envelope profile.

That a maximum height should exist demonstrates the limit of the effective mass formalism, since a literal interpretation of the effective mass equation (12) would incorrectly lead one to conclude that in the case of negative effective mass, the condensate should move against the gravitational field without bound. The incorrectness of this reasoning lies in the fact that Eq. (12) neglects all terms coupling different Bloch modes (see derivation in Appendix A). What happens in reality is that as the condensate moves under the combined influence of the external potential and of its intrinsic nonlinear potential, other Bloch modes with different quasi-momenta than those initially excited inevitably become populated. This in turn changes the effective mass of the system, or even makes this band-specific concept meaningless. In the present example, the effective mass changes its sign to positive at the maximum height as successive energy bands become excited. In the mean time, the oscillations in the density profile have a contrast less than one \[see Fig. 4(b)\] which shows that the phase of the wave function at maximum height becomes spatially continuous (otherwise, the density would vanish at each phase singularities) while the initial wave function possesses a series of $`\pi `$-phase jumps.

### IV.3 Discussion

We conclude this section with some general comments on the generation and robustness of condensate states of negative effective mass. As one can see from Eq. (19), the condensate wave function suffers a $`\pi `$-phase jump every lattice period. This immediately suggests that these states can be prepared by the phase-imprinting method that has been previously successfully implemented to generate dark solitons as in Ref. phase . However, due to edge effects and inaccuracy in controlling the intensity of the phase imprinting pulse, this method may not be fully adequate to generate a series of accurate $`\pi `$-jumps as required by Eq. (19).

A better choice may therefore consist in first preparing a condensate with uniform phase in a specific electronic state. Taking advantage of the internal degrees of freedom of the atoms, one can then drive the condensate into another internal state with an optical field of appropriate topological spatial geometry. The resulting center-of-mass wave function of the condensate having then the correct phase structure. This method has been successfully used in the past to create a series of accurate $`\pi `$-phase jumps in condensates pi1 ; pi2 ; pi3 .

We conclude this section with a brief discussion of the robustness of condensate states of negative effective mass. As we have mentioned, the state (19) is only an approximation to the true Bloch state (5). The overlap between these two states decreases as $`h^2`$ increases. Furthermore, $`|1/m^{}|`$ decreases with increasing $`h^2`$, as the system approaches the tight-binding regime. On the other hand, $`h^2`$ cannot be made too small, since the size of the gap between the first and second energy band is itself proportional to $`h^2`$, and a small gap will result in inter-band transitions to a state of positive effective mass. The optimal value od $`h^2`$ also depends on other system parameters, and as a result we found that $`h^22`$ turned out to be optimal for our purpose.

## V Conclusion

In summary, we have studied the dynamics of an atomic condensate in a periodic optical lattice using an effective mass method. We have obtained an effective equation of motion governing the time evolution of the envelope of the condensate wave function in which the periodic external potential appears in the form of an effective mass, which can be either positive or negative, and a global drift. Numerical calculation confirmed that this envelope function approach provides us with useful qualitative insight into the condensate dynamics. Our study focused on negative effective masses, in which case the condensate behaves as if it subject to an inverted nonlinearity and an inverted external potential.

We emphasize that caution must be used when applying the effective mass concept, which neglects the potentially important coupling between different Bloch modes. This is particularly so for a system with negative effective mass moving in an external potential, as illustrated in Sec. IV.2. An ordinary particle inside an external potential will move towards the potential minimum, in the mean time gaining kinetic energy. A particle with negative effective mass, however, will move towards the potential maximum and in doing so, it looses kinetic energy. This process cannot go on forever since there is only a finite amount of kinetic energy for the system to loose. In the process, the effective mass will eventually change from negative to positive.

Our study also demonstrates that, as an alternative to the established Feshbach resonance method, changing the effective mass of a Bose-Einstein condensate provides us with an additional way to control the nonlinear atom-atom interaction in a condensate. This approach should be particularly useful for systems with no Feshbach resonances at convenient magnetic field strengths, or when the presence of external magnetic fields is undesirable. Controlling atom-atom interactions via their effective mass also provides the possibility to induce anisotropic nonlinear interactions, due to the tensorial character of the effective mass in higher dimensions.

Our study here focused on the dynamics of the system. In the future, it will also be interesting to study its static properties. For doing so, one needs to provide an additional confining potential. We must keep in mind that a confining potential for a negative effective mass state is an anti-trap, instead of a trap. This suggests the possibility of trapping strong-field-seeking states using the magnetic traps implemented in current BEC experiments.

###### Acknowledgements.

This work is supported in part by the US Office of Naval Research, by the National Science Foundation, by the US Army Research Office, by NASA, and by the Joint Services Optics Program. L. O. Baksmaty wishes to thank the Horton Foundation. We would also like to thank Michael Banks for invaluable computer support.

## Appendix A Derivation of the effective mass equation

In this Appendix, we give a derivation of the effective mass equation (12). For the sake of simplicity, we restrict ourselves to a one-dimensional system and neglect the nonlinear terms, which can be added straightforwardly.

We are concerned with the solution of the Schrödinger equation

$$i\mathrm{}\frac{\psi (z,t)}{t}=\left[H_0+U(z)\right]\psi (z,t),$$
(22)
where $`H_0=\widehat{p}^2/2m+V(z)`$ is the Hamiltonian for periodic potential $`V`$ and $`U`$ is some additional potential that varies on a spatial scale much larger than the lattice period.

### A.1 The $`𝐤𝐩`$ perturbation method

Associated with each quasi-momentum $`\mathrm{}k`$ is a set of Bloch functions $`\varphi _{nk}`$ with energy $`E_{nk}`$ as defined in Eq. (9). In this section, the power expansion around the central wave vector $`k_0`$ of the energy band function $`E_{nk}`$ is investigated. In particular, we seek expressions for first and second derivatives of $`E_{nk}`$ with respect to $`k`$. These expressions will become useful in the derivation of Eq. (12).

The “$`𝐤𝐩`$” perturbation approach is a straightforward method to relate the energies and wave functions at nearby points in the quasi-momentum space. We proceed by writing the Bloch function $`\varphi _{n𝐤}`$ as

$$\varphi _{nk}(z)=e^{ikz}u_{nk}(z),$$

where $`u_{nk}`$ is the cell periodic function. From Eq. (9), we immediately have

$$\left[\frac{\widehat{p}^2}{2m}+V(z)+\frac{\mathrm{}^2k^2}{2m}+\frac{\mathrm{}k}{m}\widehat{p}\right]u_{nk}=E_{nk}u_{nk}.$$
(23)
We then perform a power expansion around $`k_0`$ up to second order, such that

$`k`$ $`=`$ $`k_0+\delta k,`$
$`E_{nk}`$ $`=`$ $`E_{nk_0}+(\delta k)\left({\displaystyle \frac{E_{nk}}{k}}\right)_0+{\displaystyle \frac{1}{2}}(\delta k)^2\left({\displaystyle \frac{^2E_{nk}}{k^2}}\right)_0,`$
$`u_{nk}`$ $`=`$ $`u_{nk_0}+(\delta k){\displaystyle \underset{mn}{}}a_m^{(1)}u_{mk_0}+(\delta k)^2{\displaystyle \underset{mn}{}}a_m^{(2)}u_{mk_0},`$
where $`(\mathrm{})_0`$ means that the corresponding derivative is evaluated at $`k=k_0`$. Next we insert these expansions into Eq. (23), and equate the coefficients of different powers of $`\delta k`$. For the zeroth order, we find

$$\left[\frac{\widehat{p}^2}{2m}+V(z)+\frac{\mathrm{}^2k_0^2}{2m}+\frac{\mathrm{}k_0}{m}\widehat{p}\right]u_{nk_0}=E_{nk_0}u_{nk_0},$$

which is simply the definition of $`u_{nk_0}`$. The first-order term yields

$`\left[\left({\displaystyle \frac{E_{nk}}{k}}\right)_0{\displaystyle \frac{\mathrm{}^2k_0}{m}}+{\displaystyle \frac{\mathrm{}}{m}}\widehat{p}\right]u_{nk_0}=`$
$`{\displaystyle \underset{mn}{}}a_m^{(1)}(E_{mk_0}E_{nk_0})u_{mk_0}.`$ (24)
Taking the dot product of Eq. (24) with $`u_{nk_0}`$ and integrating over $`z`$ in the first Brillouin zone, we find

$`{\displaystyle \frac{1}{\mathrm{}}}\left({\displaystyle \frac{E_{nk}}{k}}\right)_0={\displaystyle \frac{\mathrm{}k_0}{m}}+{\displaystyle \frac{1}{m}}p_{nn}(k_0)`$
$`={\displaystyle \frac{1}{m}}\varphi _{nk_0}|\widehat{p}|\varphi _{nk_0}=v_g,`$
where we have used the orthonormality condition of the Bloch functions and

$$p_{nl}(k)\frac{2\pi }{\mathrm{\Omega }}u_{nk}|\widehat{p}|u_{lk}.$$

One can see that the drift velocity as defined in Eq. (13) is just the gradient of the band energy at $`k_0`$.

Carrying out a similar procedure with $`u_{mk_0}`$ ($`mn`$), we get from (24) the expression for $`a_m^{(1)}`$ as

$$a_m^{(1)}=\frac{\mathrm{}}{m}\frac{p_{mn}}{E_{mk_0}E_{nk_0}}.$$

Analogously, by using the second-order equation, we can find the second derivative of $`E_{nk}`$ as

$$\left(\frac{^2E_{nk}}{k^2}\right)_0=\frac{\mathrm{}^2}{m}+\frac{\mathrm{}^2}{m^2}\underset{mn}{}\frac{2p_{mn}p_{nm}}{E_{mk_0}E_{nk_0}}.$$

Since $`1/m^{}=(^2E_{nk}/k^2)_0/\mathrm{}^2`$, we have

$$\frac{1}{m^{}}=\frac{1}{m}+\frac{1}{m^2}\underset{mn}{}\frac{2p_{mn}p_{nm}}{E_{mk_0}E_{nk_0}}.$$
(25)
### A.2 Derivation of Eq. (12)

We now proceed to derive Eq. (12, adaptingwith slight modifications from the approach of Ref.book . First, we introduce the Fourier expansion of the envelope function as

$$f_n(z,t)=A_{nk}(t)e^{i(kk_0)z}𝑑k.$$

From Eq. (11), we have

$$\psi (z,t)=\underset{n}{}A_{nk}(t)\chi _{nk}(z)e^{iE_{nk_0}t/\mathrm{}}𝑑k,$$
(26)
where

$$\chi _{nk}(z)=e^{i(kk_0)z}\varphi _{nk_0}(z)=e^{ikz}u_{nk_0}(z).$$

After some algebra, one finds the matrix elements of $`H_0`$ and $`U`$ as

$`\chi _{nk}|H_0|\chi _{lq}=`$
$`\delta (qk)\left\{\left[E_{lk_0}+{\displaystyle \frac{\mathrm{}^2(q^2k_0^2)}{2m}}\right]\delta _{nl}+{\displaystyle \frac{\mathrm{}(qk_0)}{m}}p_{nl}(k_0)\right\},`$
$`\chi _{nk}|U|\chi _{lq}=\delta _{nl}\stackrel{~}{U}(kq),`$
where

$$\stackrel{~}{U}(kq)=\frac{1}{2\pi }U(z)e^{i(qk)z}𝑑z.$$

Inserting (26) into (22), we have

$$i\mathrm{}\frac{A_{nk}(t)}{t}=\frac{\mathrm{}^2(k^2k_0^2)}{2m}A_{nk}+\frac{\mathrm{}(kk_0)}{m}\underset{l}{}p_{nl}(k_0)A_{lk}+𝑑q\stackrel{~}{U}(kq)A_{nq}.$$
(27)
The only terms that couple different bands ($`ln`$) are those involving the momentum matrix elements $`p_{nl}`$. As long as these inter-band couplings are weak enough, we can treat them as weak perturbations. To lowest order, we have

$$A_{lk}=\frac{\mathrm{}(kk_0)}{m}\frac{p_{ln}(k_0)}{E_{lk_0}E_{nk_0}}A_{nk}.$$

Inserting this expression into Eq. (27), we have

$$i\mathrm{}\frac{A_{nk}}{t}=\left[\frac{\mathrm{}^2(\delta k)^2}{2m^{}}+\mathrm{}v_g(\delta k)\right]A_{nk}+𝑑q\stackrel{~}{U}(kq)A_{nq},$$

where $`\delta k=kk_0`$. Performing an inverse Fourier transform from $`A_{nk}(t)`$ back to $`f_n(z,t)`$ and replacing $`\delta k`$ by its configuration space representation $`i/z`$, we finally obtain the one-dimensional version of Eq. (12). The nonlinear terms can be added to this equation in a straightforward way.