Topological phase transitions in superradiance lattices

Topological phases of matters are of fundamental interest and have promising applications. Fascinating topological properties of light have been unveiled in classical optical materials. However, the manifestation of topological physics in quantum optics has not been discovered. Here we study the topological phases in a two-dimensional momentum-space superradiance lattice composed of timed Dicke states (TDS) in electromagnetically induced transparency (EIT). By periodically modulating the three EIT coupling fields, we can create a Haldane model with in-situ tunable topological properties, which manifest themselves in the contrast between diffraction signals emitted by superradiant TDS. The topological superradiance lattices provide a controllable platform for simulating exotic phenomena in condensed matter physics and offer a basis of topological quantum optics and novel photonic devices.

The discovery of the quantum Hall effect (QHE) [1] reveals a new class of matter phases, topological insulators (TI's), which have been extensively studied in solid-state materials [2][3][4][5][6][7] and recently in photonic structures [8][9][10][11][12], time-periodic systems [5,[13][14][15][16] and optical lattices of cold atoms [18]. All these topological systems are lattices in real space. Our recent study shows that Scully's timed Dicke states (TDS) [19][20][21] can form a superradiance lattice (SL) in momentum space [1]. Here we report the discovery of topological phase transitions in a two-dimensional SL in electromagnetically induced transparency (EIT). By periodically modulating the three EIT coupling fields, we can create a Haldane model [2] with in-situ tunable topological properties. The Chern numbers of the energy bands and hence the topological properties of the SL manifest themselves in the contrast between diffraction signals emitted by superradiant TDS. The topological superradiance lattices (TSL) provide a controllable platform for simulating exotic phenomena in condensed matter physics and offer a basis of topological quantum optics and novel photonic devices.
The Haldane model [2] proposed in 1988 shows that the QHE is an intrinsic topological property of the energy bands. It inspired the discovery of the quantum spin Hall effect [3][4][5] and topological superconductors [7]. The model consists of a honeycomb tight-binding lattice with complex next-nearest-neighbour (NNN) hopping, which breaks the time-reversal symmetry and induces an energy gap between two bands that have opposite Chern numbers. Inversion symmetry breaking by an on-site potential can also open band gaps. The interplay of these two types of symmetry breaking leads to transitions between phases with Chern numbers 0 and ±1. Notwithstanding its foundational role in topological condensed matter physics, the Haldane model has never been realized in solid-state systems. Floquet modulation of a circular polarised light in graphene was proposed to induce the complex NNN tunnelling [13][14][15]. However, the required light is soft X-ray which would directly excite electrons and be absorbed. The Haldane model of cold atoms in optical lattices [18] and its mapping in quantum circuit [23] were recently realised with elegant but complicated experiments. Here we propose a novel quantum optics realisation of the Haldane model using two-dimensional SL's in a simple EIT configuration. Since Dicke's seminal paper in 1954 [24], superradiance has been an important topic in quantum optics. A single photon with wave vector k p can excite a spatially extended N -atom ensemble from the ground state |G = |g 1 , g 2 , . . . , g N to the TDS [19-21] where e i and g i are the excited and ground states of the ith atom at position r i , respectively.  The two Chern numbers corresponding to two topologically distinctive configurations of the three phase angles φ1, φ2 and φ3 on a unit circle.
|m with three coherent plane wave fields, we construct a honeycomb SL of TDS in momentum space [1], as shown in Fig.1. The SL Hamiltonian with rotating-wave approximation is (see Supplementary Information) where k = k p + r(k 2 − k 1 ) + s(k 3 − k 2 ) with integers r and s, Ω l is the Rabi frequency of the coupling field along k l , |m k are defined the same as in Eq.(1) but for |m states, and ∆ c = ω em − ν c is the detuning of the transition frequency ω em between |e and |m from the angular frequency ν c of the EIT coupling fields. The states |e k and |m k−k1 form two sublattices of a honeycomb structure in momentum space [see Fig. 1(c)]. Ω l causes the nearest-neighbour hopping. As we will show below, periodically modulating Ω l can introduce NNN hopping with controllable phases. On the other hand, the detuning ∆ c causes an energy offset between the two sublattices, breaking the inversion symmetry. Thus three Haldane phases with Chern numbers ±1 and 0 can be insitu realised. We first show how to induce complex NNN hopping. We consider periodically modulated Rabi frequencies Ω l (l = 1, 2, 3) of coupling fields, where Ω s and Ω d are the static and dynamic components of the Rabi frequency, ν d is the modulation frequency and φ l is the modulation phase. Here we have assumed the size of the atomic ensemble be much smaller than c/ν d (where c is the speed of light) and hence neglected the position dependence of the modulation. The periodic modulation induces Floquet sidebands [2,5] with energy separation ν d . While Ω s is the intra-sideband nearestneighbour hopping, Ω d induces the inter-sideband hopping. We choose ν d Ω s/d , ∆ c so that the Floquet sidebands are well separated in energy. Thus by second-order perturbation, the effective Hamiltonian of the NNN intrasideband transition mediated by Ω d is (see Methods) where the NNN hopping coefficient Ω lj = 2iΩ sin(φ l −φ j ) with Ω = Ω 2 d /ν d . The fact that the NNN hopping coefficient is purely imaginary is crucial to the topological phases. The loop transitions via NNN hopping accumulate nonzero phases, as shown in Fig.1 (c). H opens a band gap tr = 4 √ 3 |α|Ω where α is a dimensionless quantity determined by the summation of Ω (l+1)l , The Chern numbers of the upper and lower bands C and C are C = −C = sign(α) (see Supplementary Information). In Fig. 2 (a), we plot α with φ 1 = 0 and 0 ≤ φ 2,3 < 2π. The topological property of this SL Haldane model can be represented by the distribution of φ l on a unit circle. There are two distinct topological configurations, counter-clockwise φ 1 , φ 2 and φ 3 for C = −1, and clockwise φ 1 , φ 2 and φ 3 for C = 1, as shown in Fig.2 (b). The time reversal t → −t in Eq.(13) is equivalent to φ l → −φ l , which leads to C → −C. Unlike TI's [5,6], TSL's have no edge states in the semiclassical limit of the coupling fields (see Supplementary Information). Neither do TSL's have Fermi surfaces. Nonetheless, the TSL has its unique topological properties that are observable. The TDS have directional superradiance emission. Of all the TDS in the SL, only those |e k with c|k| − ω eg within the energy bands (ω eg is the transition frequency between |e and |g ) can satisfy both energy and momentum conservation, and have directional emission in k [19]. We call these states superradiant TDS and the other ones subradiant TDS. The topological orders lead to different light emissions from different superradiant TDS. Alternatively, we can also tune the probe field frequency to test the topological band properties at certain energy, which is analogue to tuning the Fermi surface in a fermionic system.
For the sake of simplicity, we set the wavevectors of the three EIT coupling fields to be k 1 = −k cx , k 2 =  k c (x − √ 3ŷ)/2 and k 3 = k c (x + √ 3ŷ)/2, and the probe field wavevector k p = −k 1 . In this case we have only three superradiant TDS with wavevectors k p , k + = −k 3 and k − = −k 2 , and the diffraction fields are along k ± , as shown in Fig.1 (a). We set all fields on resonance, i.e., ∆ c = 0 and ∆ p ≡ ω eg − ν p = 0. The excitation |e kp flows to |e k± and emits photons along k ± . We denote the steady state probability amplitudes of states |e k± as   c k± and define the superradiance contrast For C = 1, the excitation current flows along |e kp → |e k+ → |e k− , as shown in Fig.1 (c). Since each TDS |e k or |m k has decoherence rate γ e or γ m , respectively, the excitation is more probable in state |e k+ than in state |e k− . We therefore have η > 0. Similarly, for C = −1, η < 0. Thus the sign change of the superradiance contrast signatures the topological phase transition, as seen in Fig.3 (a). The superradiance contrast in Fig.3 (a) is consistent with Fig.2 (b) except for the two diagonal corners where the topological currents are weak and the local effect inside a unit cell dominates (see Supplementary Information). The topological phase with C = 0 can be reached by breaking the inversion symmetry. Substantially easier than in graphene electrons and cold atoms in optical lattices [18], a sublattice offset in SL's can be introduced by choosing ∆ c = 0. We set φ l = 4π(l−1)/3 and thus C = 1 when ∆ c = 0. The energy gap in absence of sublattice offset is tr = 18 Ω . Topological phase transitions occur at |∆ c | = tr [2] (see Supplementary Information). In Fig.3 (b), we plot η as a function of ∆ c / tr . For Ω = 0.01, the phase transition is smeared out by the relatively large γ e . The phase transition gets more apparent for larger Ω . At the phase transition point, η has peaks or kinks because the two bands touch at the middle of the band gap, where the probe field probes, as shown in Fig.3 (c).
The TSL's have unique features in transient light propagation under pulse probe. In Fig.4, we compare the pulse propagation in a trivial SL and in a TSL (see Supplementary Information). For a weak probe pulse, the linear susceptibility is χ (1) ∝ c kp and the linear absorption is Imχ (1) . The two multi-wave-mixing signals along k ± correspond to the nonlinear susceptibilities χ ± ∝ c k± (see Methods). We simulate the pulse propagation for the three modes along k p,± using coupled wave equations [27] (see Supplementary Information). For a trivial SL, the light propagating along k ± is symmetric, while for a TSL with C = 1, the topological currents drive the probe pulse to k + , even if the NNN hopping is two orders of magnitude smaller than the nearest-neighbour hopping.
Topological SL's can be readily realised in experiments for cold alkali atoms. Taking 85 Rb D1 line for example, we can have |g = |5 2 S 1/2 , F = 2 , |e = |5 2 P 1/2 , F = 2 and |m = |5 2 S 1/2 , F = 3 . γ e = 2.9MHz and γ m is controllable via inhomogeneous magnetic field [28]. The Rabi frequency Ω s = 3γ e = 8.6MHz (intensity 25mW/cm 2 ). The modulation frequency can be ν d = 10γ e = 28.8MHz which is large enough to separate the Floquet bands. In the µK regime, the thermal random motions have negligible Doppler shifts (∼kHz where |G = |g 1 , g 2 , ..., g N is the ground state. If the atoms are uniformly distributed and N Ω s /γ e , we can regard all points in the Brillouin zone are covered. The probe field Rabi frequency Ω p = µE p / where E p is the probe field strength and µ is the dipole moment between |e and |g . In the weak probe field limit Ω p γ e , the steady state probability amplitude (see Supplementary  Information) where the integration is over the first Brillouin zone, S is its area and is the probability amplitude of the atom at position r in state |e .
Pulse propagation. The wave-mixing susceptibility is χ ± = √ N µc k± /V 0 E p where V is the volume of the atomic ensemble and 0 is the permittivity in vacuum. We neglect the absorption of the coupling fields. The light propagation along k p,± was simulated using the coupled wave equations (see Supplementary Information) with the quasi-static approximation (steady-state values of c k are used), which is justified by the condition that the pulse duration is much longer than 1/γ e . The phase-matching condition is assumed. The incident probe pulse has Gaussian profile in both x and y directions.

EFFECTIVE HAMILTONIAN
The three-level atoms we use to construct the twodimensional (2D) superradiance lattice (SL) have a ground state |g , an excited state |e and a third state |m . The optical fields that couple |e and |m have three modes with wave vectors ). The interaction Hamiltonian with rotating-wave approximation is where κ l and a l are the vacuum coupling strength and annihilation operator of the lth mode, respectively. N is the number of atoms and r j is the position of the jth atom.
Initially the atomic ensemble is in the ground state |G = |g 1 , g 2 , ..., g N . A weak probe field with wave vector k p can prepare the atomic ensemble in the timed Dicke state (TDS) e ikp·rj |g 1 , g 2 , ..., e j , ..., g N (11) by the collective absorption of a single photon. The combinational quantum states of the atoms and three coupling modes can be written as |e kp , n 1 , n 2 , n 3 , where n l (l = 1, 2, 3) is the photon number in mode k l . This state is coupled to |m kp−k2 , n 1 , n 2 + 1, n 3 with coupling strength − κ 2 √ n 2 + 1. The state |m kp−k2 is defined the same as in equation (11) with e j replaced by m j and k p replaced by k p − k 2 . Similar coupling also exists for the other two modes. |m kp−k2 , n 1 , n 2 + 1, n 3 is in turn coupled to |e kp+k1−k2 , n 1 − 1, n 2 + 1, n 3 or |e kp+k3−k2 , n 1 , n 2 + 1, n 3 − 1 via excitation by modes k 1 or k 3 . These states form a honeycomb lattice with discrete momentum coordinates, as shown in Fig. 5, called the superradiance lattice [1]. The two sublattices of the SL correspond to TDS for |e and |m . We denote the coupling field frequency as ν c and the transition frequency between |e and |m as ω em . The energy difference between the two sublattice is ∆ c = ω em − ν c . We set the zero energy at the middle of the energies of the two sublattices. Then the energies of the |e and |m sublattices are ∆ c /2 and −∆ c /2, respectively.
Although the coupling strengths in the SL are sitedependent, if we use coherent coupling fields with large average photon numbers n l 1, the Rabi frequency of the lth mode can be approximated as the classical Rabi frequency Ω l = κ l n l . The SL with quantized photon numbers has edges when either of the three coupling field photon numbers reduces to zero. The total lattice has a triangular boundary. However, since the TDS have finite life time, which we suppose to be τ on average, once we create an excitation |e kp in the SL, the average distance this excitation can travel is in the order of Ω l τ . In this paper we have Ω l τ n l and the edges can never be reached. Interesting edge effects might exist in the fewphoton limit and it will be discussed elsewhere. Since the photon numbers are correlated to the momenta k's in the TDS, we can simplify the notations by dropping the photon numbers. The Hamiltonian of these TDS is where k = k p + r(k 2 − k 1 ) + s(k 3 − k 2 ) with integers r and s. The next-nearest-neighbour (NNN) hopping terms are induced by the periodic modulation of the three coupling field Rabi frequencies, where Ω s and Ω d are the static and dynamic components of the Rabi frequencies. ν d is the modulation frequency. φ l is the modulation phase of the lth field. k l d is the modulation wavevector. We choose the atomic ensemble much smaller than 1/|k l d | such that the position dependence of the phase −k l d · r can be neglected. We set −k l d · r = 0 for simplicity. We expand the Hamiltonian into static, positive-and negative-frequency components, where Note that H ±1 are not Hermitian themselves, but H +1 is the Hermitian conjugate of H −1 . The phase factors e ±iφ l in H ±1 play the crucial role for the complex NNN hopping terms in the Haldane model. The dynamics of the system is a Floquet problem [2,3]. Based on the Floquet theorem, the wave function can be written as where is the quasi-eigenenergy and (20) The terms with the same time evolution phase factors should be equal on both sides. We therefore have The quasi-eigenenergy can be obtained by diagonalizing the above Hamiltonian.
For the sake of simplicity, we assume the separation between the Floquet sidebands is much larger than the bandwidth, ν d Ω s,d , ∆ c , where the perturbation theory can be applied [4,5]. When the probe field is near resonance, ∆ p = ω eg − ν p Ω s where ω eg is the transition frequency between |e and |g and ν p is the probe field frequency, only the Floquet band with n = 0 in Eq.(21) is relevant. For states with eigenfrequencies near = 0, the effective Hamiltonian can be obtained by standard second-order perturbation as where with NNN Hamiltonian Explicitly, where The crucial factor i = √ −1 comes from the quantum interference between the two pathways shown in Fig.6 (a). The loop transitions via NNN hopping accumulate nonzero phases due to these complex factors, as shown in Fig.6 (b).
The effective Hamiltonian is greatly simplified in the real-space representation. We denote the real-space basis states as |e rj = |g 1 , g 2 , . . . , e j , . . . , g N , |m rj = |g 1 , g 2 , . . . , m j , . . . , g N , The effective Hamiltonian can be written as where the effective magnetic field h(r j ) = (h x , h y , h z ) with and the pseudo spin σ j = (σ x j , σ y j , σ z j ) with the Pauli matrices for the jth atom defined as σ x j = |e rj m rj | + |m rj e rj |, σ y j = −i|e rj m rj | + i|m rj e rj | and σ z j = |e rj e rj | − |m rj m rj |.

BERRY CONNECTION AND BERRY CURVATURE
The topological properties of the wavefunctions are determined by the Berry connection A and Berry curvature B, where ∇ and ∇× are the gradient and curl operators with respect to r. Similar to a spin-1/2 in a magnetic field h, the eigenenergies of the SL eigenstates of the effective Hamiltonian in equation (27) are where ± are for the upper and lower bands and h is the magnitude of h. The eigen wavefunction in the upper band can be written as This wavefunction is well defined except for the south pole of the Bloch sphere where h = (0, 0, h z ) with h z = −h. The eigen wavefunction near the south pole can be written as |ψ + (r) = e iφ(r) |ψ + (r) = 1 This way, the whole Bloch sphere is fully covered by two gauges [6].
The Berry connection in the upper band (except for the south pole of the Bloch sphere) is and near the south pole we can use The Berry curvature is where abc with a, b, c ∈ x, y, z is the Levi-Civita symbol. The Chern number is defined as the total Berry curvature in the whole first Brillouin zone which counts the winding number of the effective magnetic field h wrapping around the Bloch sphere in the whole Brillouin zone [7,8].
In Fig.7, we plot h in the Brillouin zone for a topological nontrivial SL where the band gap is opened by the NNN hopping. If we go from the K point r 9kcŷ , h moves from the north pole to the south pole of the Bloch sphere. The Chern number is one. If the band gap is opened by the on-site offset, ∆ c , h z is a constant in the Brillouin zone. h can only cover a patch on the Bloch sphere and the Chern number is zero.
The two bands of H eff have the smallest gap at the K and K points r ± , where 3 l=1 e ir±·k l = 0, and hence h x = h y = 0. At these symmetry points of the Brillouin zone the Hamiltonian is diagonal, where α = − 3 l=1 sin(φ l+1 − φ l ). Let us consider the specific case that ∆ c = 0, φ l = (l − 1)4π/3, and thus α = 3 √ 3/2. In this case H(r ± ) = ±9 Ω σ z has opposite signs at r ± , as shown in Fig.7 (a). It is obvious that at r + , the eigenstate in the upperband with eigenenergy 9 Ω is |e , while at r − , |e is the eigenstate in the lower band, as shown in Fig.8 (a).
Near r + , we can use equation (34) to describe the wavefunction. However, equation (34) cannot describe all the wavefunctions in the upper band. It has a singularity at r − where the magnetic field h points to the south pole. We can remove the singularity by tuning ∆ c . When ∆ c > 4 √ 3Ω α, the effective magnetic field h does not experience the whole Bloch sphere [see Fig.7 (b)]. In this case, h z in equation (42) has the same sign for r ± , |e is in the upper band at both r ± , as shown in Fig.8 (b). The wavefunction can be described by the same gauge. Thus the single-valued Berry connection A has no singularities on a closed surface and the Chern number is zero. Generally, the Chern number of the upper band can be written as (43) Specifically, when ∆ c = 0, C = sign(α).

DYNAMIC EVOLUTION
We can solve the dynamics of the atomic ensemble in real space with the total Hamiltonian including the probe field, withν p = ∆ c /2 − ∆ p being the probe detuning with respect to the middle of the band gap and Ω p the probe field Rabi frequency. The wavefunction in real space The dynamic equations of the probability amplitudes arė where γ e/m is the decoherence rate of |e/m states.
In the limit of weak probe field Ω p γ e , we have c G ≈ 1. In the steady state,ċ e (r j ) =ċ m (r j ) = 0. From Eq.(47), we obtain Substituting Eq.(48) in Eq.(46), we obtain In the SL coordinates, the wavefunction can be written as where the probability amplitude For uniformly distributed atoms, in the limit N Ω s /γ e , we can assume all points in the Brillouin zone are occupied by atoms, so the summation can be written as integration in the first Brillouin zone, where S is the area of the first Brillouin zone.

DIFFRACTION CONTRAST AND TOPOLOGY
In this section, we quantitatively compare the contrast between |c k+ | 2 and |c k− | 2 in topological and trivial SL's. We assumeν p = 0 and γ m = 0. Then Eq.(49) becomes We assume Ω , ∆ c Ω s and thus c e (r) is highly centred at K and K points, where we have h x = h y = 0 and According to Eq.(52), the probability amplitude c k can be approximately calculated by the integration of c e (r) in small areas near r ± , where the phase factors have been taken out of the integration since they do not change appreciably in those small areas. For a topological SL where Ω = 0 and ∆ c = 0, h z has opposite signs at r ± . We denote˜r (56) For k p = −k 1 , k + = −k 3 and k − = −k 2 , we have For a trivial SL where Ω = 0 and ∆ c = 0, h z has the same sign at r ± . We have˜r Therefore, c k+ = c k− and η = 0.

COUPLED-WAVE EQUATIONS
In this section, we calculate the propagation of a weak probe pulse. We set k p = −k 1 , k + = −k 3 and k − = −k 2 . There will be directional emission in k ± by TDS |e k± . The emitted photons in k ± interact with the atomic ensemble the same way as the probe field k p , creating TDS and resulting in directional emission in the other two modes. For example, emission at k + excites the atom ensemble and leads to emission at k − and k p . Therefore, the three optical fields in k p,± are coupled via the atoms. In the following, we derive the coupling equations.
The expectation value of the dipole moment of the jth atom at position r j is whereμ j is the dipole operator of the jth atom and µ = g j |μ j |e j . Since the atoms are homogeneously distributed, the polarization density as a function of positions is, where V is the volume of the atomic ensemble. The probe field Rabi frequency is Ω p = µE p / with E p being the probe field strength. The polarization density contains the Fourier components in k ± , P ±,p = 0 χ ± E p with χ ± = √ N µc k± /V 0 E p and 0 is the permittivity in vacuum. The notation P ±,p means the k ± Fourier component of polarization density generated by optical 1 0 1 Figure 9: The contrast of the NNN coupling strengths. The contrast (|Ω13| − |Ω12|)/(|Ω13| + |Ω12|) where Ω13 (Ω12) is the NNN hopping coefficient from |e kp to |e k + (|e k − ). The dot lines are zero points. field in k p . Once modes k ± are excited and the corresponding fields E ± are generated, they also polarize the atoms. For φ 1 = 0, φ 2 = 4π/3 and φ 3 = 2π/3, the three fields have a cyclic relation between each other. We have P p,+ = 0 χ − E + , P −,+ = 0 χ + E + , P p,− = 0 χ + E − and P +,− = 0 χ − E − .
We assume the quasi-static approximation in which the atoms are assumed in steady sates, i.e., c k is given by equation (52). This approximation is justified when the pulse duration is much longer than the decoherence time 1/γ e . We also use the slowly-varying-envelop approximation. The fields in the three relevant modes k j with j = p, + and − are denoted as E j (r)e −iνpt+ikj ·r where |k j | = ν p /c ≡ k p . The coupled-wave Maxwell equations for the three modes are wherek p/+/− are the unit vectors in the directions of k p/+/− .

LOCAL AND GLOBAL EFFECTS
In the main text, we use the chiral topological current to explain the consistence between η and the topological order. The topological current is a global effect for which we must treat the lattice as a whole. Since the three superradiant TDS |e kp,± are next nearest neighbours in the lattice, the local effect induced by the transitions inside a unit cell is also important. The direct NNN transition strength from |e kp to |e k+ or |e k− is Ω 13 or Ω 12 , respectively. If the dynamics is dominated by these direct transitions, we expect the contrast of the probability η should be similar to (|Ω 13 | − |Ω 12 |)/(|Ω 13 | + |Ω 12 |), which is plotted in Fig.9. We plot the contrast η in Fig.10 (a) to (c) for different Ω /γ m . The larger Ω is, the closer η is to Fig.9 and the more important is the local effect. For Ω γ m in Fig.10  (c), the effect of the direct NNN transitions is suppressed by the decoherence rate γ m . Thus the global effect dominates for most areas. In Fig.10 (d), η is calculated for only a unit cell. For the areas denoted with "L" in Fig.10 (c), the result is consistent with the corresponding areas in Fig.10 (d). In those areas, the topological current is weak and the local effect dominates. For the areas denoted with "G", the topological current is strong enough to overcome the local effect and the result is consistent with Fig.2 (a)