Direct observation of zitterbewegung in a Bose–Einstein condensate

Among the great discoveries of the enlightenment was the realization that physical laws are equivalent in all places, at all times, and for all scales; this remains a central tenet in contemporary science. Quantum simulation exploits this universality to study the behaviour of systems that are difficult to access or impossible to manipulate, by performing direct measurements on analogue systems composed of well-characterized and highly manipulable quantum building blocks. In this work, we used neutral rubidium atoms to simulate zitterbewegung, a trembling motion usually associated with relativistic electrons [1], and we illuminate its microscopic origins by drawing an analogy to the well-understood atomic physics of Rabi oscillations. The Dirac equation—describing the motion of free fermions—is an essential part of our current description of nature; by engineering new Dirac particles in novel settings, we expose the equation's properties by direct measurement. Simulations of the Dirac equation have been proposed for superconductors [2], semiconductors [3, 4], graphene [5], cold atoms [4, 6–14] and photonic systems [15]; and have been realized with cold atoms [16], trapped ions [17] and photons [18]. The ion and photon experiments demonstrated zitterbewegung for quantities analogous to position or time in the Dirac equation. Here, we directly observed a neutral-atom Bose–Einstein condensate (BEC) undergoing zitterbewegung in space and time.

Zitterbewegung, as observed here, is an example of a broader class of phenomena where a group of states with differing velocities are quantum mechanically coupled together and undergo Rabi-like oscillations [19–22]. As with the present case, eigenstates of the full Hamiltonian are static, but superpositions can tremble. Neutrino oscillations [23] are an example of this generalization: neutrinos are produced by the weak nuclear force in superpositions of the propagating (i.e. mass) eigenstates, each with a different mass and, therefore, velocity.

The precise control and direct measurement techniques available in systems of ultracold atoms, coupled with their accessible length and energy scales, make these systems ideal for quantum simulation. In this experiment, our quantum building blocks were Bose-condensed ⁸⁷Rb atoms. Using two counter-propagating Raman lasers (figure 1(a)) with wavelength λ = 790.1 nm, we coupled the atoms' $\left | f = 1, m_{\mathrm {F}} = \mp 1\right \rangle = \left |\uparrow \downarrow \right >$ atomic hyperfine states (comprising our effective two level system) to their external motion [24] with a four-photon Raman transition (figure 1(b)). In this environment, each atom's behaviour was governed by the one-dimensional Dirac Hamiltonian, making its motion analogous to that of a relativistic electron. The system's characteristic momentum ℏk_R = 2πℏ/λ—that of a single photon—specifies the recoil energy E_R = ℏ²k²_R/2m = h × 3.68 kHz, where m is the atomic mass. These recoil units set the scale for all physical quantities in our analogue system, such as the recoil velocity v_R = ℏk_R/m. The Raman lasers drove the four-photon $\left |\uparrow , \hbar k_x = p_x + 2\hbar k_{\mathrm { R}} \right >\leftrightarrow \left |\downarrow , \hbar k_x = p_x - 2\hbar k_{\mathrm { R}}\right >$ transition (resonant when p_x = 0), where ℏk_x is the atomic momentum along e_x and p_x will play the role of momentum in the Dirac equation. The simulated speed of light c* = 2v_R = 11.6 mm s⁻¹ was twice the atoms' recoil velocity, a factor of ≈ 10¹⁰ less than the true speed of light. The artificial rest energy m*c* ² = ℏΩ/2 ≲ 1E_R was a factor of ≈ 10¹⁷ less than the electron's rest energy (ℏΩ is the four-photon laser coupling strength). The effective Compton wavelength λ*_C = h/m*c* ≈ 1 μm, the approximate amplitude of zitterbewegung, exceeded that of an electron by a factor of ≈10⁶. These new scales enabled our direct measurement of zitterbewegung.

Zoom In Zoom Out Reset image size

The dynamics of our ultracold ⁸⁷Rb atoms were described by the one-dimensional Dirac equation

$$\begin{equation} \hat H_{\rm D}\vert\psi\rangle = (c^* \hat{p}_x \check\sigma_z + m^* {c^*}^2 \check\sigma_x)\vert\psi\rangle, \end{equation} \tag{ 1 }$$

where $\hat p_x$ is the momentum operator; $\check \sigma _{x,y,z}$ are the Pauli spin operators; and |ψ〉 is represented as a two-component spinor, whose components are defined by $\left |\uparrow \downarrow ,p_x\right >$, the m* = 0 eigenstates of $\hat H_{\mathrm { D}}$. For the massless, m* = 0 case, this equation simply describes particles (positive energy) or anti-particles (negative energy) travelling with velocity ±c*, as depicted by the dashed lines in figure 1(c). The mass term couples together these m* = 0 states, producing an avoided crossing (solid curves in figure 1(c)) with energy given by the familiar relativistic dispersion E(p_x) = ± (p²_xc*² + m*²c*⁴)^1/2, gapped at p_x = 0 by twice the rest energy. In our atomic analogue, the two massless states coupled by the effective rest energy physically corresponded to the atomic states $\left |\uparrow \downarrow \right >$ moving with velocity ±c*.

1.1. Zitterbewegung equations of motion

Zitterbewegung arises because the Pauli matrices associated with the two terms in the Dirac equation do not commute. In the Heisenberg representation of quantum mechanics the operators, not the wavefunctions, depend on time: for example $\hat v_x = \mathrm {d}\hat x/\mathrm {d}t = [\hat x, \hat H_{\mathrm { D}}]/\mathrm {i}\hbar$. In this formalism, the velocity operator obeys the differential equation

$$\begin{equation} \frac{\mathrm{d}^2\hat v_x}{\mathrm{d}t^2} + \Omega^2 \hat v_x = \frac{2\Omega {c^*}^2}{\hbar} \hat p_x \check\sigma_x. \end{equation} \tag{ 2 }$$

For an initial state $\left |\uparrow ,p_x = 0\right >$, which gives initial conditions $\langle \hat v_x \rangle = c^*$ and $\langle \mathrm {d} \hat v_x/\mathrm {d}t \rangle = 0$, the expectation values of the position and velocity observables oscillate with the zitterbewegung frequency Ω according to

$$\begin{equation} \langle \hat x (t)\rangle = x(0) + \frac{\lambda^*_{\rm C}}{4\pi} \sin(\Omega t),\quad \langle \hat v_x (t)\rangle = c^* \cos(\Omega t). \end{equation} \tag{ 3 }$$

Initial states with $\langle \hat p_x\rangle \neq 0$, or localized wave packets, follow more complex trajectories [25]. Zitterbewegung, as usually understood, refers to trembling in position; an oscillatory velocity is the obvious dual. In these experiments, we observed the out-of-phase oscillation of these conjugate quantities.

1.2. The atomic Dirac Hamiltonian

The one-dimensional Dirac Hamiltonian for a system of ⁸⁷Rb atoms can be realized by coupling different spin-momentum states. The three m_F states comprising the 5S_1/2, f = 1 electronic ground state manifold are subject to a two-photon Raman process (figure 1(a)), and the atomic dynamics along e_x are described by the three-level Hamiltonian

$$\begin{eqnarray} &&\hat H^{3\times 3} = \left(\frac{\hbar^2{\hat k_x}^2}{2m}+4E_{\rm R}\right)\check 1_3 + \frac{2\hbar^2 k_{\rm R} \hat k_x}{m} \check\sigma_{3,z} +\frac{\hbar\Omega_2}{2}\check\sigma_{3,y} +(\check 1_3- \check\sigma_{3,z} )\hbar\epsilon, \end{eqnarray} \tag{ 4 }$$

where Ω₂ is the two-photon Raman coupling strength; is the quadratic Zeeman shift that energetically displaces the m_F  = 0 state; $\check \sigma _{3,z}$ are the generalized Pauli operators for a spin-1 system; and $\check 1_3$ is the 3 × 3 identity. We concentrate on the avoided crossing at k_x = 0 between the states that adiabatically connect to $\left |m_F = -1\right >$ and $\left |m_F = +1\right >$. By adiabatically eliminating the lowest-energy eigenstate, we obtain the effective two-level Hamiltonian

$$\begin{equation} \hat H^{2\times 2} = \left[ \frac{\hbar^2 (\hat k_x^2+4k_{\rm R}^2)}{2m} +E_4 \right] \check 1 + \frac{2 \hbar^2 k_{\rm R} \hat k_x}{m} \check\sigma_{x} +\frac{\hbar\Omega }{2}\check\sigma_{z}, \end{equation} \tag{ 5 }$$

which includes a global rotation of the system $\check \sigma _x \rightarrow \check \sigma _y$, $\check \sigma _y \rightarrow \check \sigma _z$, $\check \sigma _z \rightarrow \check \sigma _x$. For k_x/k_R ≪ 1, the effective coupling is Ω = ℏΩ²₂/2(4E_L − ℏ) and E₄ = E_LΩ/Ω₂. Ignoring the uniform energy offset, we identify the parameters from the Dirac Hamiltonian (1): the effective c* = 2ℏk_R/m is twice the atomic recoil velocity, the rest energy m*c*² = ℏΩ/2 is the coupling strength and the Compton wavelength λ*_C = h/m*c* = 8πℏk_R/mΩ sets the scale for the zitterbewegung's amplitude. The equivalence of this Hamiltonian (5) and the Dirac Hamiltonian (1) provides the opportunity for our quantum simulation of relativistic electron dynamics.

To study zitterbewegung with an ultracold atomic gas, we measured the positions and velocities of atomic systems subject to the Dirac Hamiltonian for varying times after starting in an initial state with speed c*. These experiments began with N ≈ 5 × 10⁴ atom optically trapped ⁸⁷Rb BECs (f_c = 0.75(10) condensate fraction) in the |f = 1,m_F  = −1〉 ground state, subject to a uniform B₀ = 2.1 mT bias magnetic field. The atoms were confined in a harmonic trap [(ω_x,ω_y,ω_z)/2π = (38,38,130) Hz] with characteristic timescales greatly exceeding those of the zitterbewegung. We transferred these atoms (at rest) to $\left |f = 1, m_F = 0\right \rangle$ using an adiabatic rapid passage technique; a fixed frequency 15.0 MHz radiofrequency magnetic field coupled the different m_F states together as the bias magnetic field was swept through resonance. Using a pair of Raman beams counterpropagating along e_x with wavelength λ = 790.1 nm and frequency difference δω = g_F μ_BB₀ + 4E_R +  (where  = h × 32 kHz is the quadratic Zeeman shift), a 30 μs π-pulse transferred approximately 85% of the atoms from |m_F  = 0,k_x = 0〉 to |m_F  = −1,k_x = 2k_R〉 (moving with velocity v = 2ℏk_R/m = c*). Before the trap appreciably altered their velocity (200 μs), we changed the Raman lasers' frequency difference to δω = g_F μ_BB₀, bringing |m_F  = −1,k_x = 2k_R〉 and |m_F  = + 1,k_x = −2k_R〉 into four-photon resonance. We then suddenly introduced a four-photon Raman coupling between these states (figure 1(b)), and allowed the system to evolve under this new Hamiltonian for an evolution time t.

Just before transferring the BEC into $\left |m_F=-1,k_x=2k_{\mathrm { R}}\right >$, two 6.8 GHz microwave pulses spaced in time by 50 ms each out-coupled ≈10% of the atoms to the f = 2 hyperfine manifold. These atoms were separately imaged (without repumping on the f = 1–2 transition) leaving the atoms in f = 1 undisturbed. These f = 2 atoms served two purposes: (i) by setting the microwave frequency 2 kHz above (first pulse) and 2 kHz below (second pulse) resonance, we tracked shifts in the bias field that would change our four-photon Raman resonance condition. Upon analysing the data, we rejected points where the atom number difference between these two images was greater than two standard deviations from being equal; (ii) we determined the BEC's position immediately before each zitterbewegung experiment began, allowing us to cancel shot-to-shot variations in the trap position. The beginning of the three transfer pulses—two microwave outcoupling pulses, and the final four-photon Raman pulse—were each separated in time by 50 ms. As three periods of a 60 Hz cycle, this separation was chosen to reduce magnetic field background fluctuations at the power line frequency, and to facilitate rethermalization between pulses.

We measured the system either by imaging the atoms immediately following this evolution (to determine the atoms' position) or by releasing the atoms from their trap and simultaneously turning off the Raman lasers, allowing for a short time-of-flight (TOF, with duration t_TOF) before imaging (to determine the atoms' velocity). Figure 2 shows the evolution of the signal for several times of flight, and figure 3 shows in situ and after TOF (t_TOF = 550 μs) measurements at several coupling strengths; the velocity-dominated TOF images clearly show the expected cosinusoidal behaviour. For in situ measurements, the Raman and trapping beams remained on during the 40 μs absorption imaging pulses. For TOF, these were removed during TOF during which time the atoms flew ballistically for t_TOF and were subsequently absorption imaged. We used high intensity imaging, with intensity I ≈ 3I_sat (where I_sat is the saturation intensity), that reduced the effective optical depth [26] and gave better signal-to-noise in the determination of the clouds' positions.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

This simple description of zitterbewegung assumes that the range δp_x of occupied momentum states is small compared to m*c*, and only those states near the avoided-crossing structure are populated. To maintain a sufficiently narrow δp_x, the spatial size of the system must be at least δx ≫ λ*_C, which, as observed in [27], is larger than the λ*_C/4π amplitude of the zitterbewegung itself. We satisfied this requirement in our experiment by using clouds whose Thomas–Fermi radii R_x = 12(2) μm greatly exceeded the measured sub-micron zitterbewegung oscillations, and overcame the fundamental measurement challenge with good statistics. Just before initializing zitterbewegung, we measured the initial position of the BEC by out-coupling and imaging ≈ 5 × 10³ atoms. In principle, this allowed us to measure the centre of the distribution with an uncertainty estimated by $R_x / \sqrt {5\times 10^3} \approx 0.17\,\mu {\mathrm { m}}$. Our actual measurements, which include technical noise and are averages of four independent images, have a typical 0.3 μm rms uncertainty (much less than both the distribution's 12(2) μm width and our ≈ 1.75 μm imaging resolution).

From fits to data as in figure 3—with parameters joint between each in situ and TOF pair—we extracted the frequency Ω, amplitude λ*_C/4π and velocity c* of the observed zitterbewegung (shown in figure 4). The observed values are attenuated by approximately 2.5 from those predicted by (3), as explained below.

Zoom In Zoom Out Reset image size

3.1. Amplitude attenuation

For our finite-temperature system, higher momentum states are thermally occupied in the initial equilibrium system (schematically illustrated in figure 1(c) and observed in figure 5). The zitterbewegung frequency for these states is increased, and the oscillations correspondingly dephase, decreasing the observed amplitude. These finite temperature effects give rise to a non-participating fraction f_k of the atomic population drifting at c*. Indeed, figure 5 shows that the majority of the 'thermal' population surrounding the initial BEC is unaffected by the coupling. Additionally, owing to imperfect preparation of the initial | ↑,p_x = 0〉 state, a fraction f₀ = 0.15(10) remained at rest in m_F  = 0 (and therefore did not participate at all in the −1 to +1 coupling). Fluctuations in the background magnetic field also contribute to both f₀ and f_k.

Zoom In Zoom Out Reset image size

The rest fraction f₀ was determined from long TOF images (such as those in figure 5). The drifting fraction f_k was found by fitting a model

$$\begin{eqnarray} &&\fl \bar x (t,t_{\rm TOF}) = (1-f_0)\left\{(1-f_k)c^* \left[ \frac{\sin(\Omega t+\phi_0)}{\Omega} + t_{\rm TOF} \cos(\Omega t+\phi_0)\right]+f_k c^*t + \bar x_0\right\}, \end{eqnarray} \tag{ 6 }$$

where $\bar x_0$ is the initial offset position and f₀ is fixed at 0.15. Using two sets of data, one in situ and one with t_TOF = 550 μs, we performed joint fits for each laser intensity (four-photon coupling strength). In the initial analysis, we fix c* = 2ℏk_R/m and fit the data to extract the parameters f_k, Ω, ϕ₀ and $\bar x_0$. The non-participating fractions for the data shown in figures 4(a) and (b) are shown in figure 4(c). Next, we found the average of f_k as a function of Ω to use in the model. Finally, we remove the background slope due to the f_kc*t term from the same five sets of data using a simple linear fit. We refit the remaining signal to a 'fully participating' model ((6) with f₀ = f_k = 0) with fixed Ω (from the original fit) to extract the effective speed of light parameter c*. We found the zitterbewegung amplitude (of the participating atoms) c*/Ω. The model used to predict the amplitude is given by (1 − f₀)(1 − f_k)2ℏk_R/mΩ, and the uncertainty is dominated by our systematic uncertainty in f₀, which is due to magnetic field variations.

The background slope due to the f_kc*t term was subtracted from the data presented in figures 2 and 3. The curves are calculated from the values from the original fit using (6) without the f_kc*t term. The average of the extracted f_k was used in the theory curve in figure 4(c) to show the expected amplitude of in situ oscillations.

By engineering a two-level quantum system from initial states with opposite velocity, we reinterpret the 'curious' physics of zitterbewegung in analogy to the Rabi oscillations ubiquitous in atomic physics. In this language, the particles trembled because the initial state was not an eigenstate of the coupled system; once subject to the Dirac Hamiltonian, the system Rabi-oscillated between bare states of equal and opposite velocity. As the atoms' coupling was provided by resonant laser light instead of the electrons' rest energy, it is natural to think of a Rabi oscillation picture where the mass (coupling) is suddenly turned on and off. (Somewhat amusingly, the mechanism by which our laser field (a coherent state of light) generates mass is analogous to the Higgs mechanism where a Higgs condensate (a coherent matter wave) generates mass in the standard model [28].) The zitterbewegung of electrons arises because two states—particle and antiparticle states—are coupled, and the resulting eigenstates are superpositions of the two. Projections of bare electron states onto this basis result, as in the case of the atoms, in oscillations between states of opposite velocity. This straightforward analogy compels us to accept that the rest energy acts exactly as a coupling field and mixes the particle and antiparticle states into eigenstates that are superpositions of the two.

While the Dirac equation generally applies only to fermionic systems in nature, quantum simulations such as ours directly realize Dirac–boson systems in the laboratory [12, 29], permitting access to new classes of experimental systems. Though BECs near these Dirac points are short-lived [22, 30–33], strong interactions, as are present near the superfluid–Mott transition in an optical lattice, can stably populate these states [34, 35], for example leading to bosonic composite-fermion states [36, 37].

We thank V Galitski, G Juzeliūnas and R F O'Connell for useful conversations. This work was partially supported by the Office of Naval Research; by the Army Research Office with funds both from the Defense Advanced Research Projects Agency Optical Lattice Emulator program and the Atomtronics Multidisciplinary University Research Initiative; and by the NSF through the Physics Frontier Center at JQI. LJL acknowledges the Natural Sciences and Engineering Research Council of Canada, KJ-G acknowledges CONACYT (Consejo Nacional de Ciencia y Tecnología), and MCB acknowledges NIST-American Recovery and Reinvestment Act.