Superfluidity of Bose–Einstein condensate in an optical lattice: Landau–Zener tunnelling and dynamical instability

Superflow of a Bose–Einstein condensate in an optical lattice is represented by a Bloch wave, a plane wave with periodic modulation of the amplitude. We review the theoretical results of the interaction effects in the energy dispersion of the Bloch waves and in the linear stability of such waves. For sufficiently strong repulsion between the atoms, the lowest Bloch band develops a loop at the edge of the Brillouin zone, with the dramatic consequence of a finite probability of Landau–Zener tunnelling even in the limit of a vanishing external force. Superfluidity can exist in the central region of the Brillouin zone in the presence of a repulsive interaction, beyond which Landau instability takes place where the system can lower its energy by making a transition into states with smaller Bloch wavenumbers. In the outer part of the region of Landau instability, the Bloch waves are also dynamically unstable in the sense that a small initial deviation grows exponentially in time. In the inner region of Landau instability, a Bloch wave is dynamically stable in the absence of persistent external perturbations. Experimental implications of our findings will be discussed.


Introduction
As one of the most remarkable macroscopic quantum phenomena, superfluidity has attracted enormous attention since its discovery [1]- [3]. The realization of Bose-Einstein condensation in dilute atomic gases [4] has provided physicists with a new fertile ground for exploring many aspects of this fascinating phenomenon, including frictionless current [5] and vortices [6]. In particular, there has been increasing interest and effort in studying superfluidity and related phenomena in the system of a Bose-Einstein condensate (BEC) in an optical lattice, such as Landau-Zener tunnelling [7]- [10], Josephson effect [11] and dynamical instability [12]- [14]. Recently, a quantum phase transition between superfluidity and a Mott insulator [15,16] was observed in such a system. Superflow of BEC in an optical lattice is represented by a Bloch wave, which can be regarded as a plane wave modulated by the periodic potential. Bloch waves and Bloch bands are the basic language and concepts in periodic systems. They are also essential to the study of superfluidity and many of its related phenomena in this periodic BEC system.
One interesting phenomenon is the tunnelling between Bloch bands at the edge of the Brillouin zone. Previously, people studied this inter-band tunnelling experimentally by dragging cold atoms with accelerated optical lattices [17]. The cold atoms are very dilute and the interaction between them can be ignored entirely. It is then curious to know how the tunnelling will change if the cold atoms are replaced with a BEC, where the atomic density is relatively large and the interaction between atoms can no longer be neglected. Our study shows that the interaction will increase the tunnelling probability in general [8,18]. When the repulsive interaction is over a critical value, a loop appears in the lowest Bloch band, resulting in a non-zero 104.3 tunnelling probability even in the adiabatic limit [8], [18]- [21]. This breakdown of adiabaticity is the result of superfluidity and can be viewed as a hysteresis phenomenon [8,22].
The most dramatic manifestation of superfluidity is that the boson system can maintain its current speed in a very tight space, such as a narrow capillary. Landau proposed a criterion for superfluidity [2,3]: if the current moves slower than sound, it experiences no viscosity; otherwise, the system suffers an instability and loses its speed. In a homogeneous BEC, the sound speed is proportional to the square root of the BEC density. It is not clear how superfluidity can be defined and studied for a BEC in an optical lattice. Our method is to examine the energy dispersion of Bloch waves. When the Bloch waves are energy minima of the system, they represent superflows; when the Bloch waves are energy saddle points, they suffer Landau instability: external disturbances can cause the system to emit phonons, thus reducing its speed.
The BEC Bloch waves can be dynamically unstable, i.e. the system diverges away from the original Bloch state upon small instantaneous perturbation [13]. This instability does not exist either for a BEC in a free space or for a free particle in a periodic lattice; it happens only when there are both interaction and a periodic lattice. The dynamical instability has been observed experimentally [12,23] and further confirmed by other theoretical studies and experimental observations [21,24,25]. This dynamical instability is quite a general phenomenon and widely exists in other systems, such as a BEC confined in a ring [26,27]. Its general implication and relation to superfluidity will be discussed in this paper.
We will discuss these phenomena in detail and review many interesting results for the quasi-one-dimensional case, where the optical lattice is created by two counter-propagating offresonance lasers, while the lateral motion of the BEC can be either ignored or confined [7,10,28]. In section 2, a brief introduction to the mean-field theory of BEC systems is given for the sake of completeness and introduction to the notation. In section 3, we define Bloch waves and Bloch bands for BECs in optical lattices and present some general results and the numerical methods to find them. In section 4, we simplify the system to a two-level model to study the tunnelling between Bloch bands. We point out the tunnelling in the adiabatic limit is related to a general problem, the adiabatic evolution of nonlinear quantum states. In section 5, the Landau instability and dynamical instability of the Bloch waves are studied, including the tight-binding limit. Other approaches to superfluidity and experimental observation and general implications of dynamical instability are discussed. Finally in section 6, we summarize the paper.

Mean-field theory of BEC systems
Even though the BEC is quite dense compared to the cold atoms, it is still very dilute with typical densities at 10 −20 m −3 , which is thousands of times more dilute than air. Due to this diluteness, along with the extreme low temperature, the interaction between atoms can be characterized by the s-wave scattering and modelled with a δ-function. When the number of atoms is large, the BEC system is very well described by the mean-field theory given by the celebrated Gross-Pitaevskii (GP) equation [29]: where m is the atomic mass, a s is the s-wave scattering length and V (x) is the external potential imposed on the system. In our case, the potential is the optical lattice created by two laser beams and is given by where V 0 is proportional to the laser intensity and k L is the wavenumber of the laser. For simplicity, we will instead use the following dimensionless GP equation: In the above equation, all the variables are scaled to be dimensionless with the system's basic parameters: the strength of the periodic potential v is in units of 4h 2 k 2 L /m, the wavefunction ψ in units of √ n 0 , where n 0 is the averaged BEC density, x is in units of 1/2k L and t is in units of m/4hk 2 L . The coupling constant c = π n 0 a s /k 2 L . This system can also be regarded as a Hamiltonian system governed by the grand canonical Hamiltonian where µ is the chemical potential. The GP equation (2.3) can be obtained by variation of the Hamiltonian, i∂ψ/∂t = δ H/δψ * .

Bloch waves and Bloch bands
Among all possible solutions of equation (2.3), there are states which extremize the Hamiltonian (2.4) and hence satisfy the time-independent GP equation: We call these extremum states nonlinear eigenstates, which are also called nonlinear coherent modes elsewhere [30]; accordingly, µ are nonlinear eigenvalues. Bloch waves are the nonlinear eigenstates of the following form: where φ k (x) is a periodic function of period 2π and k is the Bloch wavenumber. In particular, from equation (3.1) we have the following equation for each Bloch wave state φ k : The set of eigenvalues µ(k) then forms Bloch bands.
When v = 0, the plane waves are the solutions of equation (3.1), and they represent a BEC flow with a speed of k. As is well known, when the speed is smaller than the speed of sound, k < √ c, it is a superflow; when k > √ c, it develops a Landau instability and loses the superfluidity. Bloch waves are the counterpart of these plane waves in a periodic system. Whether these Bloch waves represent superflows or not is determined by the energy dispersion of these Bloch waves, as we will show later.
In the linear case, c = 0, all eigenstates are Bloch waves [31]. The situation is very different for the periodic nonlinear system, where it is possible to have non-Bloch wave eigenstates. Furthermore, there are nonlinear Bloch waves that do not have their linear counterparts. Here is a beautiful example [14,19]. Equation (3.1) has an exact solution of a Bloch wave at the Brillouin zone edge k = 1/2: This Bloch wave state only exists when c ≥ v so that it has no linear counterpart. This shows that there are some 'extra' nonlinear Bloch waves when the nonlinearity is strong enough. Their corresponding 'extra' eigenvalues should make the nonlinear Bloch band look very different from the linear Bloch band. This is indeed the case: our numerical results show a loop structure in the lowest Bloch band when c > v, as seen in figure 1. This loop structure was first found in [8] with a two-mode approximation, then confirmed in [20,21]. In particular, in [21], a loop is also found for the second band at the Brillouin zone centre. This loop structure is a consequence of superfluidity. One clear indication is that the exact Bloch wave at the Brillouin zone edge, equation (3.4), represents a flow of non-zero speed, √ [19]. This is quite remarkable if we recall that all the Bloch waves at the zone edge carry no currents in the linear case. For a free particle, the flow e ix/2 is completely stopped due to first-order Bragg scattering by the periodic potential, leading to a zero-current Bloch wave at the edge, k = 1/2. However, when the superfluidity gets stronger, that is, when c > v, the flow e ix/2 can no longer be stopped by Bragg scattering, yielding a zone-edge Bloch wave that carries a current. The connection between superfluidity and the loop structure is discussed in depth in the context of hysteresis by Mueller [22].
These Bloch states can be prepared experimentally, at least for the ones in the lowest Bloch band, by adiabatic control. For a trapped cigar-shaped BEC, we slowly turn on an optical lattice along its longitudinal direction; we then have a BEC in the Bloch state at the centre of the Brillouin zone. Other Bloch states can be achieved by accelerating the lattice for a certain 104.6 amount of time. These experimental techniques of adiabatic turning-on and accelerating optical lattices have been demonstrated successfully with either cold atoms [17,32] or BECs [10].
There are several numerical methods to find the Bloch wave φ k . The method described in [13] is to expand φ k in Fourier series and minimize the Hamiltonian (2.4) in the space expanded by Fourier coefficients. Another one is used in [20], where equation (3.1) is first solved for different values of µ and k and then the Bloch band µ(k) is found by interpolation.
A much better method is the following. One still starts with the expansion of φ k in Fourier series: (3.5) where N is the cut-off. With the substitution of the above equation into equation (3.3), one obtains 2N + 1 equalities for the coefficients of each Fourier term e inx : Finally, the Bloch wave is found by minimizing the following sum: The results in figure 1 are obtained with this method by using N = 10.

Tunnelling between Bloch bands
Consider the scenario that we drag the BEC across the entire Brillouin zone with an accelerating lattice. If the BEC is initially in a Bloch state belonging to the lowest Bloch band and the acceleration is small enough, the system will stay in the band and undergo Bloch oscillations [7,10]. As one increases the acceleration, the BEC will have an increasing chance of tunnelling into the upper band at the edge of the Brillouin zone, where the bandgap is the smallest. Without interaction, this tunnelling is nothing other than the famous Landau-Zener tunnelling [33], a well understood phenomenon. With the interaction in our case, the tunnelling behaviour is strongly modified, in particular, by the loop structure in the Bloch band, as we will see later.
In an accelerating lattice, the system is described by where α is the acceleration. Even in the linear case, c = 0, the above equation is difficult to analyse. To reduce the mathematical complexity, we will approximate it with a two-level model without losing the essential physics.

Two-level model
The tunnelling mainly occurs around the edge of the Brillouin zone, k = 1/2, where the wavefunction is dominated by two plane wave components. We zoom in on this region and write where |a| 2 +|b| 2 = 1. Following the simple algebra described in [8] by substituting the above twomode wavefunction into equation (4.1), we obtain a two-level nonlinear Schrödinger equation: where where γ = αt and κ = |b| 2 − |a| 2 . However simplified it may appear, this two-level model captures the essence of our BEC system as it reproduces the looped Bloch band (see figure 2). In this simple model, the 'Bloch bands' are obtained by finding the eigenvalues of H (γ ): The eigenvalues µ(γ ) are also called adiabatic energy levels. How the levels change with the coupling strength c is shown in figure 2, where we again see the appearance of the loop for c > v in the lower energy level.

104.8
For this two-level model, the original tunnelling problem in the BEC system can be restated as follows: if the system starts in the lower level, as the energy bias γ changes with the sweeping rate α from the far negative end to the far positive end, what is the probability that the system ends up in the upper level?
When c = 0, this is precisely the well-known Landau-Zener model, which can be solved exactly. The tunnelling probability is [33] It is clear from this formula that, in the adiabatic limit α → 0, the tunnelling probability is zero, i.e. the system stays in the lower level throughout the entire process. This is nothing other than a special case of the quantum adiabatic theorem [34]: in the adiabatic change of a parameter the system starting in an eigenstate stays in the same instantaneous eigenstate. When c > v, this adiabaticity is broken down by the loop structure. Suppose we start with a state on the lower adiabatic level and move it up along the branch by changing γ so slowly that little tunnelling to the upper level is generated. After passing the crossing point X, the state remains in the course of moving up in energy until hitting the terminal point T, where it has no way to go any further except to jump to the upper and lower levels. This leads to a nonzero probability tunnelling into the upper level in the adiabatic limit.
This qualitative analysis is corroborated by our numerical integration of equation (4.3). The tunnelling probabilities as a function of the sweeping rate α are plotted in figure 3 for different values of c at v = 0.1. In general, the tunnelling probability is enhanced due to the interaction. However, its dependence on α changes dramatically with the increasing interaction strength c. When c < v, the tunnelling probability approaches zero in the adiabatic limit, α → 0. It appears that the approach is exponential, as supported by a recent analysis [18]. At the critical point, c = v, the tunnelling probability still turns to zero as the sweeping slows down; however, not exponentially.
When c is over the critical value, it becomes apparent in figure 3 that the tunnelling curve tends to intersect the vertical axis at a finite value, indicating the tunnelling probability is not zero in the adiabatic limit, α → 0. Of course, one can never be absolutely sure that the tunnelling probability at α → 0 is not zero by numerical calculation since there is always a limit on the smallness of α in a computer code. This uncertainty is ended in [18], where this nonzero adiabatic tunnelling probability is obtained analytically.
One can also directly integrate numerically the original Schrödinger equation (4.1) to find out the tunnelling probability and how it depends on the acceleration α. This is done in [35], where one sees a similar crossover from the exponential to non-exponential dependence of the tunnelling probability on α. We argue that this crossover behaviour can serve as indirect evidence of the existence of the loop structure in the Bloch band. We believe that this crossover should be observable in an experiment similar to the one conducted in [10].
The nonlinear Landau-Zener tunnelling studied here is likely to be as general as the linear counterpart and should exist in many physical systems. We found one example in molecular magnets [36], where the jamming effect among interacting molecular spins is the result of nonlinear Landau-Zener tunnelling.

Nonlinear quantum adiabatic theorem
The nonlinear two-level model (4.3) is a very good example of a general problem: how a nonlinear quantum system, such as a BEC system, responds to the change of a system parameter. Of particular interest is that the parameter change is adiabatic since adiabatic control is an important tool in preparing a quantum system (linear or nonlinear) in a desired quantum state [37]. For a linear quantum system, the answer is the quantum adiabatic theorem, which can be found in standard textbooks [34]. The theorem states that, if the system is initially in a non-degenerate eigenstate, it stays in the same eigenstate as the parameter changes adiabatically.
On the other hand, the answer for nonlinear quantum systems is not as straightforward. The nonlinearity proves to be a big obstacle for generalizing the linear quantum adiabatic theorem [38]. Very non-trivial modifications to the linear case are expected as it is found in the nonlinear two-level model (4.3) that tunnelling between energy levels happens even in the adiabatic limit.
In a recent paper [39], we have successfully generalized the quantum adiabatic theorem to nonlinear quantum systems by combining ideas from classical adiabatic dynamics and quantum geometric phases. It is found that the eigenstates correspond to fixed points of an equivalent classical problem. The adiabatic evolution of these eigenstates depends on the nature of the corresponding fixed points. If the fixed point is elliptic, the system can stay in the same eigenstate as in the linear case. Otherwise, when the fixed point is hyperbolic, the system will stray away from the eigenstate, leading to tunnelling between different eigenstates. This adiabatic condition 104.10 can also be specified in terms of the Bogoliubov spectrum. If the Bogoliubov spectrum is real, the adiabaticity is kept; otherwise, it is broken.

Superfluidity and Landau instability
A quantum Bose liquid or gas at very low temperatures possesses a very remarkable property known as superfluidity: its flows suffer no viscosity through capillaries or other types of tight spaces when its speed is slower than a critical value. This was discovered even before the concept of Bose-Einstein condensation [1]. Landau [2,3] gave a simple explanation for this intriguing phenomenon. He argued that a quantum current suffers friction only when the system's elementary excitation, the phonon, can lower its energy. With the Galilean transformation, he found that a quantum liquid flowing with a speed smaller than the sound speed enjoys such a frictionless transport. When it flows with a speed larger than the sound speed, the quantum liquid can lower its energy by emitting phonons, leading to the slowdown of its speed.
Such a loss of superfluidity can also be regarded as an instability: the system can no longer stay in its original state, thus keeping its flow against perturbation by the roughness of capillaries or other external disturbance. We will call this instability Landau instability. Note that this Landau instability is different from the Landau instability discussed in plasma physics [40].
In our case, the study of superfluidity or Landau instability is complicated by the presence of the optical lattice, which breaks the translational symmetry along its flow direction and modulates the BEC density. This complication is only mathematical. The physics is still the same, that is, the Landau instability is determined by whether the elementary excitation around a BEC flow lowers the system's energy or not. If it does not, the BEC flow is a local energy minimum of the system and it is a superflow. Otherwise, the flow is an energy saddle point and it suffers Landau instability (see figure 4).
Whether these Bloch waves are local energy minima is determined by their energy dispersions. To calculate the dispersion, we perturb the system around a Bloch state: Due to the periodicity of our system and the Bloch waves, these perturbations can be decomposed into different modes labelled by q: where q is also a kind of Bloch wavenumber and ranges between −1/2 and 1/2. The perturbation functions u k and v k are of a periodicity of 2π in x. By substituting the perturbed state (5.1) into the Hamiltonian (2.4)and neglecting terms of orders higher than two, we obtain a quadratic form of the energy deviation, which is block diagonal in q. Each block is given by where Schematic illustration of superfluidity, Landau instability and dynamical instability. A dynamically unstable state must be a saddle point; a state that is an energy saddle point can be dynamically stable.
If M k (q) is positive definite for all −1/2 ≤ q ≤ 1/2, the Bloch wave φ k is a local minimum and represents a superflow. Otherwise, δ E k can be negative for some q, indicating that the Bloch wave is a saddle point and describes only a normal flow. The special case v = 0, when the optical lattice is turned off, is the simple case considered by Landau [2]. In this case, the Bloch state φ k becomes a plane wave e ikx and the operator M k (q) becomes a 2 × 2 matrix: whose eigenvalues are Since λ + is always positive, the matrix M k (q) is not positive definite only when λ − ≤ 0, or equivalently, |k| ≥ q 2 /4 + c. It immediately follows that the BEC flow e ikx becomes a saddle point when the flow speed exceeds the sound speed, |k| > √ c. This recovers exactly the Landau condition for the breakdown of superfluidity [2], which has recently been confirmed experimentally on BEC [5].
The situation becomes more complicated in the case of v > 0, when the optical lattice is turned on. We resort to numerical calculations to examine whether the matrices M k (q) are positive definite and the corresponding BEC Bloch waves are local energy minima. Our focus is on the Bloch states in the lowest Bloch band, excluding the loop. The results are summarized in region expands with increasing atomic interaction c and occupies the entire Brillouin zone for sufficiently large c. On the other hand, the lattice potential strength v does not affect the superflow region very much as we see in each row. The phase boundaries for v 1 are well reproduced from the analytical expression k = q 2 /4 + c for v = 0, which is plotted as triangles in the first column.
The Landau instability of the states in the lower branch of the loop has been discussed in [21]. It is found that the area of superfluidity increases when the interaction c increases or the lattice strength v decreases.

Remarks on superfluidity
So far, we have been identifying superfluidity as a stability in the sense of Landau [2,3]. In this approach, there is an implicit assumption that the boson system has already achieved phase coherence and can be described by a complex order parameter, the macroscopic wavefunction. Whether the state described by the macroscopic wavefunction possesses superfluidity is then examined by how the system responds to perturbations. With the assumption of phase coherence, this approach is not well-equipped to study quantum phase transitions at zero temperature, such as the superfluid-Mott-insulator transition [15,16]. Nevertheless, this approach can still provide some signals for the onset of some possible phase transition as we will discuss in the context of dynamical instability in the next section.
There are various other approaches to superfluidity. One of particular interest was proposed by Fisher et al in [41], where superfluidity is defined as how strong the system is affected by an artificially imposed phase twist at boundaries. Since only the superfluid part of the system responds to the phase twist, by calculating the change of the system energy, one can find out the superfluid fraction of the Bose system. This approach has its root in the Josephson effect and Anderson's generalized concept of rigidity [42].
This phase twist approach does not assume the phase coherence; quite the opposite, it is particularly designed to calculate how much phase coherence a system has. It has been widely used in the Bose-Hubbard model to study the superfluid-Mott-insulator transition [43]. With the recognition that a BEC in an optical lattice can be described by the Bose-Hubbard model [16], this approach has now been adopted to study the superfluid-Mott-insulator transition for various systems of BECs in optical lattices [44,45].
Both approaches are of perturbative nature, i.e. they study how the Bose system responds to small perturbations. In Landau's approach, all possible disturbances are considered: one is able to determine the robustness of a superfluid as measured by the critical velocity. In contrast, only a special type of perturbation is considered in the method of phase twist, namely, an artificial phase twist imposed at the boundaries. This perturbation is sensitive to the interplay between kinetic energy and interaction. It is not clear how this phase twisting corresponds to the usual perturbations (such as disorder or roughness) which would degrade the flow of a normal fluid. In short, the relation between these two pictures is still a subject of active debate [46].

Dynamical instability
One unique feature in the system of BECs in optical lattices is dynamical instability, which does not exist in the absence of either atomic interaction, c = 0, or an optical lattice, v = 0. Many of the BEC Bloch waves can be dynamically unstable against certain perturbation modes q only 104.14 when both factors are present. By dynamical instability, we mean that small deviations from a state grow exponentially in the course of time evolution (see figure 4).
The dynamical instability can be determined from the linear stability analysis of the GP equation (2.3). Assume that the system experiences a small disturbance δφ k at a Bloch state: (5.8) where the disturbance can be similarly written as Plugging the above wavefunction into equation (2.3) and keeping only the linear terms, we arrive at a linear dynamical equation: The eigenvalues ε k (q) of the matrix σ z M k (q) determine the dynamical instability: if they are real for all −1/2 ≤ q ≤ 1/2, the Bloch state is dynamically stable; if otherwise, i.e. some of them are complex, it is dynamically unstable. One who is familiar with the quasi-particle excitation in a quantum fluid may immediately notice that σ z M k (q) also appears in the Bogoliubov approach and gives us the spectrum of the phonon excitation [2]. When all of its eigenvalues ε k (q) are real, the positive half is the phonon spectrum while the other half is deemed non-physical. However, we will refer to the modes corresponding to the negative half as anti-phonons for the sake of easy reference in the following discussion.
Again, we first look at the case v = 0, where the eigenvalues of σ z M k (q) are ε ± (q) = kq ± q 2 c + q 4 /4. (5.12) These eigenvalues are always real, implying that the BEC flows in a free space are always dynamically stable. When v = 0, the situation is totally different: the eigenvalues ε k (q) of σ z M k (q) can be complex and Bloch states can be dynamically unstable. The results are summarized in figure 5, where these ε k (q) are complex in the dark-shaded areas. The dark shaded areas lie completely inside the light shaded areas, reflecting the rigorous conclusion that only saddle-point Bloch states can have dynamical instability. Its proof is given in the appendix. The dynamical instability is the result of the resonance coupling between a phonon mode and an anti-phonon mode by first-order Bragg scattering. As shown in the appendix, the matrix σ z M k (q) is real in the momentum representation, meaning its complex eigenvalue can appear only in conjugate pairs and they must come from a pair of real eigenvalues that are degenerate prior to the coupling. Degeneracies or resonances within the phonon spectrum or within the antiphonon spectrum do not give rise to dynamical instability; they only generate gaps in the spectra. Based on this general conclusion, we have analysed dynamical instability for two extreme limits, v 1 and c v. The limit, v 1, corresponds to the cases shown in the first column of figure 5. For this limit, we can approximate the phonon spectrum and the anti-phonon spectrum with the ones 104.15 given in equation (5.12). By equating them, ε + (q − 1) = ε − (q), for the degeneracy, we find that the dynamical instability should occur on the following curves: These curves are plotted as full dots in figure 5 and they fall right in the middle of the thin darkshaded areas. To some extent, one can regard these thin dark-shaded areas as broadening of the curves (5.13). It is noted in [21] that the relation (5.13) is also the result of ε + (q −1)+ε + (−q) = 0, which involves only the physical phonons. Therefore, the physical meaning of equation (5.13) is that one can excite a pair of phonons with total energy zero and with total momentum equal to a reciprocal wavenumber of the lattice.
The other limit, c v, is the case plotted in the first row of figure 5. The open circles along the left edges of these dark-shaded areas are given by where E 1 (k) is the lowest Bloch band of In this linear periodic system, the excitation spectrum (phonon or anti-phonon) just corresponds to transitions from the Bloch states of energy E 1 (k) to other Bloch states of energy E n (k + q), or vice versa. The above equation is just the resonance condition between such excitations in the lowest band (n = 1). Alternatively, we can write the resonance condition as So, this condition may be viewed as the energy and momentum conservation for two particles interacting and decaying into two different Bloch states E 1 (k + q) and E 1 (k − q). This is the same physical picture behind equation (5.13).
One common feature of all the diagrams in figure 5 is that there is a critical Bloch wavenumber k d beyond which the Bloch waves φ k are dynamically unstable. The onset instability at k d always corresponds to q = 1/2. In other words, if we drive the Bloch state φ k from k = 0 to 1/2 the first unstable mode appearing is always q = ±1/2, which represents period doubling. Interestingly, a BEC Bloch wave solution of double periods was found recently in [47]. Only for k > k d can longer wavelength instabilities occur. The growth of these unstable modes drives the system far away from the Bloch state and spontaneously breaks the translational symmetry of the system.
The dynamical instability of the lower loop states is studied in [21]. The range of k, where dynamical instability occurs, is found to be decreasing with increasing c and decreasing v.

Tight-binding approximation
We have so far presented a theory of how to determine the Landau and dynamical instabilities for the system of BECs in optical lattices. It is valid for the full range of parameters c and v as long as the mean-field treatment is adequate. It is, nevertheless, worthwhile to consider the system in certain extreme limits, where analytical results are available. The particular limit that we now look into is the tight-binding limit, where the lattice is so strong that the BEC system can be considered as a chain of trapped BECs that are weakly linked. This limit was studied 104.16 in [24] with its focus on dynamical instability. In this case, the BEC trapped in each well of the lattice can be approximately described by (5.17) where the real wavefunction ϕ 0 (x) is very localized in the well [0, 2π ] such that The tight-binding approximation is to write Substituting it into the GP equation (2.3), we find, after neglecting a trivial energy constant [24], Based on these two results, one can easily determine the stability phase diagram of the Bloch wave solution (5.23), which is plotted in figure 6. When U/2τ < 1, there is a right borderline for dynamical instability (see the top right corner of figure 6(a)); it is given by This border disappears for U/2τ ≥ 1. For any value of U/2τ , the left border of dynamical instability is at k = 1/4 for any q while the boundary of Landau instability is described by Note that the borderline k = 1/4 is consistent with our analysis in the last subsection. It is the result of equation (5.14) when one uses the tight-binding band energy E(k) ∝ cos(2kπ). It is evident that the phase diagram figure 6 is a natural extrapolation of the results shown in figure 5, and it also follows the trend in figure 5 as one increases the interaction: the overall area of instability decreases while the dynamical instability spreads more in the area of Landau instability.  [12,25].

Experiments
This dynamical instability was observed in a recent experiment [12], whose set-up is schematically drawn in figure 7. A cigar-shaped BEC was formed in a magnetic trap superimposed with an optical lattice. This essentially prepared a BEC in a Bloch state at k = 0, the centre of the Brillouin zone. The magnetic trap was kept on to keep the quasi-one-dimensional configuration, but its centre was suddenly shifted by x along the longitudinal direction. This is equivalent to displacing the whole BEC off the centre of the harmonic trap, then releasing it.
According to the semiclassical theory [31]: where E(k) is the system energy and ω is the longitudinal frequency of the magnetic trap.
The BEC begins to oscillate in both the real space and the momentum space (Brillouin zone).
Since E(k) ∝ k 2 near the centre of the Brillouin zone, both oscillations are harmonic when the displacement x is small, as observed in the experiment. However, when x was over a certain critical value, say, x c , it was observed that the oscillations were disrupted and became dissipative. This dissipative behaviour is caused mainly by the dynamical instability in our opinion although the authors of [12] attributed it to Landau instability. When the displacement x is small, the oscillation is confined near the centre of the Brillouin zone, where Bloch states are stable as shown in figure 5. With the increase of x, the Bloch state will be driven into the dangerous zone of dynamical instability; the oscillation then can be disrupted by the growth of modes of instability. In our units, this experiment has v ∼ 0.2 and c ∼ 0.02, which is in the regime where the dynamical instability is rampant. Further experiments in [12] found that x c increases with the decrease of the lattice strength v and there was no dissipation when the density of BEC was low. These two observations basically rule out Landau instability as the main cause of dissipation. As we see in each row of figure 5, the Landau instability, as represented by the light-shaded areas, is affected very little by the change of v. Thus x c should see no changes within experimental error when v changes. Also clear in figure 5, when c is small, almost all Bloch states have Landau instability. As a result, more severe dissipation would have been observed for smaller c if the Landau instability was the main contributor.
Very recently, there was another experiment with a similar set-up, but the optical lattice used is much stronger [25]. The advantage of strong optical lattices is that dynamical instability is more severe, rendering Landau instability less relevant. Therefore, the interpretation of the experimental observation is more definite and much less controversial than the one in [12].

104.19
There are other possible ways to observe the dynamical instability. One way is to prepare a Bloch state, then monitor how its periodicity changes over with Bragg scattering of a probing laser beam by the BEC cloud [48]. Another way is indirectly, through the breakdown of Bloch oscillations [10,11,32]. Within a proper regime of parameters, the dynamical instability can be very severe and thus drive the system far away from a Bloch state within a Bloch period. This leads to the breakdown of Bloch oscillations; the observation of this breakdown is certainly an observation of dynamical instability. However, one must be careful since the nonlinear Landau-Zener tunnelling due to the loop structure, as discussed in the last section, can also lead to the destruction of Bloch oscillations. This ambiguity can be avoided by using dilute BEC such that c < v, where the loop does not exist.

Discussion on dynamical instability
As we emphasized above, dynamical instability happens only when both the optical lattice and the interaction between atoms exist. However, it does not imply that it is special to BECs in optical lattices. On the contrary, it is very general. For example, it has been identified in a BEC confined in a torus. In the following, we will demonstrate this point with a one-mode problem. To provide a slightly different angle, we will present it in the language of operator algebra and use the Bogoliubov approach [2]. A discussion of a similar two-mode problem can be found in [49].
The one-mode problem is described by the Hamiltonian written as The Bogoliubov approach is to diagonalize the Hamiltonian (5.32) with the introduction of a new set of operators, Forb andb † to be boson operators, it requires that If u and v satisfy the Hamiltonian (5.32) becomes diagonal inb andb † , up to a constant: The phonon energy ε is easily found by solving equation (5.37), and it is When A < B, this energy becomes imaginary; we again encounter dynamical instability. The simplicity of this one-mode problem allows us to look more carefully into this instability. We notice from equation (5.38) that H = 0 when ε becomes imaginary. At the same time, we find with equation (5.37) that |u| 2 − |v| 2 = 0, violating the condition (5.36) forb andb † to be boson operators. In other words, the Hamiltonian (5.32) defies the diagonalization by the Bogoliubov approach when A < B.
The physics lying behind the failure of the Bogoliubov diagonalization is the dynamical instability: the system is intrinsically unstable and, therefore, does not have any phonon excitation. We can look into this in a different representation. We writeâ andâ † in terms of space and momentum operators: As a result, the Hamiltonian assumes the standard form of a harmonic oscillator: It is evident here that, when A < B, the system is a harmonic oscillator with a negative mass, an unstable system. A quantum system with dynamical instability, such as a harmonic oscillator with a negative mass, does not exist in nature. Its appearance in the theoretical study only signals the onset of instability of the underlying unperturbed states. In our case, these unperturbed states are the BEC Bloch states. Once these Bloch states develop dynamical instability, it means that the GP equation, which predicts the existence of these Bloch states, is no longer a good description of the system. The system enters a new unknown phase and demands a new mathematical description other than the mean-field theory, the GP equation. In [25], the new phase was called 'insulator'. There may be some experimental indications to call it 'insulator'; theoretically, we are yet to develop a theory accurately depicting the new phase. Our theory, along with the one in [24,25], can only predict when or where the dynamical instability sets in and the system enters a new phase; it does not give any indication how the new phase should look like and behave.
The dynamical instability is closely related to Landau instability. As already briefly stated above, Bloch waves with dynamical instability must also have Landau instability; on the other hand, Bloch waves with Landau instability are not necessarily dynamically unstable. This comes from a rigorous result on the relation between eigenvalues of M k (q) and σ z M k (q), as proved in the appendix: when all the eigenvalues of M k (q) are all positive, the eigenvalues ε k (q) of σ z M k (q) are all real.

Summary
In summary, we have examined the superfluidity and many of its related properties of BECs in optical lattices. We studied the tunnelling between Bloch bands and found that the tunnelling 104.21 is not zero even in the adiabatic limit when the repulsive interaction c is strong enough. This dramatic tunnelling behaviour is due to the loop structure in the Bloch band, a result of strong superfluidity. We investigated Landau instability and dynamical instability of the system. With a phase diagram, we showed where Landau and dynamical instability occur and how they change with the variation of the two system parameters, c and v. We have discussed the experimental observation of dynamical instability and its implication in this system.

104.22
One can check by straightforward calculations that all the elements, A mn , B mn and C mn are real after noticing that the Bloch state φ k (x) when expanded in Fourier series has all real coefficients. This particular expression is very convenient for numerical calculation: after neglecting Fourier terms of orders higher than the cut-off, say, N , one reduces σ z M(q) to a real (4N + 2) × (4N + 2) matrix. To our experience, N = 10 is sufficient.
(ii) If ε is an eigenvalue of σ z M(q), then −ε * is an eigenvalue of σ z M(−q).
Write equation (A.2) in another form: By some simple algebraic manipulations, one immediately observes that −ε * is an eigenvalue of σ z M(−q) with the corresponding eigenvector being {v * , u * }.
(iii) When all the eigenvalues λ(q) are positive, the eigenvalues ε(q) must all be real.
We multiply both sides of equation (A.2) by σ z , then by Y † and obtain Since the matrix M(q) is Hermitian, the left-hand side of the above equation is always real. As a result, if an eigenvalue ε(q) is complex, then for its corresponding eigenvector Y , one has When all eigenvalues of M(q) are positive, the left-hand side of equation (A.8) must be positive. In this particular case, one can conclude that the eigenvalues ε(q) are all real and share the sign of Y † σ z Y .