Thank you for visiting nature.com. You are using a browser version with limited support for CSS. To obtain the best experience, we recommend you use a more up to date browser (or turn off compatibility mode in Internet Explorer). In the meantime, to ensure continued support, we are displaying the site without styles and JavaScript.

# Evidence for an atomic chiral superfluid with topological excitations

## Abstract

Topological superfluidity is an important concept in electronic materials as well as ultracold atomic gases1. However, although progress has been made by hybridizing superconductors with topological substrates, the search for a material—natural or artificial—that intrinsically exhibits topological superfluidity has been ongoing since the discovery of the superfluid 3He-A phase2. Here we report evidence for a globally chiral atomic superfluid, induced by interaction-driven time-reversal symmetry breaking in the second Bloch band of an optical lattice with hexagonal boron nitride geometry. This realizes a long-lived Bose–Einstein condensate of 87Rb atoms beyond present limits to orbitally featureless scenarios in the lowest Bloch band. Time-of-flight and band mapping measurements reveal that the local phases and orbital rotations of atoms are spontaneously ordered into a vortex array, showing evidence of the emergence of global angular momentum across the entire lattice. A phenomenological effective model is used to capture the dynamics of Bogoliubov quasi-particle excitations above the ground state, which are shown to exhibit a topological band structure. The observed bosonic phase is expected to exhibit phenomena that are conceptually distinct from, but related to, the quantum anomalous Hall effect3,4,5,6,7 in electronic condensed matter.

## Main

Quantum simulation involves the use of a comparably simple and precisely controlled synthetic quantum system to mimic poorly understood, isolated phenomena of a far more complex but less-controlled quantum system, while excluding the superimposed secondary structure that would impede a clear understanding8. Optical lattices—that is, ultracold neutral atoms trapped in laser-induced periodic potentials9,10,11—are well established systems that simulate an elementary class of many-body lattice model, the conventional s-band Hubbard model12,13. However, a large part of the physics that is relevant in electronic many-body systems is related to the coupling of electrons to magnetic fields via the Lorentz force, and hence remains inaccessible to quantum simulations in conventional optical lattices. Extensive research has been undertaken to establish a similar mechanism for neutral atoms in optical lattices, giving rise to so-called artificial magnetic fields or gauge fields14,15 using dynamical techniques—similar to what is sometimes summarized as Floquet engineering for condensed-matter physics16. This has enabled the single-particle band structures of optical lattices to be endowed with built-in local17 or even global18,19 magnetic flux and the resultant single-particle topological properties; however, this comes at the cost of notable decoherence and markedly reduced lifetimes. An alternative approach towards optical lattice simulators beyond conventional s-band Hubbard physics is the implementation of orbital degrees of freedom by making use of higher Bloch bands20,21,22,23,24. This recent direction has led to interesting new multi-orbital quantum phases, including examples that have local angular momentum. An intriguing, but as-yet unobtained, goal of this approach is to achieve a multi-orbital optical lattice quantum simulator that can, for example, capture bosonic versions of the topological superfluidity of paired electrons25, or of quantum-Hall-related physics26.

Here we take a step in this direction by demonstrating interaction-induced, spontaneous time-reversal symmetry (TRS) breaking and formation of a chiral atomic superfluid in an orbital optical lattice with global angular momentum and topological excitations and edge states. This metastable state exhibits remarkable robustness with a lifetime of hundreds of milliseconds. As an unequivocal signature of the prevalence of global angular momentum, characteristic momentum spectra are observed using time-of-flight spectroscopy. Our experimental observations directly show the spontaneous breaking of TRS in accordance with a mean-field tight-binding model and exact band calculations. The system is expected to exhibit the bosonic counterpart of the quantum anomalous Hall effect3,4,5,6,7; however, not as a result of an engineered band structure, as in the famous Haldane model27. Rather, as shown by theoretical considerations (Supplementary Information), it is the contact interaction between degenerate p orbitals that leads to a spontaneous breaking of TRS, the formation of global angular momentum, a topologically non-trivial excitation band structure and the existence of topological edge states.

Our experiments use a boron nitride optical lattice that is composed of alternating shallow and deep potential wells—denoted A and B, respectively—arranged on a hexagonal lattice. As indicated in Fig. 1a, the lattice potential is formed by three laser beams operating at a wavelength λ of 1,064 nm. The beams propagate within the x–y plane and intersect at an angle of 120°. Each beam comprises two spectral components at frequencies ω1 and ω2, with linear polarizations in the z-direction. Each set of spectral components with the same frequency forms a triangular lattice with adjustable amplitude. By tuning the frequency difference ω1 − ω2, the intensity patterns of both triangular lattices are shifted with respect to each other such that their sum produces a boron nitride lattice geometry with A and B wells. The relative potential offset of the two classes of well, denoted ΔV, can be rapidly tuned. An optical dipole trap provides weak harmonic confinement with regard to the z-direction. For details, see Methods and Supplementary information.

Initially, a Bose–Einstein condensate of 87Rb atoms was loaded into the lowest Bloch band in the centre of the first Brillouin zone, denoted Γ in Fig. 1b. The overall trap depth along the vertical z-axis—including the lattice potential, the dipole trap and gravitation—is 221 nK. As a key feature of our preparation protocol, we applied additional cooling by means of evaporation after the atoms were loaded into the lattice. The dipole trap depth decreases within 15 ms, such that the overall trap depth along the z-direction is reduced to 41 nK and the energetic atoms can escape. At this stage the atomic wavefunction is composed of s orbitals residing in the deep B wells of the lattice. In a subsequent rapid quench the atoms are excited to the second band (compare with Fig. 1b) while remaining at the Γ point, which constitutes a dynamically unstable energy maximum28. This is accomplished by tuning ΔV in approximately 100 μs according to the diagram in Fig. 2a, until the local s orbitals in the B wells come to lie between the s and p orbitals in the A wells and hence belong to the second band. This also leads to a further reduction, to 24 nK, of the overall trap depth with respect to the z-direction, and thus an additional boost of evaporative cooling. Typical momentum spectra (Fig. 2c) recorded shortly (0.5 ms) after the quench show that a large number of Bloch states of the second band are populated with no apparent coherence28. Related excitation protocols without additional evaporation have been used previously with bipartite square lattices22,23 and with hexagonal lattices29,30. For details, see Methods and Supplementary Information.

The second band gives rise to a symmetric double-well scenario in quasi-momentum space. It possesses two inequivalent global energy minima located at the corners of the first Brillouin zone—denoted KΔ and K in Fig. 1b—which exhibit perfect degeneracy, protected by TRS. A central observation of this work is that, in a re-condensation process, TRS is spontaneously broken and a condensate forms in either of the two K points with equal probability. This is shown in Fig. 2c, which displays momentum spectra recorded after variable holding times of between 0.5 ms and 755 ms. As time proceeds, sharp Bragg resonances grow at the KΔ and K points. These resonances indicate the build-up of long-range coherence, and hence the formation of condensate fractions at KΔ and K. Notably, these condensate fractions—according to the observations—are generally not equally sized, as is shown in Fig. 2c, in which the three rows show examples with either a dominant KΔ component (top), a dominant K component (bottom) or both components of similar size (middle).

We next point out that the observation of a dominant condensate fraction in either of the K points constitutes clear evidence of broken TRS. Exact band theory shows that the momentum spectra observed in the top and bottom rows of Fig. 2c are very well reproduced by the calculated momentum spectra of the Bloch functions associated with KΔ and K, respectively. For example, the calculated momentum spectrum of the Bloch function at KΔ in Fig. 2b is plotted for comparison in the top row of Fig. 2c. Hence, the observation of population in a single K point unequivocally indicates the presence of a wavefunction that is well approximated by the Bloch function associated with that K point. Note that these Bloch functions, $${\Phi }_{{{\rm{K}}}_{\Delta }}$$ and $${\Phi }_{{{\rm{K}}}_{\nabla }}$$, inherently break TRS, as stated in Fig. 1c. The first panel schematically illustrates the composition of $${\Phi }_{{{\rm{K}}}_{\Delta }}$$ in terms of local s orbitals in the shallow wells and px and py orbitals in the deep wells, which form p+ = px + ipy hybrid orbitals (p in case of K). The three distinct colours indicate local phases that differ by 2π/3, which results from an inherent vortical phase texture of $${\Phi }_{{{\rm{K}}}_{\Delta }}$$. In the second and third panels of Fig. 1c, $$|{\Phi }_{{K}_{\Delta }}|$$, which is derived from exact band theory, is plotted for two different settings of ΔV, showing that different relative populations at A and B wells can be adjusted. The fourth panel shows the vortical mass current associated with $${\Phi }_{{{\rm{K}}}_{\Delta }}$$ for the same choice of ΔV as in the third panel. Note that the vortices residing at each deep well across the entire lattice all share the same sense of rotation—an interesting property that is associated with the threefold rotation symmetry of the hexagonal lattice in connection with the twofold degeneracy of the p orbital manifold in two dimensions.

We define an experimental measure of chirality χ = (nΔ − n)/(nΔ + n) as the imbalance between $${\Phi }_{{{\rm{K}}}_{\Delta }}$$ and $${\Phi }_{{{\rm{K}}}_{\nabla }}$$ condensate populations, nΔ and n, respectively. This quantity is observed to be random in each experimental implementation with nearly zero average over many shots. However, although for short holding times—that is, in the initial phase of the condensation process—a distribution of χ with a maximum at zero is found, for long times beyond 100 ms this distribution develops pronounced maxima around ±0.5. This is shown in Fig. 3: in Fig. 3a, the average of the modulus of the chirality $$\langle |\chi |\rangle$$ is plotted against the holding time, whereas in Fig. 3b histograms of the distribution of χ are shown for holding times ranging from 35 ms to 505 ms. As is seen in Fig. 3a, $$\langle |\chi |\rangle$$ is small for short holding times of a few ms, and continuously increases until it reaches a maximum at around 0.4 for a holding time of about 265 ms, where the distribution of χ (compare with Fig. 3b) shows a pronounced bimodal structure. According to these observations, the condensation process develops in two stages: the first stage spontaneously breaks U(1) symmetry, forming a condensate that is approximately described by an equal superposition of Bloch states $${\Phi }_{{{\rm{K}}}_{\Delta }}$$ and $${\Phi }_{{{\rm{K}}}_{\nabla }}$$, thus globally preserving TRS and zero net angular momentum. Alternatively, as in the one-dimensional double-well scenario of ref. 31, domains could be formed, such that locally TRS is broken but the net angular momentum remains at zero. However, the narrow width of the observed Bragg peaks—seen in the momentum spectra with equally occupied K points in Fig. 2c—indicate well established coherence over large parts of the lattice, which is not compatible with a fine-grain domain structure. The second condensation stage breaks the global TRS by spontaneously bifurcating to form either of the two states $${\Phi }_{{{\rm{K}}}_{\Delta }}$$ or $${\Phi }_{{{\rm{K}}}_{\nabla }}$$, with equal probability.

Theoretical considerations (detailed in the Supplementary Information) suggest that the second stage of the condensation process is triggered by the ferromagnetic collisional interactions within the manifold of degenerate p orbitals in the deep lattice wells. In the presence of such collisions, the ground state according to mean-field theory is given by a condensate in either of the two K points rather than a superposition of both (see Supplementary Information). The atoms are self-organized into an order of synchronized vortices centred at each of the deep wells across the entire boron nitride lattice, as illustrated in Fig. 1c. The fraction of atoms populating the p+ orbitals give rise to a non-zero global angular momentum. To access the ground state, which requires build-up of global orbital angular momentum, angular momentum must be exchanged with the environment. Indeed, a signature of this process is observed in the experiment. As shown in Fig. 3c, during the second condensation phase, past about 10 ms holding time, we observe a stream of atoms leaving the lattice perpendicular to the lattice plane (x–y plane) along the direction of gravity (−z-direction). This supports the interpretation that kinetic energy and angular momentum are transferred out of the lattice plane into the weakly confined z-direction, and finally escape from the lattice. It is of interest to note that, for stronger collisional interactions, under the assumption of quasi-momentum conservation, a condensate at the M points (Fig. 1b) is predicted to have a lower energy (Supplementary Fig. 5). However, a phase-separated mixed state with equal contributions from both K points might nevertheless be energetically more favourable.

The plots in Fig. 2c are obtained with ΔV adjusted during step 1 in Fig. 2a, such that most of the atoms in the second band reside in local s orbitals in the shallow wells and thus do not contribute notable orbital angular momentum. This corresponds to the wavefunction illustrated in the second panel of Fig. 1c. We now discuss how, by means of a slightly modified preparation protocol, a state with large global orbital angular momentum can be experimentally formed, which possesses a wavefunction according to the third panel of Fig. 1c. By adding an additional step (step 2 in Fig. 4a) 205 ms after the preparation protocol of step 1 in Fig. 2a, one may relocate a substantial portion of the atomic population towards the vortical p orbitals in the deep wells. This step consists of adiabatically tuning ΔV (in 1 ms) until the local s orbitals in the shallow wells lie energetically above the p orbitals of the deep wells (see Methods and Supplementary Information). The absolute value of the corresponding wavefunction is shown in the third panel of Fig. 1c. In Fig. 4c we observe the resulting faster band relaxation that notably depletes the second band after only 10 ms. This observation indicates that the collision dynamics is now dominated by atoms in the p orbitals of the deep wells, which provide faster decay due to increased overlap with the s orbitals of the lowest band23. As a direct signature of dominant population of the p orbitals, their additional nodes lead to destructive interference, thus preventing the occurrence of Bragg resonances within the areas enclosed by red dashed circles in Fig. 4c. This observation follows the exact band calculations in Fig. 4b. Note the contrast with the analogous spectra in Fig. 2b and with the corresponding experimental spectra in the top row of Fig. 2c. A large global orbital angular momentum arises, to which each atom in the p+ orbitals contributes a portion ħ. A quantitative analysis of the dependence of angular momentum on ΔV is shown in Supplementary Fig. 4.

Symmetric double-well scenarios in quasi-momentum space that lead to an interaction-induced spontaneous breaking of TRS and the formation of phases with staggered chiral order have been previously reported22,31,32. The long-lived atomic orbital superfluid devised in this work exhibits global angular momentum of the ground state wavefunction, although the Hamiltonian preserves TRS. As a consequence, basic Bogoliubov-de Gennes analysis verifies the emergence of bosonic topological excitations and edge modes (Supplementary Information). This opens up the prospect of studying these signatures experimentally. A direct comparison between experiment and theory would be enabled with regards to fundamental concepts that are related to, but have no previous analogue in, condensed-matter electronic and spin materials. Previous calculations have shown that, in a checkerboard square lattice with alternating shallow and deep wells and local non-zero angular momentum order, interaction-induced gaps should open in Bogoliubov quasiparticle spectra, and edge states protected by topological invariants should occur33,34. The proposed Zeeman-like bias is, however, an experimental challenge to implement for orbital degrees of freedom. Here, the boron nitride lattice resolves this challenge by spontaneously breaking TRS globally (driven by intrinsic interaction) with an emerging global angular momentum order. This paves the way to study dynamically controlled quasiparticles as exotic as a possible counterpart of Majorana fermions in a bosonic superfluid. The latter is widely sought-after in electronic topological superconductors, and the Supplementary Information presents calculations that support this expectation.

## Methods

### Realization of optical lattice

A two-dimensional hexagonal boron nitride (BN) optical lattice potential is created by three laser beams, propagating in the xy plane and intersecting at 120° angles. Each laser beam comprises two frequency components ω1 and ω2—which are both linearly polarized along the z direction—with wavelength λ ≈ 1,064 nm and therefore with negative detuning with respect to the relevant atomic transitions of rubidium atoms at 780 nm and 795 nm. The two frequency components are derived from two independent lasers. Experimentally, we first combine the two laser beams at frequencies ω1 and ω2, respectively, before splitting them into three beams. The total electric field for all beams is then written as

$$\begin{array}{l}{\bf{E}}({\bf{r}},t)={E}_{1}{{\bf{e}}}_{z}\sum _{j}\cos ({{\bf{k}}}_{j}\cdot {\bf{r}}-{\omega }_{1}t+{\theta }_{j})\\ +{E}_{2}{{\bf{e}}}_{z}\sum _{j}\cos ({{\bf{k}}}_{j}^{{\prime} }\cdot {\bf{r}}-{\omega }_{2}t+{\theta }_{j}^{{\prime} }).\end{array}$$

Here, the wave vectors are given by $${{\bf{k}}}_{1}={k}_{{\rm{L}}}(\,-\sqrt{3}/2,1/2)$$, $${{\bf{k}}}_{2}={k}_{{\rm{L}}}(\,\sqrt{3}/2,1/2)$$, $${{\bf{k}}}_{3}={k}_{{\rm{L}}}(0,-1)$$, and $${{\bf{k}}}_{j}^{{\prime} }=({\omega }_{2}/{\omega }_{1}){{\bf{k}}}_{j}$$, where $${k}_{{\rm{L}}}={\omega }_{1}/c$$. $${\theta }_{j}={\omega }_{1}{L}_{j}/c+{\theta }_{0}$$ and $${\theta }_{j}^{{\prime} }={\omega }_{2}{L}_{j}/c+{\theta }_{0}^{{\prime} }$$, where Lj denotes the optical path length of the jth beam from the splitting point to the centre of the lattice. The corresponding laser intensity I(r) is proportional to the time averaging of the square of the electric field according to

$$\begin{array}{ll}I({\bf{r}})\,\propto \, & \,\frac{3}{2}{E}_{1}^{2}+{E}_{1}^{2}\sum _{\langle i,j\rangle }\cos [({{\bf{k}}}_{i}-{{\bf{k}}}_{j})\cdot {\bf{r}}+{\theta }_{i}-{\theta }_{j}]\\ & +\frac{3}{2}{E}_{2}^{2}+{E}_{2}^{2}\sum _{\langle i,j\rangle }\cos [({{\bf{k}}}_{i}^{{\prime} }-{{\bf{k}}}_{j}^{{\prime} })\cdot {\bf{r}}+{\theta }_{i}^{{\prime} }-{\theta }_{j}^{{\prime} }],\end{array}$$

where the summation $$\langle i,j\rangle$$ is limited to $$\langle 1,2\rangle$$, $$\langle 2,3\rangle$$, $$\langle 3,1\rangle$$. The generated optical lattice potential is proportional to the laser intensity and takes the form

$$\begin{array}{ll}{V}_{{\rm{BN}}}({\bf{r}})\,=\, & -{V}_{1}\{3+2\sum _{\langle i,j\rangle }\cos [({{\bf{k}}}_{i}-{{\bf{k}}}_{j})\cdot {\bf{r}}+({\theta }_{i}-{\theta }_{j})]\}\\ & -{V}_{2}\{3+2\sum _{\langle i,j\rangle }\cos [({{\bf{k}}}_{i}^{{\prime} }-{{\bf{k}}}_{j}^{{\prime} })\cdot {\bf{r}}+({\theta }_{i}^{{\prime} }-{\theta }_{j}^{{\prime} })]\},\end{array}$$

where V1,2 ≥ 0 for the relevant case of red detuning. Each of the two spectral components creates a triangular lattice potential, which sum up to form the total potential VBN. No interference terms arise because the frequency difference Δωω1 – ω2 is chosen in the range of a few GHz, which exceeds by far all relevant timescales of the atom dynamics. The relative position of the two triangular lattices is determined by the phase differences

$$\Delta {\theta }_{ij}=({\theta }_{i}-{\theta }_{j})-({\theta }_{i}^{{\prime} }-{\theta }_{j}^{{\prime} })=\frac{{\omega }_{1}-{\omega }_{2}}{c}({L}_{i}-{L}_{j}).$$

Note that with Δω ≈ 2π × 3 GHz, a change of the lengths Li of the order of 10 μm corresponds to irrelevant changes of Δθij of the order of 10−4 × 2π. Hence, Δθij can be readily adjusted without the need for interferometric control of the lengths Li. Experimentally, we choose convenient values for the lengths L1, L2 and L3 and lock the frequency difference between the two lasers accordingly to fine-tune the relative position of the two triangular lattices appropriately to generate the desired boron nitride optical lattice. We set $$({L}_{1}-{L}_{2},{L}_{2}-{L}_{3},{L}_{3}-{L}_{1})=(\,-6.04,3.02,3.02)$$cm and $$\Delta \omega ={\omega }_{1}-{\omega }_{2}=2{\rm{\pi }}\times 3.308$$ GHz, which leads to (Δθ12, Δθ23, Δθ31) = $$(\,-4{\rm{\pi }}/3,2{\rm{\pi }}/3,2{\rm{\pi }}/3)$$ and hence the boron nitride lattice potential

$$\begin{array}{l}{V}_{{\rm{BN}}}({\bf{r}})\,=\,-\,{V}_{1}\left\{3+2\sum _{\langle i,j\rangle }\cos \left[({{\bf{k}}}_{i}-{{\bf{k}}}_{j})\cdot {\bf{r}}-\frac{2{\rm{\pi }}}{3}\right]\right\}\\ \,-\,{V}_{2}\left\{3+2\sum _{\langle i,j\rangle }\cos \left[({{\bf{k}}}_{i}-{{\bf{k}}}_{j})\cdot {\bf{r}}+\frac{2{\rm{\pi }}}{3}\right]\right\}.\end{array}$$

The first (second) term describes the triangular lattice potential that gives rise to the A sites (B sites) in the combined boron nitride lattice shown in Fig. 1a. The potential difference between A wells and B wells can be readily adjusted on the microsecond timescale by tuning the ratio V1/V2. To avoid deformations of the lattice potential, the relative frequency difference between the two sets of triangular lattices and the laser intensities are carefully stabilized.

A Bose–Einstein condensate of typical 40,000 87Rb atoms in the state $$|F=1,{m}_{F}=-1\rangle$$ (where F is the quantum number of total angular momentum and mF = magnetic quantum number) is prepared in an optical dipole trap formed by two crossed laser beams with trapping frequencies of $$\{{\omega }_{x},{\omega }_{y},{\omega }_{z}\}=2{\rm{\pi }}\times \{26,27,71\}$$ Hz. Including the gravitational force (pointing into the −z direction) the trap depth along the −z direction is 34 nK. A bias magnetic field of 1 G is applied along the z axis. Within 120 ms, the lattice beam intensity is ramped up to V1 = 7.04ER and V2 = 8.03ER, where $${E}_{{\rm{R}}}={h}^{2}/2m{\lambda }^{2}$$. At this stage the overall trap depth along the −z direction is 221 nK. To provide additional evaporative cooling, after 5 ms, the depth of the optical dipole trap is ramped down in 15 ms, such that the overall trap depth along the −z direction is reduced to 41 nK. Excitation into the second band is obtained by swapping the depths of the A and B sites via linearly changing (V1,V2) to (7.81,7.23)ER rapidly in 0.1 ms. The trap depth along the −z direction is thereby further reduced to 24 nK, which gives rise to further evaporation.

We recorded momentum spectra via time-of-flight spectroscopy or performed band mapping that enabled us to observe the quasi-momentum distribution. These techniques were used to derive the data in Fig. 2c. To obtain the data in Fig. 4c, the experimental protocol was slightly extended. After excitation to the second band and a subsequent holding time of 205 ms, the p orbitals in the A wells are lowered by continuously increasing V1 to 8.35ER in 1 ms. Therefore, atoms are transferred from the s orbitals in the shallow B wells to the p orbitals in adjacent deep A wells. Momentum spectra of the atoms in the xy plane were obtained by switching off all potentials in less than 1 μs and, after a 20-ms-long ballistic expansion, performing absorption imaging. For band mapping, we decreased the intensity of the lattice exponentially with a time constant of 260 μs followed by 20 ms of ballistic expansion before performing absorption imaging. For the atom loss data shown in Fig. 3c, we performed in-situ absorption imaging of a plane perpendicular to the x–y plane, after the atoms were excited to the second band and held there for a variable time.

## Data availability

All data presented in the figures are available upon request from the corresponding authors.

## Code availability

The numerical simulations were performed with MATLAB R2019b. Codes are available upon request from the corresponding authors.

## References

1. 1.

Sato, M. & Ando, Y. Topological superconductors: a review. Rep. Prog. Phys. 80, 076501 (2017).

2. 2.

Wheatley, J. C. Experimental properties of superfluid 3He. Rev. Mod. Phys. 47, 415–470 (1975).

3. 3.

Yu, R. et al. Quantized anomalous Hall effect in magnetic topological insulators. Science 329, 61–64 (2010).

4. 4.

Chang, C.-Z. et al. Experimental observation of the quantum anomalous Hall effect in a magnetic topological insulator. Science 340, 167–170 (2013).

5. 5.

Liu, C.-X., Zhang, S.-C. & Qi, X.-L. The quantum anomalous hall effect: theory and experiment. Annu. Rev. Condens. Matter Phys. 7, 301–321 (2016).

6. 6.

Deng, Y. et al. Quantum anomalous Hall effect in intrinsic magnetic topological insulator MnBi2Te4. Science 367, 895–900 (2020).

7. 7.

Serlin, M. et al. Intrinsic quantized anomalous Hall effect in a moiré heterostructure. Science 367, 900–903 (2020).

8. 8.

Feynman, R. P. Simulating physics with computers. Int. J. Theor. Phys. 21, 467–488 (1982).

9. 9.

Verkerk, P. et al. Dynamics and spatial order of cold cesium atoms in a periodic optical potential. Phys. Rev. Lett. 68, 3861–3864 (1992).

10. 10.

Hemmerich, A. & Hänsch, T. W. Two-dimesional atomic crystal bound by light. Phys. Rev. Lett. 70, 410–413 (1993).

11. 11.

Grynberg, G. & Robilliard, C. Cold atoms in dissipative optical lattices. Phys. Rep. 355, 335–451 (2001).

12. 12.

Bloch, I. Ultracold quantum gases in optical lattices. Nat. Phys. 1, 23–30 (2005).

13. 13.

Lewenstein, M. et al. Ultracold atomic gases in optical lattices: mimicking condensed matter physics and beyond. Adv. Phys. 56, 243–379 (2007).

14. 14.

Dalibard, J., Gerbier, F., Juzeliūnas, G. & Öhberg, P. Colloquium: artificial gauge potentials for neutral atoms. Rev. Mod. Phys. 83, 1523–1543 (2011).

15. 15.

Galitski, V., Juzeliūnas, G. & Spielman, I. B. Artificial gauge fields with ultracold atoms. Phys. Today 72, 38–44 (2019).

16. 16.

Oka, T. & Kitamura, S. Floquet engineering of quantum materials. Annu. Rev. Condens. Matter Phys. 10, 387–408 (2019).

17. 17.

Jotzu, G. et al. Experimental realization of the topological Haldane model with ultracold fermions. Nature 515, 237–240 (2014).

18. 18.

Aidelsburger, M. et al. Realization of the Hofstadter Hamiltonian with ultracold atoms in optical lattices. Phys. Rev. Lett. 111, 185301 (2013).

19. 19.

Miyake, H., Siviloglou, G. A., Kennedy, C. J., Burton, W. C. & Ketterle, W. Realizing the Harper Hamiltonian with laser-assisted tunneling in optical lattices. Phys. Rev. Lett. 111, 185302 (2013).

20. 20.

Liu, W. V. & Wu, C. Atomic matter of nonzero-momentum Bose-Einstein condensation and orbital current order. Phys. Rev. A 74, 013607 (2006).

21. 21.

Lewenstein, M. & Liu, W. V. Orbital dance. Nat. Phys. 7, 101–103 (2011).

22. 22.

Wirth, G., Ölschläger, M. & Hemmerich, A. Evidence for orbital superfluidity in the P-band of a bipartite optical square lattice. Nat. Phys. 7, 147–153 (2011).

23. 23.

Kock, T., Hippler, C., Ewerbeck, A. & Hemmerich, A. Orbital optical lattices with bosons. J. Phys. At. Mol. Opt. Phys. 49, 042001 (2016).

24. 24.

Li, X. & Liu, W. V. Physics of higher orbital bands in optical lattices: a review. Rep. Prog. Phys. 79, 116401 (2016).

25. 25.

Qi, X. L. & Zhang, S. C. Topological insulators and superconductors. Rev. Mod. Phys. 83, 1057 (2011).

26. 26.

von Klitzing, K. et al. 40 years of the quantum Hall effect. Nat. Rev. Phys. 2, 397–401 (2020).

27. 27.

Haldane, F. D. M. Model for a quantum Hall effect without Landau levels: condensed-matter realization of the “parity anomaly”. Phys. Rev. Lett. 61, 2015–2018 (1988).

28. 28.

Sharma, V., Choudhury, S. & Mueller, E. J. Dynamics of Bose-Einstein recondensation in higher bands. Phys. Rev. A 101, 033609 (2020).

29. 29.

Weinberg, M., Staarmann, C., Ölschläger, C., Simonet, J. & Sengstock, K. Breaking inversion symmetry in a state-dependent honeycomb lattice: artificial graphene with tunable band gap. 2D Mater. 3, 024005 (2016).

30. 30.

Jin, S. et al. Evidence of Potts-Nematic superfluidity in a hexagonal sp2 optical lattice. Phys. Rev. Lett. 126, 035301 (2021).

31. 31.

Parker, C. V., Ha, L.-C. & Chin, C. Direct observation of effective ferromagnetic domains of cold atoms in a shaken optical lattice. Nat. Phys. 9, 769–774 (2013).

32. 32.

Struck, J. et al. Quantum simulation of frustrated classical magnetism in triangular optical lattices. Science 333, 996–999 (2011).

33. 33.

Xu, Z.-F., You, L., Hemmerich, A. & Liu, W. V. π-Flux dirac bosons and topological edge excitations in a bosonic chiral p-wave superfluid. Phys. Rev. Lett. 117, 085301 (2016).

34. 34.

Di Liberto, M., Hemmerich, A. & Morais Smith, C. Topological varma superfluid in optical lattices. Phys. Rev. Lett. 117, 163001 (2016).

## Acknowledgements

This work is supported by the National Key R&D Program of China (grant no. 2018YFA0307200), by the Key-Area Research and Development Program of Guangdong Province (grant no. 2019B030330001), by the NSFC (grant no. U1801661), by the high-level special funds from SUSTech (grant no. G02206401), by a fund from Guangdong province (grant no. 2019ZT08X324 (X.-Q.W., G.-Q.L., J.-Y.L. and Z.-F.X.)), by AFOSR grant no. FA9550-16-1-0006, by MURI-ARO grant no. W911NF-17-1-0323 through UC Santa Barbara, Shanghai Municipal Science and Technology Major Project (grant no. 2019SHZDZX01) (W.V.L.) and by DFG-He2334/17-1 (A.H.).

## Funding

Open access funding provided by Universität Hamburg.

## Author information

Authors

### Contributions

X.-Q.W. and J.-Y.L. carried out the experiments. G.-Q.L., Z.-F.X. and W.V.L. performed the calculations. Z.-F.X. devised the experiment with input by A.H. The lattice design was contributed by A.H. The manuscript was prepared by Z.-F.X, W.V.L. and A.H. All authors discussed the research.

### Corresponding authors

Correspondence to W. Vincent Liu, Andreas Hemmerich or Zhi-Fang Xu.

## Ethics declarations

### Competing interests

The authors declare no competing interests.

Peer review information Nature thanks the anonymous reviewers for their contribution to the peer review of this work. Peer reviewer reports are available.

Publisher’s note Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

## Supplementary information

### Supplementary Information

This file contains supplementary text, supplementary equations s1 – s31, supplementary figures s1 – s7 and supplementary references.

## Rights and permissions

Reprints and Permissions

Wang, XQ., Luo, GQ., Liu, JY. et al. Evidence for an atomic chiral superfluid with topological excitations. Nature 596, 227–231 (2021). https://doi.org/10.1038/s41586-021-03702-0

• Accepted:

• Published:

• Issue Date: